Calculation of fermionic Green functions with Grassmann higher-order tensor renormalization group

We develop calculational method for fermionic Green functions in the framework of Grassmann higher-order tensor renormalization group. The validity of the method is tested by applying it to three-dimensional free Wilson fermion system. We compare the numerical results for chiral condensate and two-point correlation functions with the exact ones obtained by analytical methods.


I. INTRODUCTION
The tensor renormalization group (TRG) is a numerical renormalization group method in the tensor network scheme, which was originally introduced by Levin and Nave [1]. This method consists of two steps. We first make a tensor network representation for the target physical system. After that, we repeat the coarse-graining of the tensor network based on the singular value decomposition (SVD). The basic concept of this local transformation is similar to the block spin transformation [2]. An advantage of the TRG method is that it does not suffer from the sign problem inherent in the Monte Carlo methods. Subsequently to the appearance of TRG, it was extended for Grassmann valued tensor networks, and the new method is called Grassmann TRG (GTRG) [3,4]. The TRG and GTRG methods have been applied to the two-dimensional relativistic quantum field theories in elementary particle physics: lattice φ 4 model [5], lattice Schwinger model [4,6], and the finite-density lattice Thirring model [7]. Meanwhile, the TRG method itself has been improved by disentangling short-range correlations in tensor network renormalization (TNR) [8] and loop-TNR [9].
Although the original idea of TRG was restricted to two-dimensional systems, Xie et al. invented a new method called higher-order TRG (HOTRG) with the use of higher-order SVD (HOSVD), which is formally applicable to any higher-dimensional system [10]. Its effectiveness is also demonstrated in Ref. [10] by applying it to two-and three-dimensional Ising models. So far, the HOTRG method has been applied to several two-and three-dimensional systems: three-dimensional Potts models [11], two-dimensional XY model [12], twodimensional O(3) model [13], and two-dimensional CP(1) model [14].
In order to use the tensor network scheme in lattice QCD, we need an algorithm which can be applied to fourdimensional relativistic fermion systems. The authors in Ref. [15] have formulated the Grassmann HOTRG (GHOTRG) based on the GTRG. It opened the way to higher-dimensional fermion systems but did not allow us to calculate physical quantities except the free energy and its differentiation. In this article we develop a method for direct calculation of fermionic Green functions in the framework of GHOTRG and test the validity of the algorithm using the relativistic free Wilson fermion system, for which the exact values of Green functions are analytically calculable.
This paper is organized as follows. We make a brief review for the GHOTRG algorithm in Sec. II. In Sec. III we explain calculational method of fermionic Green functions with the GHOTRG method. Numerical results for Green functions in the relativistic free Wilson fermion system are presented in Sec. IV. Section V is devoted to summary and outlook.

II. GRASSMANN HIGHER-ORDER TENSOR RENORMALIZATION
The detailed explanation for the GHOTRG algorithm is already given for the two-dimensional case in Ref. [15]. Although the extension to three-dimensional cases is rather straightforward, we give a brief review to fix the notations. We assume a = 1 for the lattice units in the following.

A. Tensor network representation
We first present a tensor network representation of the partition function for the three-dimensional free Wilson fermion with the Wilson parameter r = 1 which is described in detail in Ref. [15]. The action with the fermion mass m is given by where n = (n 1 , n 2 , n 3 ) is the lattice coordinate and the fermion fields ψ n ,ψ n have two spinor components. A arXiv:1711.08121v1 [hep-lat] 22 Nov 2017 tensor network representation of the partition function is written by where the indices are given by representing that ψ has two spinor components. The tensor T is defined as where the bosonic tensor T have normal number elements whose explicit form is given in Ref. [15], and the Grassmann tensor G is defined as dη i,n,s η i,n,s = dη i,n,sηi,n,s = 1.
The original fermion fields ψ,ψ have already been integrated out, while another set of independent Grassmann variables η,η have been introduced. Because of the Grassmann nature, elements of the bosonic tensor take non-trivial values only when the indices satisfy the following condition: x i,n,s + x i,n−î,s mod 2 = 0.
Otherwise the element is zero. This allows us the block diagonalization of the bosonic tensors explained in Appendix A, which contributes to save the computational and memory costs.

