Tensor Renormalization Group Centered About a Core Tensor

TL;DR

This paper introduces a modified tensor renormalization group algorithm that reduces computational cost while maintaining accuracy, enabling more efficient analysis of classical models on 2D lattices.

Contribution

A new tensor renormalization group method that coarse-grains only adjacent rows and columns, reducing computational complexity from O(χ^6) to O(χ^4).

Findings

Cost scales as O(χ^4), cheaper than TRG.

Accuracy comparable to TRG at same bond dimension.

Enables analysis of larger or more complex models.

Abstract

We propose a modified form of a tensor renormalization group algorithm for evaluating partition functions of classical statistical mechanical models on 2D lattices. This algorithm coarse-grains only the rows and columns of the lattice adjacent to a single core tensor at each step, such that the lattice size shrinks linearly with the number of coarse-graining steps as opposed to shrinking exponentially as in the usual tensor renormalization group (TRG). However, the cost of this new approach only scales as $O(\chi^4)$ in terms of the bond dimension $\chi$, significantly cheaper than the $O(\chi^6)$ cost scaling of TRG, whereas numerical benchmarking indicates that both approaches have comparable accuracy for the same bond dimension $\chi$. It follows that the new algorithm can allow for more accurate investigation of classical models on modestly sized lattices than is possible with standard TRG.

Equations

$$\epsilon =\Vert F-P_{h}F\Vert ,$$

$$\displaystyle Z=\sum_{\{s\}}\left(\prod_{\mbox{$\langle i,j\rangle$}}\exp(-\beta s_{i}s_{j})\right),\;\;s_{i}\in\{+1,-1\}$$

$$Z=\sum _{\{\sigma \}}{e}^{-H(\{\sigma \})/T},$$

$$H\left(\{\sigma \}\right)=-\sum _{⟨i,j⟩}{\sigma }_i{\sigma }_j,$$

$$A_{ijkl}=e^{\left({\sigma }_i{\sigma }_j+{\sigma }_j{\sigma }_k+{\sigma }_k{\sigma }_l+{\sigma }_l{\sigma }_i\right)/T},$$

$$Z=\sum _{ijk\cdots }{A}_{ijkl}{A}_{mnoj}{A}_{krst}{A}_{opqr}\cdots .$$

Full text

I Introduction

In the past few decades tensor networks (TN) TN1 ; TN2 ; TN3 have become a popular tool in computational and theoretical physics, both for the study of classical and quantum many-body systems. Building on the tremendous success of matrix product states (MPS) MPS1 ; MPS2 ; MPS3 and the density matrix renormalization group (DMRG) DMRG1 ; DMRG2 algorithm in studying one dimensional systems, many other classes of tensor network ansatz and algorithms have been developed. Examples include tensor network ansatz such as the multi-scale entanglement renormalization ansatz (MERA) MERA ; 2DMERA1 ; 2DMERA2 ; CritMERA and the projected entangled pair states (PEPS)PEPS1 ; PEPS2 ; PEPS3 , as well as tensor network realizations of coarse-graining algorithms such as the tensor renormalization group (TRG) TRG ; TRGplus ; TERG ; SRG ; TEFR ; TRGenv ; TRG3D ; SpinNet ; TRG3DB ; TRGcheap ; TRGani and tensor network renormalization (TNR) TNR1 ; TNR2 ; TNR3 ; TNR4 or similar methods TNRc1 ; TNRc2 ; TNRc3 ; TNRc4 ; TNRc5 . Despite ongoing efforts to improve computational performance, high computational cost is still the main obstacle for large scale calculations with many tensor network methods. Thus, even with the rapid development of computing power, there is a strong desire to reduce the computational cost through introduction of more efficient algorithms.

In this work we introduce a novel coarse-graining algorithm, which we call core-tensor renormalization group (CTRG), that can be applied to contract a $2D$ tensor network. Thus, similar the standard tensor renormalization group, this method could be used to study $2D$ classical many-body systems, where the network encodes the partition function, or $1D$ quantum systems, where the network encodes the Euclidean path integral given from a Suzuki-Trotter decomposition Suzuki1 ; Suzuki2 . However, the proposed CTRG method is demonstrated to be significantly more efficient than TRG in certain scenarios, allowing for more accurate numerical investigation of many-body systems, yet still retains the robustness and simplicity of implementation of the TRG approach.

The outline of our manuscript is as follows. We begin in Sect. II with a non-technical overview of the proposed algorithm and its benefits, before providing the full implementation details in Sect. III. Benchmark results and comparisons with TRG are provided in Sect. IV, then discussions and conclusions are presented in Sect. V.

