Finite-size scaling of O(n) systems at the upper critical dimensionality

Abstract

INTRODUCTION

Accordingly, the standard FSS formulae (2) and (3) predict that the critical susceptibility diverges as |$\chi \asymp L^{2y_h-d} = L^2$| for d > d_c. However, for the Ising model on 5D periodic hypercubes, χ was numerically observed to scale as L ≍ L^5/2 instead of L² [13–18]. The FSS for d ≥ d_c turns out to be surprisingly subtle and remains a topic of extensive controversy [13–21].

According to (8), the critical correlation function still exhibits a Gaussian-like decay, g(r, L) ≍ r^{−(d − 2)}, up to a length scale ξ₁ = L^{d/[2(d − 2)]}, and then enters an r-independent plateau whose height vanishes as L^−d/2. Since the length ξ₁ is vanishingly small compared to the linear size, ξ₁/L → 0, the plateau effectively dominates the scaling behavior of g(r, L) and the FSS of χ. The two-length scaling form (8) has been numerically confirmed for the 5D Ising model and self-avoiding random walk, with a geometric explanation based on the introduction of an unwrapped length on the torus [18]. It is also consistent with the rigorous calculations for the so-called random-length random-walk model [20]. It is noteworthy that the two-length scaling is able to explain both the FSS χ₀ ≡ χ ≍ L^5/2 for the susceptibility (the magnetic fluctuations at the zero Fourier mode) [14] and the FSS χ_k ≍ L² for the magnetic fluctuations at nonzero modes [15,17].

• Let |$\vec{\cal M} \equiv \sum _{\bf r} \vec{S}_{\bf r}$| specify the total magnetization of a spin configuration, and measure its ℓ moment as |$M_{\ell } \equiv \langle |\vec{\cal M}|^\ell \rangle$|⁠. Equation (9) predicts that |$M_{\ell } \sim L^{\ell y^*_h}+q L^{\ell y_h}$|⁠, with q a nonuniversal constant. In particular, the magnetic susceptibility |$\chi _0 \equiv L^{-d} M_2 \asymp L^{d/2} [1+\mathcal {O}(L^{(4-d)/2})]$|⁠, where the FSS from the Gaussian term |$\tilde{f}_0$| is effectively a finite-size correction, but its existence is important in analyzing numerical data [21].

• Let |$\vec{\cal M}_{\bf k} \equiv \sum _{\bf r} \vec{S}_{\bf r} e^{i {\bf k} \cdot {\bf r}}$| specify the Fourier mode of magnetization with momentum k ≠ 0, and measure its ℓ moment as |$M_{\ell , {\bf k}} \equiv \langle |\vec{\cal M}_{\bf k} |^\ell \rangle$|⁠. The magnetic fluctuations at k ≠ 0 behave as |$\chi _{\bf k} \equiv L^{-d} M_{2, {\bf k} } \sim L^{2y_h-d}=L^2$|⁠. The behaviors of χ₀ and χ_k have been confirmed for the 5D Ising model [15,17,18,20].

• The Binder cumulant |$Q \equiv \langle |\vec{\cal M} |^2 \rangle ^2/\break\langle |\vec{\cal M} |^4 \rangle$| should take the complete-graph value, as expected from the correspondence between the term with |$\tilde{f}_1$| in (9) and the complete-graph FSS. For the Ising model, the complete-graph calculations give Q = 4[Γ(3/4)/Γ(1/4)]² ≈ 0.456 947, consistent with the 5D result in [13].

Analogously, the FSS behaviors of the energy density, its higher-order fluctuations and the ℓ-moment Fourier modes at k ≠ 0 can be derived from (9).

We expect that the FSS formulae (8) and (9) are valid not only for the O(n) vector model but also for generic systems of continuous phase transitions at d > d_c. An example is given for percolation that has d_c = 6. It was observed [29] that, at criticality, the probability distributions of the largest-cluster size follow the same scaling function for 7D periodic hypercubes and on the complete graph.

In this work, we focus on the FSS for the O(n) vector model at the upper critical dimensionality d = d_c. In this marginal case, it is known that multiplicative and additive logarithmic corrections would appear in the FSS. However, exploring these logarithmic corrections turns out to be notoriously hard. The challenge comes from the lack of analytical insights, the existence of slow finite-size corrections, as well as the unavailability of very large system sizes in simulations of high-dimensional systems.

