Block-Structured AMR Software Framework
 
Loading...
Searching...
No Matches
AMReX_MLNodeLinOp_2D_K.H
Go to the documentation of this file.
1#ifndef AMREX_ML_NODE_LINOP_2D_K_H_
2#define AMREX_ML_NODE_LINOP_2D_K_H_
3#include <AMReX_Config.H>
4
5namespace amrex {
6
7template <typename T>
8void 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))) {
29 // xlo & ylo
30 if (i == dlo.x-1 && j == dlo.y-1 && (bflo[0] || bflo[1]))
31 {
32 if (bflo[0] && bflo[1])
33 {
34 a(i,j,k) = a(i+1+offset, j+1+offset, k);
35 }
36 else if (bflo[0])
37 {
38 a(i,j,k) = a(i+1+offset, j, k);
39 }
40 else if (bflo[1])
41 {
42 a(i,j,k) = a(i, j+1+offset, k);
43 }
44 }
45 // xhi & ylo
46 else if (i == dhi.x+1 && j == dlo.y-1 && (bfhi[0] || bflo[1]))
47 {
48 if (bfhi[0] && bflo[1])
49 {
50 a(i,j,k) = a(i-1-offset, j+1+offset, k);
51 }
52 else if (bfhi[0])
53 {
54 a(i,j,k) = a(i-1-offset, j, k);
55 }
56 else if (bflo[1])
57 {
58 a(i,j,k) = a(i, j+1+offset, k);
59 }
60 }
61 // xlo & yhi
62 else if (i == dlo.x-1 && j == dhi.y+1 && (bflo[0] || bfhi[1]))
63 {
64 if (bflo[0] && bfhi[1])
65 {
66 a(i,j,k) = a(i+1+offset, j-1-offset, k);
67 }
68 else if (bflo[0])
69 {
70 a(i,j,k) = a(i+1+offset, j, k);
71 }
72 else if (bfhi[1])
73 {
74 a(i,j,k) = a(i, j-1-offset, k);
75 }
76 }
77 // xhi & yhi
78 else if (i == dhi.x+1 && j == dhi.y+1 && (bfhi[0] || bfhi[1]))
79 {
80 if (bfhi[0] && bfhi[1])
81 {
82 a(i,j,k) = a(i-1-offset, j-1-offset, k);
83 }
84 else if (bfhi[0])
85 {
86 a(i,j,k) = a(i-1-offset, j, k);
87 }
88 else if (bfhi[1])
89 {
90 a(i,j,k) = a(i, j-1-offset, k);
91 }
92 }
93 else if (i == dlo.x-1 && bflo[0])
94 {
95 a(i,j,k) = a(i+1+offset, j, k);
96 }
97 else if (i == dhi.x+1 && bfhi[0])
98 {
99 a(i,j,k) = a(i-1-offset, j, k);
100 }
101 else if (j == dlo.y-1 && bflo[1])
102 {
103 a(i,j,k) = a(i, j+1+offset, k);
104 }
105 else if (j == dhi.y+1 && bfhi[1])
106 {
107 a(i,j,k) = a(i, j-1-offset, k);
108 }
109 }
110 });
111}
112
114void mlndlap_restriction (int i, int j, int k, Array4<Real> const& crse,
115 Array4<Real const> const& fine, Array4<int const> const& msk) noexcept
116{
117 int ii = i*2;
118 int jj = j*2;
119 int kk = 0;
120 if (msk(ii,jj,kk)) {
121 crse(i,j,k) = Real(0.0);
122 } else {
123 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)
124 + Real(2.)*fine(ii-1,jj ,kk) + Real(4.)*fine(ii ,jj ,kk) + Real(2.)*fine(ii+1,jj ,kk)
125 + fine(ii-1,jj+1,kk) + Real(2.)*fine(ii ,jj+1,kk) + fine(ii+1,jj+1,kk));
126 }
127}
128
129template <int rr>
131void mlndlap_restriction (int i, int j, int k, Array4<Real> const& crse,
132 Array4<Real const> const& fine, Array4<int const> const& msk,
133 Box const& fdom,
134 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
135 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
136{
137 const int ii = i*rr;
138 const int jj = j*rr;
139 if (msk(ii,jj,0)) {
140 crse(i,j,k) = Real(0.0);
141 } else {
142 const auto ndlo = amrex::lbound(fdom);
143 const auto ndhi = amrex::ubound(fdom);
144 Real tmp = Real(0.0);
145 for (int joff = -rr+1; joff <= rr-1; ++joff) {
146 Real wy = rr - std::abs(joff);
147 for (int ioff = -rr+1; ioff <= rr-1; ++ioff) {
148 Real wx = rr - std::abs(ioff);
149 int itmp = ii + ioff;
150 int jtmp = jj + joff;
151 if ((itmp < ndlo.x && (bclo[0] == LinOpBCType::Neumann ||
152 bclo[0] == LinOpBCType::inflow)) ||
153 (itmp > ndhi.x && (bchi[0] == LinOpBCType::Neumann ||
154 bchi[0] == LinOpBCType::inflow))) {
155 itmp = ii - ioff;
156 }
157 if ((jtmp < ndlo.y && (bclo[1] == LinOpBCType::Neumann ||
158 bclo[1] == LinOpBCType::inflow)) ||
159 (jtmp > ndhi.y && (bchi[1] == LinOpBCType::Neumann ||
160 bchi[1] == LinOpBCType::inflow))) {
161 jtmp = jj - joff;
162 }
163 tmp += wx*wy*fine(itmp,jtmp,0);
164 }
165 }
166 crse(i,j,k) = tmp*(Real(1.0)/Real(rr*rr*rr*rr));
167 }
168}
169
171void mlndlap_semi_restriction (int i, int j, int k, Array4<Real> const& crse,
172 Array4<Real const> const& fine, Array4<int const> const& msk, int idir) noexcept
173{
174 int kk = 0;
175 if (idir == 1) {
176 int ii = i*2;
177 int jj = j;
178 if (msk(ii,jj,kk)) {
179 crse(i,j,k) = Real(0.0);
180 } else {
181 crse(i,j,k) = Real(1./4.)*(fine(ii-1,jj,kk) + Real(2.)*fine(ii,jj,kk) + fine(ii+1,jj,kk));
182 }
183 } else if (idir == 0) {
184 int ii = i;
185 int jj = j*2;
186 if (msk(ii,jj,kk)) {
187 crse(i,j,k) = Real(0.0);
188 } else {
189 crse(i,j,k) = Real(1./4.)*(fine(ii,jj-1,kk) + Real(2.)*fine(ii,jj,kk) + fine(ii,jj+1,kk));
190 }
191 }
192}
193
194//
195// masks
196//
197
199void mlndlap_set_nodal_mask (int i, int j, int k, Array4<int> const& nmsk,
200 Array4<int const> const& cmsk) noexcept
201{
202 using namespace nodelap_detail;
203
204 int s = cmsk(i-1,j-1,k) + cmsk(i ,j-1,k)
205 + cmsk(i-1,j ,k) + cmsk(i ,j ,k);
206 if (s == 4*crse_cell) {
207 nmsk(i,j,k) = crse_node;
208 }
209 else if (s == 4*fine_cell) {
210 nmsk(i,j,k) = fine_node;
211 } else {
212 nmsk(i,j,k) = crse_fine_node;
213 }
214}
215
217void mlndlap_set_dirichlet_mask (Box const& bx, Array4<int> const& dmsk,
218 Array4<int const> const& omsk, Box const& dom,
219 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
220 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
221{
222 const auto lo = amrex::lbound(bx);
223 const auto hi = amrex::ubound(bx);
224 for (int j = lo.y; j <= hi.y; ++j) {
226 for (int i = lo.x; i <= hi.x; ++i) {
227 if (!dmsk(i,j,0)) {
228 dmsk(i,j,0) = (omsk(i-1,j-1,0) == 1 || omsk(i,j-1,0) == 1 ||
229 omsk(i-1,j ,0) == 1 || omsk(i,j ,0) == 1);
230 }
231 }}
232
233 const auto domlo = amrex::lbound(dom);
234 const auto domhi = amrex::ubound(dom);
235
236 if (bclo[0] == LinOpBCType::Dirichlet && lo.x == domlo.x) {
237 for (int j = lo.y; j <= hi.y; ++j) {
238 dmsk(lo.x,j,0) = 1;
239 }
240 }
241
242 if (bchi[0] == LinOpBCType::Dirichlet && hi.x == domhi.x) {
243 for (int j = lo.y; j <= hi.y; ++j) {
244 dmsk(hi.x,j,0) = 1;
245 }
246 }
247
248 if (bclo[1] == LinOpBCType::Dirichlet && lo.y == domlo.y) {
250 for (int i = lo.x; i <= hi.x; ++i) {
251 dmsk(i,lo.y,0) = 1;
252 }
253 }
254
255 if (bchi[1] == LinOpBCType::Dirichlet && hi.y == domhi.y) {
257 for (int i = lo.x; i <= hi.x; ++i) {
258 dmsk(i,hi.y,0) = 1;
259 }
260 }
261}
262
264void mlndlap_set_dot_mask (Box const& bx, Array4<Real> const& dmsk,
265 Array4<int const> const& omsk, Box const& dom,
266 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
267 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
268{
269 const auto lo = amrex::lbound(bx);
270 const auto hi = amrex::ubound(bx);
271 for (int j = lo.y; j <= hi.y; ++j) {
273 for (int i = lo.x; i <= hi.x; ++i) {
274 dmsk(i,j,0) = static_cast<Real>(omsk(i,j,0));
275 }}
276
277 const auto domlo = amrex::lbound(dom);
278 const auto domhi = amrex::ubound(dom);
279
280 if ((bclo[0] == LinOpBCType::Neumann || bclo[0] == LinOpBCType::inflow)
281 && lo.x == domlo.x)
282 {
283 for (int j = lo.y; j <= hi.y; ++j) {
284 dmsk(lo.x,j,0) *= Real(0.5);
285 }
286 }
287
288 if ((bchi[0] == LinOpBCType::Neumann || bchi[0] == LinOpBCType::inflow)
289 && hi.x == domhi.x)
290 {
291 for (int j = lo.y; j <= hi.y; ++j) {
292 dmsk(hi.x,j,0) *= Real(0.5);
293 }
294 }
295
296 if ((bclo[1] == LinOpBCType::Neumann || bclo[1] == LinOpBCType::inflow)
297 && lo.y == domlo.y)
298 {
300 for (int i = lo.x; i <= hi.x; ++i) {
301 dmsk(i,lo.y,0) *= Real(0.5);
302 }
303 }
304
305 if ((bchi[1] == LinOpBCType::Neumann || bchi[1] == LinOpBCType::inflow)
306 && hi.y == domhi.y)
307 {
309 for (int i = lo.x; i <= hi.x; ++i) {
310 dmsk(i,hi.y,0) *= Real(0.5);
311 }
312 }
313}
314
316void mlndlap_set_dot_mask (Box const& bx, Array4<Real> const& dmsk,
317 Array4<int const> const& omsk,
318 Array4<int const> const& fmsk, Box const& dom,
319 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bclo,
320 GpuArray<LinOpBCType, AMREX_SPACEDIM> const& bchi) noexcept
321{
322 const auto lo = amrex::lbound(bx);
323 const auto hi = amrex::ubound(bx);
324 for (int j = lo.y; j <= hi.y; ++j) {
326 for (int i = lo.x; i <= hi.x; ++i) {
327 if (fmsk(i,j,0) == 0) {
328 dmsk(i,j,0) = static_cast<Real>(omsk(i,j,0));
329 } else {
330 dmsk(i,j,0) = Real(0);
331 }
332 }}
333
334 const auto domlo = amrex::lbound(dom);
335 const auto domhi = amrex::ubound(dom);
336
337 if ((bclo[0] == LinOpBCType::Neumann || bclo[0] == LinOpBCType::inflow)
338 && lo.x == domlo.x)
339 {
340 for (int j = lo.y; j <= hi.y; ++j) {
341 dmsk(lo.x,j,0) *= Real(0.5);
342 }
343 }
344
345 if ((bchi[0] == LinOpBCType::Neumann || bchi[0] == LinOpBCType::inflow)
346 && hi.x == domhi.x)
347 {
348 for (int j = lo.y; j <= hi.y; ++j) {
349 dmsk(hi.x,j,0) *= Real(0.5);
350 }
351 }
352
353 if ((bclo[1] == LinOpBCType::Neumann || bclo[1] == LinOpBCType::inflow)
354 && lo.y == domlo.y)
355 {
357 for (int i = lo.x; i <= hi.x; ++i) {
358 dmsk(i,lo.y,0) *= Real(0.5);
359 }
360 }
361
362 if ((bchi[1] == LinOpBCType::Neumann || bchi[1] == LinOpBCType::inflow)
363 && hi.y == domhi.y)
364 {
366 for (int i = lo.x; i <= hi.x; ++i) {
367 dmsk(i,hi.y,0) *= Real(0.5);
368 }
369 }
370}
371
372}
373
374#endif
#define AMREX_PRAGMA_SIMD
Definition AMReX_Extension.H:80
#define AMREX_FORCE_INLINE
Definition AMReX_Extension.H:119
#define AMREX_HOST_DEVICE_FOR_3D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:106
#define AMREX_GPU_HOST_DEVICE
Definition AMReX_GpuQualifiers.H:20
int idir
Definition AMReX_HypreMLABecLap.cpp:1093
Array4< int const > offset
Definition AMReX_HypreMLABecLap.cpp:1089
Array4< Real > fine
Definition AMReX_InterpFaceRegister.cpp:90
Array4< Real const > crse
Definition AMReX_InterpFaceRegister.cpp:92
if(!(yy_init))
Definition amrex_iparser.lex.nolint.H:935
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_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_GPU_HOST_DEVICE bool cellCentered() const noexcept
Returns true if BoxND is cell-centered in all indexing directions.
Definition AMReX_Box.H:319
constexpr int fine_node
Definition AMReX_MLNodeLinOp_K.H:58
constexpr int crse_node
Definition AMReX_MLNodeLinOp_K.H:56
constexpr int crse_fine_node
Definition AMReX_MLNodeLinOp_K.H:57
Definition AMReX_Amr.cpp:49
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_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
BoxND< AMREX_SPACEDIM > Box
Definition AMReX_BaseFwd.H:27
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_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition AMReX_Array4.H:315
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition AMReX_Array4.H:308
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
IntVectND< AMREX_SPACEDIM > IntVect
Definition AMReX_BaseFwd.H:30
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_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_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