II Overview

In this section we provide an overview of the CTRG algorithm and discuss its cost scaling in comparison to TRG. As with TRG, the starting point of the CTRG algorithm is a $2D$ lattice of tensors. For simplicity we assume to have a homogeneous square lattice of linear dimension $L$ comprised of copies of a tensor $A_{0}$, assumed to have 4 indices of dimension $d$. Note that the method can also easily be extended to non-homogeneous networks as well as to other lattice geometries (e.g. hexagonal and triangular lattices). Through standard mappings, this network could encode either the partition of a $2D$ classical system on the square lattice, where the $A_{0}$ tensors encode the Boltzmann weights associated to plaquettes (see Appendix A), or encode the Euclidean path integral of a $1D$ quantum system, where the $A_{0}$ tensors encode the local gates given from a Suzuki-Trotter expansion Suzuki1 ; Suzuki2 . As such, contracting the $2D$ tensor network can be tantamount to calculating properties, such as expectation values, of $2D$ classical or $1D$ quantum many-body systems.

As a preliminary step of CTRG we select one tensor $A_{0}$ to act as the ‘core’ tensor, which we re-label as $A_{c}$, and similarly re-label tensors in the same column as $A_{v}$ and tensors in the same row $A_{h}$. Then an iteration of the CTRG algorithm simply consists of the following two steps as depicted in Fig. 1: (i) contract the adjacent row and column of the network into the core row/column, (ii) truncate the dimension of the core row/column, retaining at most $\chi$ dimensional indices, and obtaining new core tensors $\{A^{\prime}_{c},A^{\prime}_{v},A^{\prime}_{h}\}$. Thus starting from an $L\times L$ lattice, which is homogeneous $A_{0}$ everywhere except for a single row and column, after a single step of CTRG we obtain an $(L-1)\times(L-1)$ lattice, which is again homogeneous $A_{0}$ everywhere except for a single row and column. It follows that after performing $L$ iterations of CTRG a network of $O(1)$ tensors is obtained, which can then be exactly contracted.

As will be shown in Sect. III, the cost of a single iteration of CTRG scales as $O(\chi^{4}d^{3})$, with $d$ the index dimension of the initial network. By extension, the cost of contracting over a spatially homogeneous lattice of linear dimension $L$ using CTRG scales as $O(\chi^{4}d^{3}L)$. In comparison, the cost of TRG is known to scale as $O(\chi^{6}\log L)$ for a homogeneous $L\times L$ lattice, with $\chi$ the TRG bond dimension. However, as will be demonstrated in the benchmark results of Sect. IV, both CTRG and TRG give comparable accuracy when using the same bond dimension $\chi$. Thus one expects CTRG to be the more computationally efficient algorithm in the large bond dimension $\chi$ limit for systems where $L$ is not too large. Indeed, the benchmark results for the $2D$ classical Ising model at critical temperature show that, even for modestly large lattices of linear dimension $L=256$ Ising spins, the CTRG algorithm vastly outperforms standard TRG. In this setting the CTRG algorithm can reach a level of accuracy in less than a minute of computation time that would require many hours using standard TRG, and the disparity between computation times grows wider as the bond dimension $\chi$ is increased.

III Algorithm

In this section we describe in detail the proposed CTRG algorithm. As explained in Sect. II the starting point of the algorithm is a square lattice tensor network composed of four index bulk tensors $A_{0}$ with index dimensions $d$. We select a single row and column to serve as the core, relabeling the column tensors $A_{v}$, the row tensors $A_{h}$, and the central tensor $A_{c}$. Here we also perform an additional step of decomposing $A_{c}$ into a pair of 3-index tensors $C_{l}$ and $C_{r}$ using the singular value decomposition (SVD), and similarly decompose certain bulk tensors $A_{0}$, located on the forward sloping diagonal through the core of the lattice, into pairs of 3-index tensors as indicated in Fig. 2(a). As will be later evident, by representing the central tensor $A_{c}$ as the product of a pair of 3-index tensors $\{C_{l},C_{r}\}$ the computational cost is reduced to $O(\chi^{4})$, which would otherwise be $O(\chi^{5})$ if instead $A_{c}$ was kept 4-index tensor. We now describe how a single iteration of CTRG, which maps an $L\times L$ lattice composed of tensors $\{A_{0},A_{v},A_{h},C_{l},C_{r}\}$ into an $(L-1)\times(L-1)$ lattice composed of tensors $\{A_{0},A^{\prime}_{v},A^{\prime}_{h},C^{\prime}_{l},C^{\prime}_{r}\}$, by appropriately absorbing an adjacent bulk row/column into the core row/column. For simplicity we describe only the step where the bulk row/column from the below/left of the core is absorbed, although in practice we alternate each iteration with absorbing the bulk row/column from the above/right of the core.