For the O(n) vector model, establishing the precise FSS form at d = d_c is not only of fundamental importance in statistical mechanics and condensed-matter physics, but also of practical relevance due to the direct experimental realizations of the model, particularly in three-dimensional quantum critical systems [3–6,10,11]. For instance, to explore the stability of Anderson–Higgs excitation modes in systems with continuous symmetry breaking (n ≥ 2), a crucial theoretical question is whether or not the Gaussian r-dependent behavior g(r) ≍ r⁻² is modified by some multiplicative logarithmic corrections.

SUMMARY OF THE MAIN FINDINGS

for n ≥ 0 and n ≠ 4, where the renormalization exponents y_t = 2 and y_h = 3 are given by (4). Furthermore, the renormalization-group calculations predicted the logarithmic-correction exponents as |$\hat{y}_t={(4-n)}/{(2n+16)}$| and |$\hat{y}_h=1/4$| [32,33]. The leading FSS of χ₀ is hence given by χ₀ ≍ L²(lnL)^1/2, independent of n.

with |$\hat{p} = 2 \hat{y}_h = 1/2$|⁠. By (12), we explicitly point out that no multiplicative logarithmic correction appears in the r dependence of g(r, L) ≍ r⁻², which is still Gaussian-like. By contrast, the plateau for |$r \ge \xi _1 \sim L/({\rm ln}L)^{\hat{p}}$| is modified as |$L^{-2}({\rm ln}L)^{\hat{p}}$|⁠. In other words, along any direction of the periodic hypercube, we have |$g(r,L) \asymp r^{-2} + v L^{-2}({\rm ln}L)^{\hat{p}}$|⁠, with |$v$| a nonuniversal constant. The r⁻² decay at shorter distances in (12) is consistent with analytical calculations for the 4D weakly self-avoiding random walk and the O(n) φ⁴ model directly in the thermodynamic limit (L → ∞) [34], which predict |$g(r) \asymp r^{-2} (1 + \mathcal {O}(1/{\rm ln} r))$|⁠.

The roles of terms with |$\tilde{f}_0$| and |$\tilde{f}_1$| in (11) are analogous to those in (9). The former arises from the Gaussian fixed point, and the latter describes the ‘background’ contributions (k = 0) for the FSS of macroscopic quantities. However, note that the term with |$\tilde{f}_1$| can no longer be regarded as an exact counterpart of the FSS of the complete graph, due to the existence of multiplicative logarithmic corrections. By contrast, the exact complete-graph mechanism applies to the |$\tilde{f}_1$| term in (9), where the logarithmic correction is absent and |$\tilde{f}_1$| corresponds to the free energy of the standard complete-graph model. According to (11), the FSS of various macroscopic quantities at d = d_c can be obtained as follows.

• The magnetization density |$m \equiv L^{-d} \langle |\vec{\cal M}|\rangle \asymp L^{-1} (\ln L)^{\hat{y}_h} [1+\mathcal {O}((\ln L)^{-\hat{y}_h}) ]$|⁠.

• The magnetic susceptibility |$\chi _0 \asymp L^{2}(\ln L)^{2 \hat{y}_h} [1+\mathcal {O}((\ln L)^{-2\hat{y}_h}) ]$|⁠.

• The magnetic fluctuations at k ≠ 0 Fourier modes χ_k ≍ L².

• The Binder cumulant Q may not take the exact complete-graph value, due to the multiplicative logarithmic correction. Some evidence was observed in a recent study by Y.D. and his coworkers for the self-avoiding random walk (n = 0) on 4D periodic hypercubes, in which the maximum system size is up to L = 700.

The FSS of the energy density, its higher-order fluctuations and the ℓ-moment Fourier modes at k ≠ 0 can be obtained.

In quantities like m and χ₀, the FSS from the Gaussian fixed point effectively plays the role of finite-size corrections. Nevertheless, we note that in the analysis of numerical data, it is important to include such scaling terms.