B. Renormalization of bosonic tensor
We explain the renormalization procedure of bosonic tensors, which is mostly based on the original HOTRG [10]. Here we focus on the coarse-graining alonĝ 1-direction depicted in Fig. 1. The coarse-graining for the other directions can be done in a similar way.
We first make up a contracted tensor M by where the integrated indices are defined as for i = 2, 3. σ is a sign factor originating from anticommutation properties in the renormalization of Grassmann part, which is explained in the next subsection.
Next we apply the eigenvalue decomposition to M + and M − defined as and obtain a unitary matrix U ± and the corresponding eigenvalues λ ± : Then, while assuming that λ ± i are obtained in descending order, we define ± as with a given D cut that one can choose. If + < − then U 2 = U + is adopted as an isometory to be used for the coarse-graining, and vice versa. In a similar manner an isometory U 3 for3-direction is obtained. Restricting the new indices 1 ≤ x b 2,n , x b 3,n ≤ D cut , the coarse-grained bosonic tensor T is defined as a product of U 2 , U 3 and M : For the lattice coordinate of the index l for the tensor, here we unify n, n+1 to n and define the new unit vector as1 =1 +1 as shown in Fig. 2.

C. Renormalization of Grassmann part
We now explain the renormalization procedure of the Grassmann part focusing on the coarse-graining along1direction.

A. Chiral condensate
We first explain a method to calculate the chiral condensate defined as where V is the system volume. Thanks to the translational invariance, Eq. (33) is rewritten as Due to the presence ofψ n0 ψ n0 , the bosonic tensor on the site n 0 is modified. Now we introduce "impure tensor" 1 and define the bosonic part in it as I. Then the impure tensor with the Grassmann part is defined as Then Eq. (33) is expressed as In the following the renormalization algorithm is formulated based on this expression. First of all, the renormalization procedure for the normal tensors is explained in Sec. II. On the other hand, in case of the the impure tensor, the contracted tensor of Eq. (10) is modified as With the use of the isometries that are the same as in the normal tensors, we obtain the following new bosonic tensor: Figure 3 gives a schematic view of the coarse-graining procedure for the impure tensor. The renormalization 1 The coarse-graining for the impure tensor in the framework of TRG is discussed in Ref [16,17]. procedure for the Grassmann part is similar to that for the normal tensors. Thus the coarse-grained impure tensor is expressed as

B. Fermion two-point correlation functions
We now explain a calculational method for the fermion two-point correlation function defined as In this case, the tensors on the site n 1 and n 2 become the impure tensors, which are denoted bȳ respectively. Inclusion of odd numbers of fermions at the site n 1 and n 2 modifies the constraint described in Eq. (9) as x i,n,s + x i,n−î,s mod 2 = 1.
Then Eq. (42) is expressed as In calculation of the fermion two-point correlation functions, the renormalization procedures are classified into three types according to the positions of the impure tensors. We explain them in order.
In this case two impure tensors are merged into one impure tensor as shown in Fig. 4. The contracted tensor is obtained by where σ is the sign factor modified by Eq. (45) as follows: We use the isometries which are the same as that for the normal tensors. After the two impure tensors are merged, we obtain a tensor network with a single impure tensor described in the previous subsection.
We exemplify the case of n 1 +2 = n 2 as depicted in Fig. 5. For each impure tensor we construct the corresponding contracted tensor: M impure The unitary matrices U 2 , U 3 are obtained by the renormalization procedure for the normal tensors. In addition, we apply the eigenvalue decomposition to M + and M − defined by M impure,+ and obtain a unitary matrix U ± and the eigenvalues λ ± and ± . in the same manner of Eq. (16) and Eq. (17). If + < − then W = U + is adopted, and vice versa. The new bosonic tensorsĪ , I are obtained by using U 2 , U 3 and W : The rest of the procedure is similar to that for the normal tensors.