An iteration of the coarse-graining step is implemented using three different isometric tensors $\{Y_{h},Y_{v},Y_{c}\}$, whose output dimension is upper bounded at some specified maximal bond dimension $\chi$. The isometries $Y_{h}$ act to compress a bulk row into the core row, while isometries $Y_{h}$ act to compress a bulk column into the core column, and finally isometries $Y_{c}$ act to compress the index joining the central tensors $C_{l}$ and $C_{r}$, as depicted in Fig. 2(b-c). Let us define the small network $F$ as that composed of the central tensors $\{C_{l},C_{r}\}$ together with a single tensor $A_{h}$ and $A_{v}$ from the core column and row respectively, and a single bulk tensor $A_{0}$, as depicted in Fig. 3. The isometries are chosen as to minimize the truncation error resulting from their insertion into the network; more specifically, they should leave network $F$ invariant when applied in conjugate pairs. For instance, if we define $P_{h}\equiv Y_{h}Y_{h}^{\dagger}$ as the rank-$\chi$ projector defined from isometry $Y_{h}$ and its conjugate, then the truncation error $\varepsilon$ when acting on $F$ is given as,

[TABLE]

where network $F$ has been considered as a matrix between its left two indices and the remaining indices, as depicted in Fig. 3(a). Fortunately the solution for the optimal rank-$\chi$ projector $P_{h}\equiv Y_{h}Y_{h}^{\dagger}$ that minimizes the error $\varepsilon$ in Eq. 1 is easily obtained: isometry $Y_{h}$ should be formed from the $\chi$ dominant eigenvectors (i.e. those with largest magnitude eigenvalue) of the Hermitian matrix $FF^{\dagger}$. Similarly, the optimal isometry $Y_{v}$ to minimize the truncation error is obtained through eigenvalue decomposition (ED) of $FF^{\dagger}$ when network $F$ is treated as a matrix between its bottom two indices and the remaining indices, see Fig. 3(b). Finally, the remaining isometry $Y_{c}$ is chosen to minimize the truncation error associated to the diagonal index joining the core tensors $C_{l}$ and $C_{r}$. This is accomplished by requiring the projector $P_{c}\equiv Y_{c}Y_{c}^{\dagger}$ to minimize the error when acting on half of network $F$, split through the diagonal as depicted in Fig. 3(c), again found through the appropriate eigenvalue decomposition.

Once the set of isometries $\{Y_{h},Y_{v},Y_{c}\}$ required for the coarse-graining step have been obtained they are then inserted into the network, as depicted in Fig. 2(b), and contracted to form the coarse-grained core tensors $\{A^{\prime}_{v},A^{\prime}_{h},C^{\prime}_{l},C^{\prime}_{r}\}$ as depicted in Fig. 2(c-d). Thus the iteration is complete, and we can begin the next iteration to absorb another bulk row/column from the network into the core row/column.

Notice that during the coarse-graining step the same isometry $Y_{h}$ is used at every location within the row under consideration, likewise the same isometry $Y_{v}$ is used at every location within the column under consideration. This is despite the fact that the CTRG method appears to break the translation invariance of the initial network, where each tensor $A_{h}$ within the core row occupies a different location with respect to the intersection tensors $\{C_{l},C_{r}\}$. This could lead one to believe (erroneously) that isometries $Y_{h}$ should vary with their relative position in the row. However our proposed approach is indeed justified since the core column after $z$ coarse-graining iterations is identical, up to truncation errors, to a width $(z+1)$ strip of tensors from the original lattice. Thus, up to small errors, tensors $A_{h}$ in the core row cannot ‘sense’ the presence of the core column; it follows that all $A_{h}$ tensors possess an equivalent tensor environment regardless of their proximity to the center.