We remark that the FSS formulae (11) and (12) for d = d_c are less generic than (8) and (9) for d > d_c. For the O(n) models, a multiplicative logarithmic correction is absent in the Gaussian r dependence of g(r, L) in (12). Although the two length scales are possibly generic features of models with logarithmic finite-size corrections at upper critical dimensionality, multiplicative logarithmic corrections to the r dependence of g(r, L) require case-by-case analyses. Equation (11) can be modified in some of these models, which include the percolation and spin-glass models in six dimensions.

We proceed to verify (11) and (12) using extensive Monte Carlo (MC) simulations of the O(n) vector model. Before giving the technical details, in Fig. 1 we present complementary evidence for (11) and (12) in the case of the critical 4D XY model. In Fig. 1(a) we show the extensive data of g(r, L) for 16 ≤ L ≤ 80, of which the largest system contains about 4 × 10⁷ lattice sites. To demonstrate the multiplicative logarithmic correction in the large-distance plateau indicated by (12), we plot g(L/2, L)L² versus lnL on a log-log scale in Fig. 1(b). The excellent agreement between the MC data and the formula |$v$|₁(lnL)^1/2 + |$v$|₂ provides a first piece of evidence for the presence of the logarithmic correction with exponent |$\hat{p} = 1/2$|⁠. The second piece of evidence comes from Fig. 1(c), which suggests that the χ₀L⁻² data can be well described by the formula q₁(lnL)^1/2 + q₂. Finally, in Fig. 1(d) we plot the k ≠ 0 magnetic fluctuations χ₁ and χ₂ with k₁ = (2π/L, 0, 0, 0) and k₂ = (2π/L, 2π/L, 0, 0), respectively, which suppress the L-dependent plateau and show the r-dependent behavior of g(r, L). Indeed, the χ₁L⁻² and χ₂L⁻² data converge rapidly to constants as L increases.

Figure 1.

Evidence for conjectured formulae (11) and (12) in the example of the critical four-dimensional (4D) XY model. (a) Correlation function g(r, L) on a log-log scale. The solid line denotes r⁻² behavior. (b) Scaled correlation g(r, L)L² with r = L/2 versus lnL on a log-log scale. Thus, the horizontal axis is effectively on a double logarithmic scale of L. The solid line represents logarithmic divergence with |$\hat{p} = 1/2$|⁠. (c) Scaled magnetic susceptibility χ₀L⁻² versus lnL on a log-log scale. The solid line accounts for logarithmic divergence with |$\hat{p} = 1/2$|⁠. (d) Scaled k ≠ 0 magnetic fluctuations χ₁L⁻² and χ₂L⁻², with k₁ = (2π/L, 0, 0, 0) and k₂ = (2π/L, 2π/L, 0, 0), respectively. The horizontal lines strongly indicate the absence of logarithmic corrections in the scaling of χ_k.

Open in new tabDownload slide

NUMERICAL RESULTS AND FINITE-SIZE SCALING ANALYSES

Using a cluster MC algorithm [35], we simulate Hamiltonian (1) on 4D hypercubic lattices up to L_max = 96 (Ising, XY) and 56 (Heisenberg), and measure a variety of macroscopic quantities, including the magnetization density m, the susceptibility χ₀, the magnetic fluctuations χ₁ and χ₂ and the Binder cumulant Q. Moreover, we compute the two-point correlation function g(r, L) for the XY model up to L_max = 80 by means of a state-of-the-art worm MC algorithm [36].

Estimates of critical temperatures

By analyzing the finite-size correction Q(L, T_c) − Q_c, we find that the leading correction is nearly proportional to (lnL)^−1/2, consistent with the prediction of (11) and (12). We let Q_c be free in the fits and have Q_c = 0.45(1), close to the complete-graph result Q_c = 0.456 947. Besides, we perform simulations for the XY and Heisenberg models on the complete graph and obtain Q_c ≈ 0.635 and 0.728, respectively, also close to the fitting results of the 4D Q data. We obtain T_c(XY) = 3.314 437(6), and in Fig. 2(a) we illustrate the location of T_c by Q.

Figure 2.

