Tensor network analysis of critical coupling in two dimensional $\phi^{4}$ theory

We make a detailed analysis of the spontaneous $Z_{2}$-symmetry breaking in the two dimensional real $\phi^{4}$ theory with the tensor renormalization group approach, which allows us to take the thermodynamic limit easily and determine the physical observables without statistical uncertainties. We determine the critical coupling in the continuum limit employing the tensor network formulation for scalar field theories proposed in our previous paper. We obtain $\left[ \lambda / \mu_{\mathrm{c}}^{2} \right]_{\mathrm{cont.}} = 10.913(56)$ with the quartic coupling $\lambda$ and the renormalized critical mass $\mu_{\mathrm{c}}$. The result is compared with previous results obtained by different approaches.


Introduction
The two dimensional real φ 4 theory has been extensively studied as the simplest model with a spontaneous symmetry breaking which plays important roles in particle physics. Although the continuous symmetry as in the complex φ 4 theory cannot be broken in two dimensions [1], the real φ 4 theory has a discrete Z 2 -symmetry that can be spontaneously broken [2]. The critical coupling, which separates two phases, has been studied by various approaches with and without lattice simulations . Further improvements will be achieved in lattice simulations if one can reach the thermodynamic limit far beyond the system size the conventional Monte Carlo methods can access to.
The tensor renormalization group (TRG) [24] enables us to take the thermodynamic limit easily: the computational cost of the TRG depends on the space-time volume logarithmically while that of the Monte Carlo methods is proportional to the volume. Furthermore, the TRG method does not contain any stochastic procedure to evaluate physical quantities. Then, using the TRG, we can observe precisely the symmetry breaking effects with larger volumes avoiding statistical uncertainties.
Since the TRG was originally proposed for the two dimensional Ising model, several improved algorithms [25][26][27] and extensions to fermion systems [28,29] have been developed. Those techniques are introduced into the path-integral formulation of the field theory and many studies have been carried out [30][31][32][33][34][35][36][37][38][39][40][41][42]. It is, however, not straightforward to apply the TRG to lattice scalar theories because the indices of tensor are essentially the field variables so that their dimension is infinite for scalar fields unlike the Ising case. 1 Recently we have proposed a new tensor network (TN) formulation for constructing a finite dimensional tensor of scalar theories [43], where we directly discretize the pathintegral expression of the scalar field theory using the Gaussian quadrature formula instead of the orthonormal functions proposed in the previous work [44]. As the order of Gauss-Hermite quadrature is increased, we can properly reproduce the partition function of a two dimensional free scalar theory. The method is already well established in SUSY quantum mechanics with an interacting scalar [45], and two dimensional theories with interacting scalars now fall within the scope of application of this method.
In this paper we determine the critical coupling of the two dimensional lattice φ 4 theory employing the TRG method with our TN formulation. The critical coupling is determined from the power law behavior of the susceptibility for the expectation value of the scalar field near the critical point. We extrapolate the critical coupling to the continuum limit and compare it to other recent results obtained by different approaches [16,17,20,23].
This paper is organized as follows. We present the two dimensional continuum and lattice φ 4 theories to fix our notation in section 2. A TN formulation of the expectation value of the scalar field is given in section 3. In section 4 we first explain the strategy and the methodology to evaluate the critical coupling from the expectation value of the field. The results of susceptibility are then presented in section 4.2, and the critical couplings at several lattice spacings are listed in table 1 of section 4.3. We show the final result of the critical coupling in the continuum limit in section 4.4. Section 5 is devoted to summary and discussion. In appendix A the coarse-graining algorithm for the expectation value of the field is explained. The systematic error originated from the Gaussian quadrature is discussed in appendix B.

