Tensor Network (Part 3)

Anurag

2022/12/01 07:27

Import Modules


#
import numpy as np
from numpy import linalg as LA
# ensure 'ncon.py' in working dir
from ncon import ncon

Tree tensor networks

In this tutorial we shall focus only on tensors networks that do not possess closed loops (i.e. they are described by acyclic graphs). This class of tensor network, which we generically refer to as tree tensor networks, possess many nice properties that networks containing closed loops lack and are thus much easier to manipulate. Fig.3.1(a) presents an example of a tree tensor network. If we select a tensor to act as the center (or root node) then it is always possible to understand the tree tensor network as being composed of a set of distinct branches extending from this chosen tensor. For instance, Fig.3.1(b) depicts the four branches (including one trivial branch) extending from the order-4 tensor A from Fig.3.1(a). Importantly, connections between the different branches are not possible in networks without closed loops.

Gauge freedom

Let T be a tensor network that, under contraction of all internal indices, evaluates to some tensor D. In this tutorial we shall concern ourselves with the uniqueness of the decomposition: is there a different choice of tensors within the network that will still evaluate to the same tensor D?

Clearly the answer is yes! As shown below in Fig.3.2(a-b), on any internal index of the network one can introduce a resolution of the identity (i.e. a pair of matrices X and X^(-1)) which, by construction, does not change the final product that the network evaluates to. However absorbing one of these matrices into each adjoining tensor does change their content (while leaving the geometry of the network unchanged). Thus we conclude that there are infinitely many choices of tensors such that the network product evaluates to some fixed output tensor. We refer to this ability to introduce an arbitrary resolution of the identity on an internal index as the gauge freedom of the network

While in some respects the gauge freedom is a nuisance (as it implies tensor decompositions are never unique), it can also be exploited to simplify many types of operations on tensor networks. Indeed, most tensor network algorithms require fixing the gauge in a prescribed manner in order to function correctly. We now discuss several ways to fix the gauge degree of freedom in such a way as to create a center of orthogonality, and the utility of doing so.

Creating a center of orthogonality

Center of Orthogonality

Let T:{A,B,C,…} be a tree tensor network, then a tensor A is a center of orthogonality if, for every branch of the network attached to A, the branch forms an isometry between its open indices and the index connected to tensor A.

In the example above, the tensor A from the network T in Fig.3.3(a) is a center of orthogonality if and only if the constraints of Fig.3.3(b) are satisfied, which demand that each of the branches connected to A forms an isometry. Here, as with Tutorial 2, the conjugate B† of a tensor B denotes complex conjugation as well as opposite vertical orientation in figures.

We now discuss two different methods for changing the gauge in network T to make any tensor A into center of orthogonality, before later revealing the significance of doing so.

Setting a center of orthogonality method 1: ‘Pulling Through’

Here we describe a method for setting a tensor A within a network T as a center of orthogonality through iterative use of the QR decomposition. (Alternatively one can use the SVD to achieve the same effect, although the QR decomposition is usually preferred as it is computationally quicker). The idea behind his method is very simple: if we transform every individual tensor within a branch into a (properly oriented) isometry, then the entire branch collectively becomes an isometry and thus satisfies the requirement of Def.3.3.

Using the network of Fig.3.3(a) as an example, begin by orienting each index with an arrow that points towards the chosen center tensor A.

Then, starting from a tensor at the tip of a branch, perform a QR decomposition on the tensor (under the partition between incoming and outgoing arrows). Next redefine the tensor in question as the orthogonal ‘Q’ part of the QR decomposition and absorb the ‘R’ matrix into the tensor connected to the outgoing arrow.

Repeat this procedure, working inwards, until all tensors are isometric w.r.t. their incoming and outgoing arrows. Tensor A is now a center of orthogonality, which follows since isometric tensors necessarily satisfy the branch constraints of Fig.3.3(b).

Fig.3.3(c): through iterative use of the QR decomposition, each tensor in the network (except the tensor chosen as the center) is given as the 'Q' part of a QR decomposition and is thus isometric. The code example Ex.3.3(c) (below) illustrates how the sequence of operations in Fig.3.3(c) can be implemented numerically, and demonstrates that the initial and final networks still contract to the same tensor. Note that we use a convention such tensor indices are ordered from left-to-right along the bottom then left-to-right along the top.


# Ex.3.3(c): Creating a center of orthogonality by 'pulling through' 
# define tensors
d = 3
A = np.random.rand(d,d,d,d); B = np.random.rand(d,d,d) 
C = np.random.rand(d,d,d); D = np.random.rand(d,d,d) 
E = np.random.rand(d,d,d); F = np.random.rand(d,d,d) 
G = np.random.rand(d,d,d)
# iterate QR decomps
DQ, DR = LA.qr(D.reshape(d**2,d)); DQ = DQ.reshape(d,d,d)
EQ, ER = LA.qr(E.reshape(d**2,d)); EQ = EQ.reshape(d,d,d)
Btilda = ncon([B,DR,ER],[[1,2,-3],[-1,1],[-2,2]])
BQ, BR = LA.qr(Btilda.reshape(d**2,d)); BQ = BQ.reshape(d,d,d)
FQ, FR = LA.qr(F.reshape(d**2,d)); FQ = FQ.reshape(d,d,d)
GQ, GR = LA.qr(G.reshape(d**2,d)); GQ = GQ.reshape(d,d,d)
Ctilda = ncon([C,GR],[[1,-2,-3],[-1,1]])
CQ, CR = LA.qr(Ctilda.reshape(d**2,d)); CQ = CQ.reshape(d,d,d)
Aprime = ncon([A,BR,FR,CR],[[1,-2,2,3],[-1,1],[-3,2],[-4,3]])
# new network is formed from tensors: {Aprime,BQ,CQ,DQ,EQ,FQ,GQ}.