Locating T_c for the 4D XY model. (a) The Binder cumulant Q with finite-size corrections being subtracted, namely, Q^*(L, T) = Q(L, T) − b(lnL)^−1/2, with b ≈ 0.1069 according to a preferred least-squares fit. The shadow marks T_c and its error margin. (b) The magnetization density m rescaled by L⁻¹ versus lnL around T_c = 3.314 44 on a log-log scale.

Open in new tabDownload slide

We further examine the estimate of T_c by the FSS of other quantities, such as the magnetization density m. For the XY model, in Fig. 2(b) we give a log-log plot of the mL data versus lnL for T = T_c, as well as for T_low = 3.314 40 and T_above = 3.314 50. The significant bending-up and bending-down features clearly suggest that T_low < T_c and T_above > T_c, providing confidence for the finally quoted error margin of T_c.

The final estimates of T_c are summarized in Table 1. For n = 1, we have T_c = 6.680 300(10), which is consonant with and improves over T_c = 6.680 263(23) [37] and marginally agrees with T_c = 6.679 63(36) [38] and 6.680 339(14) [13]. For n = 2, our determination T_c = 3.314 437(6) significantly improves over T_c = 3.31 [39,40] and 3.314 [41]. For n = 3, our result T_c = 2.198 79(2) rules out T_c = 2.192(1) from a high-temperature expansion [42].

Finite-size scaling of the two-point correlation

We remark that FSS analyses for g(L/2, L) have already been performed in [16] with the formula g(L/2, L) = AL⁻²[ln(L/2 + B)]^1/2 (A and B are constants) and in [13] with a similar formula. These FSSs in the literature correspond to the first scaling term in (14). Hence, (14) serves as a forward step for complete FSS by involving the scaling term |$v$|₂L⁻², which arises from the Gaussian fixed point.

Finite-size scaling of the magnetic susceptibility

Figure 3.

The magnetic fluctuations (a) χ₀ and (b) χ₁ rescaled by L² versus lnL on a log-log scale for the critical Ising and Heisenberg models. The black lines in (a) represent the least-squares fits, and the red line in (b) denotes a constant.

Open in new tabDownload slide

We note that previous studies based on a FSS without high-order corrections produced estimates of |${\hat{y}}_h$||$(=\hat{p}/2)$|⁠, considered to be consistent with |${\hat{y}}_h=1/4$| [38,43–45]. The maximum lattice size therein was L_max = 24, four times smaller than L_max = 96 of the present study. In particular, it was reported [43] that |$2\hat{y}_h=0.45(8)$| and |$4\hat{y}_h=0.80(25)$|⁠. Nevertheless, we find that the fit |$\chi _0= q_1 L^2 ({\rm ln}L)^{2 \hat{y}_h}$| by dropping the correction term q₂L² would yield |$\hat{y}_h=0.21(1)$| (Ising), 0.20(1) (XY), and 0.19(1) (Heisenberg), which are smaller than and inconsistent with the predicted value |$\hat{y}_h=1/4$|⁠. This suggests the significance of q₂L² in the susceptibility χ₀, which arises from the r dependence of g(r, L).

Finite-size scaling of the magnetic fluctuations at nonzero Fourier modes

We consider the magnetic fluctuations χ₁ with |k₁| = 2π/L and χ₂ with |$|{\bf k}_2|= 2 \sqrt{2} \pi /L$|⁠. We have compared the FSSs of χ₀, χ₁ and χ₂ in Fig. 1(c) and (d) for the critical 4D XY model. As L increases, χ₁L⁻² and χ₂L⁻² converge rapidly, suggesting the absence of a multiplicative logarithmic correction. This is in sharp contrast to the behavior of χ₀L⁻², which diverges logarithmically. For the Ising and Heisenberg models, the FSS of the fluctuations at nonzero modes is also free of a multiplicative logarithmic correction (Fig. 3(b)).

Surprisingly, we find that the scaled fluctuations χ₁L⁻² ≈ 0.15 are equal within error bars for the Ising, XY and Heisenberg models.

Furthermore, we show in Fig. 4 χ₁ and χ₂ versus T for the 4D XY model. We observe that the magnetic fluctuations at nonzero Fourier modes reach maximum at T_c and that the χ₁L⁻² (χ₂L⁻²) data for different Ls collapse well not only at T_c but also for a wide range of |$(T-T_c)L^{y_t}$| with y_t = 2.