Two dimensional φ 4 theory
The Euclidean continuum action of the two dimensional real φ 4 theory is defined as with a real scalar field φ (x) ∈ R, the bare mass µ 0 , and the quartic coupling constant λ > 0. This theory is super renormalizable. The mass renormalization is required only at one-loop level while the renormalization of the coupling is not necessary in all orders of the perturbation theory. The suffix zero to denote the bare parameter is used only for the mass. This model has the Z 2 -symmetry (φ → −φ) classically, but it may be spontaneously broken at a quantum level. The expectation value of the scalar field φ (x) , which does not actually depend on the coordinate x from the translational invariance, is an order parameter of the symmetry breaking: φ (x) = 0 for the symmetric phase and φ (x) = 0 for the symmetry broken phase.
The dimensionful coupling λ gives a typical scale of this model while a dimensionless coupling normalized by µ 2 characterizes its physical behavior, where µ is the renormalized mass. At the tree level, two phases and the sign of κ 0 = λ/µ 2 0 are in one-to-one correspondence, and the critical point is κ 0,c = 0. The negative one corresponds to the broken phase. However, the critical point can change with the quantum corrections. The lattice simulations play a crucial role to determine the renormalized critical point κ c beyond the perturbation theory.
Let us define the lattice theory on a square lattice where a positive integer L is the lattice size. The lattice spacing a is assumed to be a = 1, and we suppress it unless necessary. We label the lattice site (n 1 , n 2 ) simply as n. The lattice scalar field φ n is defined on the site n and satisfies the periodic boundary conditions, whereρ is the unit vector of the ρ-direction. The lattice action is given by (2.5) Note that µ 0 and λ are dimensionless quantities written as a 2 µ 2 0 and a 2 λ, respectively. The mass renormalization can be calculated by the lattice perturbation theory. We employ the same renormalization condition as in refs. [2,16]. The renormalized mass µ is defined as which is the one-loop self energy. We use eq. (2.6) to determine the renormalized critical coupling in later section.
3 Tensor network formulation for two dimensional lattice φ 4 theory The expectation value of the scalar field can be expressed as a TN form. We present the detailed procedure following the formulation given in our previous work [43].
We begin with the partition function The action S h is given by where a constant external field h is introduced to investigate the phase transition taking h → 0 after L → ∞. The Boltzmann weight can be expressed as a product of local factors, The function f (φ 1 , φ 2 ) is a symmetric matrix with the continuous indices φ 1 , φ 2 ∈ R, and the continuity makes numerical treatments hard.
To define a finite dimensional tensor, we use the Gauss-Hermite quadrature formula, in which the integral of a target function g (x) with the weight function e −x 2 is approximated by a discrete sum as where y α (α = 1, 2, . . . , K) is the α-th root of the K-th Hermite polynomial H K (y) and w α is the associated weight. 2 K is the order of approximation. If g (x) is a polynomial function with the degree of 2K − 1 or less, the Gauss-Hermite quadrature is known to reproduce the exact result. For general functions, the convergence of the Gauss-Hermite quadrature is not obvious, and in this paper we numerically check the convergence for our particular case. By replacing each integral in eq. (3.2) by eq. (3.6) and using eq. (3.4), we obtain a discretized version of the partition function where {α} denotes n∈Γ L K αn=1 . This equation can be regarded as a TN representation with the bond dimension K. However, truncation of the bond dimension is practically required for the latter procedures because the computational cost of the TRG is more expensive than that of the Gauss-Hermite quadrature. Then we decompose the K × K matrix M αβ = f (y α , y β ) using the singular value decomposition (SVD): 3 where σ i are the singular values sorted as σ 1 ≥ σ 2 ≥ · · · ≥ σ K ≥ 0 and U, V are unitary matrices. The range of this summation is truncated. We will discuss this point at the end of this section. Finally we have Now let us turn to the expectation value of φ, Since Z is already expressed as a uniform TN as shown in eq. (3.9), the remaining task is to represent the numerator Z 1 as a TN. The insertion of the operator φ n in eq. (3.12) affects the form of the tensor as where an extra y α is inserted compared to eq. (3.10). We call the modified tensorT impurity tensor. Repeating the same procedures as for the partition function, we obtain the following TN form: in which one of T s in eq. (3.10) is replaced by the impurity tensorT .
The contraction of tensor indices in eq. (3.9) is exactly taken from 1 to K. In actual computations, we initially truncate the bond dimension from K to D on account of computational complexity; namely we redefine {x,t} as n∈Γ L D xn=1 D tn=1 , and then the partition function in eq. (3.9) is initially approximated. In coarse-graining steps the number of tensors in the network is reduced with truncating the bond dimensions. In this study we keep the size of tensors as D. 4 Now we have two parameters to control the accuracy of the approximations, the degree of the Gauss-Hermite quadrature K and the size of tensors D, and we will check the stability of numerical results with respect to D and K in appendix B.

