Simulation of two-dimensional quantum systems using a tree tensor network that exploits the entropic area law

This work explores the use of a tree tensor network ansatz to simulate the ground state of a local Hamiltonian on a two-dimensional lattice. By exploiting the entropic area law, the tree tensor network ansatz seems to produce quasiexact results in systems with sizes well beyond the reach of exact diagonalization techniques. We describe an algorithm to approximate the ground state of a local Hamiltonian on a $L\times L$ lattice with the topology of a torus. Accurate results are obtained for $L=\left\{4,6,8\right\}$, whereas approximate results are obtained for larger lattices. As an application of the approach, we analyze the scaling of the ground-state entanglement entropy at the quantum critical point of the model. We confirm the presence of a positive additive constant to the area law for half a torus. We also find a logarithmic additive correction to the entropic area law for a square block. The single copy entanglement for half a torus reveals similar corrections to the area law with a further term proportional to $1/L$.

Figure 1

(Color online) Example of TTN for a $2\times 2$ lattice and a $4\times 4$ lattice. Notice (right) that the TTN for a 2D lattice can always be represented as a planar graph, with the leaves or physical indices ordered on a line. The tensors labeled with $w_i$ are isometric tensors. They act locally by projecting the ground state onto its local support with dimension ${\chi }_i$ (see Sec. 2 for further explanation).

Figure 2

(Color online) (i) Diagrammatic representation of three types of isometric tensors in the TTN for a $4\times 4$ lattice in Fig. 1. (ii) Graphical representation of the constraints in Eqs. (3, 4, 5) fulfilled by the isometric tensors.

Figure 3

(Color online) The isometric TTN of Fig. 1 for a $4\times 4$ lattice ${\mathcal{L}}_0$ is associated with a coarse-graining transformation that generates a sequence of increasingly coarse-grained lattices ${\mathcal{L}}_1$, ${\mathcal{L}}_2$, and ${\mathcal{L}}_3$. Notice that in this example we have added an extra index to the top isometry $w_3$, corresponding to the single site of an extra top lattice ${\mathcal{L}}_3$, which we can use to encode in the TTN a whole subspace of ${\mathbb{V}}^{\otimes N}$ instead of a single state $|\Psi ⟩$.

Figure 4

(Color online) Interacting boundaries ${\Sigma }_{1/2}$, ${\Sigma }_{1/4}$, and ${\Sigma }_{1/8}$ corresponding to one half, one quarter and one eighth of a lattice for three different choices of boundary conditions.

Figure 5

(Color online) Isometric TTN for lattices of $6\times 6$, $8\times 8$, and $10\times 10$ sites as used in the manuscript for the purpose of benchmarking the performance of the algorithm of Sec. 4. Notice that all TTN have the same structure on the two top layers of isometries, whose manipulation dominates the computational cost of the algorithm, while they differ in the lower layers. In particular, in the $10\times 10$ lattice two lower layers of isometries are required since a single layer of isometries mapping a block of $5\times 5$ sites directly into a single effective site would have been too expensive given the present capabilities of the desktop computers used for the simulations.

Figure 6

(Color online) Computation of the expectation value $⟨\Psi |o^{\left[s\right]}|\Psi ⟩$ of a one-site operator $o^{\left[s\right]}$ acting on site $s∊\mathcal{L}$. (i) Tensor network to be contracted. (ii) Tensor network left after many of the isometries are annihilated by their Hermitian conjugate, see Fig. 2. After the steps from (iii) to (v) the expectation value is obtained.

Figure 7

(Color online) Computation of the expectation value $⟨\Psi |o^{\left[s\right]}{o}^{\left[s^{\prime }\right]}|\Psi ⟩$ corresponding to a two-site correlation function. (i) Tensor network to be contracted. (ii) Tensor network left after several isometries are annihilated by their Hermitian conjugate. After the steps from (iii) to (v) the expectation value is obtained.

Figure 8

