Block-Structured AMR Software Framework
 
Loading...
Searching...
No Matches
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
5namespace amrex {
6
7template <typename T>
8inline 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
509AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
510void 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
536template <int rr>
537AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
538void 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
588AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
589void 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
641AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
642void 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
661AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
662void 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
731AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
732void 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
808AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
809void mlndlap_set_dot_mask (Box const& bx, Array4<Real> const& dmsk,
810 Array4<int const> const& omsk,
811 Array4<int const> const& fmsk, Box const& dom,
812 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
813 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
814{
815 const auto lo = amrex::lbound(bx);
816 const auto hi = amrex::ubound(bx);
817 for (int k = lo.z; k <= hi.z; ++k) {
818 for (int j = lo.y; j <= hi.y; ++j) {
819 AMREX_PRAGMA_SIMD
820 for (int i = lo.x; i <= hi.x; ++i) {
821 if (fmsk(i,j,k) == 0) {
822 dmsk(i,j,k) = static_cast<Real>(omsk(i,j,k));
823 } else {
824 dmsk(i,j,k) = Real(0);
825 }
826 }}}
827
828 const auto domlo = amrex::lbound(dom);
829 const auto domhi = amrex::ubound(dom);
830
831 if ((bclo[0] == LinOpBCType::Neumann || bclo[0] == LinOpBCType::inflow)
832 && lo.x == domlo.x)
833 {
834 for (int k = lo.z; k <= hi.z; ++k) {
835 for (int j = lo.y; j <= hi.y; ++j) {
836 dmsk(lo.x,j,k) *= Real(0.5);
837 }}
838 }
839
840 if ((bchi[0] == LinOpBCType::Neumann || bchi[0] == LinOpBCType::inflow)
841 && hi.x == domhi.x)
842 {
843 for (int k = lo.z; k <= hi.z; ++k) {
844 for (int j = lo.y; j <= hi.y; ++j) {
845 dmsk(hi.x,j,k) *= Real(0.5);
846 }}
847 }
848
849 if ((bclo[1] == LinOpBCType::Neumann || bclo[1] == LinOpBCType::inflow)
850 && lo.y == domlo.y)
851 {
852 for (int k = lo.z; k <= hi.z; ++k) {
853 AMREX_PRAGMA_SIMD
854 for (int i = lo.x; i <= hi.x; ++i) {
855 dmsk(i,lo.y,k) *= Real(0.5);
856 }}
857 }
858
859 if ((bchi[1] == LinOpBCType::Neumann || bchi[1] == LinOpBCType::inflow)
860 && hi.y == domhi.y)
861 {
862 for (int k = lo.z; k <= hi.z; ++k) {
863 AMREX_PRAGMA_SIMD
864 for (int i = lo.x; i <= hi.x; ++i) {
865 dmsk(i,hi.y,k) *= Real(0.5);
866 }}
867 }
868
869 if ((bclo[2] == LinOpBCType::Neumann || bclo[2] == LinOpBCType::inflow)
870 && lo.z == domlo.z)
871 {
872 for (int j = lo.y; j <= hi.y; ++j) {
873 AMREX_PRAGMA_SIMD
874 for (int i = lo.x; i <= hi.x; ++i) {
875 dmsk(i,j,lo.z) *= Real(0.5);
876 }}
877 }
878
879 if ((bchi[2] == LinOpBCType::Neumann || bchi[2] == LinOpBCType::inflow)
880 && hi.z == domhi.z)
881 {
882 for (int j = lo.y; j <= hi.y; ++j) {
883 AMREX_PRAGMA_SIMD
884 for (int i = lo.x; i <= hi.x; ++i) {
885 dmsk(i,j,hi.z) *= Real(0.5);
886 }}
887 }
888}
889
890}
891
892#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_GpuLaunchMacrosC.nolint.H:106
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::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::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::cellCentered
AMREX_GPU_HOST_DEVICE bool cellCentered() const noexcept
Returns true if BoxND is cell-centered in all indexing directions.
Definition AMReX_Box.H:319
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
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::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::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::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