Other cases
This is the case that two impure tensors are renormalized independently. The contracted tensors for impure tensors are given bȳ The isometries are obtained by the renomalization procedure for the normal tensors. The rest of the procedure is similar to that for the normal tensors.

IV. NUMERICAL RESULTS
In this section, we present the numerical results. The anti-periodic boundary condition is imposed for one of the directions. We employ it to test the calculational method of the fermionic Green functions on a 256 3 lattice.

A. Free energy
We first present the results for the free energy of threedimensional free Wilson fermion system. Figure 6 compares the numerical results with the exact values for the free energy as a function of the fermion mass m. The latter is obtained by an analytical calculation in the momentum spaces. It is hard to detect the deviation between the measured values with D cut ≥ 14 and the exact ones. In Fig. 7 we directly plot the magnitude of the error defined by which clarify the difference between the numerical results with the exact values. We observe that the error decreases monotonically in the range of 0 ≤ m ≤ 2 as D cut increases from 14 to 22. It is also noted that the error diminishes monotonically as the fermion mass becomes larger. The maximal error is about 0.2% at m = 0 and the minimal one is about 0.02% at m = 2.

B. Chiral condensate
We summarize the results for the chiral condensate in Figs. 8 and 9. The former plots the numerical results as a function of the fermion mass m together with the The error reduces monotonically as the fermion mass increases. This is a similar tendency as observed for the error of the free energy. The maximal error is about 5% at m = 0 and the minimal one is about 0.1% at m = 2.
It should be noted that the error does not change so much against the increment of D cut . The reason may be that we use the common unitary matrices, which are optimized for the normal tensors but not for the impure tensors, in the coarse-graining of both tensors.

C. Fermion two-point correlation functions
We consider C 11 (n 1 , n 2 ) and C 12 (n 1 , n 2 ) as representative cases. Figures 10, 11, 12 and 13 show real and imaginary parts of C 11 (n 1 , n 2 ) and C 12 (n 1 , n 2 ), respectively, as a function of R, which is defined by n 2 = n 1 + R ·1. In Figs. 10 and 12 we find that the accuracy of real part becomes worse on m = 0 as observed for the free energy and the chiral condensate. Another important finding is that the errors become relatively large at larger R for all the choices of m. Here we point out that the magnitude of the imaginary part may give us an estimate of the magnitude of the error for the real part, since the imaginary part is expected to vanish theoretically. Figures 11 and 13, in fact, show that the magnitude of the real part at R = 8 are comparable with that of the imaginary part. This indicates that the signal-to-noise ratio of the real part at the large distance (around R = 8)   may be of order of one. In our case, where the correlation length is small, it is difficult to accurately compute the correlation function at large separation. For the case of the larger correlation length, however, one may measure the long separated correlation function, whose signal-tonoise ratio is expected to be better.

V. SUMMARY AND OUTLOOK
We have developed a method to calculate fermionic Green functions in the framework of GHOTRG using the impure tensors. The validity of the method is tested by calculating the chiral condensate and the fermion twopoint correlation functions in the three-dimensional free Wilson fermion system on a 256 3 lattice, where we employ 0 ≤ m ≤ 2 for the fermion mass parameter and 14 ≤ D cut ≤ 22 for the tensor indices. In most cases we find good consistency between the numerical results and the exact values: their deviation is at most a few % level. Only exception is the fermion two-point correlation functions at a long distance. It may be due to the difficulty of evaluating tiny correlations left at a long distance.
The method developed in this paper can be applied to any physical systems which contain fermions. We plan to apply it to more complicated systems including gauge interactions. It is also important to improve the accuracy of the method. It could be effective to use an idea of the higher-order second renormalization group (HOSRG) [10,18] for the formulation of impure tensors.