2D Fermi-Hubbard Model Overview

• The 2D Fermi-Hubbard model is a fundamental square-lattice framework describing spin-½ fermions with nearest-neighbor hopping and on-site interactions, leading to diverse phases such as antiferromagnetism, Mott insulation, and superconductivity.
• Advanced numerical and quantum simulation techniques like DCA, tensor-network methods, and dual fermion approaches are used to investigate its rich finite-temperature magnetism, pseudogap behavior, and transport anomalies.
• Research on this model reveals complex regimes including phase separation, pairing crossovers, and spectral-topology changes that challenge conventional many-body theories and guide experimental quantum simulations.

The two-dimensional Fermi-Hubbard model is the single-band lattice model of spin-$\frac12$ fermions on a square lattice with nearest-neighbor hopping and on-site interaction, conventionally written

$$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$

In the repulsive case $U>0$, it organizes antiferromagnetic, Mott, pseudogap, and doped metallic or pairing regimes near half-filling; in the attractive case $U<0$, it provides a controlled setting for the BCS–BEC crossover and for correlation-driven changes in spectral topology (IV. et al., 2017, Carlström, 2023). Because the model is simple enough to define microscopically yet difficult enough to resist a complete solution in two dimensions, it has become a standard reference system for thermodynamics, transport, fluctuation diagnostics, quantum-gas-microscope measurements, tensor-network studies, diagrammatic approaches, and quantum-simulation algorithms (LeBlanc et al., 2013, Chen et al., 2020).

1. Hamiltonian, lattice structure, and standard conventions

On the square lattice, $t$ sets the energy and temperature scale, and the nearest-neighbor tight-binding dispersion is

$\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$

For the repulsive model, the Brillouin zone is the square with components $q_x,q_y\in(-\pi,\pi]$, and the commensurate antiferromagnetic wavevector is $\mathbf{Q}=(\pi,\pi)$ (IV. et al., 2017). The particle density per site is $n=\langle n_i\rangle$, while doping away from half-filling is commonly written as $x=1-n$ or $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$0 (IV. et al., 2017, LeBlanc et al., 2013).

Half-filling means $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$1. In the particle-hole-symmetric convention often used for the square-lattice repulsive model, $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$2 corresponds to half-filling when chemical potential is referenced to $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$3 (LeBlanc et al., 2013). Many analyses adopt the idealized case $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$4, strictly local interactions, a single band, and no phonons. In that limit the square lattice has Fermi-surface nesting and particle-hole symmetry at $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$5, while the Mermin–Wagner theorem forbids antiferromagnetic long-range order at any nonzero temperature in pure two dimensions (IV. et al., 2017).

Several extensions preserve the same basic Hubbard structure while modifying symmetry or control parameters. The attractive model changes only the sign of $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$6 and is used to study pairing, pseudogap formation, and the BCS–BEC crossover (Carlström, 2023, Shi et al., 2023). A tilted version adds a linear potential $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$7 and couples mass transport to local heating through energy conservation (Guardado-Sanchez et al., 2019). SU($$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$8) generalizations replace the two spin states by $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$9 flavors and retain the same nearest-neighbor kinetic term and on-site interaction structure (Pasqualetti et al., 2023).

2. Ground-state organization near half-filling

For the square-lattice repulsive model, the zero-temperature phase structure near half-filling is highly structured but also model-specific. At $U>0$0, the ground state is an antiferromagnetic insulator for any $U>0$1 (IV. et al., 2017). At moderate interactions, $U>0$2, the ground state for $U>0$3 is a $U>0$4-wave superfluid, and for small doping $U>0$5 the model is regarded as most likely unstable to phase separation (IV. et al., 2017). The same source emphasizes that the phase-separation dome likely lies within the finite-temperature magnetic fluctuation regime rather than outside it.

These statements must be read together with the restrictions of the idealized model. The quoted phase diagram assumes a strictly two-dimensional, single-band, nearest-neighbor model with $U>0$6 and no phonons, so quantitative ranges depend on $U>0$7 and are sensitive to finite-temperature and finite-cluster biases in numerical methods (IV. et al., 2017). A recurrent issue is that large correlation lengths and phase-separation tendencies can bias finite-size calculations toward homogeneous states.

