Block-Structured AMR Software Framework
AMReX_MLEBABecLap_K.H
Go to the documentation of this file.
1 #ifndef AMREX_MLEBABECLAP_K_H_
2 #define AMREX_MLEBABECLAP_K_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_MLLinOp_K.H>
6 
7 #include <AMReX_EBCellFlag.H>
8 
9 #if (AMREX_SPACEDIM == 2)
10 #include <AMReX_MLEBABecLap_2D_K.H>
11 #else
12 #include <AMReX_MLEBABecLap_3D_K.H>
13 #endif
14 
15 namespace amrex {
16 
17 // note that the mask in these functions is different from masks in bndry registers
18 // 1 means valid data, 0 means invalid data
19 
20 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
21 void mlebabeclap_apply_bc_x (int side, Box const& box, int blen,
22  Array4<Real> const& phi,
23  Array4<int const> const& mask,
24  Array4<Real const> const& area,
25  BoundCond bct, Real bcl,
26  Array4<Real const> const& bcval,
27  int maxorder, Real dxinv, int inhomog, int icomp) noexcept
28 {
29  const auto lo = amrex::lbound(box);
30  const auto hi = amrex::ubound(box);
31  const int i = lo.x; // boundary cell
32  const int s = 1-2*side; // +1 for lo and -1 for hi
33  switch (bct) {
34  case AMREX_LO_NEUMANN:
35  {
36  for (int k = lo.z; k <= hi.z; ++k) {
37  for (int j = lo.y; j <= hi.y; ++j) {
38  if (mask(i,j,k) == 0 && mask(i+s,j,k) == 1) {
39  phi(i,j,k,icomp) = phi(i+s,j,k,icomp);
40  }
41  }
42  }
43  break;
44  }
45  case AMREX_LO_REFLECT_ODD:
46  {
47  for (int k = lo.z; k <= hi.z; ++k) {
48  for (int j = lo.y; j <= hi.y; ++j) {
49  if (mask(i,j,k) == 0 && mask(i+s,j,k) == 1) {
50  phi(i,j,k,icomp) = -phi(i+s,j,k,icomp);
51  }
52  }
53  }
54  break;
55  }
56  case AMREX_LO_DIRICHLET:
57  {
58  const int NX = amrex::min(blen+1, maxorder);
59  GpuArray<Real,4> x{-bcl * dxinv, Real(0.5), Real(1.5), Real(2.5)};
60  Array2D<Real, 0, 3, 0, 2> coef{};
61  for (int r = 0; r <= maxorder-2; ++r) {
62  poly_interp_coeff(-Real(0.5), x.data(), r+2, &(coef(0,r)));
63  }
64  for (int k = lo.z; k <= hi.z; ++k) {
65  for (int j = lo.y; j <= hi.y; ++j) {
66  if (mask(i,j,k) == 0 && mask(i+s,j,k) == 1) {
67  int order = 1;
68  bool has_cutfaces = false;
69  for (int r = 0; r <= NX-2; ++r) {
70  Real a = area(i+(1-side)+s*r,j,k);
71  if (a > Real(0.0)) {
72  ++order;
73  if (a < Real(1.0)) {
74  has_cutfaces = true;
75  }
76  } else {
77  break;
78  }
79  }
80  if (has_cutfaces) { order = amrex::min(2,order); }
81  if (order == 1) {
82  if (inhomog) {
83  phi(i,j,k,icomp) = bcval(i,j,k,icomp);
84  } else {
85  phi(i,j,k,icomp) = Real(0.0);
86  }
87  } else {
88  Real tmp = Real(0.0);
89  for (int m = 1; m < order; ++m) {
90  tmp += phi(i+m*s,j,k,icomp) * coef(m,order-2);
91  }
92  phi(i,j,k,icomp) = tmp;
93  if (inhomog) {
94  phi(i,j,k,icomp) += bcval(i,j,k,icomp)*coef(0,order-2);
95  }
96  }
97  }
98  }
99  }
100  break;
101  }
102  default: {}
103  }
104 }
105 
106 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
107 void mlebabeclap_apply_bc_y (int side, Box const& box, int blen,
108  Array4<Real> const& phi,
109  Array4<int const> const& mask,
110  Array4<Real const> const& area,
111  BoundCond bct, Real bcl,
112  Array4<Real const> const& bcval,
113  int maxorder, Real dyinv, int inhomog, int icomp) noexcept
114 {
115  const auto lo = amrex::lbound(box);
116  const auto hi = amrex::ubound(box);
117  const int j = lo.y; // boundary cell
118  const int s = 1-2*side; // +1 for lo and -1 for hi
119  switch (bct) {
120  case AMREX_LO_NEUMANN:
121  {
122  for (int k = lo.z; k <= hi.z; ++k) {
123  for (int i = lo.x; i <= hi.x; ++i) {
124  if (mask(i,j,k) == 0 && mask(i,j+s,k) == 1) {
125  phi(i,j,k,icomp) = phi(i,j+s,k,icomp);
126  }
127  }
128  }
129  break;
130  }
131  case AMREX_LO_REFLECT_ODD:
132  {
133  for (int k = lo.z; k <= hi.z; ++k) {
134  for (int i = lo.x; i <= hi.x; ++i) {
135  if (mask(i,j,k) == 0 && mask(i,j+s,k) == 1) {
136  phi(i,j,k,icomp) = -phi(i,j+s,k,icomp);
137  }
138  }
139  }
140  break;
141  }
142  case AMREX_LO_DIRICHLET:
143  {
144  const int NX = amrex::min(blen+1, maxorder);
145  GpuArray<Real,4> x{-bcl * dyinv, Real(0.5), Real(1.5), Real(2.5)};
146  Array2D<Real, 0, 3, 0, 2> coef{};
147  for (int r = 0; r <= maxorder-2; ++r) {
148  poly_interp_coeff(-Real(0.5), x.data(), r+2, &(coef(0,r)));
149  }
150  for (int k = lo.z; k <= hi.z; ++k) {
151  for (int i = lo.x; i <= hi.x; ++i) {
152  if (mask(i,j,k) == 0 && mask(i,j+s,k) == 1) {
153  int order = 1;
154  bool has_cutfaces = false;
155  for (int r = 0; r <= NX-2; ++r) {
156  Real a = area(i,j+(1-side)+s*r,k);
157  if (a > Real(0.0)) {
158  ++order;
159  if (a < Real(1.0)) {
160  has_cutfaces = true;
161  }
162  } else {
163  break;
164  }
165  }
166  if (has_cutfaces) { order = amrex::min(2,order); }
167  if (order == 1) {
168  if (inhomog) {
169  phi(i,j,k,icomp) = bcval(i,j,k,icomp);
170  } else {
171  phi(i,j,k,icomp) = Real(0.0);
172  }
173  } else {
174  Real tmp = Real(0.0);
175  for (int m = 1; m < order; ++m) {
176  tmp += phi(i,j+m*s,k,icomp) * coef(m,order-2);
177  }
178  phi(i,j,k,icomp) = tmp;
179  if (inhomog) {
180  phi(i,j,k,icomp) += bcval(i,j,k,icomp)*coef(0,order-2);
181  }
182  }
183  }
184  }
185  }
186  break;
187  }
188  default: {}
189  }
190 }
191 
192 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
193 void mlebabeclap_apply_bc_z (int side, Box const& box, int blen,
194  Array4<Real> const& phi,
195  Array4<int const> const& mask,
196  Array4<Real const> const& area,
197  BoundCond bct, Real bcl,
198  Array4<Real const> const& bcval,
199  int maxorder, Real dzinv, int inhomog, int icomp) noexcept
200 {
201  const auto lo = amrex::lbound(box);
202  const auto hi = amrex::ubound(box);
203  const int k = lo.z; // boundary cell
204  const int s = 1-2*side; // +1 for lo and -1 for hi
205  switch (bct) {
206  case AMREX_LO_NEUMANN:
207  {
208  for (int j = lo.y; j <= hi.y; ++j) {
209  for (int i = lo.x; i <= hi.x; ++i) {
210  if (mask(i,j,k) == 0 && mask(i,j,k+s) == 1) {
211  phi(i,j,k,icomp) = phi(i,j,k+s,icomp);
212  }
213  }
214  }
215  break;
216  }
217  case AMREX_LO_REFLECT_ODD:
218  {
219  for (int j = lo.y; j <= hi.y; ++j) {
220  for (int i = lo.x; i <= hi.x; ++i) {
221  if (mask(i,j,k) == 0 && mask(i,j,k+s) == 1) {
222  phi(i,j,k,icomp) = -phi(i,j,k+s,icomp);
223  }
224  }
225  }
226  break;
227  }
228  case AMREX_LO_DIRICHLET:
229  {
230  const int NX = amrex::min(blen+1, maxorder);
231  GpuArray<Real,4> x{-bcl * dzinv, Real(0.5), Real(1.5), Real(2.5)};
232  Array2D<Real, 0, 3, 0, 2> coef{};
233  for (int r = 0; r <= maxorder-2; ++r) {
234  poly_interp_coeff(-Real(0.5), x.data(), r+2, &(coef(0,r)));
235  }
236  for (int j = lo.y; j <= hi.y; ++j) {
237  for (int i = lo.x; i <= hi.x; ++i) {
238  if (mask(i,j,k) == 0 && mask(i,j,k+s) == 1) {
239  int order = 1;
240  bool has_cutfaces = false;
241  for (int r = 0; r <= NX-2; ++r) {
242  Real a = area(i,j,k+(1-side)+s*r);
243  if (a > Real(0.0)) {
244  ++order;
245  if (a < Real(1.0)) {
246  has_cutfaces = true;
247  }
248  } else {
249  break;
250  }
251  }
252  if (has_cutfaces) { order = amrex::min(2,order); }
253  if (order == 1) {
254  if (inhomog) {
255  phi(i,j,k,icomp) = bcval(i,j,k,icomp);
256  } else {
257  phi(i,j,k,icomp) = Real(0.0);
258  }
259  } else {
260  Real tmp = Real(0.0);
261  for (int m = 1; m < order; ++m) {
262  tmp += phi(i,j,k+m*s,icomp) * coef(m,order-2);
263  }
264  phi(i,j,k,icomp) = tmp;
265  if (inhomog) {
266  phi(i,j,k,icomp) += bcval(i,j,k,icomp)*coef(0,order-2);
267  }
268  }
269  }
270  }
271  }
272  break;
273  }
274  default: {}
275  }
276 }
277 
278 }
279 
280 #endif
AMReX_EBCellFlag.H
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
AMREX_GPU_HOST_DEVICE
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
mask
Array4< int const > mask
Definition: AMReX_InterpFaceRegister.cpp:93
AMREX_LO_NEUMANN
#define AMREX_LO_NEUMANN
Definition: AMReX_LO_BCTYPES.H:6
AMREX_LO_DIRICHLET
#define AMREX_LO_DIRICHLET
Definition: AMReX_LO_BCTYPES.H:5
AMREX_LO_REFLECT_ODD
#define AMREX_LO_REFLECT_ODD
Definition: AMReX_LO_BCTYPES.H:7
AMReX_MLEBABecLap_2D_K.H
AMReX_MLEBABecLap_3D_K.H
AMReX_MLLinOp_K.H
amrex::BoundCond
Maintain an identifier for boundary condition types.
Definition: AMReX_BoundCond.H:20
amrex
Definition: AMReX_Amr.cpp:49
amrex::mlebabeclap_apply_bc_z
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_apply_bc_z(int side, Box const &box, int blen, Array4< Real > const &phi, Array4< int const > const &mask, Array4< Real const > const &area, BoundCond bct, Real bcl, Array4< Real const > const &bcval, int maxorder, Real dzinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLEBABecLap_K.H:193
amrex::mlebabeclap_apply_bc_y
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_apply_bc_y(int side, Box const &box, int blen, Array4< Real > const &phi, Array4< int const > const &mask, Array4< Real const > const &area, BoundCond bct, Real bcl, Array4< Real const > const &bcval, int maxorder, Real dyinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLEBABecLap_K.H:107
amrex::min
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & min(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:21
amrex::poly_interp_coeff
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void poly_interp_coeff(T xInt, T const *AMREX_RESTRICT x, int N, T *AMREX_RESTRICT c) noexcept
Definition: AMReX_LOUtil_K.H:24
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::mlebabeclap_apply_bc_x
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_apply_bc_x(int side, Box const &box, int blen, Array4< Real > const &phi, Array4< int const > const &mask, Array4< Real const > const &area, BoundCond bct, Real bcl, Array4< Real const > const &bcval, int maxorder, Real dxinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLEBABecLap_K.H:21
amrex::Array2D
Definition: AMReX_Array.H:282
amrex::Array4
Definition: AMReX_Array4.H:61
amrex::GpuArray
Definition: AMReX_Array.H:34