(Color online) Computation of the overlap or fidelity $⟨{\Psi }_1|{\Psi }_2⟩$ between two states $|{\Psi }_1⟩$ and $|{\Psi }_2⟩$ each represented with a TTN. (i) Tensor network corresponding to $⟨{\Psi }_1|{\Psi }_2⟩$. Notice that in this case no isometry is annihilated, since the isometries of the two TTNs are not the same. After the steps from (ii) to (iv) the overlap or fidelity is obtained.

Figure 9

(Color online) Computation of the spectrum $\left\{p_{\alpha }\right\}$ of the reduced density matrix ${\rho }^A$ for a block $A$ that corresponds to one of the coarse-grained sites. (i) Tensor network corresponding to ${\rho }^A$ where $A$ is half of the lattice. (ii) Tensor network left after several isometries are annihilated with their Hermitian conjugate. (iii) since the spectrum of ${\rho }^{\left[A\right]}$ is not changed by the isometries acting on $A$, we can also eliminate those isometries and we are left with a network consisting of only two tensors, which can now be contracted together.

Figure 10

(Color online) Tensor network corresponding to the cost function $E\left(\left\{w_i\right\}\right)=⟨{\Psi }_{\left\{w_i\right\}}|H|{\Psi }_{\left\{w_i\right\}}⟩$ to be minimized.

Figure 11

(Color online) Tensor network representation for the cost function $E\left(w\right)=F\left(w\right)+G$ in Eq. (25) depending only on one isometry $w$.

Figure 12

(Color online) Examples of the three different types of two-site terms $E^{s{s}^{\prime }}$ contributing to $E\left(w\right)$: in (i) both $s$ and $s^{\prime }$ are contained within the block $A$ associated to $w$; in (ii) only one of the sites, say $s$, belongs to $A$; finally in (iii) both sites $s$ and $s^{\prime }$ are outside $A$. The terms (i) and (ii) contribute to $F^{AA}\left(w\right)$ and $F^{AB}\left(w\right)$ in Eq. (26), respectively, whereas the term (iii) contributes to the constant $G$ in 25.

Figure 13

(Color online) Examples of the two different types of two-site terms that contribute to the environment ${\rm Y}$ for the isometry $w$: (i) both $s$ and $s^{\prime }$ are contained within the block $A$ associated to $w$, and therefore this term contributes to ${{\rm Y}}^{AA}$; (ii) only one of the sites, say, $s$, belongs to $A$ and therefore this term contributes to ${{\rm Y}}^{AB}$.

Figure 14

(Color online) Expectation value $⟨{\sigma }_x{\sigma }_x⟩$ as a function of the transverse magnetic field $\lambda$ and for lattices of $4\times 4$, $6\times 6$, and $8\times 8$ sites. Notice that, as the lattice size grows, $⟨{\sigma }_x{\sigma }_x⟩$ becomes steeper and less smooth around $\lambda \approx 3$, consistent with the existence of a critical point at ${\lambda }_c\approx 3.044$ in the thermodynamic limit (Refs. 63, 64).

Figure 15

(Color online) Expectation value $⟨{\sigma }_z⟩$ for the transverse magnetization as a function of the transverse magnetic field $\lambda$ and for lattices of $4\times 4$, $6\times 6$, and $8\times 8$ sites. Again, as the lattice size grows $⟨{\sigma }_z⟩$ becomes steeper and less smooth around $\lambda \approx 3$, consistent with the existence of a critical point at ${\lambda }_c\approx 3.044$ in the thermodynamic limit (Refs. 63, 64).

Figure 16

(Color online) Error of the ground state energy per site $\Delta e\equiv e\left(\chi \right)-e_{exact}$ [cf. Eq. (33)] for a $6\times 6$ lattice plotted as a function of $1/\chi$ at $\lambda =3.05266$. $\chi$ varies in the range from $\chi =100$ to $\chi =550$. $e_{exact}$ is extracted from Ref. 1. In the plot we see that $\Delta e$ ranges from $\Delta e\sim \mathcal{O}\left({10}^{-4}\right)$ in the case of $\chi =100$ to $\Delta e\sim \mathcal{O}\left({10}^{-7}\right)$ in the case of $\chi =550$. All the TTN simulations used here require much smaller computational resources than the full exact diagonalization calculation as we explain in detail in the main text.

