A new technique to incorporate multiple fermion flavors in tensor renormalization group method for lattice gauge theories

We propose a new technique to incorporate multiple fermion flavors in the tensor renormalization group method for lattice gauge theories, where fermions are treated by the Grassmann tensor network formalism. The basic idea is to separate the site tensor into multiple layers associated with each flavor and to introduce the gauge field in each layer as replicas, which are all identified later. This formulation, after introducing an appropriate compression scheme in the network, enables us to reduce the size of the initial tensor with high efficiency compared with a naive implementation. The usefulness of this formulation is demonstrated by investigating the chiral phase transition and the Silver Blaze phenomenon in 2D Abelian gauge theories with Nf flavors of Wilson fermions up to Nf = 4.


Introduction
Nonperturbative computation in fermionic systems has always been challenging due to the anti-commuting nature of Grassmann variables.In Monte Carlo methods, the Grassmann variables have to be integrated out first, yielding the fermion determinant det M , which makes the computation very time-consuming since the matrix M has a size proportional to the system size V .While the computational cost can be made O(V ) by using the pseudofermion technique with an appropriate Hybrid Monte Carlo algorithm, the calculation is still typically a few orders of magnitude more time-consuming than corresponding bosonic systems.Moreover, in many interesting fermionic systems such as finite density systems, strongly-correlated electron systems, and theories with chiral fermions, the fermion determinant becomes complex, which causes the notorious sign problem in conventional Monte Carlo methods.In order to overcome this problem, various methods such as the complex Langevin method [1,2,3,4,5], the Lefschetz thimble method [6,7,8,9,10,11,12], and the density of state method [13,14,15] have been developed.However, each method has its pros and cons, and many models still remain out of reach.
All these problems associated with fermionic systems can be solved beautifully in the tensor network method [16,17,18,19,20,21,22,23,24,25,26,27,28,29], which is not a statistical approach based on important sampling.This method was first introduced to handle many-body systems in condensed matter physics with the main application to the calculation of the ground state based on the variational principle [30,31,32,33,34,35].However, it can also be used to directly compute the partition function with some procedures based on coarse-graining, which is similar in spirits to the real-space renormalization group and hence the name, "the tensor renormalization group (TRG) method".Notably, it enables the computation of the partition function with a computational cost that grows only logarithmically with the system size.Although the original TRG method was proposed for a two-dimensional bosonic system [17], improved versions were subsequently developed [18,19,20,29], and it has also been generalized to higher dimensional systems [21,22,23,24] and to fermionic systems, where Grassmann variables are treated directly [25,26,27,28,24] unlike in Monte Carlo methods.Using this "Grassmann tensor network", one does not have to deal with the fermion determinant, and furthermore, the sign problem does not exist in the method from the outset because it is not a statistical approach.
A recent achievement of the TRG method is its application to gauge theory, in particular in the parameter regions that are not accessible to Monte Carlo methods due to the sign problem.Notable examples include the 2D gauge theories with a θ term [36,37,38], 2D SU(2) gauge-Higgs model [39], one-flavor Schwinger model [27,40,41], 2D QCD [42], 3D SU(2) gauge theory [43], 4D Z K gauge-Higgs models [44,45] and so on.Among these applications, gauge theories with matter fields are of particular importance since typically they are not exactly solvable.However, the TRG method has so far been applied only to the case with one fermion flavor.When there are many flavors of fermions on a single lattice site, one encounters a problem that the size of the local Hilbert space, and thus the size of the initial tensor, grows exponentially with the number of flavors.This prevents us from studying theories with multiple flavors including QCD, which is a non-Abelian gauge theory with two (or three) flavors of light quarks.
In this paper, we propose a new technique that makes it possible to incorporate multiple flavors of fermions in gauge theory within the Grassmann tensor network formalism.The main idea is to separate the initial tensor into multiple layers associated with each flavor.Since the fermions with different flavors are interacting with the same gauge field, the interaction in the flavor direction becomes non-local after integrating out the gauge field.In order to avoid this problem, we introduce the gauge field in each layer as replicas and identify them all later.Once the system can be described by a tensor network, which is one dimension higher than the original theory due to the flavor direction, one can use the standard coarse-graining technique to compute the partition function.Our method is expected to be useful also in applying the TRG to the domain-wall formalism [46,47] for chiral fermions regarding the flavor direction in our method as the extra space-time dimension.We also introduce an efficient compression scheme that further reduces the size of the initial tensor drastically, especially at large K.For N f = 1, the performance of our method is found to be as good as in the previous calculations for the Schwinger model [27,40,41].
The usefulness of this formulation is demonstrated by applying it to Abelian gauge theories in two dimensions with N f flavors of Wilson fermions.First, we investigate the chiral phase transition in the N f = 2 case, which shows that the result obtained in the Z K gauge theory converges to the U(1) result obtained by using the Monte Carlo method [48] as we increase K. Next, we investigate the Silver Blaze phenomenon in the case of finite fermion density up to N f = 4.
While this paper was being completed, we encountered a paper [49] which addresses the same issue of treating multiple fermion flavors in the TRG method.There the initial tensor for all the flavors is separated into multiple layers by using the matrix product decomposition, which requires the memory of order O(e cN f ).In contrast, the memory cost of our method is of order O(1) since the layers are separated analytically and all the layers are identical.This memory cost reduction enabled the investigation of gauge theories, which was not possible in Ref. [49].
The rest of this paper is organized as follows.In section 2, we explain our basic idea to implement multiple flavors in the TRG method.In particular, we derive the initial tensor for 2D Abelian gauge theories with Wilson fermions as an example.In section 3, we describe how we perform the procedures of the TRG method using the initial tensor.In particular, we discuss how we compress the initial tensor efficiently by inserting isometries and explain how we perform coarse-graining in the flavor space.In section 4, we present the numerical results obtained by our method.After showing some results of the performance tests concerning the initial tensor compression and the coarse-graining procedure in the flavor direction, we demonstrate the usefulness of our method by investigating the chiral phase transition and the Silver Blaze phenomenon in 2D Abelian gauge theories.Section 5 is devoted to a summary and discussions.In Appendix A, we give a brief review of the Grassmann tensor network.In Appendix B, we describe the coarse-graining algorithm in detail.In Appendix C, we discuss the generalization of our technique to a model with local multi-flavor interactions, which is important, in particular, in applying our method to the domain-wall formalism.

