Bridge++  Version 1.4.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
fopr_CloverTerm_General.cpp
Go to the documentation of this file.
1 
15 
17 
18 namespace Imp {
19 #if defined USE_GROUP_SU3
20 #include "fopr_Wilson_impl_SU3.inc"
21 #elif defined USE_GROUP_SU2
22 #include "fopr_Wilson_impl_SU2.inc"
23 #elif defined USE_GROUP_SU_N
24 #include "fopr_Wilson_impl_SU_N.inc"
25 #endif
26 
27 //====================================================================
28 
29  const std::string Fopr_CloverTerm_General::class_name = "Imp::Fopr_CloverTerm_General";
30 
31 //====================================================================
33  {
34  const std::string str_vlevel = params.get_string("verbose_level");
35 
36  m_vl = vout.set_verbose_level(str_vlevel);
37 
38  //- fetch and check input parameters
39  double kappa_s, kappa_t, cSW_s, cSW_t;
40  std::vector<int> bc;
41 
42  int err = 0;
43  err += params.fetch_double("hopping_parameter_spatial", kappa_s);
44  err += params.fetch_double("hopping_parameter_temporal", kappa_t);
45  err += params.fetch_double("clover_coefficient_spatial", cSW_s);
46  err += params.fetch_double("clover_coefficient_temporal", cSW_t);
47  err += params.fetch_int_vector("boundary_condition", bc);
48 
49  if (err) {
50  vout.crucial(m_vl, "Error at %s: input parameter not found.\n", class_name.c_str());
51  exit(EXIT_FAILURE);
52  }
53 
54 
55  set_parameters(kappa_s, kappa_t, cSW_s, cSW_t, bc);
56  }
57 
58 
59 //====================================================================
60  void Fopr_CloverTerm_General::set_parameters(double kappa_s, double kappa_t,
61  double cSW_s, double cSW_t,
62  std::vector<int> bc)
63  {
64  //- print input parameters
65  vout.general(m_vl, "%s:\n", class_name.c_str());
66  vout.general(m_vl, " kappa_s = %12.8f\n", kappa_s);
67  vout.general(m_vl, " kappa_t = %12.8f\n", kappa_t);
68  vout.general(m_vl, " cSW_s = %12.8f\n", cSW_s);
69  vout.general(m_vl, " cSW_t = %12.8f\n", cSW_t);
70  for (int mu = 0; mu < m_Ndim; ++mu) {
71  vout.general(m_vl, " boundary[%d] = %2d\n", mu, bc[mu]);
72  }
73 
74  //- range check
75  // NB. kappa,cSW == 0 is allowed.
76  assert(bc.size() == m_Ndim);
77 
78  //- store values
79  m_kappa_s = kappa_s;
80  m_kappa_t = kappa_t;
81  m_cSW_s = cSW_s;
82  m_cSW_t = cSW_t;
83 
84  // m_boundary.resize(m_Ndim); // already resized in init.
85  for (int mu = 0; mu < m_Ndim; ++mu) {
86  m_boundary[mu] = bc[mu];
87  }
88  }
89 
90 
91 //====================================================================
93  {
94  m_U = (Field_G *)U;
95  set_csw();
96  }
97 
98 
99 //====================================================================
100  void Fopr_CloverTerm_General::init(std::string repr)
101  {
106  m_NinF = 2 * m_Nc * m_Nd;
107 
108  m_U = 0;
109 
110  m_repr = repr;
111 
112  m_boundary.resize(m_Ndim);
113  m_SG.resize(m_Ndim * m_Ndim);
114 
115  unique_ptr<GammaMatrixSet> gmset(GammaMatrixSet::New(m_repr));
116 
117  m_GM5 = gmset->get_GM(gmset->GAMMA5);
118 
119  m_SG[sg_index(0, 1)] = gmset->get_GM(gmset->SIGMA12);
120  m_SG[sg_index(1, 2)] = gmset->get_GM(gmset->SIGMA23);
121  m_SG[sg_index(2, 0)] = gmset->get_GM(gmset->SIGMA31);
122  m_SG[sg_index(3, 0)] = gmset->get_GM(gmset->SIGMA41);
123  m_SG[sg_index(3, 1)] = gmset->get_GM(gmset->SIGMA42);
124  m_SG[sg_index(3, 2)] = gmset->get_GM(gmset->SIGMA43);
125 
126  m_SG[sg_index(1, 0)] = m_SG[sg_index(0, 1)].mult(-1);
127  m_SG[sg_index(2, 1)] = m_SG[sg_index(1, 2)].mult(-1);
128  m_SG[sg_index(0, 2)] = m_SG[sg_index(2, 0)].mult(-1);
129  m_SG[sg_index(0, 3)] = m_SG[sg_index(3, 0)].mult(-1);
130  m_SG[sg_index(1, 3)] = m_SG[sg_index(3, 1)].mult(-1);
131  m_SG[sg_index(2, 3)] = m_SG[sg_index(3, 2)].mult(-1);
132 
133  m_SG[sg_index(0, 0)] = gmset->get_GM(gmset->UNITY);
134  m_SG[sg_index(1, 1)] = gmset->get_GM(gmset->UNITY);
135  m_SG[sg_index(2, 2)] = gmset->get_GM(gmset->UNITY);
136  m_SG[sg_index(3, 3)] = gmset->get_GM(gmset->UNITY);
137  // these 4 gamma matrices are actually not used.
138 
139  // m_fopr_w = new Fopr_Wilson(repr);
140  if (m_repr == "Dirac") {
143  } else if (m_repr == "Chiral") {
146  }
147 
148  // m_w2.reset(m_NinF,m_Nvol,1);
149  }
150 
151 
152 //====================================================================
154  {
155  // nothing to do.
156  }
157 
158 
159 //====================================================================
161  {
162  (this->*m_gm5)(v, f);
163  }
164 
165 
166 //====================================================================
168  {
169  const int Nvc = 2 * CommonParameters::Nc();
170  const int Nd = CommonParameters::Nd();
171 
172  const double *v1 = f.ptr(0);
173  double *v2 = w.ptr(0);
174 
175  int id1 = 0;
176  int id2 = Nvc;
177  int id3 = Nvc * 2;
178  int id4 = Nvc * 3;
179 
180  // threadding applied.
182  int i_thread = ThreadManager_OpenMP::get_thread_id();
183 
184  int is = m_Nvol * i_thread / Nthread;
185  int ns = m_Nvol * (i_thread + 1) / Nthread - is;
186 
187  for (int site = is; site < is + ns; ++site) {
188  // for (int site = 0; site < m_Nvol; ++site) {
189  for (int icc = 0; icc < Nvc; icc++) {
190  int in = Nvc * Nd * site;
191  v2[icc + id1 + in] = v1[icc + id3 + in];
192  v2[icc + id2 + in] = v1[icc + id4 + in];
193  v2[icc + id3 + in] = v1[icc + id1 + in];
194  v2[icc + id4 + in] = v1[icc + id2 + in];
195  }
196  }
197  }
198 
199 
200 //====================================================================
202  {
203  const int Nvc = 2 * CommonParameters::Nc();
204  const int Nd = CommonParameters::Nd();
205 
206  const double *v1 = f.ptr(0);
207  double *v2 = w.ptr(0);
208 
209  int id1 = 0;
210  int id2 = Nvc;
211  int id3 = Nvc * 2;
212  int id4 = Nvc * 3;
213 
214  // threadding applied.
216  int i_thread = ThreadManager_OpenMP::get_thread_id();
217 
218  int is = m_Nvol * i_thread / Nthread;
219  int ns = m_Nvol * (i_thread + 1) / Nthread - is;
220 
221  for (int site = is; site < is + ns; ++site) {
222  // for (int site = 0; site < m_Nvol; ++site) {
223  for (int icc = 0; icc < Nvc; icc++) {
224  int in = Nvc * Nd * site;
225  v2[icc + id1 + in] = v1[icc + id1 + in];
226  v2[icc + id2 + in] = v1[icc + id2 + in];
227  v2[icc + id3 + in] = -v1[icc + id3 + in];
228  v2[icc + id4 + in] = -v1[icc + id4 + in];
229  }
230  }
231  }
232 
233 
234 //====================================================================
236  const int mu, const int nu)
237  {
238  assert(mu != nu);
239  mult_iGM(v, m_SG[sg_index(mu, nu)], w);
240  }
241 
242 
243 //====================================================================
245  {
246  mult_csw(v, f);
247  }
248 
249 
250 //====================================================================
252  {
253  (this->*m_csw)(v, w);
254  }
255 
256 
257 //====================================================================
259  {
260  assert(w.nex() == 1);
261 
262  const int Nc = CommonParameters::Nc();
263  const int Nvc = 2 * Nc;
264  const int Ndf = 2 * Nc * Nc;
265  const int Nd = CommonParameters::Nd();
266  const int Nvol = w.nvol();
267 
268  const int id1 = 0;
269  const int id2 = Nvc;
270  const int id3 = Nvc * 2;
271  const int id4 = Nvc * 3;
272 
273  const double kappa_cSW_s = m_kappa_s * m_cSW_s;
274  const double kappa_cSW_t = m_kappa_t * m_cSW_t;
275 
276  const double *w2 = w.ptr(0);
277  double *v2 = v.ptr(0);
278 
279  double *Bx = m_Bx.ptr(0);
280  double *By = m_By.ptr(0);
281  double *Bz = m_Bz.ptr(0);
282  double *Ex = m_Ex.ptr(0);
283  double *Ey = m_Ey.ptr(0);
284  double *Ez = m_Ez.ptr(0);
285 
286  // threadding applied.
288  int i_thread = ThreadManager_OpenMP::get_thread_id();
289 
290  int is = m_Nvol * i_thread / Nthread;
291  int ns = m_Nvol * (i_thread + 1) / Nthread - is;
292 
293  for (int site = is; site < is + ns; ++site) {
294  int iv = Nvc * Nd * site;
295  int ig = Ndf * site;
296 
297  for (int ic = 0; ic < Nc; ++ic) {
298  int icr = 2 * ic;
299  int ici = icr + 1;
300  int icg = ic * Nvc + ig;
301 
302  v2[icr + id1 + iv] = 0.0;
303  v2[ici + id1 + iv] = 0.0;
304  v2[icr + id2 + iv] = 0.0;
305  v2[ici + id2 + iv] = 0.0;
306 
307  v2[icr + id3 + iv] = 0.0;
308  v2[ici + id3 + iv] = 0.0;
309  v2[icr + id4 + iv] = 0.0;
310  v2[ici + id4 + iv] = 0.0;
311 
312  // isigma_23 * Bx
313  v2[icr + id1 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id2 + iv], Nc);
314  v2[ici + id1 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id2 + iv], Nc);
315  v2[icr + id2 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id1 + iv], Nc);
316  v2[ici + id2 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id1 + iv], Nc);
317 
318  v2[icr + id3 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id4 + iv], Nc);
319  v2[ici + id3 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id4 + iv], Nc);
320  v2[icr + id4 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id3 + iv], Nc);
321  v2[ici + id4 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id3 + iv], Nc);
322 
323  // isigma_31 * By
324  v2[icr + id1 + iv] += kappa_cSW_s * mult_uv_r(&By[icg], &w2[id2 + iv], Nc);
325  v2[ici + id1 + iv] += kappa_cSW_s * mult_uv_i(&By[icg], &w2[id2 + iv], Nc);
326  v2[icr + id2 + iv] -= kappa_cSW_s * mult_uv_r(&By[icg], &w2[id1 + iv], Nc);
327  v2[ici + id2 + iv] -= kappa_cSW_s * mult_uv_i(&By[icg], &w2[id1 + iv], Nc);
328 
329  v2[icr + id3 + iv] += kappa_cSW_s * mult_uv_r(&By[icg], &w2[id4 + iv], Nc);
330  v2[ici + id3 + iv] += kappa_cSW_s * mult_uv_i(&By[icg], &w2[id4 + iv], Nc);
331  v2[icr + id4 + iv] -= kappa_cSW_s * mult_uv_r(&By[icg], &w2[id3 + iv], Nc);
332  v2[ici + id4 + iv] -= kappa_cSW_s * mult_uv_i(&By[icg], &w2[id3 + iv], Nc);
333 
334  // isigma_12 * Bz
335  v2[icr + id1 + iv] -= kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id1 + iv], Nc);
336  v2[ici + id1 + iv] += kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id1 + iv], Nc);
337  v2[icr + id2 + iv] += kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id2 + iv], Nc);
338  v2[ici + id2 + iv] -= kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id2 + iv], Nc);
339 
340  v2[icr + id3 + iv] -= kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id3 + iv], Nc);
341  v2[ici + id3 + iv] += kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id3 + iv], Nc);
342  v2[icr + id4 + iv] += kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id4 + iv], Nc);
343  v2[ici + id4 + iv] -= kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id4 + iv], Nc);
344 
345  // isigma_41 * Ex
346  v2[icr + id1 + iv] += kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id2 + iv], Nc);
347  v2[ici + id1 + iv] -= kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id2 + iv], Nc);
348  v2[icr + id2 + iv] += kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id1 + iv], Nc);
349  v2[ici + id2 + iv] -= kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id1 + iv], Nc);
350 
351  v2[icr + id3 + iv] -= kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id4 + iv], Nc);
352  v2[ici + id3 + iv] += kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id4 + iv], Nc);
353  v2[icr + id4 + iv] -= kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id3 + iv], Nc);
354  v2[ici + id4 + iv] += kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id3 + iv], Nc);
355 
356  // isigma_42 * Ey
357  v2[icr + id1 + iv] -= kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id2 + iv], Nc);
358  v2[ici + id1 + iv] -= kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id2 + iv], Nc);
359  v2[icr + id2 + iv] += kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id1 + iv], Nc);
360  v2[ici + id2 + iv] += kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id1 + iv], Nc);
361 
362  v2[icr + id3 + iv] += kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id4 + iv], Nc);
363  v2[ici + id3 + iv] += kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id4 + iv], Nc);
364  v2[icr + id4 + iv] -= kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id3 + iv], Nc);
365  v2[ici + id4 + iv] -= kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id3 + iv], Nc);
366 
367  // isigma_43 * Ez
368  v2[icr + id1 + iv] += kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id1 + iv], Nc);
369  v2[ici + id1 + iv] -= kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id1 + iv], Nc);
370  v2[icr + id2 + iv] -= kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id2 + iv], Nc);
371  v2[ici + id2 + iv] += kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id2 + iv], Nc);
372 
373  v2[icr + id3 + iv] -= kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id3 + iv], Nc);
374  v2[ici + id3 + iv] += kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id3 + iv], Nc);
375  v2[icr + id4 + iv] += kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id4 + iv], Nc);
376  v2[ici + id4 + iv] -= kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id4 + iv], Nc);
377  }
378  }
379 #pragma omp barrier
380  }
381 
382 
383 //====================================================================
385  {
386  assert(w.nex() == 1);
387 
388  const int Nc = CommonParameters::Nc();
389  const int Nvc = 2 * Nc;
390  const int Ndf = 2 * Nc * Nc;
391  const int Nd = CommonParameters::Nd();
392  const int Nvol = w.nvol();
393 
394  const int id1 = 0;
395  const int id2 = Nvc;
396  const int id3 = Nvc * 2;
397  const int id4 = Nvc * 3;
398 
399  const double kappa_cSW_s = m_kappa_s * m_cSW_s;
400  const double kappa_cSW_t = m_kappa_t * m_cSW_t;
401 
402  const double *w2 = w.ptr(0);
403  double *v2 = v.ptr(0);
404 
405  double *Bx = m_Bx.ptr(0);
406  double *By = m_By.ptr(0);
407  double *Bz = m_Bz.ptr(0);
408  double *Ex = m_Ex.ptr(0);
409  double *Ey = m_Ey.ptr(0);
410  double *Ez = m_Ez.ptr(0);
411 
412  // threadding applied.
414  int i_thread = ThreadManager_OpenMP::get_thread_id();
415 
416  int is = m_Nvol * i_thread / Nthread;
417  int ns = m_Nvol * (i_thread + 1) / Nthread - is;
418 
419  for (int site = is; site < is + ns; ++site) {
420  int iv = Nvc * Nd * site;
421  int ig = Ndf * site;
422 
423  for (int ic = 0; ic < Nc; ++ic) {
424  int icr = 2 * ic;
425  int ici = icr + 1;
426  int icg = ic * Nvc + ig;
427 
428  v2[icr + id1 + iv] = 0.0;
429  v2[ici + id1 + iv] = 0.0;
430  v2[icr + id2 + iv] = 0.0;
431  v2[ici + id2 + iv] = 0.0;
432 
433  v2[icr + id3 + iv] = 0.0;
434  v2[ici + id3 + iv] = 0.0;
435  v2[icr + id4 + iv] = 0.0;
436  v2[ici + id4 + iv] = 0.0;
437 
438  // isigma_23 * Bx
439  v2[icr + id1 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id2 + iv], Nc);
440  v2[ici + id1 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id2 + iv], Nc);
441  v2[icr + id2 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id1 + iv], Nc);
442  v2[ici + id2 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id1 + iv], Nc);
443 
444  v2[icr + id3 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id4 + iv], Nc);
445  v2[ici + id3 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id4 + iv], Nc);
446  v2[icr + id4 + iv] -= kappa_cSW_s * mult_uv_i(&Bx[icg], &w2[id3 + iv], Nc);
447  v2[ici + id4 + iv] += kappa_cSW_s * mult_uv_r(&Bx[icg], &w2[id3 + iv], Nc);
448 
449  // isigma_31 * By
450  v2[icr + id1 + iv] += kappa_cSW_s * mult_uv_r(&By[icg], &w2[id2 + iv], Nc);
451  v2[ici + id1 + iv] += kappa_cSW_s * mult_uv_i(&By[icg], &w2[id2 + iv], Nc);
452  v2[icr + id2 + iv] -= kappa_cSW_s * mult_uv_r(&By[icg], &w2[id1 + iv], Nc);
453  v2[ici + id2 + iv] -= kappa_cSW_s * mult_uv_i(&By[icg], &w2[id1 + iv], Nc);
454 
455  v2[icr + id3 + iv] += kappa_cSW_s * mult_uv_r(&By[icg], &w2[id4 + iv], Nc);
456  v2[ici + id3 + iv] += kappa_cSW_s * mult_uv_i(&By[icg], &w2[id4 + iv], Nc);
457  v2[icr + id4 + iv] -= kappa_cSW_s * mult_uv_r(&By[icg], &w2[id3 + iv], Nc);
458  v2[ici + id4 + iv] -= kappa_cSW_s * mult_uv_i(&By[icg], &w2[id3 + iv], Nc);
459 
460  // isigma_12 * Bz
461  v2[icr + id1 + iv] -= kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id1 + iv], Nc);
462  v2[ici + id1 + iv] += kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id1 + iv], Nc);
463  v2[icr + id2 + iv] += kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id2 + iv], Nc);
464  v2[ici + id2 + iv] -= kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id2 + iv], Nc);
465 
466  v2[icr + id3 + iv] -= kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id3 + iv], Nc);
467  v2[ici + id3 + iv] += kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id3 + iv], Nc);
468  v2[icr + id4 + iv] += kappa_cSW_s * mult_uv_i(&Bz[icg], &w2[id4 + iv], Nc);
469  v2[ici + id4 + iv] -= kappa_cSW_s * mult_uv_r(&Bz[icg], &w2[id4 + iv], Nc);
470 
471  // isigma_41 * Ex
472  v2[icr + id1 + iv] += kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id4 + iv], Nc);
473  v2[ici + id1 + iv] -= kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id4 + iv], Nc);
474  v2[icr + id2 + iv] += kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id3 + iv], Nc);
475  v2[ici + id2 + iv] -= kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id3 + iv], Nc);
476 
477  v2[icr + id3 + iv] += kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id2 + iv], Nc);
478  v2[ici + id3 + iv] -= kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id2 + iv], Nc);
479  v2[icr + id4 + iv] += kappa_cSW_t * mult_uv_i(&Ex[icg], &w2[id1 + iv], Nc);
480  v2[ici + id4 + iv] -= kappa_cSW_t * mult_uv_r(&Ex[icg], &w2[id1 + iv], Nc);
481 
482  // isigma_42 * Ey
483  v2[icr + id1 + iv] -= kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id4 + iv], Nc);
484  v2[ici + id1 + iv] -= kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id4 + iv], Nc);
485  v2[icr + id2 + iv] += kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id3 + iv], Nc);
486  v2[ici + id2 + iv] += kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id3 + iv], Nc);
487 
488  v2[icr + id3 + iv] -= kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id2 + iv], Nc);
489  v2[ici + id3 + iv] -= kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id2 + iv], Nc);
490  v2[icr + id4 + iv] += kappa_cSW_t * mult_uv_r(&Ey[icg], &w2[id1 + iv], Nc);
491  v2[ici + id4 + iv] += kappa_cSW_t * mult_uv_i(&Ey[icg], &w2[id1 + iv], Nc);
492 
493  // isigma_43 * Ez
494  v2[icr + id1 + iv] += kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id3 + iv], Nc);
495  v2[ici + id1 + iv] -= kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id3 + iv], Nc);
496  v2[icr + id2 + iv] -= kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id4 + iv], Nc);
497  v2[ici + id2 + iv] += kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id4 + iv], Nc);
498 
499  v2[icr + id3 + iv] += kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id1 + iv], Nc);
500  v2[ici + id3 + iv] -= kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id1 + iv], Nc);
501  v2[icr + id4 + iv] -= kappa_cSW_t * mult_uv_i(&Ez[icg], &w2[id2 + iv], Nc);
502  v2[ici + id4 + iv] += kappa_cSW_t * mult_uv_r(&Ez[icg], &w2[id2 + iv], Nc);
503  }
504  }
505 #pragma omp barrier
506  }
507 
508 
509 //====================================================================
511  {
512  set_fieldstrength(m_Bx, 1, 2);
513  set_fieldstrength(m_By, 2, 0);
514  set_fieldstrength(m_Bz, 0, 1);
515  set_fieldstrength(m_Ex, 3, 0);
516  set_fieldstrength(m_Ey, 3, 1);
517  set_fieldstrength(m_Ez, 3, 2);
518  }
519 
520 
521 //====================================================================
523  const int mu, const int nu)
524  {
526 
527  assert(Nthread == 1);
528  // this function must be called in single thread region.
529 
530  m_staple.upper(m_Cup, *m_U, mu, nu); // these staple constructions
531  m_staple.lower(m_Cdn, *m_U, mu, nu); // are multi-threaded.
532 
533  //#pragma omp parallel
534  {
535  mult_Field_Gnd(Fst, 0, *m_U, mu, m_Cup, 0);
536  multadd_Field_Gnd(Fst, 0, *m_U, mu, m_Cdn, 0, -1.0);
537 
538  mult_Field_Gdn(m_v1, 0, m_Cup, 0, *m_U, mu);
539  multadd_Field_Gdn(m_v1, 0, m_Cdn, 0, *m_U, mu, -1.0);
540  }
541 
542  m_shift.forward(m_v2, m_v1, mu);
543 
544  //#pragma omp parallel
545  {
546  axpy(Fst, 1.0, m_v2);
547  ah_Field_G(Fst, 0);
548  scal(Fst, 0.25);
549  }
550  }
551 
552 
553 //====================================================================
555  {
556  // The following counting explicitly depends on the implementation
557  // and to be recalculated when the code is modified.
558  // Present counting is based on rev.1107. [24 Aug 2014 H.Matsufuru]
559 
560  const int Lvol = CommonParameters::Lvol();
561 
562  double flop_site
563  = static_cast<double>(m_Nc * m_Nd * (2 + 48 * m_Nc));
564  double flop = flop_site * static_cast<double>(Lvol);
565 
566  return flop;
567  }
568 
569 
570 //====================================================================
571 }
572 //============================================================END=====
void scal(Field &x, const double a)
scal(x, a): x = a * x
Definition: field.cpp:282
BridgeIO vout
Definition: bridgeIO.cpp:495
void ah_Field_G(Field_G &W, const int ex)
static int get_num_threads()
returns available number of threads.
const double * ptr(const int jin, const int site, const int jex) const
Definition: field.h:142
void mult_isigma(Field_F &, const Field_F &, const int mu, const int nu)
void general(const char *format,...)
Definition: bridgeIO.cpp:195
GammaMatrix get_GM(GMspecies spec)
void set_parameters(const Parameters &params)
void multadd_Field_Gnd(Field_G &W, const int ex, const Field_G &U1, const int ex1, const Field_G &U2, const int ex2, const double ff)
Container of Field-type object.
Definition: field.h:39
int fetch_double(const string &key, double &value) const
Definition: parameters.cpp:211
int nvol() const
Definition: field.h:116
void mult_Field_Gdn(Field_G &W, const int ex, const Field_G &U1, const int ex1, const Field_G &U2, const int ex2)
void set_fieldstrength(Field_G &, const int, const int)
Class for parameters.
Definition: parameters.h:46
static int Lvol()
void lower(Field_G &, const Field_G &, const int mu, const int nu)
constructs lower staple in mu-nu plane.
Definition: staple_lex.cpp:180
static int get_thread_id()
returns thread id.
Wilson-type fermion field.
Definition: field_F.h:37
void gm5_dirac(Field &, const Field &)
const Field_G * m_U
pointer to gauge configuration.
SU(N) gauge field.
Definition: field_G.h:38
void mult_csw_dirac(Field &, const Field &)
void gm5_chiral(Field &, const Field &)
void mult_iGM(Field_F &y, const GammaMatrix &gm, const Field_F &x)
gamma matrix multiplication (i is multiplied)
Bridge::VerboseLevel m_vl
Definition: fopr.h:128
int nex() const
Definition: field.h:117
void upper(Field_G &, const Field_G &, const int mu, const int nu)
constructs upper staple in mu-nu plane.
Definition: staple_lex.cpp:156
void set_config(Field *U)
setting pointer to the gauge configuration.
void(Fopr_CloverTerm_General::* m_gm5)(Field &, const Field &)
std::vector< GammaMatrix > m_SG
void mult_csw(Field &, const Field &)
void axpy(Field &y, const double a, const Field &x)
axpy(y, a, x): y := a * x + y
Definition: field.cpp:168
void crucial(const char *format,...)
Definition: bridgeIO.cpp:178
void mult_sigmaF(Field &, const Field &)
Field_G m_v2
for calculation of field strength.
void mult_gm5(Field &v, const Field &w)
gamma_5 multiplication. [31 Mar 2017 H.Matsufuru]
void mult_csw_chiral(Field &, const Field &)
double flop_count()
this returns the number of floating point operations.
string get_string(const string &key) const
Definition: parameters.cpp:116
int fetch_int_vector(const string &key, vector< int > &value) const
Definition: parameters.cpp:294
static const std::string class_name
static VerboseLevel set_verbose_level(const std::string &str)
Definition: bridgeIO.cpp:131
void multadd_Field_Gdn(Field_G &W, const int ex, const Field_G &U1, const int ex1, const Field_G &U2, const int ex2, const double ff)
void(Fopr_CloverTerm_General::* m_csw)(Field &, const Field &)
void forward(Field &, const Field &, const int mu)
void mult_Field_Gnd(Field_G &W, const int ex, const Field_G &U1, const int ex1, const Field_G &U2, const int ex2)