Numerical results
In this section the numerical results are shown step by step. First, the technical details related to the Gauss-Hermite quadrature and the TRG is summarized in section 4.1. In section 4.2 we extract the susceptibility of the expectation value of φ with effectively taking the thermodynamic limit. Using the susceptibility we determine the critical parameters in section 4.3, and in section 4.4 we take the limit λ → 0 to obtain the critical coupling in the continuum limit.

Methods
The expectation value φ is obtained by the ratio of Z (K) and Z 1 (K), which are individually evaluated by the TRG method with D ∈ [16,64]. One can use the standard technique for Z (K) and needs a little ingenuity for Z 1 (K). The coarse-graining procedure of the tensor network with an impurity tensor is shown in refs. [34,44,47,48], and the detailed technique is also described in appendix A.
We remark that the systematic errors can be estimated by investigating the K-and D-dependences of the numerical results. There are three sources of systematic errors in the present case: (a) discretization of the scalar fields with the K-point Gauss-Hermite quadrature, (b) truncation of the bond dimension of the initial tensor discussed in the last paragraph in section 3, and (c) truncation in coarse-graining steps. The effect of (a) is appeared as K-dependence of the results while the errors from (b) and (c) are observed through the D-dependence of the results. In figure 1 we plot the singular values of the SVD in eq. (3.8). Clear hierarchy of the singular values assures that the modification of contraction range that is explained at the end of the last section effectively work. We fix K = 256 that is large enough so that the systematic errors associated with the choice of K are to be smaller than those of D as discussed in appendix B. In section 4.3 we investigate the D-dependence of the critical coupling λ/µ 2 c and estimate the systematic errors by measuring small fluctuations which occur when D changes.

Susceptibility
The critical point is given by the peak position of the susceptibility for φ defined by where φ h,L denotes the expectation value measured for given constant external field h and lattice size L. 5 Equation (4.1) is reduced to χ = lim h→0 lim L→∞ φ h,L /h since φ 0,L = 0. It is rather complicated procedure to take the double limits with respect to h and L numerically. We first find a constant behavior of φ h,L /h as L increases, which effectively represents φ h,∞ /h. Then, the reduction of h allows us to extract a converged value corresponding to the susceptibility χ.
In figure 2 we show φ h,L /h as a function of L for several values of h. One can see that it becomes a constant for larger volumes with L ≥ 10 6 , which effectively corresponds to the thermodynamic limit within the systematic errors.  Figure 3 shows the h-dependence of φ h,∞ /h with the same parameter set as in figure 2. It seems to converge to a constant for sufficiently small h ≤ 10 −11 . We can extract the limit of φ h,∞ /h to the vanishing h, which corresponds to the susceptibility χ defined in eq. (4.1) with φ 0,L = 0. Figure 4 shows the zoomed version of figure 3, and there the values of φ h,∞ /h lie on a quadratic function. Since φ h,∞ is an odd function of h, it behaves as φ h,∞ = c 1 h + c 3 h 3 + · · · with coefficients c i (i = 1, 3, . . .), and one can find that the ratio does φ h,∞ /h = c 1 + c 3 h 2 + · · · . This is the reason why the data in figure 4 show the quadratic behavior. Then we simply fit them using a quadratic function to determine χ, and its error is defined as the difference between the obtained susceptibility and φ h,∞ /h at h = 10 −11 . 7 Similar behaviors are observed for other parameter sets employed in this work (see table 1). We obtain the susceptibility for all the parameter sets taking the same procedure as here.