The basic idea to implement multiple flavors
In this section, we explain our basic idea to implement multiple flavors in the TRG method.First, we split the system into multiple layers by introducing replicas of gauge fields for each flavor, and then we construct the initial tensor based on the Grassmann tensor network formalism including the flavor direction.

Splitting the system into multiple layers
For simplicity, we will describe our idea in the case of Abelian gauge theory G ⊆ U(1) on a two-dimensional square lattice Λ 2 with the lattice spacing a although it can be readily applied to higher dimensions and non-Abelian cases.The gauge field A x,µ is represented on the lattice by the link variable U x,µ = exp(iaA x,µ ) ≡ exp(iφ x,µ ) ∈ G. Let us then consider the lattice action with N f flavors of Wilson fermions given by where β = 1/(ga) 2 is the inverse gauge coupling and and γ ν are the 2D gamma matrices given, for instance, by the Pauli matrices γ 1 = σ 1 and γ 2 = σ 2 .The fermion fields are represented by the 2-component Dirac spinors ψ for each flavor α with charge q α .We also introduce dimensionless chemical potential μα and mass mα for each flavor α.
Since fermions have 4N f internal degrees of freedom, the local Hilbert space at a given lattice site has dimension D = 2 4N f .Consequently, the size of the initial tensor grows as D 4 = 2 16N f .It is thus beneficial to separate different flavors from each other to avoid the exponential growth of the tensor size with N f .To that end, we split the link variables into N f replicas and define the partition function as where the new action S (α) is given as This enables us to split the Boltzmann weight into several layers corresponding to each replica, which are linked by the delta function.This decomposition into multiple layers has recently been considered in Ref. [49] but with the matrix product decomposition instead of using the delta function.Let us emphasize, however, that separating the layers by hand from the beginning as we do here is crucial in avoiding completely the singular value decomposition, which can be both memory-consuming and computationally expensive in lattice gauge theories.While we have not introduced local interactions among different flavors, our tensor construction can be generalized to such cases as described in Appendix C.

