Phase transition of four-dimensional Ising model with higher-order tensor renormalization group

We apply the higher-order tensor renormalization group (HOTRG) to the four-dimensional ferromagnetic Ising model, which has been attracting interests in the context of the triviality of the scalar $\phi^4_{d=4}$ theory. We investigate the phase transition of this model with HOTRG enlarging the lattice size up to $1024^4$ with parallel computation. The results for the internal energy and the magnetization are consistent with the weak first-order phase transition.


Introduction
It is well known that the critical behavior of the Ising model on the higher-dimensional hypercubic lattice is well explained with the mean-field theory. In dimensions larger than four, the effect of the background fluctuations becomes negligible and the model in the critical region exactly obeys the mean-field exponents [1,2]. At the upper critical dimension, however, multiplicative logarithmic corrections are added to the leading scaling behavior of the mean-field theory. Some of these corrections were derived by the perturbative calculation with the renormalization group method [3]. Since the Ising model is specified by the infinite coupling limit of the single-component scalar φ 4 4 theory, the model in four dimensions has been attracting the interest of particle physicists for a long time in the context of the triviality of the O(4) scalar φ 4 4 theory, which corresponds to the scalar sector of the standard model describing the generation of gauge boson and fermion mass through the Higgs mechanism [4][5][6][7][8][9][10]. There is also a recent study to discuss the triviality of the O(N ) φ 4 4 theory with the higher-loop beta function [11][12][13]. A numerical study of the Ising model on a hypercubic lattice serves as a non-perturbative test of the triviality [14,15]; if the leading scaling behavior is the mean-field type and it is modified only by the multiplicative logarithmic factor, one obtains a supporting evidence for the triviality. In fact, the numerical investigation based on the Monte Carlo simulation has successfully caught the mean-field exponents [16][17][18][19][20], but there remains some controversy over the appearance of the logarithmic corrections [19][20][21][22][23]. Actually no Monte Carlo study has confirmed the logarithmic correction in the scaling behavior of the specific heat, which is (ln |t|) 1/3 with t the reduced temperature expected from the perturbative renormalization group analysis. This is mainly because the cubic root of logarithmic divergence is too weak to detect by the finite size scaling analysis, or the specific heat may be actually bounded [19]. Indeed, the finite volume effect of the four-dimensional Ising model had been investigated from various viewpoints [9]. A detailed Monte Carlo study has found serious finite-volume effect due to non-trivial boundary effects in the four-dimensional Ising model [20]. From the viewpoint of numerical calculation, it could be possible that there remain some unrevealing aspects in the phase transition of this model and it should be worth trying different approaches other than the Monte Carlo method.
For this purpose we employ the tensor network scheme to investigate the four-dimensional classical Ising model. This scheme has various types of numerical algorithms [24], which can be divided into two streams: Hamiltonian approach and Lagrangian one. The latter enables us to evaluate the partition functions directly via tensor network representation. A typical algorithm is the tensor renormalization group (TRG) [25], which was originally proposed by Levin and Nave for the two-dimensional Ising model. The TRG method has been successfully applied to the two-dimensional field theories with the path-integral formulation in the particle physics [26][27][28][29][30][31][32][33][34][35][36][37][38][39]. The higher-order TRG (HOTRG) [40] is an improvement of TRG with the extension to higher dimensions. One of attractive features in TRG and HOTRG is that we are allowed to directly study the thermodynamic properties; we can systematically increase the system size by repeating the coarse-graining steps in the algorithms. Although earlier studies with HOTRG are restricted to two-and three-dimensional systems [41][42][43][44][45][46][47][48][49][50][51][52][53][54][55], including the three-dimensional classical Ising model [40], the algorithm itself is readily extended to a four-dimensional lattice. In this paper, we employ the HOTRG method to investigate the phase transition of the classical Ising model on the four-dimensional hypercube. The accuracy of HOTRG is controlled by the bond dimension D cut , which is varied up to 14 in this study. In order to investigate the phase transition of the model, we measure the internal energy and the magnetization through the evaluation of the tensor network with some impurity located at the center of hypercube. This paper is organized as follows. In Sec. 2 we briefly review the HOTRG method and explain an approach to evaluate the internal energy and the magnetization of the Ising model. We present numerical results in Sec. 3 and discuss the properties of the phase transition. Sec. 4 is devoted to summary and outlook.

