Tensor Network Theory Library  Beta release 1.2.1
A library of routines for performing TNT-based operations
 All Data Structures Functions Variables Groups Pages
tntMpsMpo.c
1 /*
2 Authors: Sarah Al-Assam, Stephen Clark and Dieter Jaksch
3 Date: $LastChangedDate$
4 (c) University of Oxford 2014
5 */
6 
11 /* Include the internal header for the TNT library */
12 #include "tntMpsInternal.h"
13 
14 
30 void tntMpsMpoConnect(tntNetwork mps,
31  tntNetwork mpo)
32 {
33 
34  tntNode A, O; /* the current MPS node and MPO node */
35  unsigned L; /* Length of system */
36 
37  /* First check that the MPS and the MPO have the same number of terms */
38  L = tntMpsLength(mps);
39 
40  if (tntMpsLength(mpo) != L) {
41  tntErrorPrint("Cannot connect the MPS to the MPO as there are not the same number of nodes in the two networks|" /* NO_COVERAGE */
42  "MPS contains %d nodes, MPO contains %d nodes",L,tntMpsLength(mpo)); /* NO_COVERAGE */
43  } /* NO_COVERAGE */
44 
45  /* Make a copy of the MPO network */
46  mpo = tntNetworkCopy(mpo);
47 
48  /* Find the MPO node for site 0 */
49  O = tntNodeFindFirst(mpo);
50 
51  /* Convert the network to a node group, leaving terminating connections intact */
52  tntNetworkToNodeGroup(&mpo,0);
53 
54  /* Find the MPS node corresponding to site 0 */
55  A = tntNodeFindFirst(mps);
56 
57  /* Now loop through all sites, attaching the copy of the MPO term to the MPS term */
58  while (1) {
59 
60  /* Join a copy of the MPO to the MPS */
61  tntNodeJoin(O,"U",A,"D");
62 
63  /* Find the next group of sites if they exist, or exit the loop if the last site has been reached */
64  if (A == tntNodeFindLast(mps)) {
65  break;
66  } else {
67  /* Find the next nodes along */
68  A = tntNodeFindConn(A, "R");
69  O = tntNodeFindConn(O, "R");
70  }
71  }
72 }
73 
74 /* To add once refs working * The routine performs the zipping algorithm described in \cite Stoudenmire2010 . */
91 double tntMpsMpoContract(tntNetwork mpsmpo,
92  int chi)
93 {
94 
95  tntNode A, O; /* The current node of the MPS and the MPO */
96  tntNode C, Cc; /* The node that is the result of contracting the MPS and its MPO node, and a copy of it taken to SVD */
97  tntNode U, Udag; /* First unitary that is a result of SVDing tnCc, and its conjugate */
98  unsigned L, j; /* The length of the MPS network, site counter */
99  double err_tot=0.0, err; /* total truncation error, and truncation error per SVD. */
100 
101  /* Find the length of the network */
102  L = tntMpsLength(mpsmpo);
103 
104  /* ------ contracting from the left to the right, truncating to truncation tolerace, to result in a single MPS ---- */
105 
106  /* Find the nodes corresponding to site 0 */
107  A = tntNodeFindFirst(mpsmpo);
108  O = tntNodeFindConn(A, "D");
109 
110  /* Contract these two nodes */
111  C = tntNodeContract(A, O, NULL, "LR=MS");
112 
113  /* Fuse the extra leg coming from the MPO with the MPS leg */
114  tntNodeFuse(C,"LM","L","");
115 
116  for (j = 0; j < L-1; j++) {
117 
118  /* Take a copy of tnC - this copy will not be connected to any other nodes in the network */
119  Cc = tntNodeCopy(C);
120 
121  /* Take an untruncated SVD of Cc, assigning the result back to the copy. */
122  Cc = tntNodeSVD(Cc, "LD", "R", "L", "R", "L", -1, &err);
123  err_tot += err;
124 
125  /* Make an unconnected copy of the first unitary - this will become the new MPS node to replace tnA when inserted back in the network */
126  U = tntNodeCopy(Cc);
127 
128  /* the rest of the nodes are not needed -> discard them */
129  tntNodeGroupFree(&Cc);
130 
131  /* Legs L and D of U define one dimension (lets call rows) of the Unitary matrix, leg R of U defines the other dimension (lets call columns) of the unitary matrix */
132 
133  /* Take a conjugate copy of tnU */
134  Udag = tntNodeCopy(U,1);
135 
136  /* Now map the legs of tnUdag so that leg2 -> leg 1 and leg 1 defines rows, and leg 1 -> leg 2 and leg 2 and 3 define columns.
137  * tnUdag is now hermitian conjuate of tnU, and when leg 2 of tnU contracts with leg 1 of tnUdag the result will be the identity */
138  tntNodeMapLegs(Udag, "LR=RL");
139 
140  /* Now insert --tnU--tnUdag-- = the identity (when the physical legs are fused with the outer legs), back into the MPS network */
141  /* Take a different action depending on whether we are at the start of the network, or in the middle of the network */
142  if (0 == j) {
143  /* if we inserting are at the beginning of the network, we can put insert the nodes directly, without temporarily (and illegally) connecting legs with non-matching dimensions */
144 
145  /* Insert the node at the start, which connects leg 2 of tnUdag to leg 1 of tnW */
146  tntNodeInsertAtStart(Udag, "L", "R", mpsmpo);
147 
148  /* Note that we defined legs 1 and 3 of tnC as rows, and legs 2 and 3 of tnUdag as columns - for a full insertion we also need to contract their physical legs */
149  tntNodeJoin(Udag, "D", C, "D");
150 
151  /* Now insert U at the start of the network. U becomes the new MPS node for this set (the previous MPS node is part of C) */
152  tntNodeInsertAtStart(U, "L", "R", mpsmpo);
153  } else {
154  /* if we are inserting in the middle of the network, we can't simply insert the nodes in, as this would cause a mismatch between dimensions of connecting legs */
155  /* This operation is thus forbidden and would cause an error in the program */
156 
157  /* Instead, split the connections between tnC and the previous MPS node, and then join up nodes one by one */
158  /* first identify the previous node */
159  A = tntNodeFindConn(C, "L");
160 
161  /* Split connections between tnA and tnC */
162  tntNodeSplit(A, C);
163 
164  /* connect tnA to tnU */
165  tntNodeJoin(A, "R", U, "L");
166 
167  /* connect tnU to tnUdag */
168  tntNodeJoin(U, "R", Udag, "L");
169 
170  /* Connect tnUdag to tnC, remembering they are joined on 2 legs */
171  tntNodeJoin(Udag, "R", C, "L");
172  tntNodeJoin(Udag, "D", C, "D");
173 
174  /* the network is now complete again */
175  }
176 
177  /* find the next MPS node along, and it's corresponding MPO node */
178  A = tntNodeFindConn(C, "R");
179  O = tntNodeFindConn(A, "D");
180 
181  /* Now contract all these nodes */
182  C = tntNodeListContract("LRDS",Udag,C,A,O);
183 
184  /* If tnA was the last node, tnC simply forms the final MPS node.
185  * Otherwise continue the loop to form further MPS nodes */
186  }
187 
188  /* Fuse the extra leg coming from the MPO with the MPS leg */
189  tntNodeFuse(C,"RS","R","");
190 
191  /* Now sweep back right to left, this time truncating to the maximum internal dimension chi */
192  /* Note setting desired orthogonality centre to zero enures a full right to left sweep */
193  err = tntMpsTruncate(mpsmpo,0,chi);
194  err_tot += err;
195 
196  /* return the total truncation error */
197  return err_tot;
198 
199 }
200 
201 
225 double tntMpsMpoProduct(tntNetwork mps,
226  tntNetwork mpo,
227  int chi)
228 {
229 
230  /* Create the full network */
231  tntMpsMpoConnect(mps,mpo);
232 
233  /* Contract the full network and return the error */
234  return tntMpsMpoContract(mps,chi);
235 
236 }
237 
void tntNodeFuse(tntNode A, tntLegLabel fuseLegs, tntLegLabel newLabelA, tntLegLabel newLabelB)
Definition: tntNodeConn.c:90
tntNode tntNodeCopy(tntNode A)
Definition: tntNodeUtil.c:304
void tntNodeGroupFree(tntNode *A)
Definition: tntNodeUtil.c:288
void tntNodeJoin(tntNode A, tntLegLabel legA, tntNode B, tntLegLabel legB)
Definition: tntNodeConn.c:52
tntNode tntNodeFindLast(tntNetwork nw)
Definition: tntNodeInfo.c:68
double tntMpsTruncate(tntNetwork mps, unsigned orth_centre, int chi)
void tntNodeMapLegs(tntNode A, tntLegLabel legmap)
Definition: tntNodeConn.c:125
tntNode tntNodeSVD(tntNode A, tntLegLabel Ulegs, tntLegLabel legUS, tntLegLabel legSU, tntLegLabel legSV, tntLegLabel legVS, int chi, double *err, tntLegLabel legmap)
unsigned tntMpsLength(tntNetwork wf)
Definition: tntMpsUtil.c:17
double tntMpsMpoProduct(tntNetwork mps, tntNetwork mpo, int chi)
Definition: tntMpsMpo.c:225
void tntMpsMpoConnect(tntNetwork mps, tntNetwork mpo)
Definition: tntMpsMpo.c:30
tntNode tntNodeFindConn(tntNode A, tntLegLabel legA)
Definition: tntNodeInfo.c:29
tntNode tntNodeContract(tntNode A, tntNode B, tntLegLabel legMapAC, tntLegLabel legMapBC)
tntNode tntNodeFindFirst(tntNetwork nw)
Definition: tntNodeInfo.c:50
void tntNodeInsertAtStart(tntNode I, tntLegLabel legInwstart, tntLegLabel legIfirst, tntNetwork nw)
Definition: tntNodeConn.c:213
tntNetwork tntNetworkCopy(tntNetwork nw)
Definition: tntNetwork.c:56
void tntNetworkToNodeGroup(tntNetwork *nwp, int strip_connections)
Definition: tntNetwork.c:268
double tntMpsMpoContract(tntNetwork mpsmpo, int chi)
Definition: tntMpsMpo.c:91
void tntNodeSplit(tntNode A, tntNode B)
Definition: tntNodeConn.c:161