Constructing the tensor network
We treat the fermion fields by the Grassmann tensor network [28] (See Appendix A for the details.).Let us rewrite the fermion action in the following form x,+ν ψ gauge fields and fermions indexing Figure 1: a) The two-dimensional tiling of the site tensor, which is composed of four subtensors.The shaded region represents a plaquette.b) The three-dimensional tiling of the site tensor for N f = 3.Each layer corresponds to some flavor α. c) The site tensor (2.14) with the gauge fields and fermionic variables.d) The site tensor with the index for each bond.In the diagrams, fermionic legs with an arrow pointing away from the tensor are non-conjugated fermions, while those with an arrow pointing into the tensor are conjugated fermions.
Thus we arrive at where we have defined the tensor x,1 )w(φ x,1 +φ (2.16) 1 This follows from the identity e −h θθ = dη dη e −ηη− θη+hηθ for a Grassmann-even constant h and one-component Grassmann-odd numbers θ, θ, η, η, which can be generalized to the multi-component case in a straightforward manner.
and introduced a short-hand notation for the Grassmann integral Note that the tensors P , L and S are associated with the plaquettes, links, and sites, respectively.The connection of these tensors is shown in Fig. 1-c).
Performing the integral (2.17) symbolically2 , we obtain the tensor S in the form of a polynomial of link fermions given as x,1 , φ x,2 , φ where the coefficient depends on the gauge link variables.Here we have introduced where θ 1 and θ 2 are the two components of η and k 1 , k 2 ∈ {0, 1} represent the occupation number of the components.Since some of these link fermions connect the same pair of sites, it is convenient to combine the fermion indices as (I a , J a ) → K a with the prescription ζK 4 x,4 = (−) p(J 4 ) ηI 4 x− 2,+ 2η where p(J) is the Grassmann parity of the Grassmannn number η J defined by with j a being the fermion occupation number of the a-th component.The sign factor (−) p(J) is introduced for the consistency of Grassmann tensor contraction (See (A.17)-(A.18).).
Using the Grassmann index notation (A.1), the tensor S can be expanded as (omitting the site index x to avoid redundancy) with the sign factors given in (2.24) and (2.25).Here the index i in φ (α) i refers to the index of the quadrature node in (2.11).The site and orientation indices of the link variables φ x,µ are omitted.
The tensor P can be rewritten in terms of the quadrature indices as ) . (2.29) The tensor L (2.16), which depends on two fields, has actually five legs because it is connected to two plaquettes, two sites, and to the global gauge field φ x,µ .Therefore, it can be written in terms of the quadrature indices as To summarize, the coefficient of the site tensor (2.14) is given by T (α) In the above expression, the indices with subscripts 1, 2, 3, and 4 are associated with the legs pointing in the direction + 1, + 2, − 1 and − 2, respectively, whereas m and n are associated with φ x,2 and φ x,1 , respectively.The schematic representation of the site tensor is given in Fig. 1-d).
3 The procedures of the TRG In this section, we describe how we perform the procedures of the TRG method using the initial tensor derived in the previous section.In particular, we discuss how we compress the initial tensor efficiently by inserting isometries and explain how we perform coarse-graining in the flavor space.