# check both networks evaluate to the same tensor
connectlist = [[3,-5,4,5],[1,2,3],[6,-10,5],[-1,-2,1],[-3,-4,2],[-6,-7,4],[-8,-9,6]]
H0 = ncon([A,B,C,D,E,F,G],connectlist)
H1 = ncon([Aprime,BQ,CQ,DQ,EQ,FQ,GQ],connectlist)
dH = LA.norm(H0-H1) / LA.norm(H0)

Setting a center of orthogonality method 2: ‘Direct Orthogonalization’

Here we describe a method for setting a tensor A within a network T as a center of orthogonality directly using a single eigen-decomposition for each branch, again using the network of Fig.3.3(a) as an example.

Begin by computing the positive-definite 'density matrix' ρ associated to each index about the chosen center A; this is given by contracting the open indices from a branch with the corresponding open indices from the conjugate of the branch as seen below.

Then computes the principle square root X of each of the density matrices ρ.

3) Finally, we make a change of gauge on each of the indices of tensor A using the appropriate X matrix and its corresponding inverse, as depicted below in Fig.3.3(e). Tensor A is now a center of orthogonality as the constraints of Fig.3.3(b) have been satisfied by construction.


# Ex.3.3(c): Creating a center of orthogonality with 'direct orthogonalization' 
# define tensors
d = 3
A = np.random.rand(d,d,d,d); B = np.random.rand(d,d,d) 
C = np.random.rand(d,d,d); D = np.random.rand(d,d,d) 
E = np.random.rand(d,d,d); F = np.random.rand(d,d,d) 
G = np.random.rand(d,d,d)
# compute density matrices and their principle square roots
rho1 = ncon([B,D,E,B,D,E],[[5,6,-2],[1,2,5],[3,4,6],[7,8,-1],[1,2,7],[3,4,8]])
rho2 = ncon([F,F],[[1,2,-2],[1,2,-1]])
rho3 = ncon([C,G,C,G],[[3,5,-2],[1,2,3],[4,5,-1],[1,2,4]])
d1, u1 = LA.eigh(rho1); sq_d1 = np.sqrt(abs(d1))
d2, u2 = LA.eigh(rho2); sq_d2 = np.sqrt(abs(d2))
d3, u3 = LA.eigh(rho3); sq_d3 = np.sqrt(abs(d3))
X1 = u1 @ np.diag(sq_d1) @ u1.T; X1inv = u1 @ np.diag(1/sq_d1) @ u1.T
X2 = u2 @ np.diag(sq_d2) @ u2.T; X2inv = u2 @ np.diag(1/sq_d2) @ u2.T
X3 = u3 @ np.diag(sq_d3) @ u3.T; X3inv = u3 @ np.diag(1/sq_d3) @ u3.T
# execute gauge changes
Aprime = ncon([A,X1,X2,X3],[[1,-2,2,3],[-1,1],[-3,2],[-4,3]])
Bprime = ncon([B,X1inv],[[-1,-2,1],[1,-3]])
Fprime = ncon([F,X2inv],[[-1,-2,1],[1,-3]])
Cprime = ncon([C,X3inv],[[-1,-2,1],[1,-3]])
# new network is formed from tensors: {Aprime,Bprime,Cprime,D,E,Fprime,G}

# check both networks evaluate to the same tensor
connectlist = [[3,-5,4,5],[1,2,3],[6,-10,5],[-1,-2,1],[-3,-4,2],[-6,-7,4],[-8,-9,6]]
H0 = ncon([A,B,C,D,E,F,G],connectlist)
H1 = ncon([Aprime,Bprime,Cprime,D,E,Fprime,G],connectlist)
dH = LA.norm(H0 - H1) / LA.norm(H0)

Comparison: both of the two methods discussed to create a center of orthogonality have their own advantages, and the preferred method may depend on the specific application in mind.

In practice 'direct orthogonalization' is typically computation cheaper and easier to execute. In addition this method only requires changing the gauge on the indices connected to the center, whereas the 'pulling through' method involves changing the gauge on all indices of the network. However there are some applications where it is desired to make every tensor into an isometry, as is achieved with 'pulling through'. In addition 'pulling through' can be advantageous if high precision is desired as the errors due to floating-point arithmetic are lesser (especially so if the condition number of the branch density matrices ρ is bad).

Tensor decompositions within networks

In the previous tutorial we described how the SVD can be applied to optimally decompose a tensor into a product with some restricted rank (i.e. as to minimize the Frobenius norm between the original tensor and the decomposition). Here we take this concept further and describe how, by creating a center of orthogonality, a tensor within a network can be optimally decomposed as to minimize the global error from the entire network.

Let us consider a network {A,B,C,D,E,F,G} that evaluates to tensor H, as depicted in Fig.3.4(a). Then, under replacement of A with some new tensor A', we call the new product H' as depicted in Fig.3.4(b).