The key steps of the CTRG algorithm include (i) the contraction of the $FF^{\dagger}$ networks as depicted in Fig. 3, (ii) the eigenvalue decomposition of the $FF^{\dagger}$ networks, and (iii) the contraction for new core tensors $\{A^{\prime}_{v},A^{\prime}_{h},C^{\prime}_{l},C^{\prime}_{r}\}$ as depicted in Fig. 2(d). Each of these steps can be accomplished with leading order computational cost $O(\chi^{4}d^{3})$ or less, with $d$ as the dimension of the original lattice and assuming that the indices of the core row/column have initial and final dimension $\chi$. In comparison, the computational cost of traditional TRG TRG is known to scale as $O(\chi^{6})$. The reduction in the cost as a power of $\chi$ when using CTRG can be attributed to the fact that only a single row/column of the bulk network is absorbed into the core at each step: the isometries only act to combine a $d$-dimensional and a $\chi$-dimensional index into a new $\chi$-dimensional index at each step. In contrast, the tensors used in TRG must act to combine two $\chi$-dimensional indices at each step, which is inherently more expensive as a function of $\chi$.

IV Benchmark results

We test the efficiency of the proposed approach by benchmarking with the square-lattice classical Ising model, whose partition function is defined,

[TABLE]

where $\langle i,j\rangle$ represents nearest-neighbor sites and $\beta$ is the inverse temperature, $\beta=1/T$. We encode this partition function as a square-lattice tensor network using standard methodology TNR4 , see also Appendix A.

As a first test, we explore the accuracy of CTRG in terms of the free energy density and internal energy density, both as a function of temperature $T$ and of bond dimension $\chi$. In each calculation, we repeat the RG iteration until the energy densities are converged (i.e. such that the finite size effects are smaller than the accuracy threshold of the method) and then compare with the exact known Ising model energy densities in the thermodynamic limit. This required at most $3000$ RG steps to achieve convergence near the critical temperature, $T_{c}=2/{\log{(1+\sqrt{2}})}$, which equates to lattices of linear size $L=3000$ Ising spins. The results, as presented in Fig. 4, are consistent with behavior of other tensor RG schemes TRG3D : the accuracy is lowest when near $T_{c}$, but improves systematically everywhere in the phase diagram as the bond dimension $\chi$ is increased.

Next, in Fig. 5, we compare the accuracy of the free energy given from CTRG against that given from a standard implementation of TRG TRG for the Ising model on an infinite strip of width $L=128$. When using either method, the strip is mapped to an width $L=1$ chain of infinite length after a finite number of RG steps, which can then be solved using standard transfer matrix techniques. In lieu of an exact value for the free energy density on this strip, we use the free energy given from a $\chi=640$ CTRG calculation as a proxy for the exact value (noting that the values plotted in Fig. 5 are converged with respect this reference energy). In Fig. 5 (a) we see that, for any fixed $\chi$, both TRG and CTRG yield very similar accuracies for the free energy at critical temperature $T_{c}$, with TRG only slightly more accurate for all $\chi$. Notable is that both methods converge polynomially to the exact free energy as a function of $\chi$, with the error in the free energy density $\Delta F$ scaling approximately as $\Delta F\propto\chi^{-3.4}$. Fig. 5(b) compares the accuracy of TRG and CTRG at fixed bond dimension $\chi=30$ as a function of temperature. Interestingly, it is seen that CTRG gives greatly improved accuracy over TRG when away from $T_{c}$, closely matching the performance of more sophisticated versions of TRG that make use of a larger tensor environment in order to improve accuracy SRG ; TRGenv .

Finally in Fig. 6 we compare the accuracy of TRG and CTRG as a function of computation time, with all calculations performed separately on the same hardware (3.4Ghz six-core desktop CPU), and for infinite length strips of width $L=128$ and $L=256$ Ising spins at critical temperature $T_{c}$. These results match the expectations that arise from consideration of the computational cost scaling of TRG, $O(\chi^{6}\log L)$, versus that of CTRG, $O(\chi^{4}L)$. Namely that for small $\chi$ calculations the $L$ component of the cost scaling dominates such that TRG is the more efficient approach, while in the large $\chi$ limit CTRG becomes the more efficient approach. However, we find it remarkable how early this crossover in efficiency occurs; on the width $L=256$ strip of Ising spins, CTRG is already the more efficient approach for calculations requiring more than a few seconds. Recently a version of TRG with a reduced computational cost scaling of $O(\chi^{5}\log L)$ has been proposed in Ref. TRGcheap, , which uses an iterative approach to find optimal projectors, rather than the SVD. Although this new variant of TRG would come closer to CTRG in terms of efficiency, it is still expected that it would be less efficient than CTRG in the large-$\chi$ limit. Note that all benchmark calculations were performed using code written in Julia (version 1.1.0) and made use of the ‘TensorOperations.jl’ package.