Compressing the initial tensor
As one can see from the expression (2.31), the initial tensor for lattice gauge theories has typically a large dimension due to the existence of many legs.It is therefore important to compress its size first before we perform the coarse-graining procedure.Here we use the compressing procedure based on the higher-order SVD, which is frequently used in HOTRG-type algorithms [50,21] (See also Appendix B.1.).The first step of the compressing procedure is to "squeeze" the legs of the S tensor (ζ a , i a ) into a smaller leg ξ a using the hybrid isometries that merge a fermionic leg and a bosonic leg into one fermionic leg as where the Hermitian conjugate is defined in (A.25) and the contraction ζν ζν of the bond (ζ ν , ζν ) is defined in (A.20).Here we have inserted four isometries; namely V 1 and V 2 are inserted in the inner links between S and the Kronecker delta nodes L, whereas V 3 and V 4 are inserted in the outer links (See Fig. 2-b).In order to obtain these isometries, we first define the tensors (repeating indices are not summed; See Fig. 3.) and construct the Grassmann Hermitian matrices ) Figure 3: The Q tensors used in the computation of the four isometries.Q 1 and Q 2 are used to compute the inner isometries, while Q 3 and Q 4 are used to compute the outer isometries in Fig. 2-b).The five-legged L-tensor in Fig. 2-a) are replaced by the three-legged nodes since both of these diagrams result in the same M matrices (3.6)-(3.13).Note that the position of S and the Kronecker delta nodes are swapped in the inner and outer cases for a given axis.
where the indices in the parenthesis are combined into a single index with the prescription (A.17)-(A.18),and the Hermitian conjugate of a two-legged Grassmann tensor is defined in (A.25).Note also that the coefficient tensor with reordered indices should have appropriate sign factors due to the permutation of Grassmann-odd variables.
Then we diagonalize the Hermitian matrices M ±a , which gives the unitary matrices V ±a .By comparing the singular value spectra of M +a and M −a , we define the isometry in (3.1) by the unitary matrix V ±a that corresponds to the one with the fastest falling spectrum.
Corresponding to (3.1), we have to attach the same isometries on L as Despite having many indices, these tensors are actually sparse due to the Kronecker deltas, which makes it more efficient to use sparse array algorithms to perform the calculation.
The schematic representation of the construction of S ′ , L x and L y is given in Fig. 2-b.Now that we have rewritten the site tensors in terms of four sparse subtensors, we can proceed to perform the final compression, which further reduces the size of the site tensor.This can be done in exactly the same way as we have done in compressing S. We first contract all the tensors together as where the repeated indices i 3 and i 4 are not summed over.Then we construct the Hermitian matrices By diagonalizing these matrices, we obtain the hybrid isometries (U ±µ )ξ i;ϕ that further combine the bosonic and fermionic indices together.Just like (V ±a )ζ j;ξ , the combined index ϕ is also truncated with the bond dimension χ c .Thus the final compressed tensor becomes (See Fig. 2-c and -d.) Here ζ's are the four legs pointing in the 2D space, while m and n are the bosonic legs joining different flavor layers.

Some comments on the computational cost
In this subsection, we make some comments on the computational cost.The most computationally expensive parts of the process are the contraction of M ±a (3.6)-(3.13)and M±a The decomposition of the site tensor in the modified HOTRG.
(3.19)-(3.22),which contain 11 and 12 loops, respectively.However, the number of loops in the computation of M ±a can be reduced since the index m is always the same as one of the other indices due to the Kronecker delta, which can effectively reduce the depth to 10 loops.In terms of complexity, the cost of the computation of M ±a and M±a are 16 5 K 5 and D 5 ξ K 7 , respectively, where D ξ is the bond dimension of the fermionic legs ξ.At the technical level, we can do a few more things to further reduce the computational cost.First, we implement another compression on the fermionic legs of the original S tensor (2.27) to reduce the bond dimension of the leg ζ a , which is of dimension D ζ = 16.In addition, we improve the speed of compression significantly by storing the tensors as sparse arrays and performing contractions with a sparse matrix-based algorithm [51].

The coarse-graining procedure
In order to perform the coarse-graining procedure in the flavor direction, we implement a modified version of HOTRG.Note first that the site tensor T ′ η 1 η 2 η3 η4 ;mn has six indices.Four of them are the fermionic indices in the 2D space directions, while m and n are the indices for the gauge link variables.Since the site tensors in the flavor direction are connected with the Kronecker deltas, the bosonic bond degrees of freedom are maximally entangled in the flavor direction, which implies that we cannot insert isometries to compress these bonds further.For this reason, we have to perform the coarse-graining procedure in the flavor direction before the coarse-graining in the two-dimensional space.
Since the size of the compressed initial tensor (3.23) grows like K 2 with K being the bond dimension of the bosonic legs, the coarse-graining can become costly with the traditional HOTRG algorithm [21] as we increase K. (See Appendix B for a review of the HOTRG and related methods.)We therefore modify the HOTRG method by first performing the decomposition such that we separate the legs into three groups based on their axes as where the legs along the x axis are in the X tensor, the legs along the y axis are in the Y tensor, and the legs along the flavor axis are in the P and Q tensors.With this decomposition, the isometries along the x and y axis are computed using the tensor X and Y, respectively.In order to obtain the partition function, we first perform coarse-graining procedure in the z (flavor) direction using the Grassmann higher-order TRG (gHOTRG) algorithm described in Appendix B.1, which gives us the N f -flavor tensor T ′′ ψ 1 ψ 2 ψ3 ψ4 ;mn .Then we take the trace of the bosonic indices which is described schematically in Fig. 5. Next, we perform the coarse-graining in the 2D plane using the Grassmann TRG (gTRG) which is described in detail in Appendix B.2.The bond dimensions for the gHOTRG in the flavor direction is χ f , while the bond dimension for gTRG in the 2D plane is χ xy .After performing the coarse-graining procedure sufficiently many times, we take the trace of the final tensor with anti-periodic boundary conditions in the imaginary time direction to obtain the partition function as where σ I is the sign factor defined in (A.11).
All of the computations in this paper are done using a Python package grassmanntn [52]

