The second-order reduced density matrix method (the RDM method) has performed well in determining energies and properties of atomic and molecular systems, achieving coupled-cluster singles and doubles with perturbative triples (CC SD(T)) accuracy without using the wave-function. One question that arises is how well does the RDM method perform with the same conditions that result in CCSD(T) accuracy in the strong correlation limit. The simplest and a theoretically important model for strongly correlated electronic systems is the Hubbard model. In this paper, we establish the utility of the RDM method when employing the $P$, $Q$, $G$, $T1$ and $T2^prime$ conditions in the two-dimension al Hubbard model case and we conduct a thorough study applying the $4times 4$ Hubbard model employing a coefficients. Within the Hubbard Hamilt onian we found that even in the intermediate setting, where $U/t$ is between 4 and 10, the $P$, $Q$, $G$, $T1$ and $T2^prime$ conditions re produced good ground state energies.