Tensor Network Approaches for Higher-Dimensional Quantum Lattice Models

In this section, we will show several representative TN approaches for simulating the quantum lattice models in ( d > 1) dimensions. We will mainly use the language of TN contractions. One may refer to several existing reviews 2017) for more exhaustive understanding on the TN simulations for quantum problems. We will focus on the algorithms based on PEPS, and show the key roles that the 2D TN contraction algorithms presented in Sect. 3 play in the higher-dimensional cases.


Variational Approaches of Projected-Entangled Pair State
Without losing generality, we consider a 2D quantum system with nearest-neighbor coupling on an infinite square lattice as an example. There are in general two ways to do the minimization. One way proposed firstly by Verstraete and Cirac [1] is considering the elements in the tensors as variational parameters. The tensors in the TN are updated one by one. In a similar spirit as the boundary-state methods (see Sect. 3.4), the key of this approach is to transform the global minimization to local ones, where one tensor (say P [ H eff is an "effective" Hamiltonian by computing ψ|Ĥ |ψ but after taking P [i] † in ψ| and P [i] in |ψ out. Figure 4.1 depictsĤ eff whereĤ is written as an infinite PEPO (iPEPO, also see Sect. 2.2.6 for PEPO) for a better illustration. Similarly, N eff is defined by computing ψ|ψ but after taking P [i] † and P [i] out. Obviously, the computations of bothĤ eff andN eff are in fact to contract the corresponding 2D TN's where the 2D TN contraction algorithms are needed. In [2], Corboz used CTMRG (see [3] or Sect. 3.3) to compute the contractions. In [4], Vanderstraeten et al. further developed this idea to a gradient method, where the gradient is calculated by implementing similar 2D TN contractions. The gradient is given as By imposing the normalization condition ψ|ψ = 1 and shifting the ground-state energy to zero byĤ ←Ĥ − ψ|Ĥ |ψ , the gradient is simplified as Note that the infinite PEPO (iPEPO) representation is not enforced to definê H eff . In fact, it is not easy to obtain the iPEPO of an arbitrary 2D (or 3D) Hamiltonian. The usual way is to start from the summation form of the Hamiltonian H = Ĥ ij , and compute the contribution toĤ eff from eachĤ ij separately [2].
Each term is computed by contracting a 2D TN, where one can reuse the results to improve the efficiency. Following the same clue (minimizing E), algorithms were proposed to combine TN with the QMC methods [5][6][7][8][9]. Still let us focus on those based on PEPS. One may transform Eq. (4.1) as where S = (s 1 , s 2 , · · · ) goes through all spin configurations and W (S) = s 1 s 2 · · · |ψ is the coefficient of the iPEPS for the given configuration. QMC sampling can be implemented by defining the weight function as W (S) 2 and the estimator E(S) as so that the energy becomes It is easy to see that the normalization condition of the weights S W (S) 2 = 1 is satisfied.
The task becomes to compute W (S) and S |Ĥ |S with different configurations. The computation of S |Ĥ |S is relatively easy since |S and |S are just two product states. The computation of W (S) is more tricky. When |ψ is a PEPS on a square lattice, W (S) is a 2D scalar TN by fixing all the physical indexes of the PEPS as where P [n] s n is a fourth-order tensor that only has the geometrical index. 1 The n-th physical index is taken as s n . Considering that most of the configurations are not translationally invariant, such QMC-TN methods are usually applied to finite-size models. One may use the finite TN version of the algorithms reviewed in Sect. 3.

Imaginary-Time Evolution Methods
Another way to compute the ground-state iPEPS is to do imaginary-time evolution, analog to the MPS methods presented in Sect. 3.1.3. For a d-dimensional quantum model, its ground-state simulation can be considered as computing the contraction of a (d + 1)-dimensional TN.
Firstly, let us show how the evolution operator for an infinitesimal imaginarytime step τ can be written as an iPEPO, which is in fact one layer of the 3D TN ( Fig. 4.2). The evolution of the iPEPS is to put the iPEPS at the bottom and to contract the TN layer by layer to it.
To proceed, we divide the local Hamiltonians on the square lattice into four group:Ĥ e,e = even i,jĤ One can see that each two terms in one group commute to each other. The evolution operator for an infinitesimal imaginary-time step (τ → 0) then can be written aŝ Let us assume translational invariance to the Hamiltonian, i.e.,Ĥ [i,j ] =Ĥ [two] . The element of two-body evolution operator is a fourth-order tensor U s i s j s i s j = s i s j | exp(−τĤ [two] )|s i s j . Implement SVD or QR decomposition on U (4.2) as (4.9) Then the two tensors T [L] and T [R] that form the iPEPO ofÛ is obtained as While the TN for the imaginary-time evolution with the iPEPO is a cubic TN, one may directly use the tensor U , which also gives a 3D but not cubic TN. Without losing generality, we in the following will use the iPEPO to present the algorithms for contraction a cubic TN. The algorithm can be readily applied to deal with the statistic models on cubic lattice or other problems that can be written as the contraction of a cubic TN.
The evolutionÛ |ψ is to contract the iPEPO (one layer of the tensors) to the iPEPS. In accordance to the translational invariance of the iPEPO, the iPEPS is also formed by two inequivalent tensors (denoted by P [L] and P [R] ). Locally, the tensors in the evolved iPEPS are given as with the composite indexesα x = (a x , α x ) (x = 1, 2, 3, 4). Obviously, the bond dimensions of the new tensors are increased by dim(a x ) times. It is necessary to preset a dimension cut-off χ : when the bond dimensions become larger than χ , approximations will be introduced to reduce the dimensions back to χ . One then can iterate the evolution of the iPEPS with bounded computational cost. After the iPEPS converges, it is considered that the ground state is reached. Therefore, one key step in the imaginary-time schemes (as well as the similar contraction schemes of 3D TN's) is to find the optimal truncations of the enlarged bonds. In the following, we will concentrate on the truncation of bond dimensions, and present three kinds of scheme known as full, simple, and cluster updates according to which environment the truncations are optimized [10]. 2

Full, Simple, and Cluster Update Schemes
For truncating the dimensions of the geometrical bonds of an iPEPS, the task is to minimize the distance between the iPEPSs before and after the truncation, i.e., ε = ||ψ − |ψ |. With the normalization condition of the iPEPSs, the problem can be reduced to the maximization of the fidelity Z = ψ |ψ . (4.14) As discussed in Sect. 3.1.2, Z is in fact a scalar TN.
Full Update Among the three kinds of update schemes, full update seems to be the most natural and reasonable, in which the truncation is optimized referring to the whole iPEPS [10][11][12][13][14][15][16]. Let us consider a translationally invariant iPEPS. For square lattice, the iPEPS is formed by the infinite copies of two tensors P [L] and P [R] located on the two sub-lattices, respectively. Their evolution is given by Eq. (4.10). We useP [L] andP [R] to denote the tensors with enlarged bond dimensions. Below, we follow Ref. [13] to explain the truncation process. with the dimension of the shared bond dim(α 4 ) = χ . We shall stress that Eq. (4.17) is not the SVD of M; the decomposition and truncation are optimized by the SVD of M e , hence is a non-local optimization.
With the formula given above, the task is to compute the environment tensor M e by the contraction algorithms of 2D TN's. In Ref. [13], the authors developed the SRG, where M e is computed by a modified version of TRG algorithm [17]. Other options include iTEBD [15], CTMRG [12], etc. Note that how to define the environment as well as how to truncate by the environment may have subtle differences in different works. The spirit is the same, which is to minimize the fidelity in Eq. (4.14) referring to the whole iPEPS.
Simple Update A much more efficient way known as the simple update was proposed by Jiang et al. [18]; it uses local environment to determine the truncations, providing an extremely efficient algorithm to simulate the 2D ground states. As shown in Fig. 2.8c, the iPEPS used in the simple update is formed by the tensors on the site and the spectra on the bonds: two tensors P [L] and P [R] located on the two sub-lattices, and λ [1] , λ [2] , λ [3] , and λ [4] on the four inequivalent geometrical bonds of each tensor. The evolution of the tensors in such an iPEPS is given by Eq. (4.10). λ [i] should be simultaneously evolved asλ [i] (a i ,α i ) = I a i λ α i with I a i = 1.
To truncate the fourth geometrical bond of P [L] (and P [R] ), for example, we construct a new tensor by contracting P [L] and P [R] and the adjacent spectra as  where one takes only the χ -largest singular values and the basis. P [L] and P [R] are updated as (4.20) The spectrumλ 4 is updated by λ in the SVD. The above procedure truncates dim(α 4 ) to the dimension cut-off χ , which can be readily applied to truncate any other bonds. According to the discussion about SVD in Sect. 2.2.1, the environment is the two tensors and the adjacent spectra λs in M, where the λs play the role of an "effective" environment that approximate the true environment (M e in the full update). From this viewpoint, the simple update uses local environment. Later by borrowing from the idea of the orthogonal form of the iPEPS on Bethe lattices [19][20][21][22][23][24][25][26], it was realized that the environment of the simple update is the iPEPS on the infinite trees [27][28][29], not just several tensors. We will talk about this in detail in the next chapter from the perspective of the multi-linear algebra.
Cluster Update By keeping the same dimension cut-off, the simple update is much more efficient than the full update. On the other hand, obviously, the full update possesses higher accuracy than the simple update by considering better the environment. The cluster update is between the simple and full updates, which is more flexible to balance between the efficiency and accuracy [10,27,30].
One way is to choose a finite cluster of the infinite TN and define the environment tensor by contracting the finite TN after taking a pair ofP [L] andP [R] out. One can consider to firstly use the simple update to obtain the spectra and put them on the boundary of the cluster [30]. This is equivalent to using a new boundary condition [27,29], different from the open or periodic boundary conditions of a finite cluster. Surely, the bigger the cluster becomes, more accurate but more consuming the computation will be. One may also consider an infinite-size cluster, which is formed by a certain number of rows of the tensors in the TN [10]. Again, both the accuracy and computational cost will in general increase with the number of rows. With infinite rows, such a cluster update naturally becomes the full update. Despite the progresses, there are still many open questions, for example, how to best balance the efficiency and accuracy in the cluster update.

Summary of the Tensor Network Algorithms in Higher Dimensions
In this section, we mainly focused on the iPEPS algorithm that simulates the ground states of 2D lattice models. The key step is to compute the environment tensor, which is to contract the corresponding TN. For several special cases such as trees and fractal lattices, the environment tensor corresponds to an exactly contractible TN, and thus can be computed efficiently (see Sect. 2.3.6). For the regular lattices such as square lattice, the environment tensor is computed by the TN contraction algorithms, which is normally the most consuming step in the iPEPS approaches.
The key concepts and ideas, such as environment, (simple, cluster, and full) update schemes, and the use of SVD, can be similarly applied to finite-size cases [31,32], the finite-temperature simulations [27,28,[33][34][35][36][37][38][39], and real-time simulations [31,40] in two dimensions. The computational cost of the TN approaches is quite sensitive to the spatial dimensions of the system. The simulations of 3D quantum systems are much more consuming than the 2D cases, where the task becomes to contract the 4D TN. The 4D TN contraction is extremely consuming, one may consider to generalize the simple update [29,41], or to construct finite-size effective Hamiltonians that mimic the infinite 3D quantum models [29,42] Many technical details of the approaches can be flexibly modified according to the problems under consideration. For example, the iPEPO formulation is very useful when computing a 3D statistic model, which is to contract the corresponding 3D TN. As for the imaginary-time evolution, it is usually more efficient to use the two-body evolution operators (see, e.g., [12,18]) rather than the iPEPO. The environment is not necessarily defined by the tensors; it can be defined by contracting everything of the TN except for the aimed geometrical bond [28,33]. The contraction order also significantly affects the efficiency and accuracy. One may consider to use the "single-layer" picture [10,31], or an "intersected" optimized contraction scheme [43].
Open Access This chapter is licensed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence and indicate if changes were made.
The images or other third party material in this chapter are included in the chapter's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the chapter's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.