We compute the ground-state phase diagram of the Hubbard and frustrated Hubbard models on the square lattice with density matrix embedding theory using clusters of up to 16 sites. We provide an error model to estimate the reliability of the computations and complexity of the physics at different points in the diagram. We find superconductivity in the ground state as well as competition between inhomogeneous charge, spin, and pairing states at low doping. The estimated errors in the study are below T c in the cuprates and on the scale of contributions in real materials that are neglected in the Hubbard model. DOI: 10.1103/PhysRevB.93.035126 The Hubbard model [1][2][3] is one of the simplest quantum lattice models of correlated electron materials. Its one-band realization on the square lattice plays a central role in understanding the essential physics of high-temperature superconductivity [4,5]. Rigorous, near-exact results are available in certain limits [6]: at high temperatures from series expansions [7][8][9][10], in infinite dimensions from converged dynamical mean-field theory [11][12][13][14], and at weak coupling from perturbation theory [15] and renormalization group analysis [16,17]. Further, at half filling, the model has no fermion sign problem, and unbiased determinantal quantum Monte Carlo simulations can be converged [18]. Away from these limits, however, approximations are necessary. Many numerical methods have been applied to the model at both finite and zero temperatures, including fixed-node, constrained path, determinantal, and variational quantum Monte Carlo (QMC) [19][20][21][22][23][24][25][26][27][28][29], density matrix renormalization group (DMRG) [30][31][32], dynamical cluster (DCA) [33,34], (cluster) dynamical mean-field theories (CDMFT) [35,36], and variational cluster approximations (VCA) [37,38]. (We refer to DCA/CDMFT/VCA collectively as Green's function cluster theories.) These pioneering works have suggested rich phenomenology in the phase diagram including metallic, antiferromagnetic, and d-wave (and other kinds of) superconducting phases, a pseudogap regime, inhomogeneous orders such as stripes, and charge, spin, and pairdensity waves, as well as phase separation [6,19,20,24,25,[27][28][29]32,35,[39][40][41][42][43][44][45][46][47][48][49][50][51][52][53][54][55][56][57][58]. However, as different numerical methods have yielded different pictures of the ground-state phase diagram, a precise quantitative picture of the ground-state phase diagram has yet to emerge.It is the goal of this paper to produce such a quantitative picture as best as possible across the full Hubbard model phase diagram below U = 8. Our method of choice is density matrix embedding theory (DMET), which is very accurate in this regime [59][60][61][62][63][64][65][66], employed together with clusters of up to 16 sites and thermodynamic extrapolation. We carefully calibrate errors in our calculations, giving error bars to quantify the remaining uncertainty in our phase diagram. These error bars also serve,...