Figure 17

(Color online) Approximate ground state energy per site $e$ [cf. Eq. (33)] for $\lambda =3.05$ plotted as a function of $1/\chi$, for lattices of $6\times 6$ and $8\times 8$ sites. Notice that in both cases the results seem to have converged for large $\chi$ up to several digits of accuracy. The insets show the difference $\Delta e\equiv e\left(\chi \right)-e\left({\chi }_{max}\right)$, as a function of $1/\chi$.

Figure 18

(Color online) Spectum $\left\{p_{\alpha }\right\}$ for the reduced density matrix $\rho$ of one half of the lattice. The results the ground state of $H_{\text{Ising}}$ for $\lambda =3.05$ in a $8\times 8$ lattice. Notice the relatively fast decay of the spectrum, with, e.g., $p_{\alpha }<{10}^{-4}$ for $\alpha >50$. Also, calculations with $\chi =200$ and $\chi =500$ produce spectra that are very similar for small $\alpha$. This is an indication that the largest eigenvalues ($p_{\alpha }$ for small $\alpha$) are already very close to their exact value.

Figure 19

(Color online) Spectrum of the reduced density matrix for one half the $8\times 8$ lattice for different values of $\lambda$. The calculations, conducted with $\chi =100$, show that the spectrum decays slowest for $\lambda$ near ${\lambda }_c$. It also shows that for magnetic fields smaller than ${\lambda }_c$, the spectrum develops a clear structure of plateaux.

Figure 20

(Color online) Two-point correlation functions $C\left(x,y\right)$ of Eq. (34) for the ground state of $H_{\text{Ising}}$ with transverse magnetic field $\lambda =3.05$. The correlation function between distant points in the torus remain large, as one would expect of a system that becomes critical in the thermodynamic limit. The results, obtained with a TTN with $\chi =100$, show that the invariance of the system under $90°$ rotations is preserved, in spite of the fact that the TTN manifestly breaks it at its top layers. Indeed, one can hardly distinguish $C\left(r,0\right)$ from $C\left(0,r\right)$.

Figure 21

(Color online) Expectation value $⟨{\sigma }_x{\sigma }_x⟩$ as a function of the transverse magnetic field $\lambda$ and for a $10\times 10$ lattice. The inset shows results obtained with $\chi =100$ and $\chi =200$ for values of the transverse magnetic field $\lambda$ close to ${\lambda }_c$. In this approximate regime, the TTN algorithm produces results that are not converged with respect to $\chi$ near the quantum critical point.

Figure 22

(Color online) Expectation value $⟨{\sigma }_z⟩$ as a function of the transverse magnetic field $\lambda$ and for a $10\times 10$ lattice. The inset shows results obtained with $\chi =100$ and $\chi =200$ for values of the transverse magnetic field $\lambda$ close to ${\lambda }_c$.

Figure 23

(Color online) Spectrum $\left\{p_{\alpha }\right\}$ of the reduced density matrix for one half of a $L\times L$ lattice in the ground state of $H_{\text{Ising}}$ for $\lambda =2.4$. These results, obtained with only $\chi =100$, show the presence of a plateau of exactly $2L$ eigenvalues $p_{\alpha }$, separated by two or more orders of magnitude from those of the next plateau. The structure of plateaux can be understood as a perturbative version of the entropic area law and explains why a TTN with relatively small $\chi$ can still produce converged results away from the critical point for large lattices $L\approx 10–30$. For a related discussion see also the Sec. IIIB in Ref. 45.

Figure 24

