Block-Structured AMR Software Framework
AMReX_MLEBABecLap_3D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_MLEBABECLAP_3D_K_H_
2 #define AMREX_MLEBABECLAP_3D_K_H_
3 #include <AMReX_Config.H>
4 #include <AMReX_REAL.H>
5 
6 #include <AMReX_EB_LeastSquares_3D_K.H>
7 
8 namespace amrex {
9 
10 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
11 void mlebabeclap_adotx_centroid (Box const& box, Array4<Real> const& y,
12  Array4<Real const> const& x, Array4<Real const> const& a,
13  Array4<Real const> const& bX, Array4<Real const> const& bY,
14  Array4<Real const> const& bZ,
15  Array4<EBCellFlag const> const& flag,
16  Array4<Real const> const& vfrc, Array4<Real const> const& apx,
17  Array4<Real const> const& apy, Array4<Real const> const& apz,
18  Array4<Real const> const& fcx, Array4<Real const> const& fcy,
19  Array4<Real const> const& fcz,
20  Array4<Real const> const& ccent, Array4<Real const> const& ba,
21  Array4<Real const> const& bcent, Array4<Real const> const& beb,
22  Array4<Real const> const& phieb,
23  const int& domlo_x, const int& domlo_y, const int& domlo_z,
24  const int& domhi_x, const int& domhi_y, const int& domhi_z,
25  const bool& on_x_face, const bool& on_y_face, const bool& on_z_face,
26  bool is_eb_dirichlet, bool is_eb_inhomog,
27  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
28  Real alpha, Real beta, int ncomp) noexcept
29 {
30  Real dhx = beta*dxinv[0]*dxinv[0];
31  Real dhy = beta*dxinv[1]*dxinv[1];
32  Real dhz = beta*dxinv[2]*dxinv[2];
33 
34  amrex::Loop(box, ncomp, [=] (int i, int j, int k, int n) noexcept
35  {
36  if (flag(i,j,k).isCovered())
37  {
38  y(i,j,k,n) = Real(0.0);
39  }
40  else if (flag(i,j,k).isRegular() &&
41  ((flag(i-1,j ,k ).isRegular() && flag(i+1,j ,k ).isRegular() &&
42  flag(i ,j-1,k ).isRegular() && flag(i ,j+1,k ).isRegular() &&
43  flag(i ,j ,k-1).isRegular() && flag(i ,j ,k+1).isRegular()) ))
44  {
45  y(i,j,k,n) = alpha*a(i,j,k)*x(i,j,k,n)
46  - dhx * (bX(i+1,j,k,n)*(x(i+1,j,k,n) - x(i ,j,k,n))
47  -bX(i ,j,k,n)*(x(i ,j,k,n) - x(i-1,j,k,n)))
48  - dhy * (bY(i,j+1,k,n)*(x(i,j+1,k,n) - x(i,j ,k,n))
49  -bY(i,j ,k,n)*(x(i,j ,k,n) - x(i,j-1,k,n)))
50  - dhz * (bZ(i,j,k+1,n)*(x(i,j,k+1,n) - x(i,j,k ,n))
51  -bZ(i,j,k ,n)*(x(i,j,k ,n) - x(i,j,k-1,n)));
52  }
53  else
54  {
55  Real kappa = vfrc(i,j,k);
56  Real apxm = apx(i,j,k);
57  Real apxp = apx(i+1,j,k);
58  Real apym = apy(i,j,k);
59  Real apyp = apy(i,j+1,k);
60  Real apzm = apz(i,j,k);
61  Real apzp = apz(i,j,k+1);
62 
63  // First get EB-aware slope that doesn't know about extdir
64  bool needs_bdry_stencil = (i <= domlo_x) || (i >= domhi_x) ||
65  (j <= domlo_y) || (j >= domhi_y) ||
66  (k <= domlo_z) || (k >= domhi_z);
67 
68  Real fxm = bX(i,j,k,n)*(x(i,j,k,n) - x(i-1,j,k,n));
69  if ( (apxm != Real(0.0)) && (vfrc(i,j,k) != Real(1.0) || vfrc(i-1,j,k) != Real(1.0) || vfrc(i+1,j,k) != Real(1.0)) ) {
70  Real yloc_on_xface = fcx(i,j,k,0);
71  Real zloc_on_xface = fcx(i,j,k,1);
72 
73  if(needs_bdry_stencil) {
74  fxm = grad_x_of_phi_on_centroids_extdir(i,j,k,n,x,phieb,flag,ccent,bcent,vfrc,
75  yloc_on_xface,zloc_on_xface,
76  is_eb_dirichlet,is_eb_inhomog,
77  on_x_face,domlo_x,domhi_x,
78  on_y_face,domlo_y,domhi_y,
79  on_z_face,domlo_z,domhi_z);
80  } else {
81  fxm = grad_x_of_phi_on_centroids(i,j,k,n,x,phieb,flag,ccent,bcent,
82  yloc_on_xface,zloc_on_xface,is_eb_dirichlet,is_eb_inhomog);
83  }
84 
85  fxm *= bX(i,j,k,n);
86  }
87 
88  Real fxp = bX(i+1,j,k,n)*(x(i+1,j,k,n) - x(i,j,k,n));
89  if ( (apxp != Real(0.0)) && (vfrc(i,j,k) != Real(1.0) || vfrc(i+1,j,k) != Real(1.0) || vfrc(i-1,j,k) != Real(1.0)) ) {
90  Real yloc_on_xface = fcx(i+1,j,k,0);
91  Real zloc_on_xface = fcx(i+1,j,k,1);
92 
93  if(needs_bdry_stencil) {
94  fxp = grad_x_of_phi_on_centroids_extdir(i+1,j,k,n,x,phieb,flag,ccent,bcent,vfrc,
95  yloc_on_xface,zloc_on_xface,
96  is_eb_dirichlet,is_eb_inhomog,
97  on_x_face,domlo_x,domhi_x,
98  on_y_face,domlo_y,domhi_y,
99  on_z_face,domlo_z,domhi_z);
100  } else {
101  fxp = grad_x_of_phi_on_centroids(i+1,j,k,n,x,phieb,flag,ccent,bcent,
102  yloc_on_xface,zloc_on_xface,is_eb_dirichlet,is_eb_inhomog);
103  }
104 
105  fxp *= bX(i+1,j,k,n);
106  }
107 
108  Real fym = bY(i,j,k,n)*(x(i,j,k,n) - x(i,j-1,k,n));
109  if ( (apym != Real(0.0)) && (vfrc(i,j,k) != Real(1.0) || vfrc(i,j-1,k) != Real(1.0) || vfrc(i,j+1,k) != Real(1.0)) ) {
110  Real xloc_on_yface = fcy(i,j,k,0);
111  Real zloc_on_yface = fcy(i,j,k,1);
112 
113  if(needs_bdry_stencil) {
114  fym = grad_y_of_phi_on_centroids_extdir(i,j,k,n,x,phieb,flag,ccent,bcent,vfrc,
115  xloc_on_yface,zloc_on_yface,
116  is_eb_dirichlet,is_eb_inhomog,
117  on_x_face,domlo_x,domhi_x,
118  on_y_face,domlo_y,domhi_y,
119  on_z_face,domlo_z,domhi_z);
120  } else {
121  fym = grad_y_of_phi_on_centroids(i,j,k,n,x,phieb,flag,ccent,bcent,
122  xloc_on_yface,zloc_on_yface,is_eb_dirichlet,is_eb_inhomog);
123  }
124 
125  fym *= bY(i,j,k,n);
126  }
127 
128  Real fyp = bY(i,j+1,k,n)*(x(i,j+1,k,n) - x(i,j,k,n));
129  if ( (apyp != Real(0.0)) && (vfrc(i,j,k) != Real(1.0) || vfrc(i,j+1,k) != Real(1.0) || vfrc(i,j-1,k) != Real(1.0)) ) {
130  Real xloc_on_yface = fcy(i,j+1,k,0);
131  Real zloc_on_yface = fcy(i,j+1,k,1);
132 
133  if(needs_bdry_stencil) {
134  fyp = grad_y_of_phi_on_centroids_extdir(i,j+1,k,n,x,phieb,flag,ccent,bcent,vfrc,
135  xloc_on_yface,zloc_on_yface,
136  is_eb_dirichlet,is_eb_inhomog,
137  on_x_face,domlo_x,domhi_x,
138  on_y_face,domlo_y,domhi_y,
139  on_z_face,domlo_z,domhi_z);
140  } else {
141  fyp = grad_y_of_phi_on_centroids(i,j+1,k,n,x,phieb,flag,ccent,bcent,
142  xloc_on_yface,zloc_on_yface,is_eb_dirichlet,is_eb_inhomog);
143  }
144 
145  fyp *= bY(i,j+1,k,n);
146  }
147 
148  Real fzm = bZ(i,j,k,n)*(x(i,j,k,n) - x(i,j,k-1,n));
149  if ( (apzm != Real(0.0)) && (vfrc(i,j,k) != Real(1.0) || vfrc(i,j,k-1) != Real(1.0) || vfrc(i,j,k+1) != Real(1.0)) ) {
150  Real xloc_on_zface = fcz(i,j,k,0);
151  Real yloc_on_zface = fcz(i,j,k,1);
152 
153  if(needs_bdry_stencil) {
154  fzm = grad_z_of_phi_on_centroids_extdir(i,j,k,n,x,phieb,flag,ccent,bcent,vfrc,
155  xloc_on_zface,yloc_on_zface,
156  is_eb_dirichlet,is_eb_inhomog,
157  on_x_face,domlo_x,domhi_x,
158  on_y_face,domlo_y,domhi_y,
159  on_z_face,domlo_z,domhi_z);
160  } else {
161  fzm = grad_z_of_phi_on_centroids(i,j,k,n,x,phieb,flag,ccent,bcent,
162  xloc_on_zface,yloc_on_zface,is_eb_dirichlet,is_eb_inhomog);
163  }
164 
165  fzm *= bZ(i,j,k,n);
166  }
167 
168  Real fzp = bZ(i,j,k+1,n)*(x(i,j,k+1,n) - x(i,j,k,n));
169  if ( (apzp != Real(0.0)) && (vfrc(i,j,k) != Real(1.0) || vfrc(i,j,k+1) != Real(1.0) || vfrc(i,j,k-1) != Real(1.0)) ) {
170  Real xloc_on_zface = fcz(i,j,k+1,0);
171  Real yloc_on_zface = fcz(i,j,k+1,1);
172 
173  if(needs_bdry_stencil) {
174  fzp = grad_z_of_phi_on_centroids_extdir(i,j,k+1,n,x,phieb,flag,ccent,bcent,vfrc,
175  xloc_on_zface,yloc_on_zface,
176  is_eb_dirichlet,is_eb_inhomog,
177  on_x_face,domlo_x,domhi_x,
178  on_y_face,domlo_y,domhi_y,
179  on_z_face,domlo_z,domhi_z);
180  } else {
181  fzp = grad_z_of_phi_on_centroids(i,j,k+1,n,x,phieb,flag,ccent,bcent,
182  xloc_on_zface,yloc_on_zface,is_eb_dirichlet,is_eb_inhomog);
183  }
184 
185  fzp *= bZ(i,j,k+1,n);
186  }
187 
188  Real feb = Real(0.0);
189  if (is_eb_dirichlet && flag(i,j,k).isSingleValued()) {
190  Real dapx = apxm-apxp;
191  Real dapy = apym-apyp;
192  Real dapz = apzm-apzp;
193  Real anorm = std::sqrt(dapx*dapx+dapy*dapy+dapz*dapz);
194  Real anorminv = Real(1.0)/anorm;
195  Real anrmx = dapx * anorminv;
196  Real anrmy = dapy * anorminv;
197  Real anrmz = dapz * anorminv;
198 
199  feb = grad_eb_of_phi_on_centroids_extdir(i,j,k,n,x,phieb,flag,ccent,bcent,vfrc,
200  anrmx,anrmy,anrmz,is_eb_inhomog,
201  on_x_face,domlo_x,domhi_x,
202  on_y_face,domlo_y,domhi_y,
203  on_z_face,domlo_z,domhi_z);
204  feb *= ba(i,j,k) * beb(i,j,k,n);
205  }
206 
207  y(i,j,k,n) = alpha*a(i,j,k)*x(i,j,k,n) + (Real(1.0)/kappa) *
208  (dhx*(apxm*fxm - apxp*fxp) +
209  dhy*(apym*fym - apyp*fyp) +
210  dhz*(apzm*fzm - apzp*fzp) - dhx*feb);
211  }
212  });
213 }
214 
215 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
216 void mlebabeclap_adotx (Box const& box, Array4<Real> const& y,
217  Array4<Real const> const& x, Array4<Real const> const& a,
218  Array4<Real const> const& bX, Array4<Real const> const& bY,
219  Array4<Real const> const& bZ, Array4<const int> const& ccm,
220  Array4<EBCellFlag const> const& flag,
221  Array4<Real const> const& vfrc, Array4<Real const> const& apx,
222  Array4<Real const> const& apy, Array4<Real const> const& apz,
223  Array4<Real const> const& fcx, Array4<Real const> const& fcy,
224  Array4<Real const> const& fcz, Array4<Real const> const& ba,
225  Array4<Real const> const& bc, Array4<Real const> const& beb,
226  bool is_dirichlet, Array4<Real const> const& phieb,
227  bool is_inhomog, GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
228  Real alpha, Real beta, int ncomp,
229  bool beta_on_centroid, bool phi_on_centroid) noexcept
230 {
231  Real dhx = beta*dxinv[0]*dxinv[0];
232  Real dhy = beta*dxinv[1]*dxinv[1];
233  Real dhz = beta*dxinv[2]*dxinv[2];
234 
235  bool beta_on_center = !(beta_on_centroid);
236  bool phi_on_center = !( phi_on_centroid);
237 
238  amrex::Loop(box, ncomp, [=] (int i, int j, int k, int n) noexcept
239  {
240  if (flag(i,j,k).isCovered())
241  {
242  y(i,j,k,n) = Real(0.0);
243  }
244  else if (flag(i,j,k).isRegular())
245  {
246  y(i,j,k,n) = alpha*a(i,j,k)*x(i,j,k,n)
247  - dhx * (bX(i+1,j,k,n)*(x(i+1,j,k,n) - x(i ,j,k,n))
248  -bX(i ,j,k,n)*(x(i ,j,k,n) - x(i-1,j,k,n)))
249  - dhy * (bY(i,j+1,k,n)*(x(i,j+1,k,n) - x(i,j ,k,n))
250  -bY(i,j ,k,n)*(x(i,j ,k,n) - x(i,j-1,k,n)))
251  - dhz * (bZ(i,j,k+1,n)*(x(i,j,k+1,n) - x(i,j,k ,n))
252  -bZ(i,j,k ,n)*(x(i,j,k ,n) - x(i,j,k-1,n)));
253  }
254  else
255  {
256  Real kappa = vfrc(i,j,k);
257  Real apxm = apx(i,j,k);
258  Real apxp = apx(i+1,j,k);
259  Real apym = apy(i,j,k);
260  Real apyp = apy(i,j+1,k);
261  Real apzm = apz(i,j,k);
262  Real apzp = apz(i,j,k+1);
263 
264  Real fxm = bX(i,j,k,n)*(x(i,j,k,n) - x(i-1,j,k,n));
265  if (apxm != Real(0.0) && apxm != Real(1.0)) {
266  int jj = j + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,0)));
267  int kk = k + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,1)));
268  Real fracy = (ccm(i-1,jj,k) || ccm(i,jj,k)) ? std::abs(fcx(i,j,k,0)) : Real(0.0);
269  Real fracz = (ccm(i-1,j,kk) || ccm(i,j,kk)) ? std::abs(fcx(i,j,k,1)) : Real(0.0);
270  if (beta_on_center && phi_on_center)
271  {
272  fxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*fxm +
273  fracy*(Real(1.0)-fracz)*bX(i,jj,k ,n)*(x(i,jj,k ,n)-x(i-1,jj,k ,n)) +
274  fracz*(Real(1.0)-fracy)*bX(i,j ,kk,n)*(x(i,j ,kk,n)-x(i-1,j ,kk,n)) +
275  fracy* fracz *bX(i,jj,kk,n)*(x(i,jj,kk,n)-x(i-1,jj,kk,n));
276  }
277  else if (beta_on_centroid && phi_on_center)
278  {
279  fxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*(x(i, j, k,n)-x(i-1, j, k,n)) +
280  fracy *(Real(1.0)-fracz)*(x(i,jj, k,n)-x(i-1,jj, k,n)) +
281  fracz *(Real(1.0)-fracy)*(x(i, j,kk,n)-x(i-1, j,kk,n)) +
282  fracy * fracz *(x(i,jj,kk,n)-x(i-1,jj,kk,n));
283  fxm *= bX(i,j,k,n);
284  }
285  }
286 
287  Real fxp = bX(i+1,j,k,n)*(x(i+1,j,k,n) - x(i,j,k,n));
288  if (apxp != Real(0.0) && apxp != Real(1.0)) {
289  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcx(i+1,j,k,0)));
290  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcx(i+1,j,k,1)));
291  Real fracy = (ccm(i,jj,k) || ccm(i+1,jj,k)) ? std::abs(fcx(i+1,j,k,0)) : Real(0.0);
292  Real fracz = (ccm(i,j,kk) || ccm(i+1,j,kk)) ? std::abs(fcx(i+1,j,k,1)) : Real(0.0);
293  if (beta_on_center && phi_on_center)
294  {
295  fxp = (Real(1.0)-fracy)*(Real(1.0)-fracz)*fxp +
296  fracy*(Real(1.0)-fracz)*bX(i+1,jj,k ,n)*(x(i+1,jj,k ,n)-x(i,jj,k ,n)) +
297  fracz*(Real(1.0)-fracy)*bX(i+1,j ,kk,n)*(x(i+1,j ,kk,n)-x(i,j ,kk,n)) +
298  fracy* fracz *bX(i+1,jj,kk,n)*(x(i+1,jj,kk,n)-x(i,jj,kk,n));
299  }
300  else if (beta_on_centroid && phi_on_center)
301  {
302  fxp = (Real(1.0)-fracy)*(Real(1.0)-fracz)*(x(i+1, j, k,n)-x(i, j, k,n)) +
303  fracy *(Real(1.0)-fracz)*(x(i+1,jj, k,n)-x(i,jj, k,n)) +
304  fracz *(Real(1.0)-fracy)*(x(i+1, j,kk,n)-x(i, j,kk,n)) +
305  fracy * fracz *(x(i+1,jj,kk,n)-x(i,jj,kk,n));
306  fxp *= bX(i+1,j,k,n);
307 
308  }
309  }
310 
311  Real fym = bY(i,j,k,n)*(x(i,j,k,n) - x(i,j-1,k,n));
312  if (apym != Real(0.0) && apym != Real(1.0)) {
313  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,0)));
314  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,1)));
315  Real fracx = (ccm(ii,j-1,k) || ccm(ii,j,k)) ? std::abs(fcy(i,j,k,0)) : Real(0.0);
316  Real fracz = (ccm(i,j-1,kk) || ccm(i,j,kk)) ? std::abs(fcy(i,j,k,1)) : Real(0.0);
317  if (beta_on_center && phi_on_center)
318  {
319  fym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*fym +
320  fracx*(Real(1.0)-fracz)*bY(ii,j,k ,n)*(x(ii,j,k ,n)-x(ii,j-1,k ,n)) +
321  fracz*(Real(1.0)-fracx)*bY(i ,j,kk,n)*(x(i ,j,kk,n)-x(i ,j-1,kk,n)) +
322  fracx* fracz *bY(ii,j,kk,n)*(x(ii,j,kk,n)-x(ii,j-1,kk,n));
323  }
324  else if (beta_on_centroid && phi_on_center)
325  {
326  fym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*(x( i,j, k,n)-x( i,j-1, k,n)) +
327  fracx *(Real(1.0)-fracz)*(x(ii,j, k,n)-x(ii,j-1, k,n)) +
328  fracz *(Real(1.0)-fracx)*(x(i ,j,kk,n)-x( i,j-1,kk,n)) +
329  fracx * fracz *(x(ii,j,kk,n)-x(ii,j-1,kk,n));
330  fym *= bY(i,j,k,n);
331 
332  }
333  }
334 
335  Real fyp = bY(i,j+1,k,n)*(x(i,j+1,k,n) - x(i,j,k,n));
336  if (apyp != Real(0.0) && apyp != Real(1.0)) {
337  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j+1,k,0)));
338  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j+1,k,1)));
339  Real fracx = (ccm(ii,j,k) || ccm(ii,j+1,k)) ? std::abs(fcy(i,j+1,k,0)) : Real(0.0);
340  Real fracz = (ccm(i,j,kk) || ccm(i,j+1,kk)) ? std::abs(fcy(i,j+1,k,1)) : Real(0.0);
341  if (beta_on_center && phi_on_center)
342  {
343  fyp = (Real(1.0)-fracx)*(Real(1.0)-fracz)*fyp +
344  fracx*(Real(1.0)-fracz)*bY(ii,j+1,k ,n)*(x(ii,j+1,k ,n)-x(ii,j,k ,n)) +
345  fracz*(Real(1.0)-fracx)*bY(i ,j+1,kk,n)*(x(i ,j+1,kk,n)-x(i ,j,kk,n)) +
346  fracx* fracz *bY(ii,j+1,kk,n)*(x(ii,j+1,kk,n)-x(ii,j,kk,n));
347  }
348  else if (beta_on_centroid && phi_on_center)
349  {
350  fyp = (Real(1.0)-fracx)*(Real(1.0)-fracz)*(x( i,j+1, k,n)-x( i,j, k,n)) +
351  fracx *(Real(1.0)-fracz)*(x(ii,j+1, k,n)-x(ii,j, k,n)) +
352  fracz *(Real(1.0)-fracx)*(x( i,j+1,kk,n)-x( i,j,kk,n)) +
353  fracx * fracz *(x(ii,j+1,kk,n)-x(ii,j,kk,n));
354  fyp *= bY(i,j+1,k,n);
355 
356  }
357  }
358 
359  Real fzm = bZ(i,j,k,n)*(x(i,j,k,n) - x(i,j,k-1,n));
360  if (apzm != Real(0.0) && apzm != Real(1.0)) {
361  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,0)));
362  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,1)));
363  Real fracx = (ccm(ii,j,k-1) || ccm(ii,j,k)) ? std::abs(fcz(i,j,k,0)) : Real(0.0);
364  Real fracy = (ccm(i,jj,k-1) || ccm(i,jj,k)) ? std::abs(fcz(i,j,k,1)) : Real(0.0);
365  if (beta_on_center && phi_on_center)
366  {
367  fzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*fzm +
368  fracx*(Real(1.0)-fracy)*bZ(ii,j ,k,n)*(x(ii,j ,k,n)-x(ii,j ,k-1,n)) +
369  fracy*(Real(1.0)-fracx)*bZ(i ,jj,k,n)*(x(i ,jj,k,n)-x(i ,jj,k-1,n)) +
370  fracx* fracy *bZ(ii,jj,k,n)*(x(ii,jj,k,n)-x(ii,jj,k-1,n));
371  }
372  else if (beta_on_centroid && phi_on_center)
373  {
374  fzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*(x( i, j,k,n)-x( i, j,k-1,n)) +
375  fracx *(Real(1.0)-fracy)*(x(ii, j,k,n)-x(ii, j,k-1,n)) +
376  fracy *(Real(1.0)-fracx)*(x( i,jj,k,n)-x( i,jj,k-1,n)) +
377  fracx * fracy *(x(ii,jj,k,n)-x(ii,jj,k-1,n));
378  fzm *= bZ(i,j,k,n);
379 
380  }
381  }
382 
383  Real fzp = bZ(i,j,k+1,n)*(x(i,j,k+1,n) - x(i,j,k,n));
384  if (apzp != Real(0.0) && apzp != Real(1.0)) {
385  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k+1,0)));
386  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k+1,1)));
387  Real fracx = (ccm(ii,j,k) || ccm(ii,j,k+1)) ? std::abs(fcz(i,j,k+1,0)) : Real(0.0);
388  Real fracy = (ccm(i,jj,k) || ccm(i,jj,k+1)) ? std::abs(fcz(i,j,k+1,1)) : Real(0.0);
389  if (beta_on_center && phi_on_center)
390  {
391  fzp = (Real(1.0)-fracx)*(Real(1.0)-fracy)*fzp +
392  fracx*(Real(1.0)-fracy)*bZ(ii,j ,k+1,n)*(x(ii,j ,k+1,n)-x(ii,j ,k,n)) +
393  fracy*(Real(1.0)-fracx)*bZ(i ,jj,k+1,n)*(x(i ,jj,k+1,n)-x(i ,jj,k,n)) +
394  fracx* fracy *bZ(ii,jj,k+1,n)*(x(ii,jj,k+1,n)-x(ii,jj,k,n));
395  }
396  else if (beta_on_centroid && phi_on_center)
397  {
398  fzp = (Real(1.0)-fracx)*(Real(1.0)-fracy)*(x( i, j,k+1,n)-x( i, j,k,n)) +
399  fracx *(Real(1.0)-fracy)*(x(ii, j,k+1,n)-x(ii, j,k,n)) +
400  fracy *(Real(1.0)-fracx)*(x( i,jj,k+1,n)-x( i,jj,k,n)) +
401  fracx * fracy *(x(ii,jj,k+1,n)-x(ii,jj,k,n));
402  fzp *= bZ(i,j,k+1,n);
403 
404  }
405  }
406 
407  Real feb = Real(0.0);
408  if (is_dirichlet) {
409  Real dapx = apxm-apxp;
410  Real dapy = apym-apyp;
411  Real dapz = apzm-apzp;
412  Real anorm = std::sqrt(dapx*dapx+dapy*dapy+dapz*dapz);
413  Real anorminv = Real(1.0)/anorm;
414  Real anrmx = dapx * anorminv;
415  Real anrmy = dapy * anorminv;
416  Real anrmz = dapz * anorminv;
417 
418  Real phib = is_inhomog ? phieb(i,j,k,n) : Real(0.0);
419 
420  Real bctx = bc(i,j,k,0);
421  Real bcty = bc(i,j,k,1);
422  Real bctz = bc(i,j,k,2);
423  Real dx_eb = get_dx_eb(kappa);
424 
425  Real dg = dx_eb / amrex::max(std::abs(anrmx), std::abs(anrmy),
426  std::abs(anrmz));
427  Real gx = bctx - dg*anrmx;
428  Real gy = bcty - dg*anrmy;
429  Real gz = bctz - dg*anrmz;
430  Real sx = std::copysign(Real(1.0),anrmx);
431  Real sy = std::copysign(Real(1.0),anrmy);
432  Real sz = std::copysign(Real(1.0),anrmz);
433  int ii = i - static_cast<int>(sx);
434  int jj = j - static_cast<int>(sy);
435  int kk = k - static_cast<int>(sz);
436 
437  gx = sx*gx;
438  gy = sy*gy;
439  gz = sz*gz;
440  Real gxy = gx*gy;
441  Real gxz = gx*gz;
442  Real gyz = gy*gz;
443  Real gxyz = gx*gy*gz;
444  Real phig = (Real(1.0)+gx+gy+gz+gxy+gxz+gyz+gxyz) * x(i ,j ,k ,n)
445  + (-gz - gxz - gyz - gxyz) * x(i ,j ,kk,n)
446  + (-gy - gxy - gyz - gxyz) * x(i ,jj,k ,n)
447  + (gyz + gxyz) * x(i ,jj,kk,n)
448  + (-gx - gxy - gxz - gxyz) * x(ii,j ,k ,n)
449  + (gxz + gxyz) * x(ii,j ,kk,n)
450  + (gxy + gxyz) * x(ii,jj,k ,n)
451  + (-gxyz) * x(ii,jj,kk,n);
452 
453  Real dphidn = (phib-phig)/dg;
454 
455  feb = dphidn * ba(i,j,k) * beb(i,j,k,n);
456  }
457 
458  y(i,j,k,n) = alpha*a(i,j,k)*x(i,j,k,n) + (Real(1.0)/kappa) *
459  (dhx*(apxm*fxm - apxp*fxp) +
460  dhy*(apym*fym - apyp*fyp) +
461  dhz*(apzm*fzm - apzp*fzp) - dhx*feb);
462  }
463  });
464 }
465 
466 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
467 void mlebabeclap_ebflux (int i, int j, int k, int n,
468  Array4<Real> const& feb,
469  Array4<Real const> const& x,
470  Array4<EBCellFlag const> const& flag,
471  Array4<Real const> const& vfrc,
472  Array4<Real const> const& apx,
473  Array4<Real const> const& apy,
474  Array4<Real const> const& apz,
475  Array4<Real const> const& bc,
476  Array4<Real const> const& beb,
477  Array4<Real const> const& phieb,
478  bool is_inhomog,
479  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
480 {
481  Real dhx = dxinv[0];
482 
483  if (!flag(i,j,k).isSingleValued())
484  {
485  feb(i,j,k,n) = Real(0.0);
486  }
487  else
488  {
489  Real kappa = vfrc(i,j,k);
490  Real apxm = apx(i,j,k);
491  Real apxp = apx(i+1,j,k);
492  Real apym = apy(i,j,k);
493  Real apyp = apy(i,j+1,k);
494  Real apzm = apz(i,j,k);
495  Real apzp = apz(i,j,k+1);
496 
497  Real dapx = apxm-apxp;
498  Real dapy = apym-apyp;
499  Real dapz = apzm-apzp;
500  Real anorm = std::sqrt(dapx*dapx+dapy*dapy+dapz*dapz);
501  Real anorminv = Real(1.0)/anorm;
502  Real anrmx = dapx * anorminv;
503  Real anrmy = dapy * anorminv;
504  Real anrmz = dapz * anorminv;
505 
506  Real phib = is_inhomog ? phieb(i,j,k,n) : Real(0.0);
507 
508  Real bctx = bc(i,j,k,0);
509  Real bcty = bc(i,j,k,1);
510  Real bctz = bc(i,j,k,2);
511  Real dx_eb = get_dx_eb(kappa);
512 
513  Real dg = dx_eb / amrex::max(std::abs(anrmx), std::abs(anrmy), std::abs(anrmz));
514  Real gx = bctx - dg*anrmx;
515  Real gy = bcty - dg*anrmy;
516  Real gz = bctz - dg*anrmz;
517  Real sx = std::copysign(Real(1.0),anrmx);
518  Real sy = std::copysign(Real(1.0),anrmy);
519  Real sz = std::copysign(Real(1.0),anrmz);
520  int ii = i - static_cast<int>(sx);
521  int jj = j - static_cast<int>(sy);
522  int kk = k - static_cast<int>(sz);
523 
524  gx = sx*gx;
525  gy = sy*gy;
526  gz = sz*gz;
527  Real gxy = gx*gy;
528  Real gxz = gx*gz;
529  Real gyz = gy*gz;
530  Real gxyz = gx*gy*gz;
531  Real phig = (Real(1.0)+gx+gy+gz+gxy+gxz+gyz+gxyz) * x(i ,j ,k ,n)
532  + (-gz - gxz - gyz - gxyz) * x(i ,j ,kk,n)
533  + (-gy - gxy - gyz - gxyz) * x(i ,jj,k ,n)
534  + (gyz + gxyz) * x(i ,jj,kk,n)
535  + (-gx - gxy - gxz - gxyz) * x(ii,j ,k ,n)
536  + (gxz + gxyz) * x(ii,j ,kk,n)
537  + (gxy + gxyz) * x(ii,jj,k ,n)
538  + (-gxyz) * x(ii,jj,kk,n);
539 
540  Real dphidn = dhx*(phib-phig)/dg;
541  feb(i,j,k,n) = -beb(i,j,k,n) * dphidn;
542  }
543 }
544 
545 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
546 void mlebabeclap_gsrb (Box const& box,
547  Array4<Real> const& phi, Array4<Real const> const& rhs,
548  Real alpha, Array4<Real const> const& a,
549  Real dhx, Real dhy, Real dhz,
550  Array4<Real const> const& bX, Array4<Real const> const& bY,
551  Array4<Real const> const& bZ,
552  Array4<int const> const& m0, Array4<int const> const& m2,
553  Array4<int const> const& m4,
554  Array4<int const> const& m1, Array4<int const> const& m3,
555  Array4<int const> const& m5,
556  Array4<Real const> const& f0, Array4<Real const> const& f2,
557  Array4<Real const> const& f4,
558  Array4<Real const> const& f1, Array4<Real const> const& f3,
559  Array4<Real const> const& f5,
560  Array4<const int> const& ccm, Array4<EBCellFlag const> const& flag,
561  Array4<Real const> const& vfrc,
562  Array4<Real const> const& apx, Array4<Real const> const& apy,
563  Array4<Real const> const& apz,
564  Array4<Real const> const& fcx, Array4<Real const> const& fcy,
565  Array4<Real const> const& fcz,
566  Array4<Real const> const& ba, Array4<Real const> const& bc,
567  Array4<Real const> const& beb,
568  bool is_dirichlet, bool beta_on_centroid, bool phi_on_centroid,
569  Box const& vbox, int redblack, int ncomp) noexcept
570 {
571  constexpr Real omega = 1.15;
572 
573  const auto vlo = amrex::lbound(vbox);
574  const auto vhi = amrex::ubound(vbox);
575 
576 // amrex::Loop(box, ncomp, [=] (int i, int j, int k, int n) noexcept
577  // amrex::Loop here causes gcc 8 to crash.
578  const auto lo = amrex::lbound(box);
579  const auto hi = amrex::ubound(box);
580  for (int n = 0; n < ncomp; ++n) {
581  for (int k = lo.z; k <= hi.z; ++k) {
582  for (int j = lo.y; j <= hi.y; ++j) {
583  for (int i = lo.x; i <= hi.x; ++i)
584  {
585  if ((i+j+k+redblack) % 2 == 0)
586  {
587  if (flag(i,j,k).isCovered())
588  {
589  phi(i,j,k,n) = Real(0.0);
590  }
591  else
592  {
593  Real cf0 = (i == vlo.x && m0(vlo.x-1,j,k) > 0)
594  ? f0(vlo.x,j,k,n) : Real(0.0);
595  Real cf1 = (j == vlo.y && m1(i,vlo.y-1,k) > 0)
596  ? f1(i,vlo.y,k,n) : Real(0.0);
597  Real cf2 = (k == vlo.z && m2(i,j,vlo.z-1) > 0)
598  ? f2(i,j,vlo.z,n) : Real(0.0);
599  Real cf3 = (i == vhi.x && m3(vhi.x+1,j,k) > 0)
600  ? f3(vhi.x,j,k,n) : Real(0.0);
601  Real cf4 = (j == vhi.y && m4(i,vhi.y+1,k) > 0)
602  ? f4(i,vhi.y,k,n) : Real(0.0);
603  Real cf5 = (k == vhi.z && m5(i,j,vhi.z+1) > 0)
604  ? f5(i,j,vhi.z,n) : Real(0.0);
605 
606  if (flag(i,j,k).isRegular())
607  {
608  Real gamma = alpha*a(i,j,k)
609  + dhx*(bX(i+1,j,k,n) + bX(i,j,k,n))
610  + dhy*(bY(i,j+1,k,n) + bY(i,j,k,n))
611  + dhz*(bZ(i,j,k+1,n) + bZ(i,j,k,n));
612 
613  Real rho = dhx*(bX(i+1,j ,k ,n)*phi(i+1,j ,k ,n) +
614  bX(i ,j ,k ,n)*phi(i-1,j ,k ,n))
615  + dhy*(bY(i ,j+1,k ,n)*phi(i ,j+1,k ,n) +
616  bY(i ,j ,k ,n)*phi(i ,j-1,k ,n))
617  + dhz*(bZ(i ,j ,k+1,n)*phi(i ,j ,k+1,n) +
618  bZ(i ,j ,k ,n)*phi(i ,j ,k-1,n));
619 
620  Real delta = dhx*(bX(i,j,k,n)*cf0 + bX(i+1,j,k,n)*cf3)
621  + dhy*(bY(i,j,k,n)*cf1 + bY(i,j+1,k,n)*cf4)
622  + dhz*(bZ(i,j,k,n)*cf2 + bZ(i,j,k+1,n)*cf5);
623 
624  Real res = rhs(i,j,k,n) - (gamma*phi(i,j,k,n) - rho);
625  phi(i,j,k,n) += omega*res/(gamma-delta);
626  }
627  else
628  {
629  Real kappa = vfrc(i,j,k);
630  Real apxm = apx(i,j,k);
631  Real apxp = apx(i+1,j,k);
632  Real apym = apy(i,j,k);
633  Real apyp = apy(i,j+1,k);
634  Real apzm = apz(i,j,k);
635  Real apzp = apz(i,j,k+1);
636 
637  Real fxm = -bX(i,j,k,n)*phi(i-1,j,k,n);
638  Real oxm = -bX(i,j,k,n)*cf0;
639  Real sxm = bX(i,j,k,n);
640  if (apxm != Real(0.0) && apxm != Real(1.0)) {
641  int jj = j + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,0)));
642  int kk = k + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,1)));
643  Real fracy = (ccm(i-1,jj,k) || ccm(i,jj,k))
644  ? std::abs(fcx(i,j,k,0)) : Real(0.0);
645  Real fracz = (ccm(i-1,j,kk) || ccm(i,j,kk))
646  ? std::abs(fcx(i,j,k,1)) : Real(0.0);
647  if (!beta_on_centroid && !phi_on_centroid)
648  {
649  fxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*fxm
650  + fracy *(Real(1.0)-fracz)*bX(i,jj,k ,n)*(phi(i,jj,k ,n)-phi(i-1,jj,k ,n))
651  +(Real(1.0)-fracy)* fracz *bX(i,j ,kk,n)*(phi(i,j ,kk,n)-phi(i-1,j ,kk,n))
652  + fracy * fracz *bX(i,jj,kk,n)*(phi(i,jj,kk,n)-phi(i-1,jj,kk,n));
653  }
654  else if (beta_on_centroid && !phi_on_centroid)
655  {
656  fxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*( -phi(i-1, j, k,n))
657  + fracy *(Real(1.0)-fracz)*(phi(i,jj,k ,n)-phi(i-1,jj, k,n))
658  +(Real(1.0)-fracy)* fracz *(phi(i,j ,kk,n)-phi(i-1, j,kk,n))
659  + fracy * fracz *(phi(i,jj,kk,n)-phi(i-1,jj,kk,n));
660  fxm *= bX(i,j,k,n);
661 
662  }
663  oxm = Real(0.0);
664  sxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*sxm;
665  }
666 
667  Real fxp = bX(i+1,j,k,n)*phi(i+1,j,k,n);
668  Real oxp = bX(i+1,j,k,n)*cf3;
669  Real sxp = -bX(i+1,j,k,n);
670  if (apxp != Real(0.0) && apxp != Real(1.0)) {
671  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcx(i+1,j,k,0)));
672  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcx(i+1,j,k,1)));
673  Real fracy = (ccm(i,jj,k) || ccm(i+1,jj,k))
674  ? std::abs(fcx(i+1,j,k,0)) : Real(0.0);
675  Real fracz = (ccm(i,j,kk) || ccm(i+1,j,kk))
676  ? std::abs(fcx(i+1,j,k,1)) : Real(0.0);
677  if (!beta_on_centroid && !phi_on_centroid)
678  {
679  fxp = (Real(1.0)-fracy)*(Real(1.0)-fracz)*fxp
680  + fracy *(Real(1.0)-fracz)*bX(i+1,jj,k ,n)*(phi(i+1,jj,k ,n)-phi(i,jj,k ,n))
681  +(Real(1.0)-fracy)* fracz *bX(i+1,j ,kk,n)*(phi(i+1,j ,kk,n)-phi(i,j ,kk,n))
682  + fracy * fracz *bX(i+1,jj,kk,n)*(phi(i+1,jj,kk,n)-phi(i,jj,kk,n));
683  }
684  else if (beta_on_centroid && !phi_on_centroid)
685  {
686  fxp = (Real(1.0)-fracy)*(Real(1.0)-fracz)*(phi(i+1, j, k,n) ) +
687  fracy *(Real(1.0)-fracz)*(phi(i+1,jj, k,n)-phi(i,jj, k,n)) +
688  fracz *(Real(1.0)-fracy)*(phi(i+1, j,kk,n)-phi(i, j,kk,n)) +
689  fracy * fracz *(phi(i+1,jj,kk,n)-phi(i,jj,kk,n));
690  fxp *= bX(i+1,j,k,n);
691 
692  }
693 
694  oxp = Real(0.0);
695  sxp = (Real(1.0)-fracy)*(Real(1.0)-fracz)*sxp;
696  }
697 
698  Real fym = -bY(i,j,k,n)*phi(i,j-1,k,n);
699  Real oym = -bY(i,j,k,n)*cf1;
700  Real sym = bY(i,j,k,n);
701  if (apym != Real(0.0) && apym != Real(1.0)) {
702  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,0)));
703  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,1)));
704  Real fracx = (ccm(ii,j-1,k) || ccm(ii,j,k))
705  ? std::abs(fcy(i,j,k,0)) : Real(0.0);
706  Real fracz = (ccm(i,j-1,kk) || ccm(i,j,kk))
707  ? std::abs(fcy(i,j,k,1)) : Real(0.0);
708  if (!beta_on_centroid && !phi_on_centroid)
709  {
710  fym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*fym
711  + fracx *(Real(1.0)-fracz)*bY(ii,j,k ,n)*(phi(ii,j,k ,n)-phi(ii,j-1,k ,n))
712  + (Real(1.0)-fracx)* fracz *bY(i ,j,kk,n)*(phi(i ,j,kk,n)-phi(i ,j-1,kk,n))
713  + fracx * fracz *bY(ii,j,kk,n)*(phi(ii,j,kk,n)-phi(ii,j-1,kk,n));
714  }
715  else if (beta_on_centroid && !phi_on_centroid)
716  {
717  fym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*( -phi( i,j-1, k,n))
718  + fracx *(Real(1.0)-fracz)*(phi(ii,j,k ,n)-phi(ii,j-1, k,n))
719  + (Real(1.0)-fracx)* fracz *(phi(i ,j,kk,n)-phi( i,j-1,kk,n))
720  + fracx * fracz *(phi(ii,j,kk,n)-phi(ii,j-1,kk,n));
721  fym *= bY(i,j,k,n);
722 
723  }
724  oym = Real(0.0);
725  sym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*sym;
726  }
727 
728  Real fyp = bY(i,j+1,k,n)*phi(i,j+1,k,n);
729  Real oyp = bY(i,j+1,k,n)*cf4;
730  Real syp = -bY(i,j+1,k,n);
731  if (apyp != Real(0.0) && apyp != Real(1.0)) {
732  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j+1,k,0)));
733  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j+1,k,1)));
734  Real fracx = (ccm(ii,j,k) || ccm(ii,j+1,k))
735  ? std::abs(fcy(i,j+1,k,0)) : Real(0.0);
736  Real fracz = (ccm(i,j,kk) || ccm(i,j+1,kk))
737  ? std::abs(fcy(i,j+1,k,1)) : Real(0.0);
738  if (!beta_on_centroid && !phi_on_centroid)
739  {
740  fyp = (Real(1.0)-fracx)*(Real(1.0)-fracz)*fyp
741  + fracx *(Real(1.0)-fracz)*bY(ii,j+1,k ,n)*(phi(ii,j+1,k ,n)-phi(ii,j,k ,n))
742  + (Real(1.0)-fracx)* fracz *bY(i ,j+1,kk,n)*(phi(i ,j+1,kk,n)-phi(i ,j,kk,n))
743  + fracx * fracz *bY(ii,j+1,kk,n)*(phi(ii,j+1,kk,n)-phi(ii,j,kk,n));
744  }
745  else if (beta_on_centroid && !phi_on_centroid)
746  {
747  fyp = (Real(1.0)-fracx)*(Real(1.0)-fracz)*(phi( i,j+1, k,n) )
748  + fracx *(Real(1.0)-fracz)*(phi(ii,j+1, k,n)-phi(ii,j, k,n))
749  + (Real(1.0)-fracx)* fracz *(phi( i,j+1,kk,n)-phi( i,j,kk,n))
750  + fracx * fracz *(phi(ii,j+1,kk,n)-phi(ii,j,kk,n));
751  fyp *= bY(i,j+1,k,n);
752 
753  }
754  oyp = Real(0.0);
755  syp = (Real(1.0)-fracx)*(Real(1.0)-fracz)*syp;
756  }
757 
758  Real fzm = -bZ(i,j,k,n)*phi(i,j,k-1,n);
759  Real ozm = -bZ(i,j,k,n)*cf2;
760  Real szm = bZ(i,j,k,n);
761  if (apzm != Real(0.0) && apzm != Real(1.0)) {
762  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,0)));
763  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,1)));
764  Real fracx = (ccm(ii,j,k-1) || ccm(ii,j,k))
765  ? std::abs(fcz(i,j,k,0)) : Real(0.0);
766  Real fracy = (ccm(i,jj,k-1) || ccm(i,jj,k))
767  ? std::abs(fcz(i,j,k,1)) : Real(0.0);
768  if (!beta_on_centroid && !phi_on_centroid)
769  {
770  fzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*fzm
771  + fracx *(Real(1.0)-fracy)*bZ(ii, j,k,n)*(phi(ii, j,k,n)-phi(ii, j,k-1,n))
772  +(Real(1.0)-fracx)* fracy *bZ( i,jj,k,n)*(phi( i,jj,k,n)-phi( i,jj,k-1,n))
773  + fracx * fracy *bZ(ii,jj,k,n)*(phi(ii,jj,k,n)-phi(ii,jj,k-1,n));
774  }
775  else if (beta_on_centroid && !phi_on_centroid)
776  {
777  fzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*( -phi( i, j,k-1,n))
778  + fracx *(Real(1.0)-fracy)*(phi(ii, j,k,n)-phi(ii, j,k-1,n))
779  + (Real(1.0)-fracx)* fracy *(phi( i,jj,k,n)-phi(i ,jj,k-1,n))
780  + fracx * fracy *(phi(ii,jj,k,n)-phi(ii,jj,k-1,n));
781  fzm *= bZ(i,j,k,n);
782 
783  }
784  ozm = Real(0.0);
785  szm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*szm;
786  }
787 
788  Real fzp = bZ(i,j,k+1,n)*phi(i,j,k+1,n);
789  Real ozp = bZ(i,j,k+1,n)*cf5;
790  Real szp = -bZ(i,j,k+1,n);
791  if (apzp != Real(0.0) && apzp != Real(1.0)) {
792  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k+1,0)));
793  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k+1,1)));
794  Real fracx = (ccm(ii,j,k) || ccm(ii,j,k+1))
795  ? std::abs(fcz(i,j,k+1,0)) : Real(0.0);
796  Real fracy = (ccm(i,jj,k) || ccm(i,jj,k+1))
797  ? std::abs(fcz(i,j,k+1,1)) : Real(0.0);
798  if (!beta_on_centroid && !phi_on_centroid)
799  {
800  fzp = (Real(1.0)-fracx)*(Real(1.0)-fracy)*fzp
801  + fracx *(Real(1.0)-fracy)*bZ(ii,j ,k+1,n)*(phi(ii,j ,k+1,n)-phi(ii,j ,k,n))
802  + (Real(1.0)-fracx)* fracy *bZ(i ,jj,k+1,n)*(phi(i ,jj,k+1,n)-phi(i ,jj,k,n))
803  + fracx * fracy *bZ(ii,jj,k+1,n)*(phi(ii,jj,k+1,n)-phi(ii,jj,k,n));
804  }
805  else if (beta_on_centroid && !phi_on_centroid)
806  {
807  fzp = (Real(1.0)-fracx)*(Real(1.0)-fracy)*(phi( i, j,k+1,n) )
808  + fracx *(Real(1.0)-fracy)*(phi(ii, j,k+1,n)-phi(ii, j,k,n))
809  + (Real(1.0)-fracx)* fracy *(phi( i,jj,k+1,n)-phi( i,jj,k,n))
810  + fracx * fracy *(phi(ii,jj,k+1,n)-phi(ii,jj,k,n));
811  fzp *= bZ(i,j,k+1,n);
812 
813  }
814  ozp = Real(0.0);
815  szp = (Real(1.0)-fracx)*(Real(1.0)-fracy)*szp;
816  }
817 
818  Real vfrcinv = Real(1.0)/kappa;
819  Real gamma = alpha*a(i,j,k) + vfrcinv *
820  (dhx*(apxm*sxm-apxp*sxp) +
821  dhy*(apym*sym-apyp*syp) +
822  dhz*(apzm*szm-apzp*szp));
823 
824  Real rho = -vfrcinv *
825  (dhx*(apxm*fxm-apxp*fxp) +
826  dhy*(apym*fym-apyp*fyp) +
827  dhz*(apzm*fzm-apzp*fzp));
828 
829  Real delta = -vfrcinv *
830  (dhx*(apxm*oxm-apxp*oxp) +
831  dhy*(apym*oym-apyp*oyp) +
832  dhz*(apzm*ozm-apzp*ozp));
833 
834  if (is_dirichlet) {
835  Real dapx = apxm-apxp;
836  Real dapy = apym-apyp;
837  Real dapz = apzm-apzp;
838  Real anorm = std::sqrt(dapx*dapx+dapy*dapy+dapz*dapz);
839  Real anorminv = Real(1.0)/anorm;
840  Real anrmx = dapx * anorminv;
841  Real anrmy = dapy * anorminv;
842  Real anrmz = dapz * anorminv;
843  Real bctx = bc(i,j,k,0);
844  Real bcty = bc(i,j,k,1);
845  Real bctz = bc(i,j,k,2);
846  Real dx_eb = get_dx_eb(vfrc(i,j,k));
847 
848  Real dg = dx_eb / amrex::max(std::abs(anrmx),std::abs(anrmy),
849  std::abs(anrmz));
850 
851  Real gx = bctx - dg*anrmx;
852  Real gy = bcty - dg*anrmy;
853  Real gz = bctz - dg*anrmz;
854  Real sx = std::copysign(Real(1.0),anrmx);
855  Real sy = std::copysign(Real(1.0),anrmy);
856  Real sz = std::copysign(Real(1.0),anrmz);
857  int ii = i - static_cast<int>(sx);
858  int jj = j - static_cast<int>(sy);
859  int kk = k - static_cast<int>(sz);
860 
861  gx *= sx;
862  gy *= sy;
863  gz *= sz;
864  Real gxy = gx*gy;
865  Real gxz = gx*gz;
866  Real gyz = gy*gz;
867  Real gxyz = gx*gy*gz;
868  Real phig_gamma = (Real(1.0)+gx+gy+gz+gxy+gxz+gyz+gxyz);
869  Real phig = (-gz - gxz - gyz - gxyz) * phi(i,j,kk,n)
870  + (-gy - gxy - gyz - gxyz) * phi(i,jj,k,n)
871  + (gyz + gxyz) * phi(i,jj,kk,n)
872  + (-gx - gxy - gxz - gxyz) * phi(ii,j,k,n)
873  + (gxz + gxyz) * phi(ii,j,kk,n)
874  + (gxy + gxyz) * phi(ii,jj,k,n)
875  + (-gxyz) * phi(ii,jj,kk,n);
876 
877  Real dphidn = ( -phig)/dg;
878  Real feb_gamma = -phig_gamma/dg * ba(i,j,k) * beb(i,j,k,n);
879  gamma += vfrcinv*(-dhx)*feb_gamma;
880  Real feb = dphidn * ba(i,j,k) * beb(i,j,k,n);
881  rho += -vfrcinv*(-dhx)*feb;
882  }
883 
884  Real res = rhs(i,j,k,n) - (gamma*phi(i,j,k,n) - rho);
885  phi(i,j,k,n) += omega*res/(gamma-delta);
886  }
887  }
888  }
889  }}}}
890 // });
891 }
892 
893 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
894 void mlebabeclap_flux_x (Box const& box, Array4<Real> const& fx, Array4<Real const> const& apx,
895  Array4<Real const> const& fcx, Array4<Real const> const& sol,
896  Array4<Real const> const& bX, Array4<int const> const& ccm,
897  Real dhx, int face_only, int ncomp, Box const& xbox,
898  bool beta_on_centroid, bool phi_on_centroid) noexcept
899 {
900  int lof = xbox.smallEnd(0);
901  int hif = xbox.bigEnd(0);
902  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
903  {
904  if (!face_only || lof == i || hif == i) {
905  if (apx(i,j,k) == Real(0.0)) {
906  fx(i,j,k,n) = Real(0.0);
907  } else if (apx(i,j,k) == Real(1.0)) {
908  fx(i,j,k,n) = -dhx*bX(i,j,k,n)*(sol(i,j,k,n)-sol(i-1,j,k,n));
909  } else {
910  Real fxm = bX(i,j,k,n)*(sol(i,j,k,n) - sol(i-1,j,k,n));
911  int jj = j + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,0)));
912  int kk = k + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,1)));
913  Real fracy = (ccm(i-1,jj,k) || ccm(i,jj,k)) ? std::abs(fcx(i,j,k,0)) : Real(0.0);
914  Real fracz = (ccm(i-1,j,kk) || ccm(i,j,kk)) ? std::abs(fcx(i,j,k,1)) : Real(0.0);
915  if (!beta_on_centroid && !phi_on_centroid)
916  {
917  fxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*fxm +
918  fracy *(Real(1.0)-fracz)*bX(i,jj,k ,n)*(sol(i,jj,k ,n)-sol(i-1,jj,k ,n)) +
919  fracz *(Real(1.0)-fracy)*bX(i,j ,kk,n)*(sol(i,j ,kk,n)-sol(i-1,j ,kk,n)) +
920  fracy* fracz *bX(i,jj,kk,n)*(sol(i,jj,kk,n)-sol(i-1,jj,kk,n));
921  }
922  else if (beta_on_centroid && !phi_on_centroid)
923  {
924  fxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*(sol(i, j, k,n)-sol(i-1, j, k,n)) +
925  fracy *(Real(1.0)-fracz)*(sol(i,jj, k,n)-sol(i-1,jj, k,n)) +
926  fracz *(Real(1.0)-fracy)*(sol(i, j,kk,n)-sol(i-1, j,kk,n)) +
927  fracy* fracz *(sol(i,jj,kk,n)-sol(i-1,jj,kk,n));
928  fxm *= bX(i,j,k,n);
929 
930  }
931  fx(i,j,k,n) = -dhx*fxm;
932  }
933  }
934  });
935 }
936 
937 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
938 void mlebabeclap_flux_y (Box const& box, Array4<Real> const& fy, Array4<Real const> const& apy,
939  Array4<Real const> const& fcy, Array4<Real const> const& sol,
940  Array4<Real const> const& bY, Array4<int const> const& ccm,
941  Real dhy, int face_only, int ncomp, Box const& ybox,
942  bool beta_on_centroid, bool phi_on_centroid) noexcept
943 {
944  int lof = ybox.smallEnd(1);
945  int hif = ybox.bigEnd(1);
946  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
947  {
948  if (!face_only || lof == j || hif == j) {
949  if (apy(i,j,k) == Real(0.0)) {
950  fy(i,j,k,n) = Real(0.0);
951  } else if (apy(i,j,k) == Real(1.0)) {
952  fy(i,j,k,n) = -dhy*bY(i,j,k,n)*(sol(i,j,k,n)-sol(i,j-1,k,n));
953  } else {
954  Real fym = bY(i,j,k,n)*(sol(i,j,k,n) - sol(i,j-1,k,n));
955  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,0)));
956  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,1)));
957  Real fracx = (ccm(ii,j-1,k) || ccm(ii,j,k)) ? std::abs(fcy(i,j,k,0)) : Real(0.0);
958  Real fracz = (ccm(i,j-1,kk) || ccm(i,j,kk)) ? std::abs(fcy(i,j,k,1)) : Real(0.0);
959  if (!beta_on_centroid && !phi_on_centroid)
960  {
961  fym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*fym +
962  fracx *(Real(1.0)-fracz)*bY(ii,j,k ,n)*(sol(ii,j,k ,n)-sol(ii,j-1,k ,n)) +
963  fracz *(Real(1.0)-fracx)*bY(i ,j,kk,n)*(sol(i ,j,kk,n)-sol(i ,j-1,kk,n)) +
964  fracx * fracz *bY(ii,j,kk,n)*(sol(ii,j,kk,n)-sol(ii,j-1,kk,n));
965  }
966  else if (beta_on_centroid && !phi_on_centroid)
967  {
968  fym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*(sol( i,j, k,n)-sol( i,j-1, k,n)) +
969  fracx *(Real(1.0)-fracz)*(sol(ii,j, k,n)-sol(ii,j-1, k,n)) +
970  fracz *(Real(1.0)-fracx)*(sol( i,j,kk,n)-sol( i,j-1,kk,n)) +
971  fracx * fracz *(sol(ii,j,kk,n)-sol(ii,j-1,kk,n));
972  fym *= bY(i,j,k,n);
973  }
974 
975  fy(i,j,k,n) = -dhy*fym;
976  }
977  }
978  });
979 }
980 
981 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
982 void mlebabeclap_flux_z (Box const& box, Array4<Real> const& fz, Array4<Real const> const& apz,
983  Array4<Real const> const& fcz, Array4<Real const> const& sol,
984  Array4<Real const> const& bZ, Array4<int const> const& ccm,
985  Real dhz, int face_only, int ncomp, Box const& zbox,
986  bool beta_on_centroid, bool phi_on_centroid) noexcept
987 {
988  int lof = zbox.smallEnd(2);
989  int hif = zbox.bigEnd(2);
990  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
991  {
992  if (!face_only || lof == k || hif == k) {
993  if (apz(i,j,k) == Real(0.0)) {
994  fz(i,j,k,n) = Real(0.0);
995  } else if (apz(i,j,k) == Real(1.0)) {
996  fz(i,j,k,n) = -dhz*bZ(i,j,k,n)*(sol(i,j,k,n)-sol(i,j,k-1,n));
997  } else {
998  Real fzm = bZ(i,j,k,n)*(sol(i,j,k,n) - sol(i,j,k-1,n));
999  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,0)));
1000  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,1)));
1001  Real fracx = (ccm(ii,j,k-1) || ccm(ii,j,k)) ? std::abs(fcz(i,j,k,0)) : Real(0.0);
1002  Real fracy = (ccm(i,jj,k-1) || ccm(i,jj,k)) ? std::abs(fcz(i,j,k,1)) : Real(0.0);
1003  if (!beta_on_centroid && !phi_on_centroid)
1004  {
1005  fzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*fzm +
1006  fracx*(Real(1.0)-fracy)*bZ(ii,j ,k,n)*(sol(ii,j ,k,n)-sol(ii,j ,k-1,n)) +
1007  fracy*(Real(1.0)-fracx)*bZ(i ,jj,k,n)*(sol(i ,jj,k,n)-sol(i ,jj,k-1,n)) +
1008  fracx* fracy *bZ(ii,jj,k,n)*(sol(ii,jj,k,n)-sol(ii,jj,k-1,n));
1009  }
1010  else if (beta_on_centroid && !phi_on_centroid)
1011  {
1012  fzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*(sol( i, j,k,n)-sol( i, j,k-1,n)) +
1013  fracx *(Real(1.0)-fracy)*(sol(ii, j,k,n)-sol(ii, j,k-1,n)) +
1014  fracy *(Real(1.0)-fracx)*(sol( i,jj,k,n)-sol( i,jj,k-1,n)) +
1015  fracx * fracy *(sol(ii,jj,k,n)-sol(ii,jj,k-1,n));
1016  fzm *= bZ(i,j,k,n);
1017 
1018  }
1019 
1020  fz(i,j,k,n) = -dhz*fzm;
1021  }
1022  }
1023  });
1024 }
1025 
1026 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1027 void mlebabeclap_flux_x_0 (Box const& box, Array4<Real> const& fx, Array4<Real const> const& apx,
1028  Array4<Real const> const& sol, Array4<Real const> const& bX,
1029  Real dhx, int face_only, int ncomp, Box const& xbox) noexcept
1030 {
1031  int lof = xbox.smallEnd(0);
1032  int hif = xbox.bigEnd(0);
1033  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1034  {
1035  if (!face_only || lof == i || hif == i) {
1036  if (apx(i,j,k) == Real(0.0)) {
1037  fx(i,j,k,n) = Real(0.0);
1038  } else {
1039  fx(i,j,k,n) = -dhx*bX(i,j,k,n)*(sol(i,j,k,n)-sol(i-1,j,k,n));
1040  }
1041  }
1042  });
1043 }
1044 
1045 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1046 void mlebabeclap_flux_y_0 (Box const& box, Array4<Real> const& fy, Array4<Real const> const& apy,
1047  Array4<Real const> const& sol, Array4<Real const> const& bY,
1048  Real dhy, int face_only, int ncomp, Box const& ybox) noexcept
1049 {
1050  int lof = ybox.smallEnd(1);
1051  int hif = ybox.bigEnd(1);
1052  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1053  {
1054  if (!face_only || lof == j || hif == j) {
1055  if (apy(i,j,k) == Real(0.0)) {
1056  fy(i,j,k,n) = Real(0.0);
1057  } else {
1058  fy(i,j,k,n) = -dhy*bY(i,j,k,n)*(sol(i,j,k,n)-sol(i,j-1,k,n));
1059  }
1060  }
1061  });
1062 }
1063 
1064 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1065 void mlebabeclap_flux_z_0 (Box const& box, Array4<Real> const& fz, Array4<Real const> const& apz,
1066  Array4<Real const> const& sol, Array4<Real const> const& bZ,
1067  Real dhz, int face_only, int ncomp, Box const& zbox) noexcept
1068 {
1069  int lof = zbox.smallEnd(2);
1070  int hif = zbox.bigEnd(2);
1071  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1072  {
1073  if (!face_only || lof == k || hif == k) {
1074  if (apz(i,j,k) == Real(0.0)) {
1075  fz(i,j,k,n) = Real(0.0);
1076  } else {
1077  fz(i,j,k,n) = -dhz*bZ(i,j,k,n)*(sol(i,j,k,n)-sol(i,j,k-1,n));
1078  }
1079  }
1080  });
1081 }
1082 
1083 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1084 void mlebabeclap_grad_x (Box const& box, Array4<Real> const& gx, Array4<Real const> const& sol,
1085  Array4<Real const> const& apx, Array4<Real const> const& fcx,
1086  Array4<int const> const& ccm,
1087  Real dxi, int ncomp, bool phi_on_centroid) noexcept
1088 {
1089  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1090  {
1091  if (apx(i,j,k) == Real(0.0)) {
1092  gx(i,j,k,n) = Real(0.0);
1093  } else if (apx(i,j,k) == Real(1.0)) {
1094  gx(i,j,k,n) = dxi*(sol(i,j,k,n)-sol(i-1,j,k,n));
1095  } else {
1096  Real gxm = (sol(i,j,k,n) - sol(i-1,j,k,n));
1097  int jj = j + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,0)));
1098  int kk = k + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,1)));
1099  Real fracy = (ccm(i-1,jj,k) || ccm(i,jj,k)) ? std::abs(fcx(i,j,k,0)) : Real(0.0);
1100  Real fracz = (ccm(i-1,j,kk) || ccm(i,j,kk)) ? std::abs(fcx(i,j,k,1)) : Real(0.0);
1101  if (!phi_on_centroid)
1102  {
1103  gxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*gxm +
1104  fracy*(Real(1.0)-fracz)*(sol(i,jj,k ,n)-sol(i-1,jj,k ,n)) +
1105  fracz*(Real(1.0)-fracy)*(sol(i,j ,kk,n)-sol(i-1,j ,kk,n)) +
1106  fracy* fracz *(sol(i,jj,kk,n)-sol(i-1,jj,kk,n));
1107  }
1108  gx(i,j,k,n) = dxi*gxm;
1109  }
1110  });
1111 }
1112 
1113 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1114 void mlebabeclap_grad_y (Box const& box, Array4<Real> const& gy, Array4<Real const> const& sol,
1115  Array4<Real const> const& apy, Array4<Real const> const& fcy,
1116  Array4<int const> const& ccm,
1117  Real dyi, int ncomp, bool phi_on_centroid) noexcept
1118 {
1119  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1120  {
1121  if (apy(i,j,k) == Real(0.0)) {
1122  gy(i,j,k,n) = Real(0.0);
1123  } else if (apy(i,j,k) == Real(1.0)) {
1124  gy(i,j,k,n) = dyi*(sol(i,j,k,n)-sol(i,j-1,k,n));
1125  } else {
1126  Real gym = (sol(i,j,k,n) - sol(i,j-1,k,n));
1127  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,0)));
1128  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,1)));
1129  Real fracx = (ccm(ii,j-1,k) || ccm(ii,j,k)) ? std::abs(fcy(i,j,k,0)) : Real(0.0);
1130  Real fracz = (ccm(i,j-1,kk) || ccm(i,j,kk)) ? std::abs(fcy(i,j,k,1)) : Real(0.0);
1131  if (!phi_on_centroid)
1132  {
1133  gym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*gym +
1134  fracx*(Real(1.0)-fracz)*(sol(ii,j,k ,n)-sol(ii,j-1,k ,n)) +
1135  fracz*(Real(1.0)-fracx)*(sol(i ,j,kk,n)-sol(i ,j-1,kk,n)) +
1136  fracx* fracz *(sol(ii,j,kk,n)-sol(ii,j-1,kk,n));
1137  }
1138  gy(i,j,k,n) = dyi*gym;
1139  }
1140  });
1141 }
1142 
1143 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1144 void mlebabeclap_grad_z (Box const& box, Array4<Real> const& gz, Array4<Real const> const& sol,
1145  Array4<Real const> const& apz, Array4<Real const> const& fcz,
1146  Array4<int const> const& ccm,
1147  Real dzi, int ncomp, bool phi_on_centroid) noexcept
1148 {
1149  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1150  {
1151  if (apz(i,j,k) == Real(0.0)) {
1152  gz(i,j,k,n) = Real(0.0);
1153  } else if (apz(i,j,k) == Real(1.0)) {
1154  gz(i,j,k,n) = dzi*(sol(i,j,k,n)-sol(i,j,k-1,n));
1155  } else {
1156  Real gzm = (sol(i,j,k,n) - sol(i,j,k-1,n));
1157  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,0)));
1158  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,1)));
1159  Real fracx = (ccm(ii,j,k-1) || ccm(ii,j,k)) ? std::abs(fcz(i,j,k,0)) : Real(0.0);
1160  Real fracy = (ccm(i,jj,k-1) || ccm(i,jj,k)) ? std::abs(fcz(i,j,k,1)) : Real(0.0);
1161  if (!phi_on_centroid)
1162  {
1163  gzm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*gzm +
1164  fracx*(Real(1.0)-fracy)*(sol(ii,j ,k,n)-sol(ii,j ,k-1,n)) +
1165  fracy*(Real(1.0)-fracx)*(sol(i ,jj,k,n)-sol(i ,jj,k-1,n)) +
1166  fracx* fracy *(sol(ii,jj,k,n)-sol(ii,jj,k-1,n));
1167  }
1168  gz(i,j,k,n) = dzi*gzm;
1169  }
1170  });
1171 }
1172 
1173 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1174 void mlebabeclap_grad_x_0 (Box const& box, Array4<Real> const& gx, Array4<Real const> const& sol,
1175  Array4<Real const> const& apx, Real dxi, int ncomp) noexcept
1176 {
1177  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1178  {
1179  if (apx(i,j,k) == Real(0.0)) {
1180  gx(i,j,k,n) = Real(0.0);
1181  } else {
1182  gx(i,j,k,n) = dxi*(sol(i,j,k,n)-sol(i-1,j,k,n));
1183  }
1184  });
1185 }
1186 
1187 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1188 void mlebabeclap_grad_y_0 (Box const& box, Array4<Real> const& gy, Array4<Real const> const& sol,
1189  Array4<Real const> const& apy, Real dyi, int ncomp) noexcept
1190 {
1191  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1192  {
1193  if (apy(i,j,k) == Real(0.0)) {
1194  gy(i,j,k,n) = Real(0.0);
1195  } else {
1196  gy(i,j,k,n) = dyi*(sol(i,j,k,n)-sol(i,j-1,k,n));
1197  }
1198  });
1199 }
1200 
1201 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1202 void mlebabeclap_grad_z_0 (Box const& box, Array4<Real> const& gz, Array4<Real const> const& sol,
1203  Array4<Real const> const& apz, Real dzi, int ncomp) noexcept
1204 {
1205  amrex::LoopConcurrent(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1206  {
1207  if (apz(i,j,k) == Real(0.0)) {
1208  gz(i,j,k,n) = Real(0.0);
1209  } else {
1210  gz(i,j,k,n) = dzi*(sol(i,j,k,n)-sol(i,j,k-1,n));
1211  }
1212  });
1213 }
1214 
1215 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1216 void mlebabeclap_normalize (Box const& box, Array4<Real> const& phi,
1217  Real alpha, Array4<Real const> const& a,
1218  Real dhx, Real dhy, Real dhz,
1219  Array4<Real const> const& bX, Array4<Real const> const& bY,
1220  Array4<Real const> const& bZ,
1221  Array4<const int> const& ccm, Array4<EBCellFlag const> const& flag,
1222  Array4<Real const> const& vfrc,
1223  Array4<Real const> const& apx, Array4<Real const> const& apy,
1224  Array4<Real const> const& apz,
1225  Array4<Real const> const& fcx, Array4<Real const> const& fcy,
1226  Array4<Real const> const& fcz,
1227  Array4<Real const> const& ba, Array4<Real const> const& bc,
1228  Array4<Real const> const& beb,
1229  bool is_dirichlet, bool beta_on_centroid, int ncomp) noexcept
1230 {
1231  amrex::Loop(box, ncomp, [=] (int i, int j, int k, int n) noexcept
1232  {
1233  if (flag(i,j,k).isRegular())
1234  {
1235  phi(i,j,k,n) /= alpha*a(i,j,k) + dhx*(bX(i,j,k,n) + bX(i+1,j,k,n))
1236  + dhy*(bY(i,j,k,n) + bY(i,j+1,k,n))
1237  + dhz*(bZ(i,j,k,n) + bZ(i,j,k+1,n));
1238  }
1239  else if (flag(i,j,k).isSingleValued())
1240  {
1241  Real kappa = vfrc(i,j,k);
1242  Real apxm = apx(i,j,k);
1243  Real apxp = apx(i+1,j,k);
1244  Real apym = apy(i,j,k);
1245  Real apyp = apy(i,j+1,k);
1246  Real apzm = apz(i,j,k);
1247  Real apzp = apz(i,j,k+1);
1248 
1249  Real sxm = bX(i,j,k,n);
1250  if (apxm != Real(0.0) && apxm != Real(1.0) && !beta_on_centroid) {
1251  int jj = j + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,0)));
1252  int kk = k + static_cast<int>(std::copysign(Real(1.0), fcx(i,j,k,1)));
1253  Real fracy = (ccm(i-1,jj,k) || ccm(i,jj,k))
1254  ? std::abs(fcx(i,j,k,0)) : Real(0.0);
1255  Real fracz = (ccm(i-1,j,kk) || ccm(i,j,kk))
1256  ? std::abs(fcx(i,j,k,1)) : Real(0.0);
1257  sxm = (Real(1.0)-fracy)*(Real(1.0)-fracz)*sxm;
1258  }
1259 
1260  Real sxp = -bX(i+1,j,k,n);
1261  if (apxp != Real(0.0) && apxp != Real(1.0) && !beta_on_centroid) {
1262  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcx(i+1,j,k,0)));
1263  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcx(i+1,j,k,1)));
1264  Real fracy = (ccm(i,jj,k) || ccm(i+1,jj,k))
1265  ? std::abs(fcx(i+1,j,k,0)) : Real(0.0);
1266  Real fracz = (ccm(i,j,kk) || ccm(i+1,j,kk))
1267  ? std::abs(fcx(i+1,j,k,1)) : Real(0.0);
1268  sxp = (Real(1.0)-fracy)*(Real(1.0)-fracz)*sxp;
1269  }
1270 
1271  Real sym = bY(i,j,k,n);
1272  if (apym != Real(0.0) && apym != Real(1.0) && !beta_on_centroid) {
1273  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,0)));
1274  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j,k,1)));
1275  Real fracx = (ccm(ii,j-1,k) || ccm(ii,j,k))
1276  ? std::abs(fcy(i,j,k,0)) : Real(0.0);
1277  Real fracz = (ccm(i,j-1,kk) || ccm(i,j,kk))
1278  ? std::abs(fcy(i,j,k,1)) : Real(0.0);
1279  sym = (Real(1.0)-fracx)*(Real(1.0)-fracz)*sym;
1280  }
1281 
1282  Real syp = -bY(i,j+1,k,n);
1283  if (apyp != Real(0.0) && apyp != Real(1.0) && !beta_on_centroid) {
1284  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcy(i,j+1,k,0)));
1285  int kk = k + static_cast<int>(std::copysign(Real(1.0),fcy(i,j+1,k,1)));
1286  Real fracx = (ccm(ii,j,k) || ccm(ii,j+1,k))
1287  ? std::abs(fcy(i,j+1,k,0)) : Real(0.0);
1288  Real fracz = (ccm(i,j,kk) || ccm(i,j+1,kk))
1289  ? std::abs(fcy(i,j+1,k,1)) : Real(0.0);
1290  syp = (Real(1.0)-fracx)*(Real(1.0)-fracz)*syp;
1291  }
1292 
1293  Real szm = bZ(i,j,k,n);
1294  if (apzm != Real(0.0) && apzm != Real(1.0) && !beta_on_centroid) {
1295  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,0)));
1296  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k,1)));
1297  Real fracx = (ccm(ii,j,k-1) || ccm(ii,j,k))
1298  ? std::abs(fcz(i,j,k,0)) : Real(0.0);
1299  Real fracy = (ccm(i,jj,k-1) || ccm(i,jj,k))
1300  ? std::abs(fcz(i,j,k,1)) : Real(0.0);
1301  szm = (Real(1.0)-fracx)*(Real(1.0)-fracy)*szm;
1302  }
1303 
1304  Real szp = -bZ(i,j,k+1,n);
1305  if (apzp != Real(0.0) && apzp != Real(1.0) && !beta_on_centroid) {
1306  int ii = i + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k+1,0)));
1307  int jj = j + static_cast<int>(std::copysign(Real(1.0),fcz(i,j,k+1,1)));
1308  Real fracx = (ccm(ii,j,k) || ccm(ii,j,k+1))
1309  ? std::abs(fcz(i,j,k+1,0)) : Real(0.0);
1310  Real fracy = (ccm(i,jj,k) || ccm(i,jj,k+1))
1311  ? std::abs(fcz(i,j,k+1,1)) : Real(0.0);
1312  szp = (Real(1.0)-fracx)*(Real(1.0)-fracy)*szp;
1313  }
1314 
1315  Real vfrcinv = Real(1.0)/kappa;
1316  Real gamma = alpha*a(i,j,k) + vfrcinv *
1317  (dhx*(apxm*sxm-apxp*sxp) +
1318  dhy*(apym*sym-apyp*syp) +
1319  dhz*(apzm*szm-apzp*szp));
1320 
1321  if (is_dirichlet) {
1322  Real dapx = apxm-apxp;
1323  Real dapy = apym-apyp;
1324  Real dapz = apzm-apzp;
1325  Real anorm = std::sqrt(dapx*dapx+dapy*dapy+dapz*dapz);
1326  Real anorminv = Real(1.0)/anorm;
1327  Real anrmx = dapx * anorminv;
1328  Real anrmy = dapy * anorminv;
1329  Real anrmz = dapz * anorminv;
1330  Real bctx = bc(i,j,k,0);
1331  Real bcty = bc(i,j,k,1);
1332  Real bctz = bc(i,j,k,2);
1333  Real dx_eb = get_dx_eb(vfrc(i,j,k));
1334 
1335  Real dg = dx_eb / amrex::max(std::abs(anrmx),std::abs(anrmy),
1336  std::abs(anrmz));
1337 
1338  Real gx = bctx - dg*anrmx;
1339  Real gy = bcty - dg*anrmy;
1340  Real gz = bctz - dg*anrmz;
1341  Real sx = std::copysign(Real(1.0),anrmx);
1342  Real sy = std::copysign(Real(1.0),anrmy);
1343  Real sz = std::copysign(Real(1.0),anrmz);
1344 
1345  gx *= sx;
1346  gy *= sy;
1347  gz *= sz;
1348  Real gxy = gx*gy;
1349  Real gxz = gx*gz;
1350  Real gyz = gy*gz;
1351  Real gxyz = gx*gy*gz;
1352  Real phig_gamma = (Real(1.0)+gx+gy+gz+gxy+gxz+gyz+gxyz);
1353  Real feb_gamma = -phig_gamma/dg * ba(i,j,k) * beb(i,j,k,n);
1354  gamma += vfrcinv*(-dhx)*feb_gamma;
1355  }
1356 
1357  phi(i,j,k,n) /= gamma;
1358  }
1359  });
1360 }
1361 
1362 }
1363 
1364 #endif
AMReX_EB_LeastSquares_3D_K.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
AMReX_REAL.H
abs
#define abs(x)
Definition: complex-type.h:85
amrex
Definition: AMReX_Amr.cpp:49
amrex::mlebabeclap_grad_x_0
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_grad_x_0(Box const &box, Array4< Real > const &gx, Array4< Real const > const &sol, Array4< Real const > const &apx, Real dxi, int ncomp) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:740
amrex::grad_y_of_phi_on_centroids
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_y_of_phi_on_centroids(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Real &xloc_on_yface, bool is_eb_dirichlet, bool is_eb_inhomog)
Definition: AMReX_EB_LeastSquares_2D_K.H:331
amrex::mlebabeclap_grad_y_0
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_grad_y_0(Box const &box, Array4< Real > const &gy, Array4< Real const > const &sol, Array4< Real const > const &apy, Real dyi, int ncomp) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:754
amrex::Box
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
amrex::max
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & max(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:35
amrex::mlebabeclap_adotx_centroid
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_adotx_centroid(Box const &box, Array4< Real > const &y, Array4< Real const > const &x, Array4< Real const > const &a, Array4< Real const > const &bX, Array4< Real const > const &bY, Array4< EBCellFlag const > const &flag, Array4< Real const > const &vfrc, Array4< Real const > const &apx, Array4< Real const > const &apy, Array4< Real const > const &fcx, Array4< Real const > const &fcy, Array4< Real const > const &ccent, Array4< Real const > const &ba, Array4< Real const > const &bcent, Array4< Real const > const &beb, Array4< Real const > const &phieb, const int &domlo_x, const int &domlo_y, const int &domhi_x, const int &domhi_y, const bool &on_x_face, const bool &on_y_face, bool is_eb_dirichlet, bool is_eb_inhomog, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Real alpha, Real beta, int ncomp) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:11
amrex::mlebabeclap_adotx
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_adotx(Box const &box, Array4< Real > const &y, Array4< Real const > const &x, Array4< Real const > const &a, Array4< Real const > const &bX, Array4< Real const > const &bY, Array4< const int > const &ccm, Array4< EBCellFlag const > const &flag, Array4< Real const > const &vfrc, Array4< Real const > const &apx, Array4< Real const > const &apy, Array4< Real const > const &fcx, Array4< Real const > const &fcy, Array4< Real const > const &ba, Array4< Real const > const &bc, Array4< Real const > const &beb, bool is_dirichlet, Array4< Real const > const &phieb, bool is_inhomog, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Real alpha, Real beta, int ncomp, bool beta_on_centroid, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:166
amrex::mlebabeclap_flux_x_0
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_flux_x_0(Box const &box, Array4< Real > const &fx, Array4< Real const > const &apx, Array4< Real const > const &sol, Array4< Real const > const &bX, Real dhx, int face_only, int ncomp, Box const &xbox) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:654
amrex::mlebabeclap_grad_z
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_grad_z(Box const &box, Array4< Real > const &gz, Array4< Real const > const &sol, Array4< Real const > const &apz, Array4< Real const > const &fcz, Array4< int const > const &ccm, Real dzi, int ncomp, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_3D_K.H:1144
amrex::ubound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:315
amrex::mlebabeclap_grad_x
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_grad_x(Box const &box, Array4< Real > const &gx, Array4< Real const > const &sol, Array4< Real const > const &apx, Array4< Real const > const &fcx, Array4< int const > const &ccm, Real dxi, int ncomp, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:692
amrex::lbound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:308
amrex::grad_x_of_phi_on_centroids_extdir
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_x_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &yloc_on_xface, bool is_eb_dirichlet, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y)
Definition: AMReX_EB_LeastSquares_2D_K.H:495
amrex::grad_z_of_phi_on_centroids
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_z_of_phi_on_centroids(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Real &xloc_on_zface, Real &yloc_on_zface, bool is_eb_dirichlet, bool is_eb_inhomog)
Definition: AMReX_EB_LeastSquares_3D_K.H:513
amrex::mlebabeclap_flux_z_0
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_flux_z_0(Box const &box, Array4< Real > const &fz, Array4< Real const > const &apz, Array4< Real const > const &sol, Array4< Real const > const &bZ, Real dhz, int face_only, int ncomp, Box const &zbox) noexcept
Definition: AMReX_MLEBABecLap_3D_K.H:1065
amrex::mlebabeclap_flux_y_0
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_flux_y_0(Box const &box, Array4< Real > const &fy, Array4< Real const > const &apy, Array4< Real const > const &sol, Array4< Real const > const &bY, Real dhy, int face_only, int ncomp, Box const &ybox) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:673
amrex::mlebabeclap_grad_z_0
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_grad_z_0(Box const &box, Array4< Real > const &gz, Array4< Real const > const &sol, Array4< Real const > const &apz, Real dzi, int ncomp) noexcept
Definition: AMReX_MLEBABecLap_3D_K.H:1202
amrex::grad_eb_of_phi_on_centroids_extdir
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_eb_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &nrmx, Real &nrmy, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y)
Definition: AMReX_EB_LeastSquares_2D_K.H:765
amrex::get_dx_eb
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real get_dx_eb(Real kappa) noexcept
Definition: AMReX_MLEBABecLap_K.H:11
amrex::grad_z_of_phi_on_centroids_extdir
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_z_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &xloc_on_zface, Real &yloc_on_zface, bool is_eb_dirichlet, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y, const bool on_z_face, const int domlo_z, const int domhi_z)
Definition: AMReX_EB_LeastSquares_3D_K.H:1010
amrex::mlebabeclap_flux_x
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_flux_x(Box const &box, Array4< Real > const &fx, Array4< Real const > const &apx, Array4< Real const > const &fcx, Array4< Real const > const &sol, Array4< Real const > const &bX, Array4< int const > const &ccm, Real dhx, int face_only, int ncomp, Box const &xbox, bool beta_on_centroid, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:588
amrex::mlebabeclap_gsrb
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_gsrb(Box const &box, Array4< Real > const &phi, Array4< Real const > const &rhs, Real alpha, Array4< Real const > const &a, Real dhx, Real dhy, Real dh, GpuArray< Real, AMREX_SPACEDIM > const &dx, Array4< Real const > const &bX, Array4< Real const > const &bY, Array4< int const > const &m0, Array4< int const > const &m2, Array4< int const > const &m1, Array4< int const > const &m3, Array4< Real const > const &f0, Array4< Real const > const &f2, Array4< Real const > const &f1, Array4< Real const > const &f3, Array4< const int > const &ccm, Array4< EBCellFlag const > const &flag, Array4< Real const > const &vfrc, Array4< Real const > const &apx, Array4< Real const > const &apy, Array4< Real const > const &fcx, Array4< Real const > const &fcy, Array4< Real const > const &ba, Array4< Real const > const &bc, Array4< Real const > const &beb, bool is_dirichlet, bool beta_on_centroid, bool phi_on_centroid, Box const &vbox, int redblack, int ncomp) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:371
amrex::grad_y_of_phi_on_centroids_extdir
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_y_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &xloc_on_yface, bool is_eb_dirichlet, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y)
Definition: AMReX_EB_LeastSquares_2D_K.H:632
amrex::LoopConcurrent
AMREX_GPU_HOST_DEVICE AMREX_ATTRIBUTE_FLATTEN_FOR void LoopConcurrent(Dim3 lo, Dim3 hi, F const &f) noexcept
Definition: AMReX_Loop.H:149
amrex::mlebabeclap_ebflux
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_ebflux(int i, int j, int k, int n, Array4< Real > const &feb, Array4< Real const > const &x, Array4< EBCellFlag const > const &flag, Array4< Real const > const &vfrc, Array4< Real const > const &apx, Array4< Real const > const &apy, Array4< Real const > const &bc, Array4< Real const > const &beb, Array4< Real const > const &phieb, bool is_inhomog, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:304
amrex::mlebabeclap_normalize
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_normalize(Box const &box, Array4< Real > const &phi, Real alpha, Array4< Real const > const &a, Real dhx, Real dhy, Real dh, const amrex::GpuArray< Real, AMREX_SPACEDIM > &dx, Array4< Real const > const &bX, Array4< Real const > const &bY, Array4< const int > const &ccm, Array4< EBCellFlag const > const &flag, Array4< Real const > const &vfrc, Array4< Real const > const &apx, Array4< Real const > const &apy, Array4< Real const > const &fcx, Array4< Real const > const &fcy, Array4< Real const > const &ba, Array4< Real const > const &bc, Array4< Real const > const &beb, bool is_dirichlet, bool beta_on_centroid, int ncomp) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:768
amrex::mlebabeclap_flux_z
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_flux_z(Box const &box, Array4< Real > const &fz, Array4< Real const > const &apz, Array4< Real const > const &fcz, Array4< Real const > const &sol, Array4< Real const > const &bZ, Array4< int const > const &ccm, Real dhz, int face_only, int ncomp, Box const &zbox, bool beta_on_centroid, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_3D_K.H:982
amrex::mlebabeclap_grad_y
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_grad_y(Box const &box, Array4< Real > const &gy, Array4< Real const > const &sol, Array4< Real const > const &apy, Array4< Real const > const &fcy, Array4< int const > const &ccm, Real dyi, int ncomp, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:716
amrex::grad_x_of_phi_on_centroids
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_x_of_phi_on_centroids(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Real &yloc_on_xface, bool is_eb_dirichlet, bool is_eb_inhomog)
Definition: AMReX_EB_LeastSquares_2D_K.H:246
amrex::mlebabeclap_flux_y
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebabeclap_flux_y(Box const &box, Array4< Real > const &fy, Array4< Real const > const &apy, Array4< Real const > const &fcy, Array4< Real const > const &sol, Array4< Real const > const &bY, Array4< int const > const &ccm, Real dhy, int face_only, int ncomp, Box const &ybox, bool beta_on_centroid, bool phi_on_centroid) noexcept
Definition: AMReX_MLEBABecLap_2D_K.H:621
amrex::sqrt
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE GpuComplex< T > sqrt(const GpuComplex< T > &a_z) noexcept
Return the square root of a complex number.
Definition: AMReX_GpuComplex.H:373
amrex::Loop
AMREX_GPU_HOST_DEVICE AMREX_ATTRIBUTE_FLATTEN_FOR void Loop(Dim3 lo, Dim3 hi, F const &f) noexcept
Definition: AMReX_Loop.H:124
amrex::Array4
Definition: AMReX_Array4.H:61