For the attractive model, the ground-state language shifts from antiferromagnetism and Mottness to pairing and crossover physics. In two dimensions with $U>0$8, pairing evolves continuously from a BCS regime of large, overlapping Cooper pairs to a BEC regime of tightly bound on-site pairs, and the Gaussian pairing state is reported to be a good approximation to the ground states of the attractive Hubbard model, in particular in the strong and weak coupling limits (Shi et al., 2023). This suggests that the sign of $U>0$9 does not merely reverse an interaction term; it reorganizes the relevant low-energy degrees of freedom.

3. Finite-temperature magnetism, pseudogap behavior, and real-space correlations

A central finite-temperature object is the static spin susceptibility

$U<0$0

with $U<0$1 (IV. et al., 2017). In the repulsive model, a pronounced and narrow peak in $U<0$2 indicates a large magnetic correlation length. An operational magnetic crossover temperature $U<0$3 was defined by the condition $U<0$4, marking the onset of a regime with $U<0$5 (IV. et al., 2017). Within the self-consistent skeleton-diagrammatic study of that regime, $U<0$6 can be as high as $U<0$7 at half-filling, vanishes at both small and large $U<0$8, and is maximal around $U<0$9, with the strongest finite-temperature magnetic fluctuations for $t$0–6 (IV. et al., 2017). Away from half-filling, the dominant magnetic peak becomes incommensurate and can split around $t$1, producing a diagonally oriented “carpet” of antiferromagnetic domains separated by diagonal domain walls; whether that texture survives to $t$2 remains open (IV. et al., 2017).

The pseudogap-like and non-Fermi-liquid regime has been resolved more microscopically through fluctuation diagnostics. In the ladder dual-fermion analysis, the low-frequency self-energy crossover is governed by relatively sharp spin fluctuations at $t$3, even at weak coupling (Arzhang et al., 2019). At $t$4, $t$5, $t$6, and $t$7, the antinodal spectral weight at $t$8 is about $t$9 of the nodal value, while the sign change of the standard Matsubara metric $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$0 tracks the crossover from Fermi-liquid to non-Fermi-liquid behavior (Arzhang et al., 2019). In the particle-hole-symmetric weak-coupling case $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$1, the crossover temperatures were found near $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$2 and $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$3 for nodal and antinodal momenta, with correlation length $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$4 near the crossover (Arzhang et al., 2019).

Quantum-gas-microscope experiments and finite-temperature simulations have made the same regime visible in real space. At $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$5 and temperatures $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$6–$\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$7, nearest-neighbor spin correlations are maximal at half-filling and weaken monotonically upon doping; at the lowest reported temperature, the nearest-neighbor spin correlator reaches about $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$8 (Cheuk et al., 2016). The nearest-neighbor moment correlator changes sign near a doping $\varepsilon_{\mathbf{k}}=-2t(\cos k_x+\cos k_y).$9, reflecting a crossover from anti-bunching of singly occupied sites to strong doublon-hole bunching, and the corresponding nonlocal moment correlations imply a nearest-neighbor contribution to potential-energy fluctuations of about $q_x,q_y\in(-\pi,\pi]$0 at half-filling (Cheuk et al., 2016). XTRG and DQMC calculations on open clusters at $q_x,q_y\in(-\pi,\pi]$1 further show that, upon doping, the diagonal spin correlation $q_x,q_y\in(-\pi,\pi]$2 changes sign around $q_x,q_y\in(-\pi,\pi]$3 and the next-nearest correlation $q_x,q_y\in(-\pi,\pi]$4 around $q_x,q_y\in(-\pi,\pi]$5, a pattern interpreted in terms of magnetic polarons and geometric-string distortions of the antiferromagnetic background (Chen et al., 2020).

4. Equation of state, thermodynamics, and transport