(Color online) Spectrum $\left\{p_{\alpha }\right\}$ of one half of a $10\times 10$ lattice in the ground state of $H_{\text{Ising}}$ for $\lambda =2.4$. Results obtained with $\chi =32$, 64 and 100 do not differ significantly in the first 21 eigenvalues. This shows that the presence and composition of the first plateau of $2L$ eigenvalues of Fig. 23 (with $L=10$ in this case) is robust with respect to $\chi$.

Figure 25

(Color online) Plot of $\Delta e$ defined in Eq. (35) for torus of sizes $L=\left\{4,6,8,10\right\}$ and various $\chi$. We present several choices of $\epsilon$ in the range ${10}^{-4}\le \epsilon \le 1.3\times {10}^{-3}$ represented by horizontal lines of different colors. The arrows identify the ${\chi }_{\epsilon }$ defined in the text for the particular choice of $\epsilon ={10}^{-3}$.

Figure 26

(Color online) Plot of the logarithm of ${\chi }_{\epsilon }$ as a function of $L$ for several choices of $\epsilon$ in the range ${10}^{-4}\le \epsilon \le 1.3\times {10}^{-3}$. In each case, for $L$ large enough, the data lie on a straight line confirming that the cost of the simulation increases exponentially with the size $L$ of the system.

Figure 27

(Color online) Analysis of $-S_{1/2}$ (top plot), minus the entanglement entropy of half a torus, for the critical Ising model on a $4\times 4$ lattice, as a function of $1/\chi$ and the ground state energy (bottom plot). The exact result requires $\chi =256$. If we are only able to compute both observables for $\chi \le 32$ (the region in the plots on the right of the vertical black line) the values of both observables at $\chi =32$ provide an upper bound to their exact results. A lower bound is obtained by linearly extrapolating the two observables from the two best values at $\chi =16$ and $\chi =32$ to $1/\chi =0$. This is represented in the plots by the lower red line on the plots where the extrapolated values are denoted by $S_{low}$ and $E_{low}$. In this way it is clear that $-S_{low}\le -S_{exact}\le -S\left(1/{\chi }_{max}\right)$ and $E_{low}\le E_{exact}\le E\left(1/{\chi }_{max}\right)$ as expressed in Eq. (36) in the main text. The differences $\left|S\left(1/{\chi }_{max}\right)-S_{low}\right|$ and $\left|E\left(1/{\chi }_{max}\right)-E_{low}\right|$ thus provide upper bounds to the error for both $E$ and $S$ induced by considering a smaller $\chi$ than the one required by the exact solution.

Figure 28

(Color online) Entropy $S_{1/2}\left(L\right)$ of one half of the torus as a function of the linear size $L$, corresponding to the ground state of $H_{\text{Ising}}$ with $\lambda =3.044$. The results for $L=4$, 6, 8, and 10 confirm the linear growth predicted in Eq. (37), where no logarithmic correction is expected. The results of our study also seem to rule out the presence of the term proportional to $s_{-1}$. In the inset we show the results of the fit. The asymmetric errorbars, obtained through the analysis outlined in Sec. 5D, are so small that they are hardly visible in the plot.

Figure 29

(Color online) Difference between entropies for one quarter and one half of the torus as a function of the linear size $L$, corresponding to the ground state of $H_{\text{Ising}}$ with $\lambda =3.044$. The results for $L=4$, 6, 8, and 10 allow us to confirm the logarithmic dependence predicted in Eq. (44), which is attributed to the presence of corners in the boundary of our block for one quarter of the lattice.

Figure 30

(Color online) Single copy entanglement $E^{\left(1\right)}\left(L\right)$ for one half of the torus as a function of the linear size $L$, corresponding to the ground state of $H_{\text{Ising}}$ with $\lambda =3.044$. The results for $L=4$, 6, 8, and 10 confirm the linear growth predicted in Eq. (48). In this case the term $e_{-1}$, analogous to the term $s_{-1}$ for the entanglement entropy, is negative as predicted by the theory. In the inset we present the results of the fit. The asymmetric error bars, obtained through the analysis outlined in Sec. 5D, are so small that they are hardly visible on the plot.