Numerical results
In this section, we present our numerical results obtained by the method introduced in the previous sections.First, we perform performance tests concerning the initial tensor compression described in section 3.1 and the coarse-graining procedure in the flavor direction described in section 3.3.Then we demonstrate the usefulness of our method by investigating the chiral phase transition and the Silver Blaze phenomenon in 2D Abelian gauge theories.In what follows, we assume that all the flavors of fermions have the same charge q α = q, mass mα = m and chemical potential μα = μ.

Performance tests
Let us first demonstrate the efficiency of the initial tensor compression.For that, we compute log Z for various parameters and measure the errors by comparing the results obtained with and without compression.It is found that the relative error of the compression is less than 10 −15 for all the cases if we choose the compression bond dimension χ c = 64.The efficiency of the compression for this choice is summarized in Table 1.It is clear that our compression scheme is efficient, in particular for large K, where the original tensor can easily become too large to be handled by the currently available computers.
Figure 6: The singular value spectrum associated with the HOTRG isometry truncation in the free electron gas model and the Z 2 gauge theory.Here we show only the spectra of the x-axis truncation, which is identical to that of the y-axis truncation.
Next, we discuss the performance of coarse-graining in the flavor direction.Here we perform the gHOTRG up to N f = 4 for the free electron gas model (K = 1) and the Z 2 gauge theory with β = 0, q = 1, m = 1 and μ = 0.In Fig. 6, we plot the singular value spectra associated with the SVD when the isometries (U x and U y in (B.14) and (B.15)) are used during the step N f = 1 → 2 and 2 → 4 with χ f = 64.One can see that the tail of the singular value spectrum grows quickly with N f , which indicates that fermions from different layers have strong degeneracy.Note that introducing gauge interaction makes the singular value spectrum decays faster.For the calculations in the subsequent subsections, we use χ f = 64 for N f = 2 and χ f = 32 for N f = 4 for the flavor coarse-graining, and χ xy = 64 for the two-dimensional coarse-graining.

The chiral phase transition
In order to demonstrate the usefulness of our method, we first apply it to the chiral phase transition in two-flavor Z 2 , Z 4 , and U(1) gauge theories.Let us define the hopping parameter by Although the chiral symmetry is broken explicitly by the Wilson term in (2.3), it is expected to be restored at some critical hopping parameter κ c .We can easily identify κ c in the TRG   Monte Carlo [48] 0.3296983759 Table 2: The critical hopping parameter for various gauge theories with N f = 2 at β = 0. Our U(1) result is computed with the 4-nodes Gauss-Legendre quadrature.
method by the location of the peak of chiral susceptibility given as a function of the hopping parameter κ, where the derivatives can be taken numerically.The peak of the chiral susceptibility exhibits a critical behavior in the large volume limit as demonstrated in Ref. [27] with N f = 1.In Fig. 7. we observe similar behavior in various gauge theories with N f = 1 and 2 at β = 0, where we have used χ f = χ xy = 64.For the U(1) case, we use the 4-nodes Gauss-Legendre quadrature to discretize the group integral.In the large-K limit, Z K gauge theory converges to the U(1) gauge theory.For N f = 1, the critical hopping parameter κ c = 0.3806(1) obtained by our method in both Z 2 and Z 4 theories is consistent with the result κ c = 0.380665(59) obtained for the U(1) ≃ Z ∞ theory in Ref. [27], which indicates that the convergence occurs already at K = 2 for β = 0.For N f = 2, we obtain κ c for the Z 2 , Z 4 , and U(1) theories at β = 0 and compare with the Monte Carlo results for the U(1) theory [48] in table 2.

