Bridge++  Ver.2.1.3
mult_Domainwall_5din_matinv_openacc-inc.h
Go to the documentation of this file.
1 
10 #ifndef MULT_DOMAINWALL_5DIN_MATINV_ACC_INCLUDED
11 #define MULT_DOMAINWALL_5DIN_MATINV_ACC_INCLUDED
12 
13 //====================================================================
15  real_t *RESTRICT vp, real_t *RESTRICT wp,
16  int jd, int Ns,
17  real_t *RESTRICT mat_inv, int *Nsize)
18 {
19  int Nx = Nsize[0];
20  int Ny = Nsize[1];
21  int Nz = Nsize[2];
22  int Nt = Nsize[3];
23  int Nst = Nx * Ny * Nz * Nt;
24  int Nst_pad = CEIL_NWP(Nst);
25 
26  int Nin5 = NVCD * Ns;
27  int size = Nin5 * Nst_pad;
28  int mat_size = ND2 * Ns;
29  int mat_size2 = mat_size * mat_size;
30 
31 #pragma acc data present(vp[0:size], wp[0:size]) \
32  copyin(Ns, Nin5, Nst, Nst_pad, mat_inv[0:mat_size2])
33  {
34 #pragma acc parallel num_workers(NUM_WORKERS) vector_length(VECTOR_LENGTH)
35  {
36 
37 #pragma acc loop gang worker vector
38  for (int idx = 0; idx < Ns * NC * Nst_pad; ++idx) {
39  int idx2_wp = idx / NWP;
40  int idx_in = idx % NWP;
41  int ic = idx2_wp % NC;
42  int is1 = (idx2_wp/NC) % Ns;
43  int idx_out = idx2_wp/(NC * Ns);
44  int site = idx_in + NWP * idx_out;
45  if(site < Nst){
46 
47  for (int id = 0; id < ND2; ++id) {
48 
49  real_t vt1r = 0.0;
50  real_t vt1i = 0.0;
51  real_t vt2r = 0.0;
52  real_t vt2i = 0.0;
53  int idx1 = id + ND2 * is1;
54  for (int is2 = 0; is2 < Ns; ++is2) {
55  int idx2, ivcd1, ivcd2;
56  real_t mat1, mat2;
57  if(jd == 1){
58  idx2 = 0 + ND2 * is2;
59  mat1 = mat_inv[idx1 + (ND2 * Ns) * idx2];
60  idx2 = 1 + ND2 * is2;
61  mat2 = mat_inv[idx1 + (ND2 * Ns) * idx2];
62  }else{
63  idx2 = 0 + ND2 * is2;
64  mat1 = mat_inv[idx2 + (ND2 * Ns) * idx1];
65  idx2 = 1 + ND2 * is2;
66  mat2 = mat_inv[idx2 + (ND2 * Ns) * idx1];
67  }
68  ivcd1 = 0 + 2 * (ic + NC * 0);
69  ivcd2 = 0 + 2 * (ic + NC * 2);
70  vt1r += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
71  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
72 
73  ivcd1 = 1 + 2 * (ic + NC * 0);
74  ivcd2 = 1 + 2 * (ic + NC * 2);
75  vt1i += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
76  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
77 
78  ivcd1 = 0 + 2 * (ic + NC * 1);
79  ivcd2 = 0 + 2 * (ic + NC * 3);
80  vt2r += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
81  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
82 
83  ivcd1 = 1 + 2 * (ic + NC * 1);
84  ivcd2 = 1 + 2 * (ic + NC * 3);
85  vt2i += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
86  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
87  }
88  int id1 = 2 * id;
89  int id2 = 2 * id + 1;
90  vp[IDX2(Nin5, 0 + 2*(ic+NC*id1) + NVCD * is1, site)] = vt1r;
91  vp[IDX2(Nin5, 1 + 2*(ic+NC*id1) + NVCD * is1, site)] = vt1i;
92  vp[IDX2(Nin5, 0 + 2*(ic+NC*id2) + NVCD * is1, site)] = vt2r;
93  vp[IDX2(Nin5, 1 + 2*(ic+NC*id2) + NVCD * is1, site)] = vt2i;
94  }
95  }
96  }
97 
98  }
99  }
100 
101 }
102 
103 //====================================================================
105  real_t *RESTRICT vp, real_t *RESTRICT wp,
106  int jd, int Ns,
107  real_t *RESTRICT mat_inv, int *Nsize)
108 {
109  int Nx = Nsize[0];
110  int Ny = Nsize[1];
111  int Nz = Nsize[2];
112  int Nt = Nsize[3];
113  int Nst = Nx * Ny * Nz * Nt;
114  int Nst_pad = CEIL_NWP(Nst);
115 
116  int Nin5 = NVCD * Ns;
117  int size = Nin5 * Nst_pad;
118  int mat_size = ND2 * Ns;
119  int mat_size2 = mat_size * mat_size;
120 
121 #pragma acc data present(vp[0:size], wp[0:size]) \
122  copyin(Ns, Nin5, Nst, Nst_pad, mat_inv[0:mat_size2])
123  {
124 #pragma acc parallel num_workers(NUM_WORKERS) vector_length(VECTOR_LENGTH)
125  {
126 
127 #pragma acc loop gang worker vector
128  for (int idx = 0; idx < NC * Nst_pad; ++idx) {
129  int idx2_wp = idx / NWP;
130  int idx_in = idx % NWP;
131  int ic = idx2_wp % NC;
132  int idx_out = idx2_wp / NC;
133  int site = idx_in + NWP * idx_out;
134  if(site < Nst){
135 
136  for (int is1 = 0; is1 < Ns; ++is1) {
137  {
138  for (int id = 0; id < ND2; ++id) {
139 
140  real_t vt1r = 0.0;
141  real_t vt1i = 0.0;
142  real_t vt2r = 0.0;
143  real_t vt2i = 0.0;
144  int idx1 = id + ND2 * is1;
145  for (int is2 = 0; is2 < Ns; ++is2) {
146  int idx2, ivcd1, ivcd2;
147  real_t mat1, mat2;
148  if(jd == 1){
149  idx2 = 0 + ND2 * is2;
150  mat1 = mat_inv[idx1 + (ND2 * Ns) * idx2];
151  idx2 = 1 + ND2 * is2;
152  mat2 = mat_inv[idx1 + (ND2 * Ns) * idx2];
153  }else{
154  idx2 = 0 + ND2 * is2;
155  mat1 = mat_inv[idx2 + (ND2 * Ns) * idx1];
156  idx2 = 1 + ND2 * is2;
157  mat2 = mat_inv[idx2 + (ND2 * Ns) * idx1];
158  }
159  ivcd1 = 0 + 2 * (ic + NC * 0);
160  ivcd2 = 0 + 2 * (ic + NC * 2);
161  vt1r += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
162  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
163 
164  ivcd1 = 1 + 2 * (ic + NC * 0);
165  ivcd2 = 1 + 2 * (ic + NC * 2);
166  vt1i += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
167  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
168 
169  ivcd1 = 0 + 2 * (ic + NC * 1);
170  ivcd2 = 0 + 2 * (ic + NC * 3);
171  vt2r += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
172  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
173 
174  ivcd1 = 1 + 2 * (ic + NC * 1);
175  ivcd2 = 1 + 2 * (ic + NC * 3);
176  vt2i += mat1 * wp[IDX2(Nin5, (ivcd1 + NVCD * is2), site)]
177  + mat2 * wp[IDX2(Nin5, (ivcd2 + NVCD * is2), site)];
178  }
179  int id1 = 2 * id;
180  int id2 = 2 * id + 1;
181  vp[IDX2(Nin5, 0 + 2*(ic+NC*id1) + NVCD * is1, site)] = vt1r;
182  vp[IDX2(Nin5, 1 + 2*(ic+NC*id1) + NVCD * is1, site)] = vt1i;
183  vp[IDX2(Nin5, 0 + 2*(ic+NC*id2) + NVCD * is1, site)] = vt2r;
184  vp[IDX2(Nin5, 1 + 2*(ic+NC*id2) + NVCD * is1, site)] = vt2i;
185  }
186  }
187  }
188 
189  }
190  }
191 
192  }
193  }
194 
195  }
196 
197 #endif
198 //============================================================END=====
NVCD
#define NVCD
Definition: define_params_SU3.h:20
mult_domainwall_5din_ee_inv_dirac_4d
void mult_domainwall_5din_ee_inv_dirac_4d(real_t *RESTRICT vp, real_t *RESTRICT wp, int jd, int Ns, real_t *RESTRICT mat_inv, int *Nsize)
Definition: mult_Domainwall_5din_matinv_openacc-inc.h:104
CEIL_NWP
#define CEIL_NWP(nst)
Definition: define_params.h:47
NC
#define NC
Definition: field_F_imp_SU2-inc.h:15
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
real_t
double real_t
Definition: bridgeACC_AField_double.cpp:14
IDX2
#define IDX2(nin, in, ist)
Definition: define_index.h:28
NWP
#define NWP
Definition: define_params.h:32
mult_domainwall_5din_ee_inv_dirac_5d
void mult_domainwall_5din_ee_inv_dirac_5d(real_t *RESTRICT vp, real_t *RESTRICT wp, int jd, int Ns, real_t *RESTRICT mat_inv, int *Nsize)
Definition: mult_Domainwall_5din_matinv_openacc-inc.h:14
ND2
#define ND2
Definition: define_params_SU3.h:18