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
tntMpsCreateBosonOp.c
1 /*
2 Authors: Sarah Al-Assam, Stephen Clark and Dieter Jaksch
3 $LastChangedDate$
4 (c) University of Oxford 2014
5 */
6 
12 /* Include the header for the TNT library */
13 #include "tntMpsInternal.h"
14 
29 void tntMpsCreateBosonOp(unsigned n_max,
30  tntNode *b,
31  tntNode *bdag,
32  tntNode *n,
33  tntNode *os_int,
34  tntNode *eye)
35 {
36  tntComplexArray barr; /* Arrays to create the operators: only b is created from elements, and the remaining arrays created from this */
37  unsigned physdim = n_max+1; /* Physical dimension is 2 for a spin half system */
38  unsigned j; /* Row number */
39  tntNode tnA, tnB, tnC, tnD; /* Used for contracting nodes */
40  tntNode bl, bdagl; /* Local copy of nodes */
41  tntIntArray qninf; /* Quantum number information for the physical legs */
42 
43  /* Quantum number information to set for invariant operators only */
44  if (TNT_SYMM_U1 == tntSymmTypeGet()) {
45  qninf = tntIntArrayAlloc(physdim);
46  for (j = 1; j < physdim; j++) qninf.vals[j] = j;
47  }
48 
49  /* Check whether each node pointer is NULL, before creating each node
50  (apart from b and bdag, which need to be created to create all the other nodes, but will be freed and set back to NULL at the end
51  of the routine if not required). */
52 
53  /* \b = [ 0 1 0 0 0 ...]
54  [ 0 0 sqrt(2) 0 0 ...]
55  [ 0 0 0 sqrt(3) 0 ...]
56  [ 0 0 0 0 2 ...]
57  [ 0 0 0 0 0 ...]
58  [ ..... ] */
59  /* Create empty array */
60  barr = tntComplexArrayAlloc(physdim, physdim);
61 
62  /* Put elements in row j-1, col j */
63  for (j = 1; j < physdim; j++) {
64  barr.vals[j-1 + j*physdim].re = sqrt((double) j);
65  }
66 
67  /* Create node from elements */
68  bl = tntNodeCreate(&barr, "DU", physdim, physdim);
69 
70  /* Free the array - it is no longer needed */
71  tntComplexArrayFree(&barr);
72 
73  /* take transpose of b to create bdag by swapping legs */
74  bdagl = tntNodeCopy(bl);
75  tntNodeMapLegs(bdagl,"UD=DU");
76 
77  /* n = bdag*b (U(1) invariant)
78  Make node through contracting annihilation and creation operators */
79  if (NULL != n) {
80  tnA = tntNodeCopy(bdagl);
81  tnB = tntNodeCopy(bl);
82  tntNodeJoin(tnA,"U",tnB,"D");
83 
84  *n = tntNodeContract(tnA, tnB);
85 
86  if (TNT_SYMM_U1 == tntSymmTypeGet()) {
87  tntNodeSetQN(*n,"D",&qninf,TNT_QN_IN);
88  tntNodeSetQN(*n,"U",&qninf,TNT_QN_OUT);
89  }
90  }
91 
92  /* os_int = bdag*bdag*b*b (U(1) invariant)
93  Make node through contracting annihilation and creation operators */
94  if (NULL != os_int) {
95  /* Make a copy of all required operators */
96  tnA = tntNodeCopy(bl);
97  tnB = tntNodeCopy(bl);
98  tnC = tntNodeCopy(bdagl);
99  tnD = tntNodeCopy(bdagl);
100 
101  /* Join them together in sequence i.e. s.t. such that annihilation operator is closest to MPS (i.e. the top node) */
102  tntNodeJoin(tnA,"D",tnB,"U");
103  tntNodeJoin(tnB,"D",tnC,"U");
104  tntNodeJoin(tnC,"D",tnD,"U");
105 
106  *os_int = tntNodeListContract(NULL, tnA, tnB, tnC, tnD);
107 
108  if (TNT_SYMM_U1 == tntSymmTypeGet()) {
109  tntNodeSetQN(*os_int,"D",&qninf,TNT_QN_IN);
110  tntNodeSetQN(*os_int,"U",&qninf,TNT_QN_OUT);
111  }
112  }
113 
114  /* eye = [ 1 0]
115  [ 0 1] (U(1) invariant) */
116  if (NULL != eye) {
117  *eye = tntNodeCreateEyeOp(bl);
118 
119  if (TNT_SYMM_U1 == tntSymmTypeGet()) {
120  tntNodeSetQN(*eye,"D",&qninf,TNT_QN_IN);
121  tntNodeSetQN(*eye,"U",&qninf,TNT_QN_OUT);
122  }
123  }
124 
125  /* Free the qn information if it was used */
126  if (TNT_SYMM_U1 == tntSymmTypeGet()) tntIntArrayFree(&qninf);
127 
128  /* Free the nodes for b and bdag if not required, otherwise assign to pointer arguments */
129  if (NULL == b) tntNodeFree(&bl);
130  else *b = bl;
131 
132  if (NULL == bdag) tntNodeFree(&bdagl);
133  else *bdag = bdagl;
134 
135  return;
136 }
tntNode tntNodeCopy(tntNode A)
Definition: tntNodeUtil.c:304
void tntNodeSetQN(tntNode A, tntLegLabel legA, tntIntArray *qvals, int legdir)
Definition: tntNodeQN.c:37
void tntNodeJoin(tntNode A, tntLegLabel legA, tntNode B, tntLegLabel legB)
Definition: tntNodeConn.c:52
tntNode tntNodeCreateEyeOp(tntNode A)
Definition: tntNodeUtil.c:195
void tntNodeMapLegs(tntNode A, tntLegLabel legmap)
Definition: tntNodeConn.c:125
void tntMpsCreateBosonOp(unsigned n_max, tntNode *b, tntNode *bdag, tntNode *n, tntNode *os_int, tntNode *eye)
void tntNodeFree(tntNode *A)
Definition: tntNodeUtil.c:275
tntComplexArray tntComplexArrayAlloc(unsigned numrows, unsigned numcols)
Definition: tntArray.c:392
void tntIntArrayFree(tntIntArray *arr)
Definition: tntArray.c:620
tntIntArray tntIntArrayAlloc(unsigned numrows, unsigned numcols)
Definition: tntArray.c:29
tntNode tntNodeContract(tntNode A, tntNode B, tntLegLabel legMapAC, tntLegLabel legMapBC)
tntNode tntNodeCreate(tntComplexArray *nodeVals, tntLegLabel leglabels, unsigned dimleg1, unsigned dim...)
Definition: tntNodeUtil.c:354
void tntComplexArrayFree(tntComplexArray *arr)
Definition: tntArray.c:650
int tntSymmTypeGet(void)
Definition: tntSys.c:498