Silver Blaze phenomenon at finite density
Next we apply our method to the finite density case with multiple flavors.In particular, it is expected that physical observables in the thermodynamic limit and at zero temperature are independent of the chemical potential up to some threshold due to the gapped spectrum in the confined phase.This is known as the Silver Blaze phenomenon [53], which is difficult to reproduce by Monte Carlo methods due to the sign problem.We investigate this phenomenon by our method to demonstrate that the sign problem is indeed solved.
Here we consider the Z 2 gauge theory and the free electron model with β = 1, q = 1, m = 1.We calculate the pressure and the number density Figure 9: ρ(μ)/N f is plotted against μ for V = 128 × 128 in the free electron model and the Z 2 gauge theory with N f = 1, 2, and 4. The non-monotonicity may be attributed to the truncation effect, which becomes more significant for larger N f .
where the derivative in (4.4) is taken numerically.In Fig. 8 we plot the differential pressure ∆P (μ) = P (μ) − P (0) and the number density ρ(μ) against μ up to the volume V = 128 × 128 in the Z 2 gauge theory with N f = 2.We observe a clear Silver Blaze phenomenon for μ ≲ 0.8.We also find that the number density saturates to the value ρ = 2 for μ > 1.4 as expected from the number of degrees of freedom at each lattice site.
In Fig. 9 we plot ρ(μ)/N f against μ in the free electron model and the Z 2 gauge theory with N f = 1, 2, and 4. The number density saturates to the value ρ = N f for all cases.Note also that in the free electron gas model, ρ(μ)/N f are expected to be the same for all N f , which is not the case in the Z 2 gauge theory due to interactions.

Summary
In this paper, we have proposed a new technique to incorporate multiple flavors in the TRG method for lattice gauge theories.The problem of the initial tensor, which grows in size exponentially with the number of flavors N f , has been overcome by separating the initial tensor into N f layers with replicas of the gauge field for each layer, which are identified later.This effectively makes the system one dimension higher due to the flavor direction.Consequently, the tensor can still be large, in particular, due to the gauge field legs in the extra dimension.In order to overcome this problem, we use a compression scheme, which proceeds in two steps by truncating first the subtensors and then the whole tensor.We have shown that this enables us to compress the size of the initial tensor by many orders of magnitude without sacrificing the accuracy.Notably, the compression is found to be more effective for larger K in Z K gauge theories, where the original tensor can be too large to perform any calculation with currently available computers.
As another important performance test, we have studied the singular value spectrum of the flavor coarse-graining procedure and find that introducing gauge interaction makes the spectrum decays faster.In order to demonstrate the usefulness of our method, we have investigated the chiral phase transition in two-dimensional Abelian gauge theories with N f = 2 by computing the critical value of the hopping parameter, which turns out to be consistent with the known value obtained by the Monte Carlo method.We have also applied our method to the case of finite density with N f = 1, 2, and 4 in the Z 2 gauge theory.In particular, we were able to observe the Silver Blaze phenomenon, which is difficult to reproduce by Monte Carlo methods due to the sign problem.
We consider that our new technique will make the TRG method applicable to many interesting gauge theories with multiple flavors that have not been explored yet.Since the main idea can be generalized to the non-Abelian case, we hope that it will be useful also in investigating QCD, where two (or three) flavors of light quarks have to be incorporated.We also expect that, by implementing better renormalization schemes such as the bondweighting methods [20,54] or stochastic sampling approach [55,56,57], the systematic error from bond truncation can be reduced.Last but not the least, we expect that our technique is useful in applying the TRG method to the domain-wall formalism for chiral fermions, where the extra dimension can be regarded as the flavor direction in our method.In that case, we need to introduce local interactions in the flavor direction as described in Appendix C. If such interactions make the singular-value spectrum in the flavor direction decay faster, one can go to larger N f , which is crucial in the domain-wall formalism.We hope to report on this in the future publication.

B Details of the coarse-graining procedures
In this section, we present the details of the coarse-graining procedures.There are actually two kinds of them used in this paper, namely the HOTRG method, which is used when we block the tensor in the z(flavor)-direction, and the TRG method, which is used in two-dimensional blocking.

