Bridge++  Ver.2.1.3
afield-inc.h
Go to the documentation of this file.
1 
10 #ifndef ACCEL_AFIELD_INC_INCLUDED
11 #define ACCEL_AFIELD_INC_INCLUDED
12 
13 #include <cstdlib>
14 
15 #include "complexTraits.h"
16 #include "lib/ResourceManager/threadManager.h"
17 
18 #include "lib_alt_Accel/inline/afield_th-inc.h"
19 #include "lib_alt_Accel/Field/aindex_lex.h"
20 #include "lib_alt_Accel/Field/aindex_eo.h"
21 
22 #include "lib_alt_Accel/BridgeACC/bridgeACC_AField.h"
23 
24 //====================================================================
25 template <typename REALTYPE>
26 void update_host(AField<REALTYPE, ACCEL>& v)
27 {
28  v.update_host();
29 }
30 
31 //====================================================================
32 template <typename REALTYPE>
33 void update_device(AField<REALTYPE, ACCEL>& v)
34 {
35  v.update_device();
36 }
37 
38 //====================================================================
39 template <typename REALTYPE>
40 void copy(AField<REALTYPE, ACCEL>& v, const AField<REALTYPE, ACCEL> &w)
41 {
42  v.copy(w);
43 }
44 
45 //====================================================================
46 template <typename REALTYPE>
47 void copy(AField<REALTYPE, ACCEL>& v, const int ex,
48  const AField<REALTYPE, ACCEL> &w, const int ex_w)
49 {
50  v.copy(ex, w, ex_w);
51 }
52 
53 //====================================================================
54 template <typename REALTYPE>
55 void axpy(AField<REALTYPE, ACCEL>& v, const REALTYPE a,
56  const AField<REALTYPE, ACCEL> &w)
57 {
58  v.axpy(a, w);
59 }
60 
61 //====================================================================
62 template <typename REALTYPE>
63 void axpy(AField<REALTYPE, ACCEL>& v,
64  const typename ComplexTraits<REALTYPE>::complex_t a,
65  const AField<REALTYPE, ACCEL> &w)
66 {
67  v.axpy(a, w);
68 }
69 
70 //====================================================================
71 template <typename REALTYPE>
72 void axpy(AField<REALTYPE, ACCEL>& v, const int ex,
73  const REALTYPE a,
74  const AField<REALTYPE, ACCEL> &w, const int ex_w)
75 {
76  v.axpy(ex, a, w, ex_w);
77 }
78 
79 //====================================================================
80 template <typename REALTYPE>
81 void aypx(const REALTYPE a, AField<REALTYPE, ACCEL>& v,
82  const AField<REALTYPE, ACCEL>& w)
83 {
84  v.aypx(a, w);
85 }
86 
87 //====================================================================
88 template <typename REALTYPE>
89 void aypx(const typename ComplexTraits<REALTYPE>::complex_t a,
90  AField<REALTYPE, ACCEL>& v,
91  const AField<REALTYPE, ACCEL>& w)
92 {
93  v.aypx(a, w);
94 }
95 
96 //====================================================================
97 template <typename REALTYPE>
98 void scal(AField<REALTYPE, ACCEL>& v, const REALTYPE a)
99 {
100  v.scal(a);
101 }
102 
103 //====================================================================
104 template <typename REALTYPE>
105 void scal(AField<REALTYPE, ACCEL>& v,
106  const typename ComplexTraits<REALTYPE>::complex_t a)
107 {
108  v.scal(a);
109 }
110 
111 //====================================================================
112 template <typename REALTYPE>
113 REALTYPE dot(AField<REALTYPE, ACCEL>& v, AField<REALTYPE, ACCEL>& w)
114 {
115  return v.dot(w);
116 }
117 
118 //====================================================================
119 template <typename REALTYPE>
120 typename ComplexTraits<REALTYPE>::complex_t dotc(
121  const AField<REALTYPE, ACCEL>& v,
122  const AField<REALTYPE, ACCEL>& w)
123 {
124  REALTYPE vw_r, vw_i;
125  v.dotc(vw_r, vw_i, w);
126  return cmplx(vw_r, vw_i);
127 }
128 
129 //====================================================================
130 template <typename REALTYPE>
131 REALTYPE norm2(const AField<REALTYPE, ACCEL>& v)
132 {
133  return v.norm2();
134 }
135 
136 //====================================================================
137 template <typename REALTYPE>
138 void xI(AField<REALTYPE, ACCEL>& v)
139 {
140  v.xI();
141 }
142 
143 //====================================================================
144 template <typename REALTYPE>
145 void conjg(AField<REALTYPE, ACCEL>& v)
146 {
147  v.conjg();
148 }
149 
150 //====================================================================
151 template <class INDEX, class AFIELD>
152 void convert(INDEX& index, AFIELD& v, const Field& w)
153 {
154 #pragma omp barrier
155 
156  int Nin = w.nin();
157  int Nvol = w.nvol();
158  int Nex = w.nex();
159 
160  vout.paranoiac("convert to ACCEL from Field start.\n");
161  vout.paranoiac(" AFIELD = %s\n", AFIELD::class_name.c_str());
162  vout.paranoiac(" Nin = %d Nvol = %d Nex = %d\n", Nin, Nvol, Nex);
163 
164  typedef typename AFIELD::real_t real_t;
165  real_t* vp = v.ptr(0);
166 
167  int Nvol_pad = v.nvol_pad();
168  if(Nvol_pad != index.nvol_pad()){
169  vout.crucial("convert: inconsistent Nvol_pad in AField and AIndex");
170  exit(EXIT_FAILURE);
171  }
172 
173  int ith, nth, is, ns;
174  set_threadtask(ith, nth, is, ns, Nvol_pad);
175 
176  for(int ex = 0; ex < Nex; ++ex){
177  for(int site = is; site < ns; ++site){
178 
179  if(site < Nvol){
180  for(int in = 0; in < Nin; ++in){
181  int iw = in + Nin * (site + Nvol * ex);
182  int iv = index.idx(in, Nin, site, ex);
183  vp[iv] = w.cmp(iw);
184  }
185  }else{
186  for(int in = 0; in < Nin; ++in){
187  int iv = index.idx(in, Nin, site, ex);
188  vp[iv] = real_t(0.0);
189  }
190  }
191 
192  }
193  }
194 
195 #pragma omp barrier
196 
197  if(ith == 0){
198  size_t nv = Nin * Nvol_pad * Nex;
199  BridgeACC::copy_to_device(vp, nv);
200  }
201 
202  vout.paranoiac("convert to ACCEL from Field finished.\n");
203 
204 #pragma omp barrier
205 
206 }
207 
208 //====================================================================
209 template <class INDEX, class AFIELD>
210 void convert_spinor(INDEX& index, AFIELD& v, const Field& w)
211 {
212  int Nc = CommonParameters::Nc();
213  int Nd = CommonParameters::Nd();
214  assert(w.nin() == 2*Nc*Nd);
215  assert(v.nin() == 2*Nc*Nd);
216 
217  convert(index, v, w);
218 
219 }
220 
221 //====================================================================
222 template <class INDEX, class AFIELD>
223 void convert_gauge(INDEX& index, AFIELD& v, const Field& w)
224 {
225  int Nc = CommonParameters::Nc();
226  assert(w.nin() == 2*Nc*Nc);
227  assert(v.nin() == 2*Nc*Nc);
228 
229  convert(index, v, w);
230 
231 }
232 
233 //====================================================================
234 template<class INDEX2, class AFIELD2, class INDEX1, class AFIELD1>
235 void convert(INDEX2& index2, AFIELD2& v2,
236  const INDEX1& index1, const AFIELD1& v1)
237 {
238  int Nin = v1.nin();
239  int Nvol = v1.nvol();
240  int Nex = v1.nex();
241  //assert(v2.check_size(Nin, Nvol, Nex)); // this does not hold
242 
243  typename AFIELD1::real_t *v1p = const_cast<AFIELD1*>(&v1)->ptr(0);
244 
245  int ith, nth, is, ns;
246  set_threadtask(ith, nth, is, ns, Nvol);
247 
248 #pragma omp barrier
249 
250  if(ith == 0){
251  int nv = Nin * Nvol * Nex;
252  BridgeACC::copy_from_device(v1p, nv);
253  }
254 
255 #pragma omp barrier
256 
257  typename AFIELD2::real_t *v2p = v2.ptr(0);
258 
259  for(int ex = 0; ex < Nex; ++ex){
260  for(int site = is; site < ns; ++site){
261  for(int in = 0; in < Nin; ++in){
262  int iv1 = index1.idx(in, Nin, site, ex);
263  int iv2 = index2.idx(in, Nin, site, ex);
264  v2p[iv2] = v1p[iv1];
265  }
266  }
267  }
268 
269 #pragma omp barrier
270 
271  if(ith == 0){
272  int nv = Nin * Nvol * Nex;
273  BridgeACC::copy_to_device(v2p, nv);
274  }
275 
276 #pragma omp barrier
277 
278 }
279 
280 //====================================================================
281 template<>
282 void convert(AIndex_lex<double,ACCEL>& index2,
283  AField<double,ACCEL>& v2,
284  const AIndex_lex<float,ACCEL>& index1,
285  const AField<float,ACCEL>& v1)
286 {
287  int Nin = v1.nin();
288  int Nvol = v1.nvol();
289  int Nex = v1.nex();
290 
291  int Nvol_pad = CEIL_NWP(Nvol);
292 
293  double *v2p = v2.ptr(0);
294  float *v1p = const_cast<AField<float,ACCEL>* >(&v1)->ptr(0);
295 
296  int ith, nth;
297  set_thread(ith, nth);
298 
299 #pragma omp barrier
300 
301  if(ith == 0){
302  int nvx = Nvol_pad * Nex;
303  BridgeACC::copy(v2p, 0, v1p, 0, Nin, nvx);
304  }
305 
306 #pragma omp barrier
307 
308 }
309 
310 //====================================================================
311 template<>
312 void convert(AIndex_lex<float,ACCEL>& index2,
313  AField<float,ACCEL>& v2,
314  const AIndex_lex<double,ACCEL>& index1,
315  const AField<double,ACCEL>& v1)
316 {
317  int Nin = v1.nin();
318  int Nvol = v1.nvol();
319  int Nex = v1.nex();
320  //assert(v2.check_size(Nin, Nvol, Nex)); // this does not hold
321 
322  int Nvol_pad = CEIL_NWP(Nvol);
323 
324  float *v2p = v2.ptr(0);
325  double *v1p = const_cast<AField<double,ACCEL>* >(&v1)->ptr(0);
326 
327  int ith, nth;
328  set_thread(ith, nth);
329 
330 #pragma omp barrier
331 
332  if(ith == 0){
333  int nvx = Nvol_pad * Nex;
334  BridgeACC::copy(v2p, 0, v1p, 0, Nin, nvx);
335  }
336 
337 #pragma omp barrier
338 
339 }
340 
341 //====================================================================
342 template <class INDEX2, class AFIELD2, class INDEX1, class AFIELD1>
343 void convert_h(INDEX2& index2, AFIELD2& v2,
344  const INDEX1& index1, const AFIELD1& v1)
345 {
346  int Nin = v1.nin();
347  int Nvol = v1.nvol();
348  int Nex = v1.nex();
349 
350  int Nvol_pad = CEIL_NWP(Nvol);
351 
352  int ith, nth, is, ns;
353  set_threadtask(ith, nth, is, ns, Nvol);
354 
355 #pragma omp barrier
356 
357  typename AFIELD1::real_t *v1p = const_cast<AFIELD1*>(&v1)->ptr(0);
358 
359  if(ith == 0){
360  int nv = Nin * Nvol * Nex;
361  BridgeACC::copy_from_device(v1p, nv);
362  }
363 
364 #pragma omp barrier
365 
366  typename AFIELD2::real_t *v2p = v2.ptr(0);
367 
368  for(int ex = 0; ex < Nex; ++ex){
369  for(int site = is; site < ns; ++site){
370  for(int in = 0; in < Nin; ++in){
371  int iv1 = index1.idxh(in, Nin, site, ex);
372  int iv2 = index2.idxh(in, Nin, site, ex);
373  v2p[iv2] = v1p[iv1];
374  }
375  }
376  }
377 
378 #pragma omp barrier
379 
380  if(ith == 0){
381  int nv = Nin * Nvol * Nex;
382  BridgeACC::copy_to_device(v2p, nv);
383  }
384 
385 #pragma omp barrier
386 
387 }
388 
389 
390 //====================================================================
391 template<>
392 void convert_h(AIndex_eo<double,ACCEL>& index2,
393  AField<double,ACCEL>& v2,
394  const AIndex_eo<float,ACCEL>& index1,
395  const AField<float,ACCEL>& v1)
396 {
397  int Nin = v1.nin();
398  int Nvol = v1.nvol();
399  int Nex = v1.nex();
400 
401  int Nvol_pad = CEIL_NWP(Nvol);
402 
403  double *v2p = v2.ptr(0);
404  float *v1p = const_cast<AField<float,ACCEL>* >(&v1)->ptr(0);
405 
406  int ith, nth;
407  set_thread(ith, nth);
408 
409 #pragma omp barrier
410 
411  if(ith == 0){
412  int nvx = Nvol_pad * Nex;
413  BridgeACC::copy(v2p, 0, v1p, 0, Nin, nvx);
414  }
415 
416 #pragma omp barrier
417 
418 }
419 
420 //====================================================================
421 template<>
422 void convert_h(AIndex_eo<float,ACCEL>& index2,
423  AField<float,ACCEL>& v2,
424  const AIndex_eo<double,ACCEL>& index1,
425  const AField<double,ACCEL>& v1)
426 {
427  int Nin = v1.nin();
428  int Nvol = v1.nvol();
429  int Nex = v1.nex();
430 
431  int Nvol_pad = CEIL_NWP(Nvol);
432 
433  float *v2p = v2.ptr(0);
434  double *v1p = const_cast<AField<double,ACCEL>* >(&v1)->ptr(0);
435 
436  int ith, nth;
437  set_thread(ith, nth);
438 
439 #pragma omp barrier
440 
441  if(ith == 0){
442  int nvx = Nvol_pad * Nex;
443  BridgeACC::copy(v2p, 0, v1p, 0, Nin, nvx);
444  }
445 
446 #pragma omp barrier
447 
448 }
449 
450 //====================================================================
451 template <class INDEX, class AFIELD>
452 void reverse(INDEX& index, Field& v, const AFIELD& w)
453 {
454 #pragma omp barrier
455 
456  int Nin = v.nin();
457  int Nvol = v.nvol();
458  int Nex = v.nex();
459 
460  vout.paranoiac("reverse to Field from ACCEL start.\n");
461 
462  int Nvol_pad = w.nvol_pad();
463  if(Nvol_pad != index.nvol_pad()){
464  vout.crucial("reverse: inconsistent Nvol_pad in AField and AIndex");
465  exit(EXIT_FAILURE);
466  }
467 
468  int ith, nth, is, ns;
469  set_threadtask(ith, nth, is, ns, Nvol);
470 
471  typedef typename AFIELD::real_t real_t;
472  real_t *wp = const_cast<AFIELD*>(&w)->ptr(0);
473 
474  if(ith == 0){
475  int nv = Nin * Nvol_pad * Nex;
476  BridgeACC::copy_from_device(wp, nv);
477  }
478 
479 #pragma omp barrier
480 
481  for(int ex = 0; ex < Nex; ++ex){
482  for(int site = is; site < ns; ++site){
483  for(int in = 0; in < Nin; ++in){
484  int iv = in + Nin * (site + Nvol * ex);
485  int iw = index.idx(in, Nin, site, ex);
486  v.set(iv, double(wp[iw]));
487  }
488  }
489  }
490 
491  vout.paranoiac("reverse to Field from ACCEL finished.\n");
492 
493 #pragma omp barrier
494 
495 }
496 
497 //====================================================================
498 template <class AFIELD>
499 void reverse(AIndex_lex<typename AFIELD::real_t,ACCEL>& index,
500  Field& v, const AFIELD& w)
501 {
502 #pragma omp barrier
503 
504  int Nin = v.nin();
505  int Nvol = v.nvol();
506  int Nex = v.nex();
507 
508  vout.paranoiac("reverse to Field from ACCEL start.\n");
509 
510  typedef typename AFIELD::real_t real_t;
511  double *vp = v.ptr(0);
512  real_t *wp = const_cast<AFIELD*>(&w)->ptr(0);
513 
514  int Nvol_pad = w.nvol_pad();
515  if(Nvol_pad != index.nvol_pad()){
516  vout.crucial("reverse: inconsistent Nvol_pad in AField and AIndex");
517  exit(EXIT_FAILURE);
518  }
519 
520  int ith, nth;
521  set_thread(ith, nth);
522 
523  if(ith == 0){
524  for(int ex = 0; ex < Nex; ++ex){
525  double *vp2 = &vp[Nin * Nvol * ex];
526  real_t *wp2 = &wp[Nin * Nvol_pad * ex];
527  BridgeACC::reverse(vp2, wp2, Nin, Nvol, Nvol_pad);
528  }
529  }
530 
531  vout.paranoiac("reverse to Field from ACCEL finished.\n");
532 
533 #pragma omp barrier
534 
535 }
536 
537 //====================================================================
538 template <class INDEX, class AFIELD>
539 void reverse_spinor(INDEX& index, Field& v, const AFIELD& w)
540 {
541  int Nc = CommonParameters::Nc();
542  int Nd = CommonParameters::Nd();
543  assert(w.nin() == 2*Nc*Nd);
544  assert(v.nin() == 2*Nc*Nd);
545 
546  reverse(index, v, w);
547 
548 }
549 
550 //====================================================================
551 template <class INDEX, class AFIELD>
552 void reverse_gauge(INDEX& index, Field& v, const AFIELD& w)
553 {
554  int Nc = CommonParameters::Nc();
555  assert(w.nin() == 2*Nc*Nc);
556  assert(v.nin() == 2*Nc*Nc);
557 
558  reverse(index, v, w);
559 
560 }
561 
562 //============================================================END=====
563 #endif
AField< REALTYPE, ACCEL >::xI
void xI()
multiplying imaginary unit.
Definition: afield-tmpl.h:751
afield_th-inc.h
convert_gauge
void convert_gauge(INDEX &index, AFIELD &v, const Field &w)
Definition: afield-inc.h:223
AField< REALTYPE, ACCEL >::aypx
void aypx(const real_t, const AField< real_t, ACCEL > &)
Definition: afield-tmpl.h:456
norm2
REALTYPE norm2(const AField< REALTYPE, ACCEL > &v)
Definition: afield-inc.h:131
Field::set
void set(const int jin, const int site, const int jex, double v)
Definition: field.h:175
AIndex_lex
Definition: aindex_lex_base.h:17
BridgeACC::copy
void copy(double *v, double *w, int nin, int nvol)
update_device
void update_device(AField< REALTYPE, ACCEL > &v)
Definition: afield-inc.h:33
AField< REALTYPE, ACCEL >::update_host
void update_host() const
Definition: afield-tmpl.h:152
AField
Definition: afield_base.h:16
Field::nex
int nex() const
Definition: field.h:128
AField< REALTYPE, ACCEL >::dotc
void dotc(real_t &, real_t &, const AField< real_t, ACCEL > &) const
complex inner-product: real and imaginary parts in order.
Definition: afield-tmpl.h:666
aindex_eo.h
BridgeACC::copy_from_device
void copy_from_device(double *v, int nv)
update_host
void update_host(AField< REALTYPE, ACCEL > &v)
Definition: afield-inc.h:26
reverse
void reverse(INDEX &index, Field &v, const AFIELD &w)
Definition: afield-inc.h:452
Field::nin
int nin() const
Definition: field.h:126
Field::real_t
double real_t
Definition: field.h:51
CEIL_NWP
#define CEIL_NWP(nst)
Definition: define_params.h:47
AField< REALTYPE, ACCEL >::axpy
void axpy(const real_t, const AField< real_t, ACCEL > &)
Definition: afield-tmpl.h:325
Bridge::BridgeIO::paranoiac
void paranoiac(const char *format,...)
Definition: bridgeIO.cpp:300
AField< REALTYPE, ACCEL >::copy
void copy(const Field &w)
Definition: afield-tmpl.h:247
aypx
void aypx(const REALTYPE a, AField< REALTYPE, ACCEL > &v, const AField< REALTYPE, ACCEL > &w)
Definition: afield-inc.h:81
AField< REALTYPE, ACCEL >::conjg
void conjg()
taking complex conjugate.
Definition: afield-tmpl.h:772
CommonParameters::Nc
static int Nc()
Definition: commonParameters.h:115
reverse_gauge
void reverse_gauge(INDEX &index, Field &v, const AFIELD &w)
Definition: afield-inc.h:552
dotc
ComplexTraits< REALTYPE >::complex_t dotc(const AField< REALTYPE, ACCEL > &v, const AField< REALTYPE, ACCEL > &w)
Definition: afield-inc.h:120
convert
void convert(INDEX &index, AFIELD &v, const Field &w)
Definition: afield-inc.h:152
reverse_spinor
void reverse_spinor(INDEX &index, Field &v, const AFIELD &w)
Definition: afield-inc.h:539
Field::class_name
static const std::string class_name
Definition: field.h:52
threadManager.h
AField< REALTYPE, ACCEL >::norm2
real_t norm2(void) const
square norm squared (|v|^2).
Definition: afield-tmpl.h:698
copy
void copy(AField< REALTYPE, ACCEL > &v, const AField< REALTYPE, ACCEL > &w)
Definition: afield-inc.h:40
aindex_lex.h
Field::nvol
int nvol() const
Definition: field.h:127
AField< REALTYPE, ACCEL >::scal
void scal(const real_t)
Definition: afield-tmpl.h:565
ComplexTraits
Definition: complexTraits.h:16
real_t
double real_t
Definition: bridgeACC_AField_double.cpp:14
conjg
void conjg(AField< REALTYPE, ACCEL > &v)
Definition: afield-inc.h:145
Field::cmp
double cmp(const int jin, const int site, const int jex) const
Definition: field.h:143
Field::ptr
const double * ptr(const int jin, const int site, const int jex) const
Definition: field.h:153
CommonParameters::Nd
static int Nd()
Definition: commonParameters.h:116
AField< REALTYPE, ACCEL >::dot
real_t dot(const AField< real_t, ACCEL > &)
Definition: afield-tmpl.h:638
AField< REALTYPE, ACCEL >::update_device
void update_device()
Definition: afield-tmpl.h:168
dot
REALTYPE dot(AField< REALTYPE, ACCEL > &v, AField< REALTYPE, ACCEL > &w)
Definition: afield-inc.h:113
xI
void xI(AField< REALTYPE, ACCEL > &v)
Definition: afield-inc.h:138
convert_h
void convert_h(INDEX2 &index2, AFIELD2 &v2, const INDEX1 &index1, const AFIELD1 &v1)
Definition: afield-inc.h:343
Bridge::BridgeIO::crucial
void crucial(const char *format,...)
Definition: bridgeIO.cpp:242
complexTraits.h
Field
Container of Field-type object.
Definition: field.h:46
scal
void scal(AField< REALTYPE, ACCEL > &v, const REALTYPE a)
Definition: afield-inc.h:98
convert_spinor
void convert_spinor(INDEX &index, AFIELD &v, const Field &w)
Definition: afield-inc.h:210
BridgeACC::copy_to_device
void copy_to_device(double *v, int nv)
AField< REALTYPE, ACCEL >
Definition: afield.h:31
axpy
void axpy(AField< REALTYPE, ACCEL > &v, const REALTYPE a, const AField< REALTYPE, ACCEL > &w)
Definition: afield-inc.h:55
Bridge::vout
BridgeIO vout
Definition: bridgeIO.cpp:572
AIndex_eo
Definition: aindex_eo_base.h:17
BridgeACC::reverse
void reverse(double *v, double *w, int nin, int nvol, int nvol_pad)