Large-cluster DCA calculations provide a thermodynamic equation of state for the repulsive square-lattice model in the thermodynamic limit, including $q_x,q_y\in(-\pi,\pi]$6, $q_x,q_y\in(-\pi,\pi]$7, double occupancy $q_x,q_y\in(-\pi,\pi]$8, nearest-neighbor spin correlations, and $q_x,q_y\in(-\pi,\pi]$9 for $\mathbf{Q}=(\pi,\pi)$0, $\mathbf{Q}=(\pi,\pi)$1, and $\mathbf{Q}=(\pi,\pi)$2 near half-filling (LeBlanc et al., 2013). At $\mathbf{Q}=(\pi,\pi)$3, a pronounced incompressible region with $\mathbf{Q}=(\pi,\pi)$4 appears near half-filling at low temperature; for $\mathbf{Q}=(\pi,\pi)$5, $\mathbf{Q}=(\pi,\pi)$6 is nearly flat as $\mathbf{Q}=(\pi,\pi)$7 varies from roughly $\mathbf{Q}=(\pi,\pi)$8 to $\mathbf{Q}=(\pi,\pi)$9 (LeBlanc et al., 2013). The same study identifies a behavioral shift in the energy below a pseudogap crossover temperature $n=\langle n_i\rangle$0, accompanied by entropy suppression and a low-temperature increase of double occupancy in the thermodynamic-limit extrapolation (LeBlanc et al., 2013). The specific heat has a high-temperature “charge” peak near $n=\langle n_i\rangle$1 and a low-temperature “spin” feature tied to the same crossover (LeBlanc et al., 2013).

Transport in the weak-coupling regime shows a distinct but equally structured set of regularities. In quantum kinetic theory for the 2D single-band model, the electrical resistivity is $n=\langle n_i\rangle$2-linear at high temperature,

$n=\langle n_i\rangle$3

and the thermal resistivity behaves as

$n=\langle n_i\rangle$4

for generic filling (Kiely et al., 2021). At half-filling, perfect nesting and a van Hove singularity produce nearly $n=\langle n_i\rangle$5-linear electrical resistivity down to the lowest temperatures in that framework, while the low-temperature Wiedemann–Franz ratio tends to $n=\langle n_i\rangle$6 at half-filling and to $n=\langle n_i\rangle$7 away from half-filling (Kiely et al., 2021). By contrast, for $n=\langle n_i\rangle$8 the zero-temperature Fermi surface is too small to allow Umklapp scattering, and resistivity is exponentially suppressed with an Umklapp gap $n=\langle n_i\rangle$9 (Kiely et al., 2021).

A distinct hydrodynamic regime emerges when the repulsive 2D model is placed under a uniform tilt. In quantum-gas-microscope experiments on a tilted square lattice with $x=1-n$0 and average filling $x=1-n$1, the decay time of prepared density waves crosses over from diffusive scaling $x=1-n$2 at weak tilt to subdiffusive scaling $x=1-n$3 at strong tilt (Guardado-Sanchez et al., 2019). In the strong-tilt regime, the decay rate is controlled by a thermal diffusivity through

$x=1-n$4

reflecting the fact that particle currents generate local heating and are bottlenecked by heat diffusion (Guardado-Sanchez et al., 2019). This identifies a transport regime in which density relaxation is no longer set by ordinary charge diffusion alone.

5. Numerical and quantum-simulation methods

Methodologically, the 2D Fermi-Hubbard model has served as a benchmark for controlled but mutually complementary approaches. A self-consistent skeleton-diagrammatic GGGW scheme dresses the single-particle propagator, particle-particle and particle-hole pair propagators, and screened interaction simultaneously; at moderate $x=1-n$5 it reproduces thermodynamic observables at the few-percent level and the amplitudes of $x=1-n$6 and $x=1-n$7 within $x=1-n$8 of DDMC benchmarks at half-filling (IV. et al., 2017). Large-cluster DCA with CT-AUX provides thermodynamic-limit extrapolations of the equation of state and spin correlations without analytic continuation for thermodynamic observables (LeBlanc et al., 2013). Ladder dual fermions give momentum-resolved fluctuation diagnostics on fine grids and isolate the $x=1-n$9 spin channel as the dominant source of non-Fermi-liquid self-energy structure (Arzhang et al., 2019).

Tensor-network and variational approaches cover complementary regions of parameter space. XTRG, combined with DQMC, reaches temperatures down to $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$00 on open clusters up to $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$01 at $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$02, quantitatively reproducing finite-temperature spin and charge correlators measured in ultracold-atom experiments (Chen et al., 2020). For the attractive model, a fermionic Gaussian variational method on 2D square lattices with periodic boundary conditions yields total energies with maximal systematic error $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$03 in the intermediate regime around $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$04, with negligible dependence on lattice size up to $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$05, and maps directly to repulsive-case energies and double occupancies through a partial particle-hole transformation (Shi et al., 2023). A fermionic PEPS-plus-imaginary-time-TEBD construction has also been used on square lattices with periodic boundaries up to $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$06, explicitly tracking Jordan–Wigner signs for vertical hopping by a doubling construction, although its small-bond-dimension energies remain visibly above exact results (Chung, 2013).

