Bridge++  Ver.2.1.3
afield_Gauge-inc.h
Go to the documentation of this file.
1 
10 #ifndef QXS_AFIELD_GAUGE_INC_INCLUDED
11 #define QXS_AFIELD_GAUGE_INC_INCLUDED
12 
13 #include <cstdlib>
14 
15 #include "lib/ResourceManager/threadManager.h"
16 
17 #include "lib_alt_QXS/Field/afield.h"
18 #include "lib_alt_QXS/Field/aindex_lex.h"
19 #include "lib_alt_QXS/inline/afield_th-inc.h"
20 
21 #include "lib_alt_QXS/Field/afield_Gauge_index-inc.h"
22 
23 //namespace {
24 // inline int idxGr(int ic1, int ic2){ return 2*(ic1 + NC * ic2); }
25 // inline int idxGi(int ic1, int ic2){ return 1 + 2*(ic1 + NC * ic2); }
26 //}
27 
28 namespace QXS_Gauge{
29 
30 //====================================================================
31 template <typename REALTYPE>
32 void set_boundary(AField<REALTYPE, QXS>& ulex,
33  const std::vector<int>& boundary)
34 {
35  typedef REALTYPE real_t;
36 
37 #pragma omp barrier
38 
39  int Nx = CommonParameters::Nx();
40  int Ny = CommonParameters::Ny();
41  int Nz = CommonParameters::Nz();
42  int Nt = CommonParameters::Nt();
43  int Ndim = CommonParameters::Ndim();
44  int Nvol = Nx * Ny * Nz * Nt;
45 
46  int Nin = ulex.nin();
47 
48  if(!ulex.check_size(Nin, Nvol, Ndim)){
49  vout.crucial("set_boundary: wrong size of input field\n");
50  exit(EXIT_FAILURE);
51  }
52 
53  AIndex_lex<real_t, QXS> index;
54 
55  int mu = 0;
56  int ipex = Communicator::ipe(mu);
57  int npex = Communicator::npe(mu);
58 
59  if(boundary[mu] != 1 && ipex == npex-1){
60 
61  real_t bc = real_t(boundary[mu]);
62  int Nyzt = Ny * Nz * Nt;
63 
64  int ith, nth, is, ns;
65  set_threadtask(ith, nth, is, ns, Nyzt);
66 
67  for(int iyzt = is; iyzt < ns; ++iyzt){
68  int site = Nx-1 + Nx * iyzt;
69  for(int in = 0; in < Nin; ++in){
70  int idx = index.idx_G(in, site, mu);
71  real_t uv = ulex.cmp(idx);
72  uv = uv * bc;
73  ulex.set(idx, uv);
74  }
75  }
76  }
77 
78  mu = 1;
79  int ipey = Communicator::ipe(mu);
80  int npey = Communicator::npe(mu);
81 
82  if(boundary[mu] != 1 && ipey == npey-1){
83 
84  real_t bc = real_t(boundary[mu]);
85  int Nxzt = Nx * Nz * Nt;
86 
87  int ith, nth, is, ns;
88  set_threadtask(ith, nth, is, ns, Nxzt);
89 
90  for(int ixzt = is; ixzt < ns; ++ixzt){
91  int ix = ixzt % Nx;
92  int izt = ixzt / Nx;
93  int site = ix + Nx * (Ny-1 + Ny * izt);
94  for(int in = 0; in < Nin; ++in){
95  int idx = index.idx_G(in, site, mu);
96  real_t uv = ulex.cmp(idx);
97  uv = uv * bc;
98  ulex.set(idx, uv);
99  }
100  }
101  }
102 
103  mu = 2;
104  int ipez = Communicator::ipe(mu);
105  int npez = Communicator::npe(mu);
106 
107  if(boundary[mu] != 1 && ipez == npez-1){
108 
109  real_t bc = real_t(boundary[mu]);
110  int Nxy = Nx * Ny;
111  int Nxyt = Nxy * Nt;
112 
113  int ith, nth, is, ns;
114  set_threadtask(ith, nth, is, ns, Nxyt);
115 
116  for(int ixyt = is; ixyt < ns; ++ixyt){
117  int ixy = ixyt % Nxy;
118  int it = ixyt / Nxy;
119  int site = ixy + Nxy * (Nz-1 + Nz * it);
120  for(int in = 0; in < Nin; ++in){
121  int idx = index.idx_G(in, site, mu);
122  real_t uv = ulex.cmp(idx);
123  uv = uv * bc;
124  ulex.set(idx, uv);
125  }
126  }
127  }
128 
129  mu = 3;
130  int ipet = Communicator::ipe(mu);
131  int npet = Communicator::npe(mu);
132 
133  if(boundary[mu] != 1 && ipet == npet-1){
134 
135  real_t bc = real_t(boundary[mu]);
136  int Nxyz = Nx * Ny * Nz;
137 
138  int ith, nth, is, ns;
139  set_threadtask(ith, nth, is, ns, Nxyz);
140 
141  for(int ixyz = is; ixyz < ns; ++ixyz){
142  int site = ixyz + Nxyz * (Nt-1);
143  for(int in = 0; in < Nin; ++in){
144  int idx = index.idx_G(in, site, mu);
145  real_t uv = ulex.cmp(idx);
146  uv = uv * bc;
147  ulex.set(idx, uv);
148  }
149  }
150  }
151 
152 #pragma omp barrier
153 
154 }
155 
156 //====================================================================
157 template <typename REALTYPE>
158 void mult_Gnn(AField<REALTYPE,QXS>& u, const int exu,
159  const AField<REALTYPE,QXS>& v, const int exv,
160  const AField<REALTYPE,QXS>& w, const int exw)
161 {
162 #pragma omp barrier
163 
164  int Nst = w.nvol();
165  int Nstv = Nst/VLEN;
166 
167  typedef AField<REALTYPE,QXS> AFIELD;
168 
169  REALTYPE* up = u.ptr(NDF * Nst * exu);
170  REALTYPE* vp = const_cast<AFIELD*>(&v)->ptr(NDF * Nst * exv);
171  REALTYPE* wp = const_cast<AFIELD*>(&w)->ptr(NDF * Nst * exw);
172 
173  int ith, nth, is, ns;
174  set_threadtask(ith, nth, is, ns, Nstv);
175 
176  svbool_t pg = set_predicate();
177 
178  for(int site = is; site < ns; ++site){
179  REALTYPE* upt = &up[VLEN * NDF * site];
180  REALTYPE* vpt = &vp[VLEN * NDF * site];
181  REALTYPE* wpt = &wp[VLEN * NDF * site];
182 
183  svreal_t w00r, w00i, w01r, w01i, w02r, w02i;
184  load_vec(pg, w00r, &wpt[VLEN * idxGr(0,0)]);
185  load_vec(pg, w00i, &wpt[VLEN * idxGi(0,0)]);
186  load_vec(pg, w01r, &wpt[VLEN * idxGr(0,1)]);
187  load_vec(pg, w01i, &wpt[VLEN * idxGi(0,1)]);
188  load_vec(pg, w02r, &wpt[VLEN * idxGr(0,2)]);
189  load_vec(pg, w02i, &wpt[VLEN * idxGi(0,2)]);
190 
191  svreal_t w10r, w10i, w11r, w11i, w12r, w12i;
192  load_vec(pg, w10r, &wpt[VLEN * idxGr(1,0)]);
193  load_vec(pg, w10i, &wpt[VLEN * idxGi(1,0)]);
194  load_vec(pg, w11r, &wpt[VLEN * idxGr(1,1)]);
195  load_vec(pg, w11i, &wpt[VLEN * idxGi(1,1)]);
196  load_vec(pg, w12r, &wpt[VLEN * idxGr(1,2)]);
197  load_vec(pg, w12i, &wpt[VLEN * idxGi(1,2)]);
198 
199  svreal_t w20r, w20i, w21r, w21i, w22r, w22i;
200  load_vec(pg, w20r, &wpt[VLEN * idxGr(2,0)]);
201  load_vec(pg, w20i, &wpt[VLEN * idxGi(2,0)]);
202  load_vec(pg, w21r, &wpt[VLEN * idxGr(2,1)]);
203  load_vec(pg, w21i, &wpt[VLEN * idxGi(2,1)]);
204  load_vec(pg, w22r, &wpt[VLEN * idxGr(2,2)]);
205  load_vec(pg, w22i, &wpt[VLEN * idxGi(2,2)]);
206 
207  for(int ic1 = 0; ic1 < NC; ++ic1){
208 
209  svreal_t v0r, v0i, v1r, v1i, v2r, v2i;
210  load_vec(pg, v0r, &vpt[VLEN * idxGr(ic1,0)]);
211  load_vec(pg, v0i, &vpt[VLEN * idxGi(ic1,0)]);
212  load_vec(pg, v1r, &vpt[VLEN * idxGr(ic1,1)]);
213  load_vec(pg, v1i, &vpt[VLEN * idxGi(ic1,1)]);
214  load_vec(pg, v2r, &vpt[VLEN * idxGr(ic1,2)]);
215  load_vec(pg, v2i, &vpt[VLEN * idxGi(ic1,2)]);
216 
217  svreal_t ut0r, ut0i, ut1r, ut1i, ut2r, ut2i;
218  mul_vec( pg, ut0r, v0r, w00r);
219  mul_vec( pg, ut0i, v0r, w00i);
220  mul_vec( pg, ut1r, v0r, w01r);
221  mul_vec( pg, ut1i, v0r, w01i);
222  mul_vec( pg, ut2r, v0r, w02r);
223  mul_vec( pg, ut2i, v0r, w02i);
224 
225  ymax_vec(pg, ut0r, v0i, w00i);
226  axpy_vec(pg, ut0i, v0i, w00r);
227  ymax_vec(pg, ut1r, v0i, w01i);
228  axpy_vec(pg, ut1i, v0i, w01r);
229  ymax_vec(pg, ut2r, v0i, w02i);
230  axpy_vec(pg, ut2i, v0i, w02r);
231 
232  axpy_vec(pg, ut0r, v1r, w10r);
233  axpy_vec(pg, ut0i, v1r, w10i);
234  axpy_vec(pg, ut1r, v1r, w11r);
235  axpy_vec(pg, ut1i, v1r, w11i);
236  axpy_vec(pg, ut2r, v1r, w12r);
237  axpy_vec(pg, ut2i, v1r, w12i);
238 
239  ymax_vec(pg, ut0r, v1i, w10i);
240  axpy_vec(pg, ut0i, v1i, w10r);
241  ymax_vec(pg, ut1r, v1i, w11i);
242  axpy_vec(pg, ut1i, v1i, w11r);
243  ymax_vec(pg, ut2r, v1i, w12i);
244  axpy_vec(pg, ut2i, v1i, w12r);
245 
246  axpy_vec(pg, ut0r, v2r, w20r);
247  axpy_vec(pg, ut0i, v2r, w20i);
248  axpy_vec(pg, ut1r, v2r, w21r);
249  axpy_vec(pg, ut1i, v2r, w21i);
250  axpy_vec(pg, ut2r, v2r, w22r);
251  axpy_vec(pg, ut2i, v2r, w22i);
252 
253  ymax_vec(pg, ut0r, v2i, w20i);
254  axpy_vec(pg, ut0i, v2i, w20r);
255  ymax_vec(pg, ut1r, v2i, w21i);
256  axpy_vec(pg, ut1i, v2i, w21r);
257  ymax_vec(pg, ut2r, v2i, w22i);
258  axpy_vec(pg, ut2i, v2i, w22r);
259 
260  save_vec(pg, &upt[VLEN * idxGr(ic1,0)], ut0r);
261  save_vec(pg, &upt[VLEN * idxGi(ic1,0)], ut0i);
262  save_vec(pg, &upt[VLEN * idxGr(ic1,1)], ut1r);
263  save_vec(pg, &upt[VLEN * idxGi(ic1,1)], ut1i);
264  save_vec(pg, &upt[VLEN * idxGr(ic1,2)], ut2r);
265  save_vec(pg, &upt[VLEN * idxGi(ic1,2)], ut2i);
266  }
267 
268  }
269 
270 #pragma omp barrier
271 }
272 
273 //====================================================================
274 template <typename REALTYPE>
275 void mult_Gnd(AField<REALTYPE,QXS>& u, const int exu,
276  const AField<REALTYPE,QXS>& v, const int exv,
277  const AField<REALTYPE,QXS>& w, const int exw)
278 {
279 #pragma omp barrier
280 
281  int Nst = w.nvol();
282  int Nstv = Nst/VLEN;
283 
284  typedef AField<REALTYPE,QXS> AFIELD;
285 
286  REALTYPE* up = u.ptr(NDF * Nst * exu);
287  REALTYPE* vp = const_cast<AFIELD*>(&v)->ptr(NDF * Nst * exv);
288  REALTYPE* wp = const_cast<AFIELD*>(&w)->ptr(NDF * Nst * exw);
289 
290  int ith, nth, is, ns;
291  set_threadtask(ith, nth, is, ns, Nstv);
292 
293  svbool_t pg = set_predicate();
294 
295  for(int site = is; site < ns; ++site){
296  REALTYPE* upt = &up[VLEN * NDF * site];
297  REALTYPE* vpt = &vp[VLEN * NDF * site];
298  REALTYPE* wpt = &wp[VLEN * NDF * site];
299 
300  svreal_t w00r, w00i, w01r, w01i, w02r, w02i;
301  load_vec(pg, w00r, &wpt[VLEN * idxGr(0,0)]);
302  load_vec(pg, w00i, &wpt[VLEN * idxGi(0,0)]);
303  load_vec(pg, w01r, &wpt[VLEN * idxGr(0,1)]);
304  load_vec(pg, w01i, &wpt[VLEN * idxGi(0,1)]);
305  load_vec(pg, w02r, &wpt[VLEN * idxGr(0,2)]);
306  load_vec(pg, w02i, &wpt[VLEN * idxGi(0,2)]);
307 
308  svreal_t w10r, w10i, w11r, w11i, w12r, w12i;
309  load_vec(pg, w10r, &wpt[VLEN * idxGr(1,0)]);
310  load_vec(pg, w10i, &wpt[VLEN * idxGi(1,0)]);
311  load_vec(pg, w11r, &wpt[VLEN * idxGr(1,1)]);
312  load_vec(pg, w11i, &wpt[VLEN * idxGi(1,1)]);
313  load_vec(pg, w12r, &wpt[VLEN * idxGr(1,2)]);
314  load_vec(pg, w12i, &wpt[VLEN * idxGi(1,2)]);
315 
316  svreal_t w20r, w20i, w21r, w21i, w22r, w22i;
317  load_vec(pg, w20r, &wpt[VLEN * idxGr(2,0)]);
318  load_vec(pg, w20i, &wpt[VLEN * idxGi(2,0)]);
319  load_vec(pg, w21r, &wpt[VLEN * idxGr(2,1)]);
320  load_vec(pg, w21i, &wpt[VLEN * idxGi(2,1)]);
321  load_vec(pg, w22r, &wpt[VLEN * idxGr(2,2)]);
322  load_vec(pg, w22i, &wpt[VLEN * idxGi(2,2)]);
323 
324  for(int ic1 = 0; ic1 < NC; ++ic1){
325 
326  svreal_t v0r, v0i, v1r, v1i, v2r, v2i;
327  load_vec(pg, v0r, &vpt[VLEN * idxGr(ic1,0)]);
328  load_vec(pg, v0i, &vpt[VLEN * idxGi(ic1,0)]);
329  load_vec(pg, v1r, &vpt[VLEN * idxGr(ic1,1)]);
330  load_vec(pg, v1i, &vpt[VLEN * idxGi(ic1,1)]);
331  load_vec(pg, v2r, &vpt[VLEN * idxGr(ic1,2)]);
332  load_vec(pg, v2i, &vpt[VLEN * idxGi(ic1,2)]);
333 
334  svreal_t ut0r, ut0i, ut1r, ut1i, ut2r, ut2i;
335  mul_vec( pg, ut0r, v0r, w00r);
336  mul_vec( pg, ut0i, v0i, w00r);
337  mul_vec( pg, ut1r, v0r, w10r);
338  mul_vec( pg, ut1i, v0i, w10r);
339  mul_vec( pg, ut2r, v0r, w20r);
340  mul_vec( pg, ut2i, v0i, w20r);
341 
342  axpy_vec(pg, ut0r, v0i, w00i);
343  ymax_vec(pg, ut0i, v0r, w00i);
344  axpy_vec(pg, ut1r, v0i, w10i);
345  ymax_vec(pg, ut1i, v0r, w10i);
346  axpy_vec(pg, ut2r, v0i, w20i);
347  ymax_vec(pg, ut2i, v0r, w20i);
348 
349  axpy_vec(pg, ut0r, v1r, w01r);
350  axpy_vec(pg, ut0i, v1i, w01r);
351  axpy_vec(pg, ut1r, v1r, w11r);
352  axpy_vec(pg, ut1i, v1i, w11r);
353  axpy_vec(pg, ut2r, v1r, w21r);
354  axpy_vec(pg, ut2i, v1i, w21r);
355 
356  axpy_vec(pg, ut0r, v1i, w01i);
357  ymax_vec(pg, ut0i, v1r, w01i);
358  axpy_vec(pg, ut1r, v1i, w11i);
359  ymax_vec(pg, ut1i, v1r, w11i);
360  axpy_vec(pg, ut2r, v1i, w21i);
361  ymax_vec(pg, ut2i, v1r, w21i);
362 
363  axpy_vec(pg, ut0r, v2r, w02r);
364  axpy_vec(pg, ut0i, v2i, w02r);
365  axpy_vec(pg, ut1r, v2r, w12r);
366  axpy_vec(pg, ut1i, v2i, w12r);
367  axpy_vec(pg, ut2r, v2r, w22r);
368  axpy_vec(pg, ut2i, v2i, w22r);
369 
370  axpy_vec(pg, ut0r, v2i, w02i);
371  ymax_vec(pg, ut0i, v2r, w02i);
372  axpy_vec(pg, ut1r, v2i, w12i);
373  ymax_vec(pg, ut1i, v2r, w12i);
374  axpy_vec(pg, ut2r, v2i, w22i);
375  ymax_vec(pg, ut2i, v2r, w22i);
376 
377  save_vec(pg, &upt[VLEN * idxGr(ic1,0)], ut0r);
378  save_vec(pg, &upt[VLEN * idxGi(ic1,0)], ut0i);
379  save_vec(pg, &upt[VLEN * idxGr(ic1,1)], ut1r);
380  save_vec(pg, &upt[VLEN * idxGi(ic1,1)], ut1i);
381  save_vec(pg, &upt[VLEN * idxGr(ic1,2)], ut2r);
382  save_vec(pg, &upt[VLEN * idxGi(ic1,2)], ut2i);
383  }
384 
385  }
386 
387 #pragma omp barrier
388 }
389 
390 //====================================================================
391 template <typename REALTYPE>
392 void mult_Gdn(AField<REALTYPE,QXS>& u, const int exu,
393  const AField<REALTYPE,QXS>& v, const int exv,
394  const AField<REALTYPE,QXS>& w, const int exw)
395 {
396 #pragma omp barrier
397 
398  int Nst = w.nvol();
399  int Nstv = Nst/VLEN;
400 
401  typedef AField<REALTYPE,QXS> AFIELD;
402 
403  REALTYPE* up = u.ptr(NDF * Nst * exu);
404  REALTYPE* vp = const_cast<AFIELD*>(&v)->ptr(NDF * Nst * exv);
405  REALTYPE* wp = const_cast<AFIELD*>(&w)->ptr(NDF * Nst * exw);
406 
407  int ith, nth, is, ns;
408  set_threadtask(ith, nth, is, ns, Nstv);
409 
410  svbool_t pg = set_predicate();
411 
412  for(int site = is; site < ns; ++site){
413  REALTYPE* upt = &up[VLEN * NDF * site];
414  REALTYPE* vpt = &vp[VLEN * NDF * site];
415  REALTYPE* wpt = &wp[VLEN * NDF * site];
416 
417  svreal_t w00r, w00i, w01r, w01i, w02r, w02i;
418  load_vec(pg, w00r, &wpt[VLEN * idxGr(0,0)]);
419  load_vec(pg, w00i, &wpt[VLEN * idxGi(0,0)]);
420  load_vec(pg, w01r, &wpt[VLEN * idxGr(0,1)]);
421  load_vec(pg, w01i, &wpt[VLEN * idxGi(0,1)]);
422  load_vec(pg, w02r, &wpt[VLEN * idxGr(0,2)]);
423  load_vec(pg, w02i, &wpt[VLEN * idxGi(0,2)]);
424 
425  svreal_t w10r, w10i, w11r, w11i, w12r, w12i;
426  load_vec(pg, w10r, &wpt[VLEN * idxGr(1,0)]);
427  load_vec(pg, w10i, &wpt[VLEN * idxGi(1,0)]);
428  load_vec(pg, w11r, &wpt[VLEN * idxGr(1,1)]);
429  load_vec(pg, w11i, &wpt[VLEN * idxGi(1,1)]);
430  load_vec(pg, w12r, &wpt[VLEN * idxGr(1,2)]);
431  load_vec(pg, w12i, &wpt[VLEN * idxGi(1,2)]);
432 
433  svreal_t w20r, w20i, w21r, w21i, w22r, w22i;
434  load_vec(pg, w20r, &wpt[VLEN * idxGr(2,0)]);
435  load_vec(pg, w20i, &wpt[VLEN * idxGi(2,0)]);
436  load_vec(pg, w21r, &wpt[VLEN * idxGr(2,1)]);
437  load_vec(pg, w21i, &wpt[VLEN * idxGi(2,1)]);
438  load_vec(pg, w22r, &wpt[VLEN * idxGr(2,2)]);
439  load_vec(pg, w22i, &wpt[VLEN * idxGi(2,2)]);
440 
441  for(int ic1 = 0; ic1 < NC; ++ic1){
442 
443  svreal_t v0r, v0i, v1r, v1i, v2r, v2i;
444  load_vec(pg, v0r, &vpt[VLEN * idxGr(0,ic1)]);
445  load_vec(pg, v0i, &vpt[VLEN * idxGi(0,ic1)]);
446  load_vec(pg, v1r, &vpt[VLEN * idxGr(1,ic1)]);
447  load_vec(pg, v1i, &vpt[VLEN * idxGi(1,ic1)]);
448  load_vec(pg, v2r, &vpt[VLEN * idxGr(2,ic1)]);
449  load_vec(pg, v2i, &vpt[VLEN * idxGi(2,ic1)]);
450 
451  svreal_t ut0r, ut0i, ut1r, ut1i, ut2r, ut2i;
452  mul_vec( pg, ut0r, v0r, w00r);
453  mul_vec( pg, ut0i, v0r, w00i);
454  mul_vec( pg, ut1r, v0r, w01r);
455  mul_vec( pg, ut1i, v0r, w01i);
456  mul_vec( pg, ut2r, v0r, w02r);
457  mul_vec( pg, ut2i, v0r, w02i);
458 
459  axpy_vec(pg, ut0r, v0i, w00i);
460  ymax_vec(pg, ut0i, v0i, w00r);
461  axpy_vec(pg, ut1r, v0i, w01i);
462  ymax_vec(pg, ut1i, v0i, w01r);
463  axpy_vec(pg, ut2r, v0i, w02i);
464  ymax_vec(pg, ut2i, v0i, w02r);
465 
466  axpy_vec(pg, ut0r, v1r, w10r);
467  axpy_vec(pg, ut0i, v1r, w10i);
468  axpy_vec(pg, ut1r, v1r, w11r);
469  axpy_vec(pg, ut1i, v1r, w11i);
470  axpy_vec(pg, ut2r, v1r, w12r);
471  axpy_vec(pg, ut2i, v1r, w12i);
472 
473  axpy_vec(pg, ut0r, v1i, w10i);
474  ymax_vec(pg, ut0i, v1i, w10r);
475  axpy_vec(pg, ut1r, v1i, w11i);
476  ymax_vec(pg, ut1i, v1i, w11r);
477  axpy_vec(pg, ut2r, v1i, w12i);
478  ymax_vec(pg, ut2i, v1i, w12r);
479 
480  axpy_vec(pg, ut0r, v2r, w20r);
481  axpy_vec(pg, ut0i, v2r, w20i);
482  axpy_vec(pg, ut1r, v2r, w21r);
483  axpy_vec(pg, ut1i, v2r, w21i);
484  axpy_vec(pg, ut2r, v2r, w22r);
485  axpy_vec(pg, ut2i, v2r, w22i);
486 
487  axpy_vec(pg, ut0r, v2i, w20i);
488  ymax_vec(pg, ut0i, v2i, w20r);
489  axpy_vec(pg, ut1r, v2i, w21i);
490  ymax_vec(pg, ut1i, v2i, w21r);
491  axpy_vec(pg, ut2r, v2i, w22i);
492  ymax_vec(pg, ut2i, v2i, w22r);
493 
494  save_vec(pg, &upt[VLEN * idxGr(ic1,0)], ut0r);
495  save_vec(pg, &upt[VLEN * idxGi(ic1,0)], ut0i);
496  save_vec(pg, &upt[VLEN * idxGr(ic1,1)], ut1r);
497  save_vec(pg, &upt[VLEN * idxGi(ic1,1)], ut1i);
498  save_vec(pg, &upt[VLEN * idxGr(ic1,2)], ut2r);
499  save_vec(pg, &upt[VLEN * idxGi(ic1,2)], ut2i);
500  }
501 
502  }
503 
504 #pragma omp barrier
505 }
506 
507 //====================================================================
508 template <typename REALTYPE>
509 void at_G(AField<REALTYPE,QXS>& u, const int ex)
510 { // anti-Hermitian traceless
511 
512 #pragma omp barrier
513 
514  int Nst = u.nvol();
515  int Nstv = Nst/VLEN;
516 
517  typedef AField<REALTYPE,QXS> AFIELD;
518 
519  REALTYPE* up = u.ptr(NDF * Nst * ex);
520 
521  int ith, nth, is, ns;
522  set_threadtask(ith, nth, is, ns, Nstv);
523 
524  svbool_t pg = set_predicate();
525  svreal_t zero;
526  set_vec(pg, zero, 0.0);
527 
528  for(int site = is; site < ns; ++site){
529 
530  REALTYPE* upt = &up[VLEN * NDF * site];
531 
532  svreal_t tri;
533  set_vec(pg, tri, real_t(0.0));
534 
535  for(int ic1 = 0; ic1 < NC; ++ic1){
536  for(int ic2 = ic1; ic2 < NC; ++ic2){
537  if(ic1 != ic2){
538  svreal_t ur, ui, unr, uni, utr, uti;
539  load_vec(pg, unr, &upt[VLEN * idxGr(ic1,ic2)]);
540  load_vec(pg, uni, &upt[VLEN * idxGi(ic1,ic2)]);
541  load_vec(pg, utr, &upt[VLEN * idxGr(ic2,ic1)]);
542  load_vec(pg, uti, &upt[VLEN * idxGi(ic2,ic1)]);
543  sub_vec(pg, ur, unr, utr);
544  add_vec(pg, ui, uni, uti);
545  scal_vec(pg, ur, real_t(0.5));
546  scal_vec(pg, ui, real_t(0.5));
547  save_vec(pg, &upt[VLEN * idxGr(ic1,ic2)], ur);
548  save_vec(pg, &upt[VLEN * idxGi(ic1,ic2)], ui);
549  scal_vec(pg, ur, real_t(-1.0));
550  save_vec(pg, &upt[VLEN * idxGr(ic2,ic1)], ur);
551  save_vec(pg, &upt[VLEN * idxGi(ic2,ic1)], ui);
552  }else{
553  svreal_t uti;
554  load_vec(pg, uti, &upt[VLEN * idxGi(ic1,ic1)]);
555  add_vec(pg, tri, uti);
556  }
557  }
558  }
559 
560  svreal_t uti;
561  for(int ic1 = 0; ic1 < NC; ++ic1){
562  load_vec(pg, uti, &upt[VLEN * idxGi(ic1,ic1)]);
563  axpy_vec(pg, uti, -1.0/3.0, tri);
564  save_vec(pg, &upt[VLEN * idxGr(ic1,ic1)], zero);
565  save_vec(pg, &upt[VLEN * idxGi(ic1,ic1)], uti);
566  }
567 
568  }
569 
570 #pragma omp barrier
571 
572 }
573 
574 } // namespace QXS_Gauge
575 
576 //============================================================END=====
577 #endif
CommonParameters::Ny
static int Ny()
Definition: commonParameters.h:106
CommonParameters::Nz
static int Nz()
Definition: commonParameters.h:107
QXS_Gauge::at_G
void at_G(AField< REALTYPE, QXS > &u, const int ex)
Definition: afield_Gauge-inc.h:509
AField< REALTYPE, QXS >::ptr
real_t * ptr(int i)
Definition: afield.h:146
CommonParameters::Ndim
static int Ndim()
Definition: commonParameters.h:117
aindex_lex.h
AIndex_lex
Definition: aindex_lex_base.h:17
ixy
int ixy
Definition: mult_Domainwall_eo_xyz_openacc-inc.h:513
izt
int izt
Definition: mult_Domainwall_eo_xyz_openacc-inc.h:262
AField< REALTYPE, QXS >
Definition: afield.h:35
AField< REALTYPE, QXS >::nvol
int nvol() const
returning size of site d.o.f.
Definition: afield.h:112
VLEN
#define VLEN
Definition: bridgeQXS_Clover_coarse_double.cpp:12
AFIELD
Field AFIELD
Definition: eigensolver.cpp:18
NDF
#define NDF
Definition: field_F_imp_SU2-inc.h:17
Vsimd_t
Definition: vsimd_double-inc.h:13
afield.h
iyzt
int iyzt
Definition: mult_Domainwall_eo_xyz_openacc-inc.h:14
CommonParameters::Nx
static int Nx()
Definition: commonParameters.h:105
ix
int ix
Definition: mult_Wilson_xyz_openacc-inc.h:16
QXS_Gauge
Definition: afield_Gauge-inc.h:28
NC
#define NC
Definition: field_F_imp_SU2-inc.h:15
QXS_Gauge::mult_Gnd
void mult_Gnd(AField< REALTYPE, QXS > &u, const int exu, const AField< REALTYPE, QXS > &v, const int exv, const AField< REALTYPE, QXS > &w, const int exw)
Definition: afield_Gauge-inc.h:275
QXS_Gauge::set_boundary
void set_boundary(AField< REALTYPE, QXS > &ulex, const std::vector< int > &boundary)
Definition: afield_Gauge-inc.h:32
CommonParameters::Nt
static int Nt()
Definition: commonParameters.h:108
AIndex_eo_accel::idx
int idx(const int in, const int Nin, const int ist, const int leo, const int Nvol2, const int ex)
Definition: aindex_eo.h:28
Communicator::npe
static int npe(const int dir)
logical grid extent
Definition: communicator.cpp:112
threadManager.h
QXS_Gauge::mult_Gnn
void mult_Gnn(AField< REALTYPE, QXS > &u, const int exu, const AField< REALTYPE, QXS > &v, const int exv, const AField< REALTYPE, QXS > &w, const int exw)
Definition: afield_Gauge-inc.h:158
it
int it
Definition: mult_Wilson_xyz_openacc-inc.h:461
afield_Gauge_index-inc.h
real_t
double real_t
Definition: bridgeACC_AField_double.cpp:14
AField< REALTYPE, QXS >::check_size
bool check_size(const int nin, const int nvol, const int nex) const
checking size parameters.
Definition: afield.h:121
svbool_t
Definition: vsimd_double-inc.h:30
AField< REALTYPE, QXS >::set
void set(const int index, const real_t a)
Definition: afield.h:160
Communicator::ipe
static int ipe(const int dir)
logical coordinate of current proc.
Definition: communicator.cpp:105
ixyz
int ixyz
Definition: mult_Domainwall_eo_t_dirac_openacc-inc.h:13
Bridge::BridgeIO::crucial
void crucial(const char *format,...)
Definition: bridgeIO.cpp:242
Field
Container of Field-type object.
Definition: field.h:46
AField< REALTYPE, QXS >::nin
int nin() const
returning size of inner (on site) d.o.f.
Definition: afield.h:109
AField< REALTYPE, QXS >::cmp
real_t cmp(const int index) const
reference of data element
Definition: afield.h:157
QXS_Gauge::mult_Gdn
void mult_Gdn(AField< REALTYPE, QXS > &u, const int exu, const AField< REALTYPE, QXS > &v, const int exv, const AField< REALTYPE, QXS > &w, const int exw)
Definition: afield_Gauge-inc.h:392
Bridge::vout
BridgeIO vout
Definition: bridgeIO.cpp:572