B.1 Grassmann HOTRG
Let us first give a brief review of the non-Grassmann HOTRG method.The main idea is to insert a resolution of the identity on the bonds we wish to truncate.For example, given a two-dimensional tensor T ijkl , where we wish to truncate the leg i, the most naive way is to perform an SVD and insert an identity U † U on that leg as where the legs a and b are truncated.The schematic representation of this process is shown in Fig. 10-a to -c.If the legs i and k in T ijkl are contracted on a periodic lattice, which means that they are pointing in the opposite direction, the isometry U ib can be moved to contract with k instead as Since the legs b and c are truncated, the tensor T ′ becomes smaller.This process is shown in Fig. 10-c to -d.Note that the isometry is not unique.If we apply the same procedure on the leg k instead, we will arrive at a different isometry.We have to choose the one that makes the singular value spectrum decay the fastest and thus allows the smallest truncation.
A more efficient but equivalent way to do this is to define a Hermitian matrix Substituting the first equality of (B.1) in (B.3), we get In other words, we can obtain U by diagonalizing the matrix M instead of performing an SVD for a four-legged tensor T , which is much slower.We can straightforwardly apply such a coarse-graining procedure to the 3-dimensional case.Given a 6-leg tensor T i 1 i 2 i 3 j 1 j 2 j 3 with i µ and j µ pointing in the opposite direction, we first contract two T tensors in z direction, and then attach the isometry (U µ ) (ii ′ ) ĩ to merge the double bond (i µ i ′ µ ) into a truncated bond ĩµ .Note that the standard HOTRG procedure can be generalized to the Grassmann tensor network [28,42].
For the flavor coarse-graining, however, we explain a slightly modified version of the Grassmann HOTRG to improve computational efficiency in dealing with large-size initial tensors.The key point is first to decompose the tensor of each layer into small sub-tensors and then block the legs of subtensors from the layers α and α ′ in the x and y directions separately.Namely, we perform gSVD for T (α) (See Fig. 11-a.)where the square root of singular value matrices are absorbed into U (m) and V (m) .Using U (m) and V (m) , one can then obtain the blocked tensor with renormalized legs as (See For convenience, we also note the explicit form of the coefficient tensors below TJ 1 J 2 J 3 J 4 = I 3 I 1 U ′ I 3 I 1 J 1 J 2 V ′ J 3 J 4 I 1 I 3 σ I 1 σ I 3 (−) p(J 1 )+p(J 2 ) , (B.24)
Since the general structure of the tensor network is similar to that in Fig. 1, appropriate compression and coarse-graining techniques can be straightforwardly applied to this site tensor. )

Figure 2 :
Figure 2: The schematic representation of the initial tensor compression.a) The original tensor.b) Four sets of isometries are inserted between L's and S. c) Another set of isometries is applied to compress the whole tensor.d) The compressed initial tensor.The tensor V a in the diagram above is the isometry with the fastest-falling singular value between V +a and V −a .Similarly, U µ is the isometry with the fastest-falling singular value between U +µ and U −µ .

Figure 5 :
Figure 5: The summary of the flavor coarse-graining procedure.

Figure 7 :
Figure 7: Chiral susceptibility is plotted as a function of κ at β = 0 for different gauge theories and N f (labeled in each plot) with lattice volume up to V = 128 2 .The dotted lines are shown to guide the eyes.The critical hopping parameter κ c can be identified as the location of the peak in the infinite volume limit.Note that the plots for N f = 1 are almost identical for the Z 2 and Z 4 cases.

Figure 8 :
Figure 8: Differential pressure ∆P (μ) and number density ρ(μ) of the Z 2 gauge theory with N f = 2 at various volume V = L × L. The small non-monotonicity may be attributed to the finite truncation effect.

Figure 10 :
Figure 10: Schematic representation of (B.1)-(B.2).Small triangles, large triangles, and circles represent U ia , V ajkl and λ a , respectively.The dashed lines represent the truncated legs.a) A four-legged tensor T ijkl .b) An SVD is performed for one of the legs, i. c) A resolution of the identity U † U is inserted between U and λ. d) An isometry is moved to the other side.

Figure 14 :
Figure 14: Connection of the tensors with the local multi-flavor interaction.a) The site tensor.b) The connection between layers corresponding to the N f = 3 case.

Table 1 :
specialized in handling the Grassmann tensor network.Summary of the initial tensor compression for various input parameters β, μ, N f , and K.The relative error of the compression is less than 10 −15 for all the cases.The size of the original tensor is obtained by the formula 16 4 K 10 , whereas that of the compressed tensor is obtained by (D x D y K) 2 .Here D x and D y represent the bond dimension of the legs I 1 and I 3 and the bond dimension of the legs I 2 and I 4 , respectively, of the compressed coefficient tensor T (α) ′ I 1 I 2 I 3 I 4 ;mn , while K represents the bond dimension of the bosonic legs m and n.