Quantum simulation has become an additional methodological branch rather than a replacement for classical methods. On a 16-qubit superconducting processor, a low-depth symmetry-preserving variational algorithm for a $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$07 ladder observed antiferromagnetic order at half-filling and its suppression away from half-filling (Stanisic et al., 2021). A QETU implementation for the $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$08 model on a 9-qubit grid-like architecture achieved ground-state fidelities above $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$09 in noiseless simulations with a sufficiently high polynomial degree, though direct energy estimation remained shot-noise limited at about $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$10 with $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$11 shots (Müller et al., 2024). In a different direction, a hybrid quantum-classical analog simulation on neutral-atom Rydberg hardware reformulated the half-filled 2D model into free fermions coupled to auxiliary spins and reproduced a metal-to-Mott crossover on a $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$12 lattice through a self-consistent auxiliary-spin solver (Julià-Farré et al., 7 Oct 2025). The common pattern is not methodological convergence but decomposition of the problem into regimes where different approximations are controlled.

6. Extensions, generalizations, and open questions

The attractive 2D Hubbard model has become a distinct arena for studying how strong correlations alter single-particle structure. At $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$13, diagrammatic Monte Carlo showed that the spectrum can become disconnected or multiply connected in momentum space because the electron operator fractionalizes into singlon-like and doublon-like components whose spectral weights collapse in momentum-selective fashion (Carlström, 2023). In that framework, open Fermi-arc-like structures are not themselves computed in all regimes, because some cases are gapped, but the nontrivial spectral connectivity is identified as the prerequisite for arcs when $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$14 intersects the excitation set (Carlström, 2023). This suggests a route to arc phenomenology that does not rely on explicit translation-symmetry breaking.

Generalizations in symmetry and lattice structure broaden the model’s scope without removing its Hubbard core. In a 2D single-layer square optical lattice, the SU($$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$15) Fermi-Hubbard model with $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$16 has been characterized experimentally through its equation of state, compressibility, occupation probabilities, and density fluctuations for $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$17, $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$18, $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$19, and $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$20, with Mott-like incompressible regions near integer filling and model-independent fluctuation thermometry consistent with EoS-based fits (Pasqualetti et al., 2023). On decorated honeycomb and triangular lattices, a half-filled Hubbard model can map to a 120-degree compass model; in that asymptotic regime, the triangular case yields collinear stripe order corresponding to a time-reversal-invariant $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$21-form-factor charge-density wave, whereas the decorated honeycomb case yields a unique, gapped fermionic symmetry-protected topological phase protected by time-reversal and reflection symmetries (Chen et al., 2016).

Several open problems remain explicit in the current literature. For the repulsive square-lattice model, the precise low-temperature location of the phase-separation dome relative to magnetic and pairing regimes remains unsettled, and finite-size methods may bias toward homogeneous superfluid states or miss phase separation altogether when $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$22 becomes large (IV. et al., 2017). In the attractive model, mapping the precise arc locus as a function of $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$23 and $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$24 remains future work (Carlström, 2023). In an altermagnetic extension of the 2D Hubbard model with sublattice-odd next-nearest-neighbor hopping, Hartree–Fock analysis finds broad mixed-phase regions even at half-filling; for $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$25, $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$26, and $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$27, the coexistence window at the crossing chemical potential spans densities from about $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$28 to $$H=-t\sum_{\langle ij\rangle,\sigma}\left(c^\dagger_{i\sigma}c_{j\sigma}+\text{h.c.}\right)+U\sum_i n_{i\uparrow}n_{i\downarrow}-\mu\sum_i n_i.$$29, so half-filling lies inside a mixed phase rather than a homogeneous one (Langmann et al., 19 Jun 2025). A plausible implication is that inhomogeneity, rather than a single uniform broken-symmetry state, may organize parts of the extended 2D Hubbard phase diagram that are experimentally accessible in cold atoms.

The 2D Fermi-Hubbard model therefore remains less a single solved phase diagram than a unifying framework in which antiferromagnetism, pseudogap behavior, doublon-hole correlations, transport anomalies, pairing tendencies, spectral-topology changes, symmetry-enriched descendants, and quantum-simulation strategies can all be formulated with a common microscopic language. Its importance derives not from analytical closure but from the unusually sharp way in which many different many-body phenomena can be traced back to the same square-lattice Hamiltonian.