Block-Structured AMR Software Framework
AMReX_MLNodeLinOp_3D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_ML_NODE_LINOP_3D_K_H_
2 #define AMREX_ML_NODE_LINOP_3D_K_H_
3 #include <AMReX_Config.H>
4 
5 namespace amrex {
6 
7 template <typename T>
8 inline void mlndlap_bc_doit (Box const& vbx, Array4<T> const& a, Box const& domain,
9  GpuArray<bool,AMREX_SPACEDIM> const& bflo,
10  GpuArray<bool,AMREX_SPACEDIM> const& bfhi) noexcept
11 {
12  Box gdomain = domain;
13  for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
14  if (! bflo[idim]) { gdomain.growLo(idim,1); }
15  if (! bfhi[idim]) { gdomain.growHi(idim,1); }
16  }
17 
18  if (gdomain.strictly_contains(vbx)) { return; }
19 
20  const int offset = domain.cellCentered() ? 0 : 1;
21 
22  const auto dlo = amrex::lbound(domain);
23  const auto dhi = amrex::ubound(domain);
24 
25  Box const& sbox = amrex::grow(vbx,1);
26  AMREX_HOST_DEVICE_FOR_3D(sbox, i, j, k,
27  {
28  if (! gdomain.contains(IntVect(i,j,k))) {
29  // xlo & ylo & zlo
30  if (i == dlo.x-1 && j == dlo.y-1 && k == dlo.z-1 && (bflo[0] || bflo[1] || bflo[2]))
31  {
32  if (bflo[0] && bflo[1] && bflo[2])
33  {
34  a(i,j,k) = a(i+1+offset, j+1+offset, k+1+offset);
35  }
36  else if (bflo[0] && bflo[1])
37  {
38  a(i,j,k) = a(i+1+offset, j+1+offset, k);
39  }
40  else if (bflo[0] && bflo[2])
41  {
42  a(i,j,k) = a(i+1+offset, j, k+1+offset);
43  }
44  else if (bflo[1] && bflo[2])
45  {
46  a(i,j,k) = a(i, j+1+offset, k+1+offset);
47  }
48  else if (bflo[0])
49  {
50  a(i,j,k) = a(i+1+offset, j, k);
51  }
52  else if (bflo[1])
53  {
54  a(i,j,k) = a(i, j+1+offset, k);
55  }
56  else if (bflo[2])
57  {
58  a(i,j,k) = a(i, j, k+1+offset);
59  }
60  }
61  // xhi & ylo & zlo
62  else if (i == dhi.x+1 && j == dlo.y-1 && k == dlo.z-1 && (bfhi[0] || bflo[1] || bflo[2]))
63  {
64  if (bfhi[0] && bflo[1] && bflo[2])
65  {
66  a(i,j,k) = a(i-1-offset, j+1+offset, k+1+offset);
67  }
68  else if (bfhi[0] && bflo[1])
69  {
70  a(i,j,k) = a(i-1-offset, j+1+offset, k);
71  }
72  else if (bfhi[0] && bflo[2])
73  {
74  a(i,j,k) = a(i-1-offset, j, k+1+offset);
75  }
76  else if (bflo[1] && bflo[2])
77  {
78  a(i,j,k) = a(i, j+1+offset, k+1+offset);
79  }
80  else if (bfhi[0])
81  {
82  a(i,j,k) = a(i-1-offset, j, k);
83  }
84  else if (bflo[1])
85  {
86  a(i,j,k) = a(i, j+1+offset, k);
87  }
88  else if (bflo[2])
89  {
90  a(i,j,k) = a(i, j, k+1+offset);
91  }
92  }
93  // xlo & yhi & zlo
94  else if (i == dlo.x-1 && j == dhi.y+1 && k == dlo.z-1 && (bflo[0] || bfhi[1] || bflo[2]))
95  {
96  if (bflo[0] && bfhi[1] && bflo[2])
97  {
98  a(i,j,k) = a(i+1+offset, j-1-offset, k+1+offset);
99  }
100  else if (bflo[0] && bfhi[1])
101  {
102  a(i,j,k) = a(i+1+offset, j-1-offset, k);
103  }
104  else if (bflo[0] && bflo[2])
105  {
106  a(i,j,k) = a(i+1+offset, j, k+1+offset);
107  }
108  else if (bfhi[1] && bflo[2])
109  {
110  a(i,j,k) = a(i, j-1-offset, k+1+offset);
111  }
112  else if (bflo[0])
113  {
114  a(i,j,k) = a(i+1+offset, j, k);
115  }
116  else if (bfhi[1])
117  {
118  a(i,j,k) = a(i, j-1-offset, k);
119  }
120  else if (bflo[2])
121  {
122  a(i,j,k) = a(i, j, k+1+offset);
123  }
124  }
125  // xhi & yhi & zlo
126  else if (i == dhi.x+1 && j == dhi.y+1 && k == dlo.z-1 && (bfhi[0] || bfhi[1] || bflo[2]))
127  {
128  if (bfhi[0] && bfhi[1] && bflo[2])
129  {
130  a(i,j,k) = a(i-1-offset, j-1-offset, k+1+offset);
131  }
132  else if (bfhi[0] && bfhi[1])
133  {
134  a(i,j,k) = a(i-1-offset, j-1-offset, k);
135  }
136  else if (bfhi[0] && bflo[2])
137  {
138  a(i,j,k) = a(i-1-offset, j, k+1+offset);
139  }
140  else if (bfhi[1] && bflo[2])
141  {
142  a(i,j,k) = a(i, j-1-offset, k+1+offset);
143  }
144  else if (bfhi[0])
145  {
146  a(i,j,k) = a(i-1-offset, j, k);
147  }
148  else if (bfhi[1])
149  {
150  a(i,j,k) = a(i, j-1-offset, k);
151  }
152  else if (bflo[2])
153  {
154  a(i,j,k) = a(i, j, k+1+offset);
155  }
156  }
157  // xlo & ylo & zhi
158  else if (i == dlo.x-1 && j == dlo.y-1 && k == dhi.z+1 && (bflo[0] || bflo[1] || bfhi[2]))
159  {
160  if (bflo[0] && bflo[1] && bfhi[2])
161  {
162  a(i,j,k) = a(i+1+offset, j+1+offset, k-1-offset);
163  }
164  else if (bflo[0] && bflo[1])
165  {
166  a(i,j,k) = a(i+1+offset, j+1+offset, k);
167  }
168  else if (bflo[0] && bfhi[2])
169  {
170  a(i,j,k) = a(i+1+offset, j, k-1-offset);
171  }
172  else if (bflo[1] && bfhi[2])
173  {
174  a(i,j,k) = a(i, j+1+offset, k-1-offset);
175  }
176  else if (bflo[0])
177  {
178  a(i,j,k) = a(i+1+offset, j, k);
179  }
180  else if (bflo[1])
181  {
182  a(i,j,k) = a(i, j+1+offset, k);
183  }
184  else if (bfhi[2])
185  {
186  a(i,j,k) = a(i, j, k-1-offset);
187  }
188  }
189  // xhi & ylo & zhi
190  else if (i == dhi.x+1 && j == dlo.y-1 && k == dhi.z+1 && (bfhi[0] || bflo[1] || bfhi[2]))
191  {
192  if (bfhi[0] && bflo[1] && bfhi[2])
193  {
194  a(i,j,k) = a(i-1-offset, j+1+offset, k-1-offset);
195  }
196  else if (bfhi[0] && bflo[1])
197  {
198  a(i,j,k) = a(i-1-offset, j+1+offset, k);
199  }
200  else if (bfhi[0] && bfhi[2])
201  {
202  a(i,j,k) = a(i-1-offset, j, k-1-offset);
203  }
204  else if (bflo[1] && bfhi[2])
205  {
206  a(i,j,k) = a(i, j+1+offset, k-1-offset);
207  }
208  else if (bfhi[0])
209  {
210  a(i,j,k) = a(i-1-offset, j, k);
211  }
212  else if (bflo[1])
213  {
214  a(i,j,k) = a(i, j+1+offset, k);
215  }
216  else if (bfhi[2])
217  {
218  a(i,j,k) = a(i, j, k-1-offset);
219  }
220  }
221  // xlo & yhi & zhi
222  else if (i == dlo.x-1 && j == dhi.y+1 && k == dhi.z+1 && (bflo[0] || bfhi[1] || bfhi[2]))
223  {
224  if (bflo[0] && bfhi[1] && bfhi[2])
225  {
226  a(i,j,k) = a(i+1+offset, j-1-offset, k-1-offset);
227  }
228  else if (bflo[0] && bfhi[1])
229  {
230  a(i,j,k) = a(i+1+offset, j-1-offset, k);
231  }
232  else if (bflo[0] && bfhi[2])
233  {
234  a(i,j,k) = a(i+1+offset, j, k-1-offset);
235  }
236  else if (bfhi[1] && bfhi[2])
237  {
238  a(i,j,k) = a(i, j-1-offset, k-1-offset);
239  }
240  else if (bflo[0])
241  {
242  a(i,j,k) = a(i+1+offset, j, k);
243  }
244  else if (bfhi[1])
245  {
246  a(i,j,k) = a(i, j-1-offset, k);
247  }
248  else if (bfhi[2])
249  {
250  a(i,j,k) = a(i, j, k-1-offset);
251  }
252  }
253  // xhi & yhi & zhi
254  else if (i == dhi.x+1 && j == dhi.y+1 && k == dhi.z+1 && (bfhi[0] || bfhi[1] || bfhi[2]))
255  {
256  if (bfhi[0] && bfhi[1] && bfhi[2])
257  {
258  a(i,j,k) = a(i-1-offset, j-1-offset, k-1-offset);
259  }
260  else if (bfhi[0] && bfhi[1])
261  {
262  a(i,j,k) = a(i-1-offset, j-1-offset, k);
263  }
264  else if (bfhi[0] && bfhi[2])
265  {
266  a(i,j,k) = a(i-1-offset, j, k-1-offset);
267  }
268  else if (bfhi[1] && bfhi[2])
269  {
270  a(i,j,k) = a(i, j-1-offset, k-1-offset);
271  }
272  else if (bfhi[0])
273  {
274  a(i,j,k) = a(i-1-offset, j, k);
275  }
276  else if (bfhi[1])
277  {
278  a(i,j,k) = a(i, j-1-offset, k);
279  }
280  else if (bfhi[2])
281  {
282  a(i,j,k) = a(i, j, k-1-offset);
283  }
284  }
285  // xlo & ylo
286  else if (i == dlo.x-1 && j == dlo.y-1 && (bflo[0] || bflo[1]))
287  {
288  if (bflo[0] && bflo[1])
289  {
290  a(i,j,k) = a(i+1+offset, j+1+offset, k);
291  }
292  else if (bflo[0])
293  {
294  a(i,j,k) = a(i+1+offset, j, k);
295  }
296  else if (bflo[1])
297  {
298  a(i,j,k) = a(i, j+1+offset, k);
299  }
300  }
301  // xhi & ylo
302  else if (i == dhi.x+1 && j == dlo.y-1 && (bfhi[0] || bflo[1]))
303  {
304  if (bfhi[0] && bflo[1])
305  {
306  a(i,j,k) = a(i-1-offset, j+1+offset, k);
307  }
308  else if (bfhi[0])
309  {
310  a(i,j,k) = a(i-1-offset, j, k);
311  }
312  else if (bflo[1])
313  {
314  a(i,j,k) = a(i, j+1+offset, k);
315  }
316  }
317  // xlo & yhi
318  else if (i == dlo.x-1 && j == dhi.y+1 && (bflo[0] || bfhi[1]))
319  {
320  if (bflo[0] && bfhi[1])
321  {
322  a(i,j,k) = a(i+1+offset, j-1-offset, k);
323  }
324  else if (bflo[0])
325  {
326  a(i,j,k) = a(i+1+offset, j, k);
327  }
328  else if (bfhi[1])
329  {
330  a(i,j,k) = a(i, j-1-offset, k);
331  }
332  }
333  // xhi & yhi
334  else if (i == dhi.x+1 && j == dhi.y+1 && (bfhi[0] || bfhi[1]))
335  {
336  if (bfhi[0] && bfhi[1])
337  {
338  a(i,j,k) = a(i-1-offset, j-1-offset, k);
339  }
340  else if (bfhi[0])
341  {
342  a(i,j,k) = a(i-1-offset, j, k);
343  }
344  else if (bfhi[1])
345  {
346  a(i,j,k) = a(i, j-1-offset, k);
347  }
348  }
349  // xlo & zlo
350  else if (i == dlo.x-1 && k == dlo.z-1 && (bflo[0] || bflo[2]))
351  {
352  if (bflo[0] && bflo[2])
353  {
354  a(i,j,k) = a(i+1+offset, j, k+1+offset);
355  }
356  else if (bflo[0])
357  {
358  a(i,j,k) = a(i+1+offset, j, k);
359  }
360  else if (bflo[2])
361  {
362  a(i,j,k) = a(i, j, k+1+offset);
363  }
364  }
365  // xhi & zlo
366  else if (i == dhi.x+1 && k == dlo.z-1 && (bfhi[0] || bflo[2]))
367  {
368  if (bfhi[0] && bflo[2])
369  {
370  a(i,j,k) = a(i-1-offset, j, k+1+offset);
371  }
372  else if (bfhi[0])
373  {
374  a(i,j,k) = a(i-1-offset, j, k);
375  }
376  else if (bflo[2])
377  {
378  a(i,j,k) = a(i, j, k+1+offset);
379  }
380  }
381  // xlo & zhi
382  else if (i == dlo.x-1 && k == dhi.z+1 && (bflo[0] || bfhi[2]))
383  {
384  if (bflo[0] && bfhi[2])
385  {
386  a(i,j,k) = a(i+1+offset, j, k-1-offset);
387  }
388  else if (bflo[0])
389  {
390  a(i,j,k) = a(i+1+offset, j, k);
391  }
392  else if (bfhi[2])
393  {
394  a(i,j,k) = a(i, j, k-1-offset);
395  }
396  }
397  // xhi & zhi
398  else if (i == dhi.x+1 && k == dhi.z+1 && (bfhi[0] || bfhi[2]))
399  {
400  if (bfhi[0] && bfhi[2])
401  {
402  a(i,j,k) = a(i-1-offset, j, k-1-offset);
403  }
404  else if (bfhi[0])
405  {
406  a(i,j,k) = a(i-1-offset, j, k);
407  }
408  else if (bfhi[2])
409  {
410  a(i,j,k) = a(i, j, k-1-offset);
411  }
412  }
413  // ylo & zlo
414  else if (j == dlo.y-1 && k == dlo.z-1 && (bflo[1] || bflo[2]))
415  {
416  if (bflo[1] && bflo[2])
417  {
418  a(i,j,k) = a(i, j+1+offset, k+1+offset);
419  }
420  else if (bflo[1])
421  {
422  a(i,j,k) = a(i, j+1+offset, k);
423  }
424  else if (bflo[2])
425  {
426  a(i,j,k) = a(i, j, k+1+offset);
427  }
428  }
429  // yhi & zlo
430  else if (j == dhi.y+1 && k == dlo.z-1 && (bfhi[1] || bflo[2]))
431  {
432  if (bfhi[1] && bflo[2])
433  {
434  a(i,j,k) = a(i, j-1-offset, k+1+offset);
435  }
436  else if (bfhi[1])
437  {
438  a(i,j,k) = a(i, j-1-offset, k);
439  }
440  else if (bflo[2])
441  {
442  a(i,j,k) = a(i, j, k+1+offset);
443  }
444  }
445  // ylo & zhi
446  else if (j == dlo.y-1 && k == dhi.z+1 && (bflo[1] || bfhi[2]))
447  {
448  if (bflo[1] && bfhi[2])
449  {
450  a(i,j,k) = a(i, j+1+offset, k-1-offset);
451  }
452  else if (bflo[1])
453  {
454  a(i,j,k) = a(i, j+1+offset, k);
455  }
456  else if (bfhi[2])
457  {
458  a(i,j,k) = a(i, j, k-1-offset);
459  }
460  }
461  // yhi & zhi
462  else if (j == dhi.y+1 && k == dhi.z+1 && (bfhi[1] || bfhi[2]))
463  {
464  if (bfhi[1] && bfhi[2])
465  {
466  a(i,j,k) = a(i, j-1-offset, k-1-offset);
467  }
468  else if (bfhi[1])
469  {
470  a(i,j,k) = a(i, j-1-offset, k);
471  }
472  else if (bfhi[2])
473  {
474  a(i,j,k) = a(i, j, k-1-offset);
475  }
476  }
477  else if (i == dlo.x-1 && bflo[0])
478  {
479  a(i,j,k) = a(i+1+offset, j, k);
480  }
481  else if (i == dhi.x+1 && bfhi[0])
482  {
483  a(i,j,k) = a(i-1-offset, j, k);
484  }
485  else if (j == dlo.y-1 && bflo[1])
486  {
487  a(i,j,k) = a(i, j+1+offset, k);
488  }
489  else if (j == dhi.y+1 && bfhi[1])
490  {
491  a(i,j,k) = a(i, j-1-offset, k);
492  }
493  else if (k == dlo.z-1 && bflo[2])
494  {
495  a(i,j,k) = a(i, j, k+1+offset);
496  }
497  else if (k == dhi.z+1 && bfhi[2])
498  {
499  a(i,j,k) = a(i, j, k-1-offset);
500  }
501  }
502  });
503 }
504 
505 //
506 // restriction
507 //
508 
509 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
510 void mlndlap_restriction (int i, int j, int k, Array4<Real> const& crse,
511  Array4<Real const> const& fine, Array4<int const> const& msk) noexcept
512 {
513  int ii = i*2;
514  int jj = j*2;
515  int kk = k*2;
516  if (msk(ii,jj,kk)) {
517  crse(i,j,k) = Real(0.0);
518  } else {
519  crse(i,j,k) = Real(1./64.)*(fine(ii-1,jj-1,kk-1)+fine(ii+1,jj-1,kk-1)
520  +fine(ii-1,jj+1,kk-1)+fine(ii+1,jj+1,kk-1)
521  +fine(ii-1,jj-1,kk+1)+fine(ii+1,jj-1,kk+1)
522  +fine(ii-1,jj+1,kk+1)+fine(ii+1,jj+1,kk+1))
523  + Real(1./32.)*(fine(ii ,jj-1,kk-1)+fine(ii ,jj+1,kk-1)
524  +fine(ii ,jj-1,kk+1)+fine(ii ,jj+1,kk+1)
525  +fine(ii-1,jj ,kk-1)+fine(ii+1,jj ,kk-1)
526  +fine(ii-1,jj ,kk+1)+fine(ii+1,jj ,kk+1)
527  +fine(ii-1,jj-1,kk )+fine(ii+1,jj-1,kk )
528  +fine(ii-1,jj+1,kk )+fine(ii+1,jj+1,kk ))
529  + Real(1./16.)*(fine(ii-1,jj,kk)+fine(ii+1,jj,kk)
530  +fine(ii,jj-1,kk)+fine(ii,jj+1,kk)
531  +fine(ii,jj,kk-1)+fine(ii,jj,kk+1))
532  + Real(1./8.)*fine(ii,jj,kk);
533  }
534 }
535 
536 template <int rr>
537 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
538 void mlndlap_restriction (int i, int j, int k, Array4<Real> const& crse,
539  Array4<Real const> const& fine, Array4<int const> const& msk,
540  Box const& fdom,
541  GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
542  GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
543 {
544  const int ii = i*rr;
545  const int jj = j*rr;
546  const int kk = k*rr;
547  if (msk(ii,jj,kk)) {
548  crse(i,j,k) = Real(0.0);
549  } else {
550  const auto ndlo = amrex::lbound(fdom);
551  const auto ndhi = amrex::ubound(fdom);
552  Real tmp = Real(0.0);
553  for (int koff = -rr+1; koff <= rr-1; ++koff) {
554  Real wz = rr - std::abs(koff);
555  for (int joff = -rr+1; joff <= rr-1; ++joff) {
556  Real wy = rr - std::abs(joff);
557  for (int ioff = -rr+1; ioff <= rr-1; ++ioff) {
558  Real wx = rr - std::abs(ioff);
559  int itmp = ii + ioff;
560  int jtmp = jj + joff;
561  int ktmp = kk + koff;
562  if ((itmp < ndlo.x && (bclo[0] == LinOpBCType::Neumann ||
563  bclo[0] == LinOpBCType::inflow)) ||
564  (itmp > ndhi.x && (bchi[0] == LinOpBCType::Neumann ||
565  bchi[0] == LinOpBCType::inflow))) {
566  itmp = ii - ioff;
567  }
568  if ((jtmp < ndlo.y && (bclo[1] == LinOpBCType::Neumann ||
569  bclo[1] == LinOpBCType::inflow)) ||
570  (jtmp > ndhi.y && (bchi[1] == LinOpBCType::Neumann ||
571  bchi[1] == LinOpBCType::inflow))) {
572  jtmp = jj - joff;
573  }
574  if ((ktmp < ndlo.z && (bclo[2] == LinOpBCType::Neumann ||
575  bclo[2] == LinOpBCType::inflow)) ||
576  (ktmp > ndhi.z && (bchi[2] == LinOpBCType::Neumann ||
577  bchi[2] == LinOpBCType::inflow))) {
578  ktmp = kk - koff;
579  }
580  tmp += wx*wy*wz*fine(itmp,jtmp,ktmp);
581  }
582  }
583  }
584  crse(i,j,k) = tmp/Real(rr*rr*rr*rr*rr*rr);
585  }
586 }
587 
588 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
589 void mlndlap_semi_restriction (int i, int j, int k, Array4<Real> const& crse,
590  Array4<Real const> const& fine, Array4<int const> const& msk, int idir) noexcept
591 {
592  if (idir == 2)
593  {
594  int ii = i*2;
595  int jj = j*2;
596  int kk = k;
597  if (msk(ii,jj,kk)) {
598  crse(i,j,k) = Real(0.0);
599  } else { // use 2-D formula
600  crse(i,j,k) = Real(1./16.)*(fine(ii-1,jj-1,kk) + Real(2.)*fine(ii ,jj-1,kk) + fine(ii+1,jj-1,kk)
601  + Real(2.)*fine(ii-1,jj ,kk) + Real(4.)*fine(ii ,jj ,kk) + Real(2.)*fine(ii+1,jj ,kk)
602  + fine(ii-1,jj+1,kk) + Real(2.)*fine(ii ,jj+1,kk) + fine(ii+1,jj+1,kk));
603  }
604  }
605  else if (idir == 1)
606  {
607  int ii = i*2;
608  int jj = j;
609  int kk = k*2;
610  if (msk(ii,jj,kk)) {
611  crse(i,j,k) = Real(0.0);
612  } else { // use 2-D formula
613  crse(i,j,k) = Real(1./16.)*(fine(ii-1,jj,kk-1) + Real(2.)*fine(ii ,jj,kk-1) + fine(ii+1,jj,kk-1)
614  + Real(2.)*fine(ii-1,jj ,kk) + Real(4.)*fine(ii ,jj,kk ) + Real(2.)*fine(ii+1,jj,kk )
615  + fine(ii-1,jj,kk+1) + Real(2.)*fine(ii ,jj,kk+1) + fine(ii+1,jj,kk+1));
616  }
617  }
618  else if (idir == 0)
619  {
620  int ii = i;
621  int jj = j*2;
622  int kk = k*2;
623  if (msk(ii,jj,kk)) {
624  crse(i,j,k) = Real(0.0);
625  } else { // use 2-D formula
626  crse(i,j,k) = Real(1./16.)*(fine(ii,jj-1,kk-1) + Real(2.)*fine(ii ,jj,kk-1) + fine(ii,jj+1,kk-1)
627  + Real(2.)*fine(ii,jj-1 ,kk) + Real(4.)*fine(ii ,jj,kk ) + Real(2.)*fine(ii,jj+1,kk )
628  + fine(ii,jj-1,kk+1) + Real(2.)*fine(ii ,jj,kk+1) + fine(ii,jj+1,kk+1));
629  }
630  }
631  else
632  {
633  amrex::Abort("mlndlap_semi_restriction semi direction wrong semi-direction. ");
634  }
635 }
636 
637 //
638 // masks
639 //
640 
641 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
642 void mlndlap_set_nodal_mask (int i, int j, int k, Array4<int> const& nmsk,
643  Array4<int const> const& cmsk) noexcept
644 {
645  using namespace nodelap_detail;
646 
647  int s = cmsk(i-1,j-1,k-1) + cmsk(i ,j-1,k-1)
648  + cmsk(i-1,j ,k-1) + cmsk(i ,j ,k-1)
649  + cmsk(i-1,j-1,k ) + cmsk(i ,j-1,k )
650  + cmsk(i-1,j ,k ) + cmsk(i ,j ,k );
651  if (s == 8*crse_cell) {
652  nmsk(i,j,k) = crse_node;
653  }
654  else if (s == 8*fine_cell) {
655  nmsk(i,j,k) = fine_node;
656  } else {
657  nmsk(i,j,k) = crse_fine_node;
658  }
659 }
660 
661 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
662 void mlndlap_set_dirichlet_mask (Box const& bx, Array4<int> const& dmsk,
663  Array4<int const> const& omsk, Box const& dom,
664  GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
665  GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
666 {
667  const auto lo = amrex::lbound(bx);
668  const auto hi = amrex::ubound(bx);
669  for (int k = lo.z; k <= hi.z; ++k) {
670  for (int j = lo.y; j <= hi.y; ++j) {
671  AMREX_PRAGMA_SIMD
672  for (int i = lo.x; i <= hi.x; ++i) {
673  if (!dmsk(i,j,k)) {
674  dmsk(i,j,k) = (omsk(i-1,j-1,k-1) == 1 || omsk(i,j-1,k-1) == 1 ||
675  omsk(i-1,j ,k-1) == 1 || omsk(i,j ,k-1) == 1 ||
676  omsk(i-1,j-1,k ) == 1 || omsk(i,j-1,k ) == 1 ||
677  omsk(i-1,j ,k ) == 1 || omsk(i,j ,k ) == 1);
678  }
679  }}}
680 
681  const auto domlo = amrex::lbound(dom);
682  const auto domhi = amrex::ubound(dom);
683 
684  if (bclo[0] == LinOpBCType::Dirichlet && lo.x == domlo.x) {
685  for (int k = lo.z; k <= hi.z; ++k) {
686  for (int j = lo.y; j <= hi.y; ++j) {
687  dmsk(lo.x,j,k) = 1;
688  }}
689  }
690 
691  if (bchi[0] == LinOpBCType::Dirichlet && hi.x == domhi.x) {
692  for (int k = lo.z; k <= hi.z; ++k) {
693  for (int j = lo.y; j <= hi.y; ++j) {
694  dmsk(hi.x,j,k) = 1;
695  }}
696  }
697 
698  if (bclo[1] == LinOpBCType::Dirichlet && lo.y == domlo.y) {
699  for (int k = lo.z; k <= hi.z; ++k) {
700  AMREX_PRAGMA_SIMD
701  for (int i = lo.x; i <= hi.x; ++i) {
702  dmsk(i,lo.y,k) = 1;
703  }}
704  }
705 
706  if (bchi[1] == LinOpBCType::Dirichlet && hi.y == domhi.y) {
707  for (int k = lo.z; k <= hi.z; ++k) {
708  AMREX_PRAGMA_SIMD
709  for (int i = lo.x; i <= hi.x; ++i) {
710  dmsk(i,hi.y,k) = 1;
711  }}
712  }
713 
714  if (bclo[2] == LinOpBCType::Dirichlet && lo.z == domlo.z) {
715  for (int j = lo.y; j <= hi.y; ++j) {
716  AMREX_PRAGMA_SIMD
717  for (int i = lo.x; i <= hi.x; ++i) {
718  dmsk(i,j,lo.z) = 1;
719  }}
720  }
721 
722  if (bchi[2] == LinOpBCType::Dirichlet && hi.z == domhi.z) {
723  for (int j = lo.y; j <= hi.y; ++j) {
724  AMREX_PRAGMA_SIMD
725  for (int i = lo.x; i <= hi.x; ++i) {
726  dmsk(i,j,hi.z) = 1;
727  }}
728  }
729 }
730 
731 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
732 void mlndlap_set_dot_mask (Box const& bx, Array4<Real> const& dmsk,
733  Array4<int const> const& omsk, Box const& dom,
734  GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
735  GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
736 {
737  const auto lo = amrex::lbound(bx);
738  const auto hi = amrex::ubound(bx);
739  for (int k = lo.z; k <= hi.z; ++k) {
740  for (int j = lo.y; j <= hi.y; ++j) {
741  AMREX_PRAGMA_SIMD
742  for (int i = lo.x; i <= hi.x; ++i) {
743  dmsk(i,j,k) = static_cast<Real>(omsk(i,j,k));
744  }}}
745 
746  const auto domlo = amrex::lbound(dom);
747  const auto domhi = amrex::ubound(dom);
748 
749  if ((bclo[0] == LinOpBCType::Neumann || bclo[0] == LinOpBCType::inflow)
750  && lo.x == domlo.x)
751  {
752  for (int k = lo.z; k <= hi.z; ++k) {
753  for (int j = lo.y; j <= hi.y; ++j) {
754  dmsk(lo.x,j,k) *= Real(0.5);
755  }}
756  }
757 
758  if ((bchi[0] == LinOpBCType::Neumann || bchi[0] == LinOpBCType::inflow)
759  && hi.x == domhi.x)
760  {
761  for (int k = lo.z; k <= hi.z; ++k) {
762  for (int j = lo.y; j <= hi.y; ++j) {
763  dmsk(hi.x,j,k) *= Real(0.5);
764  }}
765  }
766 
767  if ((bclo[1] == LinOpBCType::Neumann || bclo[1] == LinOpBCType::inflow)
768  && lo.y == domlo.y)
769  {
770  for (int k = lo.z; k <= hi.z; ++k) {
771  AMREX_PRAGMA_SIMD
772  for (int i = lo.x; i <= hi.x; ++i) {
773  dmsk(i,lo.y,k) *= Real(0.5);
774  }}
775  }
776 
777  if ((bchi[1] == LinOpBCType::Neumann || bchi[1] == LinOpBCType::inflow)
778  && hi.y == domhi.y)
779  {
780  for (int k = lo.z; k <= hi.z; ++k) {
781  AMREX_PRAGMA_SIMD
782  for (int i = lo.x; i <= hi.x; ++i) {
783  dmsk(i,hi.y,k) *= Real(0.5);
784  }}
785  }
786 
787  if ((bclo[2] == LinOpBCType::Neumann || bclo[2] == LinOpBCType::inflow)
788  && lo.z == domlo.z)
789  {
790  for (int j = lo.y; j <= hi.y; ++j) {
791  AMREX_PRAGMA_SIMD
792  for (int i = lo.x; i <= hi.x; ++i) {
793  dmsk(i,j,lo.z) *= Real(0.5);
794  }}
795  }
796 
797  if ((bchi[2] == LinOpBCType::Neumann || bchi[2] == LinOpBCType::inflow)
798  && hi.z == domhi.z)
799  {
800  for (int j = lo.y; j <= hi.y; ++j) {
801  AMREX_PRAGMA_SIMD
802  for (int i = lo.x; i <= hi.x; ++i) {
803  dmsk(i,j,hi.z) *= Real(0.5);
804  }}
805  }
806 }
807 
808 }
809 
810 #endif
AMREX_PRAGMA_SIMD
#define AMREX_PRAGMA_SIMD
Definition: AMReX_Extension.H:80
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
AMREX_HOST_DEVICE_FOR_3D
#define AMREX_HOST_DEVICE_FOR_3D(...)
Definition: AMReX_GpuLaunch.nolint.H:50
AMREX_GPU_HOST_DEVICE
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
idir
int idir
Definition: AMReX_HypreMLABecLap.cpp:1093
offset
Array4< int const > offset
Definition: AMReX_HypreMLABecLap.cpp:1089
fine
Array4< Real > fine
Definition: AMReX_InterpFaceRegister.cpp:90
crse
Array4< Real const > crse
Definition: AMReX_InterpFaceRegister.cpp:92
if
if(!(yy_init))
Definition: amrex_iparser.lex.nolint.H:935
amrex::BoxND::growHi
AMREX_GPU_HOST_DEVICE BoxND & growHi(int idir, int n_cell=1) noexcept
Grow the BoxND on the high end by n_cell cells in direction idir. NOTE: n_cell negative shrinks the B...
Definition: AMReX_Box.H:659
amrex::BoxND::growLo
AMREX_GPU_HOST_DEVICE BoxND & growLo(int idir, int n_cell=1) noexcept
Grow the BoxND on the low end by n_cell cells in direction idir. NOTE: n_cell negative shrinks the Bo...
Definition: AMReX_Box.H:648
amrex::nodelap_detail::fine_node
constexpr int fine_node
Definition: AMReX_MLNodeLinOp_K.H:58
amrex::nodelap_detail::crse_node
constexpr int crse_node
Definition: AMReX_MLNodeLinOp_K.H:56
amrex::nodelap_detail::crse_fine_node
constexpr int crse_fine_node
Definition: AMReX_MLNodeLinOp_K.H:57
amrex_eb_flux_reg_nd_module::crse_cell
integer, parameter crse_cell
Definition: AMReX_EBFluxRegister_nd.F90:7
amrex_eb_flux_reg_nd_module::fine_cell
integer, parameter fine_cell
Definition: AMReX_EBFluxRegister_nd.F90:9
amrex
Definition: AMReX_Amr.cpp:49
amrex::mlndlap_bc_doit
void mlndlap_bc_doit(Box const &vbx, Array4< T > const &a, Box const &domain, GpuArray< bool, AMREX_SPACEDIM > const &bflo, GpuArray< bool, AMREX_SPACEDIM > const &bfhi) noexcept
Definition: AMReX_MLNodeLinOp_1D_K.H:8
amrex::mlndlap_set_nodal_mask
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlndlap_set_nodal_mask(int i, int, int, Array4< int > const &nmsk, Array4< int const > const &cmsk) noexcept
Definition: AMReX_MLNodeLinOp_1D_K.H:93
amrex::Box
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
amrex::grow
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > grow(const BoxND< dim > &b, int i) noexcept
Grow BoxND in all directions by given amount.
Definition: AMReX_Box.H:1211
amrex::abs
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE T abs(const GpuComplex< T > &a_z) noexcept
Return the absolute value of a complex number.
Definition: AMReX_GpuComplex.H:356
amrex::mlndlap_restriction
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlndlap_restriction(int i, int, int, Array4< Real > const &crse, Array4< Real const > const &fine, Array4< int const > const &msk) noexcept
Definition: AMReX_MLNodeLinOp_1D_K.H:41
amrex::ubound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:315
amrex::lbound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:308
amrex::mlndlap_semi_restriction
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlndlap_semi_restriction(int, int, int, Array4< Real > const &, Array4< Real const > const &, Array4< int const > const &, int) noexcept
Definition: AMReX_MLNodeLinOp_1D_K.H:86
amrex::IntVect
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
amrex::Abort
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
amrex::mlndlap_set_dirichlet_mask
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlndlap_set_dirichlet_mask(Box const &bx, Array4< int > const &dmsk, Array4< int const > const &omsk, Box const &dom, GpuArray< LinOpBCType, AMREX_SPACEDIM > const &bclo, GpuArray< LinOpBCType, AMREX_SPACEDIM > const &bchi) noexcept
Definition: AMReX_MLNodeLinOp_1D_K.H:109
amrex::mlndlap_set_dot_mask
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlndlap_set_dot_mask(Box const &bx, Array4< Real > const &dmsk, Array4< int const > const &omsk, Box const &dom, GpuArray< LinOpBCType, AMREX_SPACEDIM > const &bclo, GpuArray< LinOpBCType, AMREX_SPACEDIM > const &bchi) noexcept
Definition: AMReX_MLNodeLinOp_1D_K.H:136