Three-dimensional Hubbard model in the thermodynamic limit

We employ the numerical linked-cluster expansion to study finite-temperature properties of the uniform cubic lattice Hubbard model in the thermodynamic limit for a wide range of interaction strengths and densities. We carry out the expansion to the 9th order and find that the convergence of the series extends to lower temperatures as the strength of the interaction increases, giving us access to regions of the parameter space that are difficult to reach by most other numerical methods. We study the precise trends in the specific heat, the double occupancy, and magnetic correlations at temperatures as low as 0.2 of the hopping amplitude in the strong-coupling regime. We show that in this regime, accurate estimates for transition temperatures to the Néel ordered phase, in agreement with the predicted asymptotic behavior, can be deduced from the low-temperature magnetic structure factor. We also find evidence for possible instability to the magnetically ordered phase away from, but close to, half filling. Our results have important implications for parametrizing fermionic systems in optical lattice experiments and for benchmarking other numerical methods.

Figure 1

(a) Specific heat of the 3D Hubbard model at half filling vs temperature for several values of the interaction strength. Thick solid lines represent the average between the last two orders after Euler resummation of the last 6 terms of the series, and Wynn resummations with 3 and 4 cycles of improvement. The error bars mark the confidence region where all the values used in the average fall. Thin dotted lines are the last two orders of the bare NLCE sums, and thin solid lines are the 9th and 10th orders of the HTSE [41, 42, 43]. (b) Double occupancy at half filling vs temperature. $U$ increases from top to bottom with the same values as in (a). (c) Double occupancy at $U=12$ away from half filling vs temperature.

Figure 2

(a) Double occupancy at half filling as a function of the interaction strength at different temperatures. $D$ rapidly decreases by increasing the interaction strength as one approaches the Mott phase. (b) Nearest-neighbor spin correlations at half filling vs $U$. Similarly to the 2D model [34], they initially rise in the weak-coupling regime by increasing the interaction strength and eventually decrease as $1/U$ in the strong-coupling regime. Its peak around $U=10$ becomes more pronounced as the temperature is lowered.

Figure 3

Antiferromagnetic structure factor at half filling as a function of temperature for several values of the interaction strength. Wynn resummation with 3 and 4 cycles of improvement have been used. The results are shown for temperatures above where the two estimates match within a few percent. The low-temperature structure factor increases as $U$ increases to 9, then reduces upon further increasing of the strength of the interaction. The dashed vertical line indicates the location of the transition temperature for $U=9$. Inset shows the inverse of $S^{AF}$ for select values of $U$. A simple fit to $A/{(T-T_N)}^B$ for $T<0.6$ suggests $T_N=0.35$ for $U=9$, which is consistent with current best estimates. Thin dotted lines show bare NLCE results for the last two orders.

Figure 4

Antiferromagnetic structure factor at half filling as a function of the interaction strength at different temperatures. The results are shown after Wynn resummations with 4 cycles of improvement (they also match results from Wynn with 3 cycles of improvement roughly within the size of symbols). The dashed line is a fit to $A/T$ (constant $A$) using values of $S^{AF}$ for the three largest $U$'s at $T=0.44$, representing the theoretical asymptotic trend in the strong-coupling regime.

Figure 5

Antiferromagnetic structure factor with $U=9$ as a function of density at several temperatures. As the temperature is lowered, the structure factor develops a sharp peak around half filling. Inset: Inverse of the structure factor vs temperature at three different densities.

Figure 6

(a) Estimated Néel transition temperature vs the interaction strength. Filled squares are obtained from extrapolations of the NLCE $S^{AF}$ to low temperatures, and empty circles are obtained from fits to $S^{AF}$ as a function of the NLCE order (see text). Filled triangles and diamonds are data from the DCA [13] and DQMC [26], respectively. The dotted line is a guide to the eye. The dashed line is the theoretical asymptotic function. (b) The bare AF structure factor for $U=12$ vs the NLCE order at different temperatures (with a uniform grid) above and below the expected transition temperature ($T_N\sim 0.30$). Dashed lines represent fits of results in the even orders to a 2nd-degree polynomial. (c) Inverse of $S^{AF}$ as a function of temperature for the 3D Heisenberg model with $J=1$ [which sets the unit of energy in (c) and (d)]. The last two orders of the bare sums and results after Wynn resummations with 6 cycles of improvement are shown. The latter points to a divergence at $T_N\sim 0.96$. (d) The staggered structure factor for the Heisenberg model as a function of the NLCE order at different temperatures (in a uniform grid) around $T_N$. Here, dashed lines are 2nd-degree polynomial fits of data in all orders.