Critical mass for each λ
The scaling behavior of the susceptibility χ as a function of µ 0 near the critical point in the symmetric phase is given by 7 This error will propagate to the next analysis, and in the end it will become that of λ/µ 2 c for given λ. The propagated error, however, does not affect the final result because it is very small compared to a fluctuation originated from D as shown in figure 6.  with the critical bare mass µ 0,c , the critical exponent γ, and a constant A. This formula can be used to determine µ 0,c by fitting the numerical result of χ. It may be more beneficial to express the above formula in the following way: Since the φ 4 theory is expected to belong to the universality class of the two dimensional Ising model, whose critical exponent is exactly known to be γ Ising = 7/4 = 1.75 [49], the modified expression could make it easier to check the consistency between γ and γ Ising ; one can clearly see that γ/γ Ising = 1 in the right hand side of eq. (4.4) if χ −1/γ Ising behave as a straight line with respect to µ 2 0 around the critical point. Figure 5 shows the susceptibility as a function of µ 2 0 at λ = 0.05 and D = 32. It is clear that the results are on a straight line which implies γ = γ Ising . We determine the values of µ 2 0,c by fitting the susceptibility using the expression of eq. (4.4) with µ 2 0,c and A ′ as free parameters and fixed γ = γ Ising . The bare critical mass µ 2 0,c for another value of λ is obtained by repeating the procedures above. 8 Then the renormalized critical mass µ 2 c is determined by solving the self-consistent equation (2.6) using µ 2 0,c and λ as inputs. Table 1 summarizes the results of µ 2 0,c obtained for 0.005 ≤ λ ≤ 0.1 and 16 ≤ D ≤ 64. We also present the dimensionless critical coupling λ/µ 2 c . Since χ 2 /d.o.f. ≤ 1 for all λ and D, the legitimity of the linear fitting shown in figure 5 are confirmed, and we can consider that fixing γ = γ Ising is a reasonable assumption. The critical bare mass µ 2 0,c monotonically decreases as λ approaches to zero while it does not show any smooth dependence on D at each λ.
In figure 6 we plot the D-dependence of λ/µ 2 c for λ = 0.05. In the large D region λ/µ 2 c converges to a value while oscillating. This kind of behavior is also observed in the TRG computations of other models such as the Ising and the Potts models [50]. As also shown in figure 6 we determine the central value of λ/µ 2 c as an average between the maximum and minimum values in the oscillating region. The error is estimated using the half width of the oscillation. Although the critical coupling for each D has a rather small systematic error as given in table 1, its magnitude is negligible compared to the fluctuation of the critical coupling with the variation of D. Thus, we present the results of the renormalized critical coupling λ/µ 2 c only with the systematic errors associated with D.  (27) 0.37 0.005 52 −0.012841211 (17) 10.82077 (11) 0.013 0.005 56 −0.0128303006(66) 10.749790 (43) 0.28 0.005 60 −0.0128349426(48) 10.779915 (31) 0.0065 Table 1: Bare critical mass µ 2 0,c , renormalized critical coupling λ/µ 2 c , and χ 2 /d.o.f. of the linear fitting described in the main body of the text for given λ and D. Errors are originated from that of the susceptibility as explained in section 4.2. Figure 6. D-dependence of λ/µ 2 c for λ = 0.05. The error estimated from the fluctuation is also shown as the colored band.

Continuum limit of the critical coupling
We finally obtain the value of the critical coupling in the continuum limit as . (4.5) Note that the continuum limit is understood as a 2 λ → 0. Figure 7 shows the λ-dependence of the critical coupling with the systematic errors. Let us take the continuum limit (λ → 0) of the critical coupling by a linear extrapolation. The result of the linear fit for λ ≤ 0.05 shown in the figure gives a reasonable chi-squared value: χ 2 /d.o.f. ≈ 0.026. The critical coupling in the continuum limit is found to be λ µ 2 c cont.   [16,17,20,23] for small λ. We have reached the smallest lattice spacing: λ = 0.005. The error bar, however, is relatively large compared to the latest Monte Carlo result around λ ∼ 0.01. Note that around λ = 0.0312 one can directly compare four data points given by different methods. As a result, the values by Schaich and Loinaz [16], Bronzin et al. [23] and us seems roughly consistent with each other while that by Bosetti et al. [20] is sizably deviated from the three results.

Summary and outlook
We have studied the two dimensional lattice φ 4 theory using the TRG with our new formulation of making a finite dimensional tensor for scalar field theories given in ref. [43]. The TRG method, whose computational cost depends on the space-time volume logarithmically, allows us to access to large volume lattices so that the thermodynamic limit and the power law behavior of the susceptibility can be precisely studied. This work is the first study to determine the critical coupling in the continuum limit employing the TRG method.  [20], and Bronzin et al. [23]) and this work. At λ = 0, data points are horizontally shifted to ensure the visibility. Note that the results by Wozar and Wipf cannot be compared at non-zero λ since they used the SLAC derivative for scalar bosons, but, in the continuum limit (λ = 0), their results are consistent with naively discretized ones within errors.
The dimensionless critical coupling was obtained as  with sufficiently small error albeit the simplest form of the TRG is employed. The error is the systematic one coming from the finite D effects. Our result shows a reasonable consistency with the recent results obtained by different approaches. The simplest TRG algorithm suffers from the growth of the systematic errors around the criticality. Alternative coarse-graining procedures such as the tensor network renormalization (TNR) [26] and loop-TNR [27] might be useful to obtain more precise results. These methods effectively work around the critical point and, in principle, are applicable to any two dimensional model irrespective of the field contents. We could expect further improvements for the accuracy of the critical coupling.