Figure 4.

Data collapses for the magnetic fluctuations (a) χ₁ and (b) χ₂ rescaled by |$L^{2 y_h-d}$| and |$L^{y_t}$| (y_h = 3, y_t = 2, d = 4) for the 4D XY model. The insets show the scaled fluctuations versus T, and the dashed lines denote T_c.

Open in new tabDownload slide

DISCUSSIONS

We propose formulae (11) and (12) for the FSS of the O(n) universality class at the upper critical dimensionality, which are tested against extensive MC simulations with n = 1, 2, 3. From the FSS of the magnetic fluctuations at zero and nonzero Fourier modes, the two-point correlation function and the Binder cumulant, we obtain complementary and solid evidence supporting (11) and (12). As byproducts, the critical temperatures for n = 1, 2, 3 are all located up to an unprecedented precision.

An immediate application of (12) is to the massive amplitude excitation mode (often called the Anderson–Higgs boson) due to the spontaneous breaking of the continuous O(n) symmetry [46], which is at the frontier of condensed matter research. At the pressure-induced quantum critical point (QCP) in the dimerized quantum antiferromagnet TlCuCl₃, the 3D O(3) amplitude mode was probed by neutron spectroscopy and a rather narrow peak width of about 15% of the excitation energy was revealed, giving no evidence for the logarithmic reduction of the width-mass ratio [3]. This was later confirmed by a quantum MC study of a 3D model Hamiltonian of O(3) symmetry [5,6]. Indeed, (12) provides an explanation why the logarithmic-correction reduction in the Higgs resonance was not observed at the 3D QCP. In numerical studies of the Higgs excitation mode at the 3D QCP, the correlation function g(τ ≡ |τ₁ − τ₂|) is measured along the imaginary-time axis β, and numerical analytical continuation is used to deal with the g(τ) data. In practice, simulations are carried out at very low temperature β → ∞, and it is expected that g(τ) ≍ τ⁻² for a significantly wide range of τ. Furthermore, it is the τ-dependent behavior of g(τ), instead of the L dependence, that plays a decisive role in numerical analytical continuation.

In the thermodynamic limit, the two-point correlation function decays as |$g(r) \sim r^{-2} \tilde{g}(r/\xi )$|⁠, where the scaling function |$\tilde{g}(r/\xi )$| quickly drops to zero as r/ξ ≫ 1. It can be seen that no multiplicative logarithmic correction exists in the algebraic decaying behavior. On the other hand, as the criticality is approached (t → 0), the correlation length diverges as |$\xi (t) \sim t^{-1/2}|{\rm ln}t|^{\hat{\nu }}$|⁠, and |$\hat{\nu } = {(n+2)}/{2(n+8)} > 0$| implies that ξ diverges faster than t^−1/2 [30,33]. Since the susceptibility can be calculated by summing the correlation as |$\chi _0 \sim \int _0^{\xi } g(r)r^{d-1}dr \sim \xi ^2$|⁠, we have |$\chi _0(t) \sim t^{-1} |\ln t|^{\hat{\gamma }}$| with |$\hat{\gamma } = 2 \hat{\nu }$|⁠. The thermodynamic scaling of χ₀(t) can also be obtained from the FSS formula (10) or (11), which gives |$\chi _0(t,L) \sim L^{2y_h-4} (\ln L)^{2\hat{y}_h} \tilde{\chi _0} (tL^{y_t} (\ln L)^{\hat{y}_t})$|⁠. By fixing |$tL^{y_t} (\ln L)^{\hat{y}_t}$| at some constant, we obtain the relation |$L \sim t^{-1/y_t} |\ln t|^{-\hat{y}_t/y_t}$|⁠. Substituting this into the FSS of χ₀(t, L) yields |$\chi _0(t) \sim t^{\gamma } |\ln t|^{\hat{\gamma }}$| with γ = (2y_h − 4)/y_t and |$\hat{\gamma } = -\gamma \hat{y}_t+2\hat{y}_h$|⁠. With |$(y_t, y_h, \hat{y}_t, \hat{y}_h)\,=\,(2, 3,\, {(4-n)}/{(2n+16)}, \frac{1}{4})$|⁠, we have γ = 1 and |$\hat{\gamma } = {(n+2)}/{(n+8)}$|⁠. The thermodynamic scaling with logarithmic corrections has been demonstrated in [4] in terms of the magnetization m of an O(3) Hamiltonian.