HOTRG w/ and w/o impurity
The partition function of the four-dimensional ferromagnetic Ising model is given by where σ i is the two-state classical spin variable on the lattice site i, ij specifies the sum over all the nearest-neighboring spin pairs, β is the inverse temperature 1/T and h is the external magnetic field. The subscript N is the size of a system. Based on the eigenvalue decomposition T = U ΛU T , one defines the eight-rank local tensor locating on each lattice site as The indices of these tensors are called bond indices. Now, we obtain the tensor network representation of Eq. (2.1) as where we assume the periodic boundary condition and the right-hand side means all the bond indices are contracted so as to restore the model defined on the four-dimensional hypercube. One way to evaluate Eq. (2.3) is HOTRG with the use of the higher-order singular value decomposition (HOSVD) [40]. In the HOTRG procedure, nearest two local tensors along the x-, y-, z-and t-directions are mapped to the coarse-grained one sequencially. Hence the lattice size is reduced by a factor of 2 after each step of coarse-graining. After repeating n steps of coarse-graining, one obtains the partition function with the system size of N = 2 n ; that is, The right-hand side is again the sum over all the bond indices so as to restore the structure of the four-dimensional lattice model with the periodic boundary condition and this is easily done by defining the trace of the coarse-grained tensor as There are two ways to evaluate the expectation values such as the internal energy and the magnetization. One is the numerical differentiation with respect to β and h. The other is the direct evaluation of the expectation value using the corresponding tensor network representation. For instance, we can obtain the internal energy through the evaluation of the nearest-neighbor local energy term σ i σ j with the HOTRG method as follows. We first define the additional local tensor as With the use of this local tensor, the tensor network representation for the local energy is given by where S i,j represent the tensors on the lattice sites i and j, respectively, and T pure tensor. The denominator is evaluated by the plain HOTRG method. In order to coarse-grain the impure tensor network of Eq. (2.7), we assume the local energy term is fixed at the center of lattice (i = 1, j = 2 = 1 +ŷ withŷ the unit vector in y-direction) during the HOTRG calculation. At the first step, we define the coarsegrained impure tensor S (1) 1 by contracting two initial impurities. For simplicity we give the corresponding expression in the two-dimensional case; where U (1) is a block-spin transformation determined within the original algorithm of HOTRG [40]. In the following steps, S and T (n) 2 . We again show the corresponding formula in two-dimensional case for simplicity; (2.9) Finally, the local energy is approximately given by . (2.10) The meaning of the trace is the same as in Eq. (2.5). Since the original model has the translational invariance, σ i σ j × d, where d is the dimensionality, should give the absolute value of internal energy. The one-point function σ i to measure the magnetization is also evaluated in the same way. In this case, we are allowed to apply the four-dimensional counterpart of Eq. (2.9) from the first coarse-graining step because the initial expression of σ i has the form where S (0) i locates only on the lattice site i (i = 1) and k runs the rest of N − 1 sites. After sufficient iterations, σ i is evaluated by Thanks to the translational invariance, σ i directly corresponds to the spatial average of the Ising spin. Note that Eqs. (2.10) and (2.12) have the same expression, but they are evaluated by different coarse-graining procedures. In Ref. [40], computational costs and memory space requirements in two-and threedimensional HOTRG are given. Computational costs are O(D 7 cut ) and O(D 11 cut ) and memory space requirements are O(D 4 cut ) and O(D 6 cut ), respectively. In straightforward expansion of the HOTRG algorithm in Ref. [40] to four dimensions, the computational cost is O(D 15 cut ) and the memory space requirement is O(D 8 cut ). In our implementation, the computational cost in each process is O(D 13 cut ) and the memory space requirement in each process is O(D 7 cut ). Reduction of the order of computational cost is achieved by using D 2 cut processes. This implementation is basically based on an idea which will be shown in Ref. [56] .We have carried out a detailed measurement of the internal energy and the magnetization with D cut = 13 employing the fine resolution of the temperature ∆T = 6.25 × 10 −6 around the transition temperature. We have repeated the calculation with D cut = 14 to confirm the qualitative features obtained with D cut = 13. But, in this case the temperature resolution remains coarser as ∆T = 3.0 × 10 −5 due to the computational cost. In the following we focus on the results with D cut = 13.
As found in Sec. 3 all the measured physical quantities seem to lose the volume dependence beyond n ≈ 30 so that the lattice size of N = 2 40 = 1024 4 is large enough to be taken as the thermodynamic limit.