A Coarse-graining of tensor network including impurity tensors
We describe the coarse-graining algorithm for a tensor network with an impurity tensor such as Z 1 (K) in eq. (3.14), which is given with a fixed integer D for truncating the SVD.
In the first few steps in the TRG, N 2 might not be larger than D. In this case D is replaced by N 2 for N 2 < D. 10 Note that the coarse-grained tensor T new is defined by gathering S [i] (i = 1, 2, 3, 4) and its bond dimension is at most D. One can see that the tensor network of T is approximately expressed as a coarse network of T new . The number of tensors decreases at each coarse-graining step since T provides two S's while T new is made of four S's. Repeating the procedure above over and over again, we finally find that the network is expressed as a single tensor whose indices are contracted with itself. The value of the tensor network is thus obtained. Now let us turn to the coarse-graining of the tensor network with an impurity tensor. As explained below, an important point is that the number of impurity tensors does not increase beyond four (the impurity tensors are located in four corners of a plaquette) even if the coarse-graining step is repeated many times.
Suppose we consider a coarse-graining procedure for the tensor network with four impurity tensors as shown in figure 10 (left), where the impurity tensors are represented Figure 9. Coarse-graining step for uniform tensor network. Circles represent tensors, and closed indices are contracted. A local block on the network is taken into account here. In the first step, every rank-four tensor is decomposed into two rank-three tensors. By contracting the four rankthree tensors, a new rank-four tensor is obtained in the second step. Here, decomposed tensors with tilde are defined byS [1(3) as the red, orange, yellow, and green circles while the normal tensors are given by the blue circles. In the first step, the rank-four tensors are decomposed into rank-three tensors. As seen in the middle panel in figure 10, those rank-three impurity tensors form a closed loop. Thanks to this structure, after taking the contraction, one obtains a network which again contains only four impurity tensors as shown in the right panel in figure 10. In this way, one can suppress the spread of the impurity tensors. Figure 10. Coarse-graining of tensor network with four impurity tensors. Red, orange, yellow and, green circles denote the impurity tensors, and the normal tensors are shown by blue ones. Note that the impurity tensors colored in red, orange, yellow, and green are inherited in both decomposition and contraction processes. A key point is that the number of impurity tensors does not change through the coarse-graining.
Starting from a tensor network with a single impurity tensor such as eq. (3.14), one obtains a network with two neighboring impurity tensors after a course-graining step. In the next step it is coarse-grained to yield a network with three impurity tensors located on three corners of a plaquette. In this way, the number of impurity tensors stepwisely increases until it reaches four. After that, as shown above, the number of impurity tensors does not change anymore.

B Systematic error from discretization of scalar fields
As shown in section 3 we discretize the scalar fields using the Gauss-Hermite quadrature in order to create a finite dimensional tensor, which gives the approximated expressions for the partition function and expectation values. The order of approximation is characterized by the degree of the Hermite polynomial K. In this appendix we discuss the systematic error associated with the approximation, namely the finite K-effect, for the expectation value of the field. 11 For an one dimensional integral in eq. (3.6), the Gauss-Hermite quadrature gives an exact result for 2K − 1 degree (or less) polynomial target functions. In the case of multi dimensional integrals, unfortunately, the convergence is not trivial. We have to check whether or not the numerical results do not have large systematic errors from the finite K-effect with K = 256, which is employed in this study.
In figure 11 the K-dependence of φ is presented.It is clear that the results do not depend on K for K ∈ [64, 256] while much larger D-dependence is observed. Thus we can conclude that our choice K = 256 is large enough, and the main source of systematic errors is the truncation of the SVD for constructing the initial tensor and the coarse-graining steps in the TRG.