Tensor Network Theory Library  Beta release 1.2.1 A library of routines for performing TNT-based operations
tntMpsPmpoMps.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 #ifdef MAKINGPUBLICDOCS
15
58 #else
59 tntComplex tntMpsPmpoMpsProduct_(
60 #endif
61  tntNetwork mps,
62  tntNodeArray *op,
63  tntIntArray *sitenum,
65  int orth_centre,
66  tntNodeArray *betas, /* !< precontracted matrices for the left part of the network. Leave out or pass NULL if not used */
67  tntNodeArray *gammas) /* !< precontracted matrices for the left part of the network. Leave out or pass NULL if not used */
68 {
69
70  tntNetwork sandwich; /* The full MPS-MPO-MPS network */
71  int L; /* Number of sites in the MPS */
72  int site_left, site_right; /* Minimum and maximum site numbers for contracting the network */
73  unsigned loop; /* Used for looping over the elements in the arrays */
74  // unsigned usempo, allsymm, mpssymm; /* Various flags for symmetry properties */
75  tntNode beta, gamma; /* Value to use for beta and gamma */
76
77  /* Find length of the MPS */
78  L = (int) tntMpsLength(mps);
79
80  /* Check orth_centre is not too big */
81  if (orth_centre >= L) tntErrorPrint("The site for orth_centre (%d) is more than or equal to the length of the MPS (%d)", orth_centre, L); /* NO_COVERAGE */
82
83  if (op->sz != sitenum->sz) tntErrorPrint("The number of sites listed (%d) does not match the number of operators listed (%d)", op->sz, sitenum->sz); /* NO_COVERAGE */
84
85  /* Initialise values for finding the minimum and maximum site number in the sitenum array */
86  site_left = L-1;
87  site_right = 0;
88
89  /* Check that the values given for the site numbers are sensible. At the same time, find the minimum and maximum site number that has an operator on */
90  for (loop = 0; loop < op->sz; loop++) {
91  if (sitenum->vals[loop] < 0) {
92  tntErrorPrint("The site (%d) for operator %d is less than zero", sitenum->vals[loop], loop); /* NO_COVERAGE */
93  } else if (sitenum->vals[loop] >= L) {
94  tntErrorPrint("The site (%d) for operator %d is more than or equal to the length of the MPS (%d)", sitenum->vals[loop], loop, L); /* NO_COVERAGE */
95  } else {
96  /* Find the minimum and maximum site number so far encountered */
97  if (sitenum->vals[loop] < site_left) {
98  site_left = sitenum->vals[loop];
99  }
100  if (sitenum->vals[loop] > site_right) {
101  site_right = sitenum->vals[loop];
102  }
103  }
104  }
105
106  /* Now determine what parts of the network need to be kept depending on the orthormalilty centre and whether beta and gamma are provided */
107  if (orth_centre >= site_left) {
108  beta = NULL;
109  } else {
110  if (NULL == betas) {
111  beta = NULL;
112  site_left = -1;
113  } else {
114  beta = betas->vals[site_left];
115  }
116  }
117
118  if (orth_centre <= site_right && orth_centre > -1) {
119  gamma = NULL;
120  } else {
121  if (NULL == gammas) {
122  gamma = NULL;
123  site_right = -1;
124  } else {
125  gamma = gammas->vals[site_right];
126  }
127  }
128
129  /* Create a site-wide MPO to preserve QN information if the mps is quantum number conserving */
130 // if (usempo) {
131 // if (beta != NULL) {
132 // tntNodeMapLegs(beta,"MS=NT");
134 // }
135 // if (gamma != NULL) {
136 // tntNodeMapLegs(gamma,"MS=NT");
138 // }
139 //
140 // tntNetwork pmpo = tntMpsCreateProductMpo(tntMpsLength(mps), op, sitenum);
141 //
142 // return tntMpsMpoMpsPartProduct(mps, pmpo, site_left, site_right, beta, gamma);
143 //
144 // } else {
145
146  /* the full network will be altered then deleted during the network contraction, so use a copy so that the original mps is untouched. */
147  sandwich = tntNetworkCopy(mps);
148
149  /* Take the product of the MPS with the PMPO */
150  tntMpsPmpoProduct(sandwich,op,sitenum);
151
152  /* Connect the altered MPS to a flipped copy of the original MPS */
153  tntMpsMpsConnect(sandwich,mps);
154
155  /* Contract the full network */
156  return tntMpsMpsContract(sandwich, site_left, site_right, beta, gamma);
157  // }
158 }
159
160
161
170 double tntMpsNodeNodeMpsProduct(tntNetwork mps,
171  tntNode op,
172  unsigned sitenum)
173 {
174
175  double val; /* value to return */
176  tntIntArray sitenum_arr = tntIntArrayAlloc(2); /* Current site number - two operators applied to this site */
177  tntNodeArray op_arr = tntNodeArrayAlloc(2); /* Array to hold current operator and h.c. */
178
179
180  sitenum_arr.vals[0] = sitenum_arr.vals[1] = sitenum;
181  op_arr.vals[0] = tntNodeCopy(op);
182
183  /* Take the hermitian conjugate of the second node */
184  op_arr.vals[1] = tntNodeCopy(op,1); /* complex conjugate */
185  tntNodeMapLegs(op_arr.vals[1],"DU=UD"); /* transpose */
186
188  val = tntComplexToReal(tntMpsPmpoMpsProduct(mps, &op_arr, &sitenum_arr, 0));
189
190  tntIntArrayFree(&sitenum_arr);
191  tntNodeArrayFree(&op_arr);
192
193  return val;
194 }
double tntComplexToReal(tntComplex var)
Definition: tntUtil.c:108
tntNode tntNodeCopy(tntNode A)
Definition: tntNodeUtil.c:304
void tntNodeArrayFree(tntNodeArray *arr)
Definition: tntArray.c:668
tntNodeArray tntNodeArrayAlloc(unsigned num_elems)
Definition: tntArray.c:428
void tntNodeMapLegs(tntNode A, tntLegLabel legmap)
Definition: tntNodeConn.c:125
tntComplex tntMpsMpsContract(tntNetwork sandwich, int site_left, int site_right, tntNode beta, tntNode gamma)
Definition: tntMpsMps.c:150
void tntMpsMpsConnect(tntNetwork mpsA, tntNetwork mpsB)
Definition: tntMpsMps.c:38
unsigned tntMpsLength(tntNetwork wf)
Definition: tntMpsUtil.c:17
Definition: tnt.h:65
void tntIntArrayFree(tntIntArray *arr)
Definition: tntArray.c:620
tntIntArray tntIntArrayAlloc(unsigned numrows, unsigned numcols)
Definition: tntArray.c:29
tntComplex tntMpsPmpoMpsProduct(tntNetwork mps, tntNodeArray *op, tntIntArray *sitenum, int orth_centre, tntNodeArray *betas, tntNodeArray *gammas)
Definition: tntMpsPmpoMps.c:57
tntNetwork tntNetworkCopy(tntNetwork nw)
Definition: tntNetwork.c:56
void tntMpsPmpoProduct(tntNetwork mps, tntNodeArray *op, tntIntArray *sitenum)
Definition: tntMpsPmpo.c:35
void tntSysQNClearWarnOff(void)
Definition: tntUtil.c:35
double tntMpsNodeNodeMpsProduct(tntNetwork mps, tntNode op, unsigned sitenum)