V Discussion

An important observation that can be drawn from the benchmark results is that a CTRG calculation with bond dimension $\chi$ seems to perform comparably well to a TRG calculation of the same bond dimension, when applied to systems at or near criticality. This is not surprising since both methods are related to tree tensor network (TTN) ansatz TTN1 ; TTN2 . It can be understood that TRG produces an approximation to the initial square-lattice tensor network as a binary tree (an observation which is more readily apparent when considering the HOTRG TRG3D version of the method). Similarly, it can be understood that CTRG produces an approximation to the initial width-$L$ lattice as a depth-$L$ tree, see Fig. 7, which would be equivalent to an MPS under folding through the center. Since both forms of tree can only support block entanglement entropy $S$ which is upper bounded as $S\leq\log(\chi)$, it makes sense that they share the same accuracy limitations when applied to critical systems, which are known to have logarithmically diverging block entanglement entropy scaling Ent1 ; Ent2 .

A more recent advance in tensor RG schemes has been the introduction of a disentangling step, as pioneered by tensor network renormalization (TNR) and related methods TNR1 ; TNR2 ; TNR3 ; TNR4 ; TNRc1 ; TNRc2 ; TNRc3 ; TNRc4 ; TNRc5 . TNR is known to produce an approximation to a square-lattice tensor network as a multi-scale entanglement renormalization ansatz (MERA) TNR2 , thus is a natural choice for investigating critical systems as it can support the correct entanglement scaling Ent1 ; Ent2 . A practical consequence of having the proper entanglement scaling is that the free energy converges exponentially in bond dimension $\chi$ to the exact value, i.e. as $\Delta F\propto\exp(-k\chi)$, as opposed to only polynomially quickly for TRG and CTRG. Thus, in the limit of large $\chi$, disentangling based methods such as TNR are still expected to outperform CTRG when applied to critical systems. However there are still some scenarios for which CTRG could the more efficient approach, such as (i) away from criticality where impact of disentangling is not so significant, or (ii) at intermediate values of $\chi$ where computational cost associated to disentangling outweighs the associated accuracy improvements. Moreover, from a practical standpoint, CTRG has the benefit of being much easier to implement than TNR and related approaches that use disentangling. This follows as CTRG inherits much of the simplicity of TRG; the isometries that implement each RG step are each given (deterministically) through a single eigen-decomposition as seen in Fig. 3. In contrast the disentanglers required by TNR are given via an iterative method that repeats many times until the tensors are deemed sufficiently converged TNR4 . Such iterative methods have the possibility of becoming trapped in sub-optimal solutions, a possibility which is entirely avoided when using TRG or CTRG.

A potentially useful extension of CTRG would be towards networks in higher spatial dimension such as a cubic lattice tensor network, which could represent either the partition function of a $3D$ classical system or the path integral for a $2D$ quantum system. In this setting the HOTRG TRG3D has proven to be a useful numerical tool, but it is hampered by a computational cost which scales as $O(\chi^{11})$ in terms of bond dimension $\chi$, although more efficient variations have recently been developed TRGani . We hypothesize that a version of CTRG generalized for higher spatial dimensions could reproduce results of equivalent accuracy to HOTRG, but with a much lower cost scaling in bond dimension $\chi$. This possibility remains an intriguing avenue for future research.

This research was undertaken thanks in part to funding from the Canada First Research Excellence Fund.

Appendix A Tensor Network Representation of Partition Functions

In this appendix we describe how the partition function $Z$ of a $2D$ classical statistical system at temperature $T$,

[TABLE]

can be expressed as a network of tensors. As a specific example we consider the square-lattice classical Ising model, with Hamiltonian functional,

[TABLE]

where $\sigma_{i}\in\{+1,-1\}$ is an Ising spin on site $i$. The partition function can be represented as a square-lattice tensor network composed of copies of a four-index tensor $A_{ijkl}$, where tensors sit in the center of every second plaquette of Ising spins according to a checkerboard tiling as depicted in Fig. 8, such that there is one tensor $A$ for every two spins. Each tensor $A$ encodes the four Boltzmann weights $e^{\sigma_{i}\sigma_{j}/T}$ of the Ising spin interactions on the edges of the plaquette on which it sits,

[TABLE]

such that the partition function $Z$ is thus as the sum over all indices from the network of $A$ tensors,

[TABLE]

Notice that this method for constructing tensor networks representations of partition functions can be employed for any model with nearest-neighbor interactions.