Infinite projected entangled pair states algorithm improved: Fast full update and gauge fixing

The infinite projected entangled pair states (iPEPS) algorithm [J. Jordan et al., Phys. Rev. Lett. 101, 250602 (2008)] has become a useful tool in the calculation of ground-state properties of two-dimensional quantum lattice systems in the thermodynamic limit. Despite its many successful implementations, the method has some limitations in its present formulation which hinder its application to some highly entangled systems. The purpose of this paper is to unravel some of these issues, in turn enhancing the stability and efficiency of iPEPS methods. For this, we first introduce the fast full update scheme, where effective environment and iPEPS tensors are both simultaneously updated (or evolved) throughout time. As we shall show, this implies two crucial advantages: (i) dramatic computational savings and (ii) improved overall stability. In addition, we extend the application of the local gauge fixing, successfully implemented for finite-size PEPS [M. Lubasch et al., Phys. Rev. B 90, 064425 (2014)], to the iPEPS algorithm. We see that the gauge fixing not only further improves the stability of the method but also accelerates the convergence of the alternating least-squares sweeping in the (either “full” or “fast full”) tensor update scheme. The improvement in terms of computational cost and stability of the resulting “improved” iPEPS algorithm is benchmarked by studying the ground-state properties of the quantum Heisenberg and transverse-field Ising models on an infinite square lattice.

Figure 1

(a) Graphical representation of a PEPS on a $5\times 5$ square lattice. Each ball is a tensor, and lines correspond to tensor indices. Lines from one tensor to another correspond to common summed indices (or contracted indices). The free index, or open leg, in each tensor is called the physical index and corresponds to the local degrees of freedom of the local Hilbert space at every site. (b) An infinite PEPS with a two-tensor unit cell. The two tensors, $A$ and $B$, are repeated on the infinite 2D lattice.

Figure 2

(a) Left: Two-dimensional lattice of tensors formed from $a$ and $b$; right: contractions to obtain tensors $a$ and $b$. (b) Environment of a given link on the lattice (here an $r$ link). (c) Effective environment of a given link on the lattice (here an $r$ link). (d) Six-tensor representation of the effective environment around the link.

Figure 3

Tensor network diagrams showing how to compute $T,R$, and $S$ in Eq. (12).

Figure 4

(a) Subtensors for iPEPS tensors $A$ and $B$. Decompose $A=X{a}_R$ and $B=b_{L}Y$ by means of decompositions such as QR or SVD. (b) The action of the two-body gate $g$ on the iPEPS tensors $A$ and $B$ connected by an $r$ link is equivalent to its action on only the subtensors $a_R$ and $b_L$, leaving $X$ and $Y$ unaffected.

Figure 5

Tensor network diagrams showing how to compute $T,R$, and $S$ in Eq. (14).

Figure 6

Tensor network diagrams for the fast full update. Details are explained in the main text.

Figure 7

(a) Contraction producing the norm tensor. The leading computational cost of this contraction is $O\left(d^4{D}^4{\chi }^2\right)+O\left(d^2{D}^6{\chi }^2\right)+O\left(d^2{D}^3{\chi }^3\right)$ [45]. (b) Diagrammatic representation of the cost function defined in Eq. (15).

Figure 8

(a) Positive approximant for the norm tensor ${\tilde{\mathcal{N}}}_{LR}$. (b) Apply the QR decomposition for the left and the right bond indices of tensor $Z$, so that $Z=Q_{L}R=L{Q}_R$. (c) Insert the identities $I=L^{-1}L=R^{-1}R$ into the left and right bond indices of the tensor $Z$ to obtain $\tilde{Z}=L^{-1}Z{R}^{-1}$ and the new subtensors ${\tilde{a}}_R=L{a}_R,{\tilde{b}}_L=b_{L}R$. (d) A new norm tensor is obtained as ${\tilde{\mathcal{N}}}_{LR}=\tilde{Z}{\tilde{Z}}^{†}$. In order to keep the compatibility once the subtensors have been updated, we need to multiply the tensors $X$ and $Y$ by matrices $L^{-1},R^{-1}$, so that $\tilde{X}=X{L}^{-1}$ and $\tilde{Y}=R^{-1}Y$, respectively. This is done before recovering the tensors $\tilde{A}=\tilde{X}{\tilde{a}}_R$ and $\tilde{B}={\tilde{b}}_L\tilde{Y}$.

Figure 9

Mean value of the relative change ${\overline{\epsilon }}_d=|d_{u+1}-d_u|/|d_{\text{init}}|$ of the cost function $d$ and the local fidelity ${\overline{\epsilon }}_{ab}=|f_{ab}^{u+1}-f_{ab}^u|/|f_{ab}^{\text{init}}|$ in the iPEPS algorithm with gauge (solid symbols) or without gauge (open symbols) for (a) and (b) the Ising model with transverse field and (c) and (d) the Heisenberg model, respectively. For the Ising model the magnetic field is $h=3.01$ and $(D,\chi )=(3,30)$, whereas for the Heisenberg model $(D,\chi )=(5,50)$. The time step is $\delta =0.005$. The mean value is taken over all steps in imaginary time.

Figure 10

Local magnetization $m_z$ as a function of the transverse field $h$ in the quantum Ising model, with gauge (solid circles) and without gauge (stars). The inset is a zoom around criticality.

Figure 11

Two-point correlator $S_{zz}\left(x\right)$ of the quantum Ising model close to the critical point $h=3.044$, computed with gauge (solid circles) and without gauge (stars).

Figure 12

(a) Relative error of the energy per link $\overline{\epsilon }=({\epsilon }_0-\epsilon )/{\epsilon }_0$ for the Heisenberg model as a function of the iPEPS bond dimension $D$. Here ${\epsilon }_0=-0.669437\left(5\right)$ is the quantum Monte Carlo result [70]. (b) Staggered magnetization $m$ of the Heisenberg model as a function of the bond dimension $D$, compared to the Monte Carlo result $m_0=0.30703$ [70].

Figure 13

Total running time (in seconds) for five time steps in imaginary-time evolution with the Heisenberg Hamiltonian, with the FU and FFU (plus gauge fixing) updates. The speedup factor of the new update scheme with respect to the old one is shown in the inset.