  The variational determination of the two-particle density matrix is an
interesting, but not yet fully explored technique that allows to obtain
ground-state properties of a quantum many-body system without reference to an
$N$-particle wave function. The one-dimensional fermionic Hubbard model has
been studied before with this method, using standard two- and three-index
conditions on the density matrix [J. R. Hammond {\it et al.}, Phys. Rev. A 73,
062505 (2006)], while a more recent study explored so-called subsystem
constraints [N. Shenvi {\it et al.}, Phys. Rev. Lett. 105, 213003 (2010)].
These studies reported good results even with only standard two-index
conditions, but have always been limited to the half-filled lattice. In this
Letter we establish the fact that the two-index approach fails for other
fillings. In this case, a subset of three-index conditions is absolutely needed
to describe the correct physics in the strong-repulsion limit. We show that
applying lifting conditions [J.R. Hammond {\it et al.}, Phys. Rev. A 71, 062503
(2005)] is the most economical way to achieve this, while still avoiding the
computationally much heavier three-index conditions. A further extension to
spin-adapted lifting conditions leads to increased accuracy in the intermediate
repulsion regime. At the same time we establish the feasibility of such studies
to the more complicated phase diagram in two-dimensional Hubbard models.

Strongly Correlated Electrons