For the critical Ising model in five dimensions, an unwrapped distance r_u was introduced to account for the winding numbers across a finite torus [18]. The unwrapped correlation was shown to behave as |$g(r_{\rm u}) \sim r_{\rm u}^{2-d} \tilde{g} (r_{\rm u}/\xi _{\rm u})$|⁠, where the unwrapped correlation length diverges as ξ_u ∼ L^d/4. This differs from typical correlation functions that are cut off by a linear system size of approximately L. We expect that at d_c = 4 the unwrapped correlation length diverges as |$\xi _{\rm u} \sim L (\ln L)^{\hat{y}_h}$|⁠, which gives the critical susceptibility as |$\chi _0(L) \sim L^2 (\ln L)^{2\hat{y}_h}$|⁠.

Besides, (12) is useful for predicting various critical behaviors. As an instance, it was observed that an impurity immersed in a 2D O(2) quantum critical environment can evolve into a quasiparticle of fractionalized charge, as the impurity-environment interaction is tuned to a boundary critical point [47–49]. Equation (12) precludes the emergence of such a quantum-fluctuation-induced quasiparticle at the 3D O(2) QCP.

We mention an open question about the specific heat of the 4D Ising model. The FSS formula (10) predicts that the critical specific heat diverges as C ≍ (ln L)^1/3. By contrast, an MC study demonstrated that the critical specific heat is bounded [37]. The complete scaling form (11) is potentially useful for reconciling the inconsistency.

Finally, it would be possible to extend the present scheme to other systems of critical phenomena, as the existence of upper critical dimensionality is a common feature therein. These systems include the percolation and spin-glass models at their upper critical dimensionality d_c = 6. We leave this for a future study.

METHODS

Throughout the paper, the raw data for any temperature T and linear size L are obtained by means of MC simulations, for which the Wolff cluster algorithm [35] and the Prokof’ev–Svistunov worm algorithm [36] are employed complementarily. Both algorithms are state-of-the-art tools in their own territories.

The O(n) vector model (1) in its original spin representation is efficiently sampled by the Wolff cluster algorithm, which is the single-cluster version of the widely utilized nonlocal cluster algorithms. The present study uses the standard procedure of the algorithm, as in the original paper [35] where the algorithm was invented. In some situations, we also use the conventional Metropolis algorithm [50] for benchmarks. The macroscopic physical quantities of interest have been introduced in aforementioned sections for the spin representation.

The two-point correlation function for the XY model (n = 2) is sampled by means of the Prokof’ev–Svistunov worm algorithm, which was invented for a variety of classical statistical models [36]. By means of a high-temperature expansion, we perform an exact transformation for the original XY spin model to a graphic model in directed-flow representation. We then introduce two defects for enlarging the state space of directed flows. The Markov chain process of evolution is built upon biased random walks of defects, which satisfy the detailed balance condition. It is defined that the evolution hits the original directed-flow state space when the two defects meet at a site. The details for the exact transformation and a step-by-step procedure for the algorithm have been presented in [51].

Acknowledgements

Y.D. is indebted to valuable discussions with Timothy Garoni, Jens Grimm and Zongzheng Zhou.

FUNDING

This work was supported by the National Natural Science Foundation of China (11774002, 11625522 and 11975024), the National Key R&D Program of China (2016YFA0301604 and 2018YFA0306501), and the Department of Education in Anhui Province.

AUTHOR CONTRIBUTIONS

J.-P.L., K.C. and Y.D. designed the research and established the formulae for finite-size scaling. J.-P.L., W.X. and Y.S. performed the simulations. J.-P.L., W.X., Y.S. and Y. D. analyzed the results. J.-P.L. and Y.D. wrote the manuscript. All the authors participated in the revisions of the manuscript.

Conflict of interest statement. None declared.

REFERENCES