Néel temperature and thermodynamics of the half-filled three-dimensional Hubbard model by diagrammatic determinant Monte Carlo

We study the thermodynamics of the three-dimensional Hubbard model at half filling on approach to the Néel transition by means of large-scale unbiased diagrammatic determinant Monte Carlo simulations. We obtain the transition temperature in the strongly correlated regime, as well as the temperature dependence of the energy, entropy, double occupancy, and nearest-neighbor spin correlation function. Our results improve the accuracy of previous unbiased studies and present accurate benchmarks in the ongoing effort to realize the antiferromagnetic state of matter with ultracold atoms in optical lattices.

Figure 1

Finite-size scaling for $T_N$ at $U=8t$ by Binder crossing analysis. Points are Monte Carlo results, lines are linear fits. The uncertainty for the Néel temperature is estimated conservatively by varying the Monte Carlo points within their respective error bars.

Figure 2

Finite-size scaling for $T_N$ at $U=5t$ by Binder crossing analysis. Points with error bars are Monte Carlo results, solid lines are linear fits. FSS procedures based on Eqs. (25, 26, 27) result in the estimates $T_N/t=0.2175\left(44\right)$ and $T_N/t=0.2211\left(26\right)$, respectively. See also Fig. 3.

Figure 3

Scaling—according to Eqs. (25) and (26)—of estimates of the critical temperature at $U=5t$ obtained from Binder crossings between the lines for different system sizes in Fig. 2. See text for discussion. The square corresponds to the crossing point between lines for $L=4$ and $L=10$, which substantially deviates from the linear scaling exhibited by all the crossings for $L>4$ (circles), demonstrating that the $L=4$ system is too small to be consistent with the critical scaling described by Eq. (24). Correspondingly, all the other crossings with the $L=4$ line are omitted from the figure.

Figure 4

Finite-size scaling for $T_N$ at $U=4t$ by Binder crossing analysis. In this case we are not able to reliably extract the Néel temperature and can provide only an upper limit, $T_N<0.17t$. Notice that the crossing of $L=6$ and $L=8$ is clearly outside the range of applicability of either (25) or (27).

Figure 5

Comparison of estimates for $T_N$ by different unbiased approaches. Also shown is the strong-coupling limiting behavior, $T_N=3.83t^2/U$. See text for discussion.

Figure 6

Dependence of the energy per particle of noninteracting fermions at half filling [described by Eq. (1) with $U=0,\mu =0$] on the inverse of the linear system size $L$.

Figure 7

Example of the dependence of energy on the inverse of the linear system size $L$ obtained with the periodic boundary conditions (PBCs, circles) and using the twist-averaged boundary conditions (TABCs, triangles) for $U=8$, $T=0.3875$ and $U=4$, $T=0.2$. The solid line is an extrapolation using the formula $E\left(L\right)=E^{\prime }\left(\infty \right)+C[E_0\left(L\right)-E_0\left(\infty \right)]$ (see text). For $U=8$, $T=0.3875$ the parameters are $E^{\prime }\left(\infty \right)=-0.5965$, $C=0.6$, while for $U=4$, $T=0.2$, $E^{\prime }\left(\infty \right)=-1.135$, $C=0.95$. The claimed thermodynamic-limit results—$E\left(\infty \right)=-0.5960\left(16\right)$ for $U=8$, $T=0.3875$ (using the PBC data) and $E\left(\infty \right)=-1.1390\left(9\right)$ for $U=4$, $T=0.2$ (using the TABC data)—are shown by the horizontal bands. More generally, we use PBC data at $U\ge 6$ and TABC data at $U\le 5$ to obtain the thermodynamic-limit values, as explained in Sec. 4a.

Figure 8

Example of the dependence of energy on the inverse of the linear system size $L$ obtained with the periodic boundary conditions (PBCs, circles) for $U=8$, $T=0.3875$ and with the twist-averaged boundary conditions (TABCs, triangles) for $U=4$, $T=0.2$, compared to the result of a simulation based on free-particle propagators of an infinite system ($G_L^{\left(0\right)}\to G_{\infty }^{\left(0\right)}$, squares). The error bars are smaller than the symbols. The dashed line is a linear fit yielding $E^{\prime \prime }\left(\infty \right)=-0.5962\left(9\right)$ [$E^{\prime \prime }\left(\infty \right)=-1.1368\left(32\right)$], in perfect agreement with the claimed conservative estimate $E\left(\infty \right)=-0.5960\left(16\right)$ [$E\left(\infty \right)=-1.1390\left(9\right)$] for $U=8$, $T=0.3875$ ($U=4$, $T=0.2$) shown by the horizontal band. The data point for the smallest system size $L=4$ at $U=4$, $T=0.2$ obtained with $G_L^{\left(0\right)}\to G_{\infty }^{\left(0\right)}$ deviates from the linear scaling followed by larger systems and therefore is excluded from the fit.

Figure 9

Energy (extrapolated to the TD limit) versus temperature at $U=8,6,5,4$. The lines represent the results of the high-temperature series expansion series of orders $2,4,6,8,10$ labeled correspondingly.

Figure 10

Entropy (extrapolated to the TD limit) versus temperature at $U=8,6,5,4$. In each panel, the vertical line shows the position of the critical temperature $T_N$ with its width given by the error bar, while the horizontal dashed lines represent the corresponding bounds of the critical entropy $S_N$ listed in Table . The rest of the lines represent the results of the high-temperature expansion series of orders $2,4,6,8,10$ labeled correspondingly.

Figure 11

Lines of constant entropy in the $T$ vs $U$ plane.

Figure 12

Double occupancy $\langle n_{\mathbf{x}\uparrow }{n}_{\mathbf{x}\downarrow }\rangle$ (squares) and the nearest-neighbor spin correlation function $|\langle S_{\mathbf{x}}^z{S}_{\mathbf{x}+{\mathbf{e}}_i}^z\rangle |$ (circles) versus $T$ at $U=8,6,5,4$ extrapolated to the TD limit. In each panel, the vertical line shows the position of the critical temperature for a given value of $U$ listed in Table with its width corresponding to the error bar.