Bridge++  Ver. 1.1.x
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
contract_4spinor.cpp
Go to the documentation of this file.
1 #include "contract_4spinor.h"
2 #include "index_lex.h"
3 
4 #include <cassert>
5 
6 #define NC 3
7 #define C1 0
8 #define C2 1
9 #define C3 2
10 
11 //====================================================================
12 void contract_at_t(dcomplex& corr, const GammaMatrix& gm,
13  const Field_F& v1, const Field_F& v2, int time)
14 {
15  int Nc = CommonParameters::Nc();
16  int Nd = CommonParameters::Nd();
17  int Nx = CommonParameters::Nx();
18  int Ny = CommonParameters::Ny();
19  int Nz = CommonParameters::Nz();
20  int Nt = CommonParameters::Nt();
21  int Nvol = CommonParameters::Nvol();
22 
23  assert(Nvol == v1.nvol());
24  assert(Nvol == v2.nvol());
25  assert(time < Nt);
26 
27  Index_lex index;
28  valarray<int> gmindex(Nd);
29  valarray<double> corr_r(Nd), corr_i(Nd);
30 
31  for (int i = 0; i < Nd; ++i) {
32  gmindex[i] = gm.index(i);
33  corr_r[i] = 0.0;
34  corr_i[i] = 0.0;
35  }
36 
37  for (int z = 0; z < Nz; ++z) {
38  for (int y = 0; y < Ny; ++y) {
39  for (int x = 0; x < Nx; ++x) {
40  int site = index.site(x, y, z, time);
41 
42  for (int s0 = 0; s0 < Nd; ++s0) {
43  int s1 = gmindex[s0];
44 
45  for (int c1 = 0; c1 < Nc; ++c1) {
46  corr_r[s0] += v1.cmp_r(c1, s1, site) * v2.cmp_r(c1, s0, site)
47  + v1.cmp_i(c1, s1, site) * v2.cmp_i(c1, s0, site);
48 
49  corr_i[s0] += v1.cmp_r(c1, s1, site) * v2.cmp_i(c1, s0, site)
50  - v1.cmp_i(c1, s1, site) * v2.cmp_r(c1, s0, site);
51  }
52  }
53  }
54  }
55  }
56 
57  corr = dcomplex(0.0, 0.0);
58  for (int s0 = 0; s0 < Nd; ++s0) {
59  corr += gm.value(s0) * dcomplex(corr_r[s0], corr_i[s0]);
60  }
61 }
62 
63 
64 //====================================================================
65 void contract_at_t(dcomplex& corr, const GammaMatrix& gm, int d3,
66  const Field_F& v1, const Field_F& v2, const Field_F& v3,
67  int time)
68 {
69  int Nc = CommonParameters::Nc();
70  int Nd = CommonParameters::Nd();
71  int Nx = CommonParameters::Nx();
72  int Ny = CommonParameters::Ny();
73  int Nz = CommonParameters::Nz();
74  int Nt = CommonParameters::Nt();
75  int Nvol = CommonParameters::Nvol();
76 
77  assert(Nc == 3);
78  assert(Nvol == v1.nvol());
79  assert(Nvol == v2.nvol());
80  assert(Nvol == v3.nvol());
81  assert(time < Nt);
82 
83  Index_lex index;
84  valarray<int> gmindex(Nd);
85  valarray<double> cr(Nd), ci(Nd);
86 
87  for (int i = 0; i < Nd; ++i) {
88  gmindex[i] = gm.index(i);
89  cr[i] = 0.0;
90  ci[i] = 0.0;
91  }
92 
93  for (int z = 0; z < Nz; ++z) {
94  for (int y = 0; y < Ny; ++y) {
95  for (int x = 0; x < Nx; ++x) {
96  int site = index.site(x, y, z, time);
97 
98  for (int d1 = 0; d1 < Nd; ++d1) {
99  int d2 = gmindex[d1];
100 
101  cr[d1] += (v1.cmp_r(C1, d1, site) * v2.cmp_r(C2, d2, site)
102  - v1.cmp_i(C1, d1, site) * v2.cmp_i(C2, d2, site)) * v3.cmp_r(C3, d3, site)
103  - (v1.cmp_r(C1, d1, site) * v2.cmp_i(C2, d2, site)
104  + v1.cmp_i(C1, d1, site) * v2.cmp_r(C2, d2, site)) * v3.cmp_i(C3, d3, site);
105  ci[d1] += (v1.cmp_r(C1, d1, site) * v2.cmp_r(C2, d2, site)
106  - v1.cmp_i(C1, d1, site) * v2.cmp_i(C2, d2, site)) * v3.cmp_i(C3, d3, site)
107  + (v1.cmp_r(C1, d1, site) * v2.cmp_i(C2, d2, site)
108  + v1.cmp_i(C1, d1, site) * v2.cmp_r(C2, d2, site)) * v3.cmp_r(C3, d3, site);
109 
110  cr[d1] += (v1.cmp_r(C2, d1, site) * v2.cmp_r(C3, d2, site)
111  - v1.cmp_i(C2, d1, site) * v2.cmp_i(C3, d2, site)) * v3.cmp_r(C1, d3, site)
112  - (v1.cmp_r(C2, d1, site) * v2.cmp_i(C3, d2, site)
113  + v1.cmp_i(C2, d1, site) * v2.cmp_r(C3, d2, site)) * v3.cmp_i(C1, d3, site);
114  ci[d1] += (v1.cmp_r(C2, d1, site) * v2.cmp_r(C3, d2, site)
115  - v1.cmp_i(C2, d1, site) * v2.cmp_i(C3, d2, site)) * v3.cmp_i(C1, d3, site)
116  + (v1.cmp_r(C2, d1, site) * v2.cmp_i(C3, d2, site)
117  + v1.cmp_i(C2, d1, site) * v2.cmp_r(C3, d2, site)) * v3.cmp_r(C1, d3, site);
118 
119  cr[d1] += (v1.cmp_r(C3, d1, site) * v2.cmp_r(C1, d2, site)
120  - v1.cmp_i(C3, d1, site) * v2.cmp_i(C1, d2, site)) * v3.cmp_r(C2, d3, site)
121  - (v1.cmp_r(C3, d1, site) * v2.cmp_i(C1, d2, site)
122  + v1.cmp_i(C3, d1, site) * v2.cmp_r(C1, d2, site)) * v3.cmp_i(C2, d3, site);
123  ci[d1] += (v1.cmp_r(C3, d1, site) * v2.cmp_r(C1, d2, site)
124  - v1.cmp_i(C3, d1, site) * v2.cmp_i(C1, d2, site)) * v3.cmp_i(C2, d3, site)
125  + (v1.cmp_r(C3, d1, site) * v2.cmp_i(C1, d2, site)
126  + v1.cmp_i(C3, d1, site) * v2.cmp_r(C1, d2, site)) * v3.cmp_r(C2, d3, site);
127 
128  cr[d1] -= (v1.cmp_r(C3, d1, site) * v2.cmp_r(C2, d2, site)
129  - v1.cmp_i(C3, d1, site) * v2.cmp_i(C2, d2, site)) * v3.cmp_r(C1, d3, site)
130  - (v1.cmp_r(C3, d1, site) * v2.cmp_i(C2, d2, site)
131  + v1.cmp_i(C3, d1, site) * v2.cmp_r(C2, d2, site)) * v3.cmp_i(C1, d3, site);
132  ci[d1] -= (v1.cmp_r(C3, d1, site) * v2.cmp_r(C2, d2, site)
133  - v1.cmp_i(C3, d1, site) * v2.cmp_i(C2, d2, site)) * v3.cmp_i(C1, d3, site)
134  + (v1.cmp_r(C3, d1, site) * v2.cmp_i(C2, d2, site)
135  + v1.cmp_i(C3, d1, site) * v2.cmp_r(C2, d2, site)) * v3.cmp_r(C1, d3, site);
136 
137  cr[d1] -= (v1.cmp_r(C2, d1, site) * v2.cmp_r(C1, d2, site)
138  - v1.cmp_i(C2, d1, site) * v2.cmp_i(C1, d2, site)) * v3.cmp_r(C3, d3, site)
139  - (v1.cmp_r(C2, d1, site) * v2.cmp_i(C1, d2, site)
140  + v1.cmp_i(C2, d1, site) * v2.cmp_r(C1, d2, site)) * v3.cmp_i(C3, d3, site);
141  ci[d1] -= (v1.cmp_r(C2, d1, site) * v2.cmp_r(C1, d2, site)
142  - v1.cmp_i(C2, d1, site) * v2.cmp_i(C1, d2, site)) * v3.cmp_i(C3, d3, site)
143  + (v1.cmp_r(C2, d1, site) * v2.cmp_i(C1, d2, site)
144  + v1.cmp_i(C2, d1, site) * v2.cmp_r(C1, d2, site)) * v3.cmp_r(C3, d3, site);
145 
146  cr[d1] -= (v1.cmp_r(C1, d1, site) * v2.cmp_r(C3, d2, site)
147  - v1.cmp_i(C1, d1, site) * v2.cmp_i(C3, d2, site)) * v3.cmp_r(C2, d3, site)
148  - (v1.cmp_r(C1, d1, site) * v2.cmp_i(C3, d2, site)
149  + v1.cmp_i(C1, d1, site) * v2.cmp_r(C3, d2, site)) * v3.cmp_i(C2, d3, site);
150  ci[d1] -= (v1.cmp_r(C1, d1, site) * v2.cmp_r(C3, d2, site)
151  - v1.cmp_i(C1, d1, site) * v2.cmp_i(C3, d2, site)) * v3.cmp_i(C2, d3, site)
152  + (v1.cmp_r(C1, d1, site) * v2.cmp_i(C3, d2, site)
153  + v1.cmp_i(C1, d1, site) * v2.cmp_r(C3, d2, site)) * v3.cmp_r(C2, d3, site);
154  }
155  }
156  }
157  }
158 
159  corr = dcomplex(0.0, 0.0);
160  for (int s0 = 0; s0 < Nd; ++s0) {
161  corr += gm.value(s0) * dcomplex(cr[s0], ci[s0]);
162  }
163 }