Bridge++  Ver. 1.1.x
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
corr2pt_Wilson_SF.cpp
Go to the documentation of this file.
1 
15 #include "corr2pt_Wilson_SF.h"
16 #include "communicator.h"
17 
18 using std::valarray;
19 using Bridge::vout;
20 
21 //====================================================================
22 
68 //====================================================================
69 double Corr2pt_Wilson_SF::fAfP(const valarray<Field_F>& sq1,
70  const valarray<Field_F>& sq2)
71 {
72  int Lt = CommonParameters::Lt();
73 
74  valarray<dcomplex> mcorr(Lt);
75  GammaMatrix qn_src, qn_sink;
76 
77  int Lvol = CommonParameters::Lvol();
78  int m_Nc = CommonParameters::Nc();
79  double norm = 0.5 / Lvol * Lt;
80  // double norm=0.5/Lvol*Lt*(2*0.130*2*0.130);
81  // vout.general(m_vl,"norm=%lf\n",norm);
82 
83  vout.general(m_vl, "AO correlator:\n");
84 
85  qn_src = m_gmset->get_GM(m_gmset->GAMMA5);
86  qn_sink = m_gmset->get_GM(m_gmset->GAMMA54);
87 
88  meson_corr(mcorr, qn_sink, qn_src, sq1, sq1);
89 
90  for (int t = 0; t < mcorr.size(); ++t) {
91  vout.general(m_vl, "fA %4d %20.12e %20.12e\n",
92  t, -norm * real(mcorr[t]), -norm * imag(mcorr[t]));
93  }
94 
95 
96  vout.general(m_vl, "PO correlator:\n");
97 
98  qn_src = m_gmset->get_GM(m_gmset->GAMMA5);
99  qn_sink = m_gmset->get_GM(m_gmset->GAMMA5);
100  // qn_src = m_gmset->get_GM(m_gmset->UNITY);
101  // qn_sink = m_gmset->get_GM(m_gmset->UNITY);
102 
103  meson_corr(mcorr, qn_sink, qn_src, sq1, sq1);
104 
105  for (int t = 0; t < mcorr.size(); ++t) {
106  vout.general(m_vl, "fP %4d %20.12e %20.12e\n",
107  t, norm * real(mcorr[t]), norm * imag(mcorr[t]));
108  }
109 
110 
111  vout.general(m_vl, "AO' correlator:\n");
112 
113  qn_src = m_gmset->get_GM(m_gmset->GAMMA5);
114  qn_sink = m_gmset->get_GM(m_gmset->GAMMA54);
115 
116  meson_corr(mcorr, qn_sink, qn_src, sq2, sq2);
117 
118  for (int t = 0; t < mcorr.size(); ++t) {
119  vout.general(m_vl, "fA' %4d %20.12e %20.12e\n",
120  t, norm * real(mcorr[Lt - 1 - t]), norm * imag(mcorr[Lt - 1 - t]));
121  }
122 
123 
124  vout.general(m_vl, "PO' correlator:\n");
125 
126  qn_src = m_gmset->get_GM(m_gmset->GAMMA5);
127  qn_sink = m_gmset->get_GM(m_gmset->GAMMA5);
128  // qn_src = m_gmset->get_GM(m_gmset->UNITY);
129  // qn_sink = m_gmset->get_GM(m_gmset->UNITY);
130 
131  meson_corr(mcorr, qn_sink, qn_src, sq2, sq2);
132 
133  for (int t = 0; t < mcorr.size(); ++t) {
134  vout.general(m_vl, "fP' %4d %20.12e %20.12e\n",
135  t, norm * real(mcorr[Lt - 1 - t]), norm * imag(mcorr[Lt - 1 - t]));
136  }
137 
138 
139  //- NB. use meson correlator at t=1 for a non-trivial check.
140  // double result = real(mcorr[0]);
141  double result = real(mcorr[1]);
142 
143  return result;
144 }
145 
146 
147 //====================================================================
148 double Corr2pt_Wilson_SF::meson_corr(valarray<dcomplex>& meson,
149  const GammaMatrix& qn_sink,
150  const GammaMatrix& qn_src,
151  const valarray<Field_F>& sq1,
152  const valarray<Field_F>& sq2)
153 {
154  int Nc = CommonParameters::Nc();
155  int Nd = CommonParameters::Nd();
156  int Nvol = CommonParameters::Nvol();
157  int Lt = CommonParameters::Lt();
158  int Nx = CommonParameters::Nx();
159  int Ny = CommonParameters::Ny();
160  int Nz = CommonParameters::Nz();
161  int Nt = CommonParameters::Nt();
162 
163  assert(meson.size() == Lt);
164 
165  GammaMatrix gm_src, gm_sink, gm5;
166  gm5 = m_gmset->get_GM(m_gmset->GAMMA5);
167  gm_src = qn_src.mult(gm5);
168  gm_sink = gm5.mult(qn_sink);
169 
170  Index_lex index;
171  valarray<dcomplex> corr_local(Nt);
172  valarray<double> corr_r(Nd), corr_i(Nd);
173  valarray<int> s2(Nd);
174 
175  Field corrF(2, Nvol, 2);
176 
177  corr_local = dcomplex(0.0);
178 
179  for (int c0 = 0; c0 < Nc; ++c0) {
180  for (int d0 = 0; d0 < Nd; ++d0) {
181  int d1 = gm_src.index(d0);
182  // int d1 = qn_src.index(d0);
183 
184  for (int t = 0; t < Nt; ++t) {
185  corr_r = 0.0;
186  corr_i = 0.0;
187  for (int z = 0; z < Nz; ++z) {
188  for (int y = 0; y < Ny; ++y) {
189  for (int x = 0; x < Nx; ++x) {
190  int site = index.site(x, y, z, t);
191 
192  for (int s0 = 0; s0 < Nd; ++s0) {
193  int s1 = gm_sink.index(s0);
194  //int s1 = qn_sink.index(s0);
195 
196  for (int c1 = 0; c1 < Nc; ++c1) {
197  corr_r[s0] += sq1[c0 + Nc * d0].cmp_r(c1, s1, site)
198  * sq2[c0 + Nc * d1].cmp_r(c1, s0, site)
199  + sq1[c0 + Nc * d0].cmp_i(c1, s1, site)
200  * sq2[c0 + Nc * d1].cmp_i(c1, s0, site);
201 
202  corr_i[s0] += sq1[c0 + Nc * d0].cmp_r(c1, s1, site)
203  * sq2[c0 + Nc * d1].cmp_i(c1, s0, site)
204  - sq1[c0 + Nc * d0].cmp_i(c1, s1, site)
205  * sq2[c0 + Nc * d1].cmp_r(c1, s0, site);
206  }
207  }
208  }
209  }
210  }
211  for (int s0 = 0; s0 < Nd; ++s0) {
212  dcomplex gmf = gm_src.value(d0) * gm_sink.value(s0);
213  //dcomplex gmf = qn_src.value(d0)*qn_sink.value(s0);
214  dcomplex corr = cmplx(corr_r[s0], corr_i[s0]);
215  corr_local[t] += gmf * corr;
216  //corr_local[t] += corr;
217  }
218  //vout.general(m_vl,"%d %lf %lf\n",t,real(corr_local[t]),imag(corr_local[t]));
219  }
220  }
221  }
222 
223  valarray<dcomplex> corr_tmp(Lt);
224 
225  int ipet = Communicator::ipe(3);
226  for (int t = 0; t < Lt; ++t) {
227  corr_tmp[t] = cmplx(0.0, 0.0);
228  }
229  for (int t = 0; t < Nt; ++t) {
230  int t2 = t + ipet * Nt;
231  // vout.general(m_vl,"%d %d\n",t,t2);
232  corr_tmp[t2] = corr_local[t];
233  }
234  for (int t = 0; t < Lt; ++t) {
235  double crr = real(corr_tmp[t]);
236  double cri = imag(corr_tmp[t]);
237  crr = Communicator::reduce_sum(crr);
238  cri = Communicator::reduce_sum(cri);
239  meson[t] = cmplx(crr, cri);
240  // vout.general(m_vl,"%d %lf %lf\n",t,crr,cri);
241  }
242 
243 
244  //- NB. use meson correlator at t=1 for a non-trivial check.
245  // double result = real(mcorr[0]);
246  double result = real(meson[1]);
247 
248 
249  return result;
250 }
251 
252 
253 //====================================================================
254 //============================================================END=====