Numerical results
We first evaluate the free energy with the plain HOTRG method. The convergence behavior is investigated by defining the following quantity; (3.1) Figure 1 shows a typical convergence behavior of ln Z N in the vicinity of the transition temperature. We observe that δf decreases monotonically as a function of D cut . We now turn to the determination of the transition temperature. Let us assume that one has just obtained the coarse-grained tensor T (n) 1;xx yy zz tt whose coarse-graining direction was t-direction. Choosing a D cut × D cut matrix as  we calculate the following quantity; Tr A (n) 2 , which counts the number of the largest singular value of A (n) . This is an indicator of the symmetry-breaking [57]. We calculate X (n) iteratively until it converges. A typical convergence behavior of X (n) is shown in Fig. 2. Notice that we sequentially redefine A (n) corresponding to the direction of coarse-graining in the practical calculation. Figure 3 shows the transition temperature T c as a function of D cut . The error bars, provided by the temperature resolution, are all smaller than the corresponding symbols. Since T c (D cut ) is estimated by X (n) with sufficiently large n, typically beyond n = 30, there remains little finite-volume effect. In this work, we have obtained T c (D cut = 13) = 6.650365(5) on 1024 4 lattice. The recent Monte Carlo study [19] obtained β c = 0.1496947(5) corresponding to T c = 6.68026 (2), which shows a slight deviation from our result with HOTRG up to D cut = 13. Note that the value of T c in Ref. [19] was obtained by the infinite-volume extrapolation using the results on relatively small lattices with L 4 ≤ 80 4 . Let us move on to the evaluation of the internal energy, which can be obtained by numerical differentiation or the coarse-graining of the impure network of Eq. (2.7). We have compared both methods varying the temperature resolution and found that the latter successfully keeps the numerical accuracy as the resolution becomes finer. In the following, we show the results with the impure tensor method. Figure 4 traces the volume dependence of the internal energy with D cut = 13. The converging behavior toward the thermodynamic limit is clearly observed. Since the system size N is given by 2 n , a hypercubic structure is restored in the condition of n mod4 = 0. Figure 5 shows the internal energy as a function of temperature for various lattice sizes with D cut = 13. In the case of n ≥ 24 (L ≥ 64), a finite jump emerges with mutual crossings of curves between different volumes around the transition temperature. These are characteristic features of the first-order phase transition as discussed in Ref. [58]. The similar volume dependence and a finite jump at L ≥ 64 have been also confirmed in case of D cut = 14. The numerical value of the finite jump ∆E(D cut = 13) in the infinite-volume limit is ∆E(D cut = 13) = 0.0034(5), which is obtained by the linear extrapolation toward the transition temperature both from the low and high temperature regions. The resolution of the temperature at the boundary between the two phases is ∆T = 6.25 × 10 −6 .
We also investigate the spontaneous magnetization, which is an order parameter to detect the symmetry-breaking phase. Figure 6 shows a typical volume dependence of magnetization toward the thermodynamic limit. We have evaluated σ i with h = 1.0×10 −9 and 2.0 × 10 −9 at each temperature and coarse-graining step. After taking the infinite-volume limit we extrapolate the value of σ i toward the h → 0 limit. Figure 7 shows the resulting spontaneous magnetization as a function of temperature. The transition temperature is consistent with both estimates by X (n) and the internal energy. We have observed a finite jump in the magnetization, whose numerical value is obtained by the linear extrapolation toward the transition temperature both from the low and high temperature regions; ∆m(D cut = 13) = 0.037 (2).
The resolution of the temperature at the boundary between the two phases is again ∆T = 6.25 × 10 −6 . Note that we have tried several choices of the external field other than h =   O(10 −9 ) and confirmed that the behavior of the magnetization is robust against the change of the magnitude of h.

Summary and outlook
We have analized the phase transition of the four-dimensional ferromagnetic Ising model employing HOTRG on L 4 ≤ 1024 4 lattices. The transition temperature is successfully determined by measuring the degeneracy of the largest singular value of the pure tensor. We have also investigated the temperature dependence of the internal energy and magnetization with the impure tensor method. We have found a finite jump for the internal energy together with mutual crossings of curves between different volumes around the transition temperature. A finite jump is also observed in the magnetization. These are characteristic features of the first-order phase transition. The numerical results obtained by the impure tensor method are consistent with the weak first-order phase transition. The resulting estimate for the transition temperature in the thermodynamic limit shows a slight deviation from the recent Monte Carlo prediction [19] obtained from the infinite-volume extrapolation of the data on relatively small lattices up to 80 4 . In future investigation, the HOTRG calculation with D cut > 14 should allow us to achieve a direct and essential improvement of this study. Our impure tensor method can be also improved by considering all the patterns of coarse-graining for the network including some impurities [59]. Another possible approach is to develop the best optimization of the Frobenius norm of impure tensor, which would be a realistic way to improve our impure tensor algorithm from the viewpoint of the computational cost of HOTRG.