Block-Structured AMR Software Framework
AMReX_MLEBTensor_2D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_ML_EB_TENSOR_2D_K_H_
2 #define AMREX_ML_EB_TENSOR_2D_K_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_MLEBABecLap_K.H>
6 
7 namespace amrex {
8 
10 void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
11  Array4<Real const> const& vel,
12  Array4<Real const> const& etax,
13  Array4<Real const> const& kapx,
14  Array4<Real const> const& apx,
15  Array4<EBCellFlag const> const& flag,
16  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
17 {
18  const Real dyi = dxinv[1];
19  const auto lo = amrex::lbound(box);
20  const auto hi = amrex::ubound(box);
21  constexpr Real twoThirds = 2./3.;
22 
23  int k = 0;
24  for (int j = lo.y; j <= hi.y; ++j) {
26  for (int i = lo.x; i <= hi.x; ++i) {
27  if (apx(i,j,0) == Real(0.0))
28  {
29  fx(i,j,0,0) = Real(0.0);
30  fx(i,j,0,1) = Real(0.0);
31  }
32  else
33  {
34  int jhip = j + flag(i ,j,0).isConnected(0, 1,0);
35  int jhim = j - flag(i ,j,0).isConnected(0,-1,0);
36  int jlop = j + flag(i-1,j,0).isConnected(0, 1,0);
37  int jlom = j - flag(i-1,j,0).isConnected(0,-1,0);
38  Real whi = mlebtensor_weight(jhip-jhim);
39  Real wlo = mlebtensor_weight(jlop-jlom);
40  Real dudy = mlebtensor_dy_on_xface(i,j,k,0,vel,dyi,
41  whi,wlo,jhip,jhim,jlop,jlom);
42  Real dvdy = mlebtensor_dy_on_xface(i,j,k,1,vel,dyi,
43  whi,wlo,jhip,jhim,jlop,jlom);
44  Real divu = dvdy;
45  Real xif = kapx(i,j,0);
46  Real mun = Real(0.75)*(etax(i,j,0,0)-xif);// restore the original eta
47  Real mut = etax(i,j,0,1);
48  fx(i,j,0,0) = -mun*(-twoThirds*divu) - xif*divu;
49  fx(i,j,0,1) = -mut*dudy;
50  }
51  }
52  }
53 }
54 
56 void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
57  Array4<Real const> const& vel,
58  Array4<Real const> const& etay,
59  Array4<Real const> const& kapy,
60  Array4<Real const> const& apy,
61  Array4<EBCellFlag const> const& flag,
62  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
63 {
64  const Real dxi = dxinv[0];
65  const auto lo = amrex::lbound(box);
66  const auto hi = amrex::ubound(box);
67  constexpr Real twoThirds = 2./3.;
68 
69  int k = 0;
70  for (int j = lo.y; j <= hi.y; ++j) {
72  for (int i = lo.x; i <= hi.x; ++i) {
73  if (apy(i,j,0) == Real(0.0))
74  {
75  fy(i,j,0,0) = Real(0.0);
76  fy(i,j,0,1) = Real(0.0);
77  }
78  else
79  {
80  int ihip = i + flag(i,j ,0).isConnected( 1,0,0);
81  int ihim = i - flag(i,j ,0).isConnected(-1,0,0);
82  int ilop = i + flag(i,j-1,0).isConnected( 1,0,0);
83  int ilom = i - flag(i,j-1,0).isConnected(-1,0,0);
84  Real whi = mlebtensor_weight(ihip-ihim);
85  Real wlo = mlebtensor_weight(ilop-ilom);
86  Real dudx = mlebtensor_dx_on_yface(i,j,k,0,vel,dxi,
87  whi,wlo,ihip,ihim,ilop,ilom);
88  Real dvdx = mlebtensor_dx_on_yface(i,j,k,1,vel,dxi,
89  whi,wlo,ihip,ihim,ilop,ilom);
90  Real divu = dudx;
91  Real xif = kapy(i,j,0);
92  Real mun = Real(0.75)*(etay(i,j,0,1)-xif);// restore the original eta
93  Real mut = etay(i,j,0,0);
94  fy(i,j,0,0) = -mut*dvdx;
95  fy(i,j,0,1) = -mun*(-twoThirds*divu) - xif*divu;
96  }
97  }
98  }
99 }
100 
102 void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
103  Array4<Real const> const& vel,
104  Array4<Real const> const& etax,
105  Array4<Real const> const& kapx,
106  Array4<Real const> const& apx,
107  Array4<EBCellFlag const> const& flag,
108  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
109  Array4<Real const> const& bvxlo,
110  Array4<Real const> const& bvxhi,
112  0,2*AMREX_SPACEDIM,
113  0,AMREX_SPACEDIM> const& bct,
114  Dim3 const& dlo, Dim3 const& dhi) noexcept
115 {
116  const Real dyi = dxinv[1];
117  const auto lo = amrex::lbound(box);
118  const auto hi = amrex::ubound(box);
119  constexpr Real twoThirds = 2./3.;
120 
121  int k = 0;
122  for (int j = lo.y; j <= hi.y; ++j) {
124  for (int i = lo.x; i <= hi.x; ++i) {
125  if (apx(i,j,0) == Real(0.0))
126  {
127  fx(i,j,0,0) = Real(0.0);
128  fx(i,j,0,1) = Real(0.0);
129  }
130  else
131  {
132  int jhip = j + flag(i ,j,0).isConnected(0, 1,0);
133  int jhim = j - flag(i ,j,0).isConnected(0,-1,0);
134  int jlop = j + flag(i-1,j,0).isConnected(0, 1,0);
135  int jlom = j - flag(i-1,j,0).isConnected(0,-1,0);
136  Real whi = mlebtensor_weight(jhip-jhim);
137  Real wlo = mlebtensor_weight(jlop-jlom);
138  Real dudy = mlebtensor_dy_on_xface(i,j,k,0,vel,dyi,
139  bvxlo,bvxhi,bct,dlo,dhi,
140  whi,wlo,jhip,jhim,jlop,jlom);
141  Real dvdy = mlebtensor_dy_on_xface(i,j,k,1,vel,dyi,
142  bvxlo,bvxhi,bct,dlo,dhi,
143  whi,wlo,jhip,jhim,jlop,jlom);
144  Real divu = dvdy;
145  Real xif = kapx(i,j,0);
146  Real mun = Real(0.75)*(etax(i,j,0,0)-xif);// restore the original eta
147  Real mut = etax(i,j,0,1);
148  fx(i,j,0,0) = -mun*(-twoThirds*divu) - xif*divu;
149  fx(i,j,0,1) = -mut*dudy;
150  }
151  }
152  }
153 }
154 
156 void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
157  Array4<Real const> const& vel,
158  Array4<Real const> const& etay,
159  Array4<Real const> const& kapy,
160  Array4<Real const> const& apy,
161  Array4<EBCellFlag const> const& flag,
162  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
163  Array4<Real const> const& bvylo,
164  Array4<Real const> const& bvyhi,
166  0,2*AMREX_SPACEDIM,
167  0,AMREX_SPACEDIM> const& bct,
168  Dim3 const& dlo, Dim3 const& dhi) noexcept
169 {
170  const Real dxi = dxinv[0];
171  const auto lo = amrex::lbound(box);
172  const auto hi = amrex::ubound(box);
173  constexpr Real twoThirds = 2./3.;
174 
175  int k = 0;
176  for (int j = lo.y; j <= hi.y; ++j) {
178  for (int i = lo.x; i <= hi.x; ++i) {
179  if (apy(i,j,0) == Real(0.0))
180  {
181  fy(i,j,0,0) = Real(0.0);
182  fy(i,j,0,1) = Real(0.0);
183  }
184  else
185  {
186  int ihip = i + flag(i,j ,0).isConnected( 1,0,0);
187  int ihim = i - flag(i,j ,0).isConnected(-1,0,0);
188  int ilop = i + flag(i,j-1,0).isConnected( 1,0,0);
189  int ilom = i - flag(i,j-1,0).isConnected(-1,0,0);
190  Real whi = mlebtensor_weight(ihip-ihim);
191  Real wlo = mlebtensor_weight(ilop-ilom);
192  Real dudx = mlebtensor_dx_on_yface(i,j,k,0,vel,dxi,
193  bvylo,bvyhi,bct,dlo,dhi,
194  whi,wlo,ihip,ihim,ilop,ilom);
195  Real dvdx = mlebtensor_dx_on_yface(i,j,k,1,vel,dxi,
196  bvylo,bvyhi,bct,dlo,dhi,
197  whi,wlo,ihip,ihim,ilop,ilom);
198  Real divu = dudx;
199  Real xif = kapy(i,j,0);
200  Real mun = Real(0.75)*(etay(i,j,0,1)-xif);// restore the original eta
201  Real mut = etay(i,j,0,0);
202  fy(i,j,0,0) = -mut*dvdx;
203  fy(i,j,0,1) = -mun*(-twoThirds*divu) - xif*divu;
204  }
205  }
206  }
207 }
208 
210 void mlebtensor_cross_terms (Box const& box, Array4<Real> const& Ax,
211  Array4<Real const> const& fx,
212  Array4<Real const> const& fy,
213  Array4<Real const> const& vel,
214  Array4<Real const> const& velb,
215  Array4<Real const> const& etab,
216  Array4<Real const> const& kapb,
217  Array4<int const> const& ccm,
218  Array4<EBCellFlag const> const& flag,
219  Array4<Real const> const& vol,
220  Array4<Real const> const& apx,
221  Array4<Real const> const& apy,
222  Array4<Real const> const& fcx,
223  Array4<Real const> const& fcy,
224  Array4<Real const> const& bc,
225  bool is_dirichlet, bool is_inhomog,
226  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
227  Real bscalar) noexcept
228 {
229  const Real dxi = dxinv[0];
230  const Real dyi = dxinv[1];
231  const auto lo = amrex::lbound(box);
232  const auto hi = amrex::ubound(box);
233 
234  for (int j = lo.y; j <= hi.y; ++j) {
236  for (int i = lo.x; i <= hi.x; ++i) {
237  if (flag(i,j,0).isRegular())
238  {
239  Ax(i,j,0,0) += bscalar*(dxi*(fx(i+1,j ,0,0) - fx(i,j,0,0))
240  + dyi*(fy(i ,j+1,0,0) - fy(i,j,0,0)));
241  Ax(i,j,0,1) += bscalar*(dxi*(fx(i+1,j ,0,1) - fx(i,j,0,1))
242  + dyi*(fy(i ,j+1,0,1) - fy(i,j,0,1)));
243  }
244  else if (flag(i,j,0).isSingleValued())
245  {
246  Real fxm_0 = fx(i,j,0,0);
247  Real fxm_1 = fx(i,j,0,1);
248  if (apx(i,j,0) > Real(0.0) && apx(i,j,0) < Real(1.0)) {
249  int jj = j + (fcx(i,j,0,0) >= Real(0.0) ? 1 : -1);
250  Real fracy = (ccm(i-1,jj,0) || ccm(i,jj,0)) ? std::abs(fcx(i,j,0,0)) : Real(0.0);
251  fxm_0 = (Real(1.0)-fracy)*fxm_0 + fracy*fx(i,jj,0,0);
252  fxm_1 = (Real(1.0)-fracy)*fxm_1 + fracy*fx(i,jj,0,1);
253  }
254 
255  Real fxp_0 = fx(i+1,j,0,0);
256  Real fxp_1 = fx(i+1,j,0,1);
257  if (apx(i+1,j,0) > Real(0.0) && apx(i+1,j,0) < Real(1.0)) {
258  int jj = j + (fcx(i+1,j,0,0) >= Real(0.0) ? 1 : -1);
259  Real fracy = (ccm(i,jj,0) || ccm(i+1,jj,0)) ? std::abs(fcx(i+1,j,0,0)) : Real(0.0);
260  fxp_0 = (Real(1.0)-fracy)*fxp_0 + fracy*fx(i+1,jj,0,0);
261  fxp_1 = (Real(1.0)-fracy)*fxp_1 + fracy*fx(i+1,jj,0,1);
262  }
263 
264  Real fym_0 = fy(i,j,0,0);
265  Real fym_1 = fy(i,j,0,1);
266  if (apy(i,j,0) > Real(0.0) && apy(i,j,0) < Real(1.0)) {
267  int ii = i + (fcy(i,j,0,0) >= Real(0.0) ? 1 : -1);
268  Real fracx = (ccm(ii,j-1,0) || ccm(ii,j,0)) ? std::abs(fcy(i,j,0,0)) : Real(0.0);
269  fym_0 = (Real(1.0)-fracx)*fym_0 + fracx*fy(ii,j,0,0);
270  fym_1 = (Real(1.0)-fracx)*fym_1 + fracx*fy(ii,j,0,1);
271  }
272 
273  Real fyp_0 = fy(i,j+1,0,0);
274  Real fyp_1 = fy(i,j+1,0,1);
275  if (apy(i,j+1,0) > Real(0.0) && apy(i,j+1,0) < Real(1.0)) {
276  int ii = i + (fcy(i,j+1,0,0) >= Real(0.0) ? 1 : -1);
277  Real fracx = (ccm(ii,j,0) || ccm(ii,j+1,0)) ? std::abs(fcy(i,j+1,0,0)) : Real(0.0);
278  fyp_0 = (Real(1.0)-fracx)*fyp_0 + fracx*fy(ii,j+1,0,0);
279  fyp_1 = (Real(1.0)-fracx)*fyp_1 + fracx*fy(ii,j+1,0,1);
280  }
281 
282  Real kappa = vol(i,j,0);
283  Real feb_0 = Real(0.0);
284  Real feb_1 = Real(0.0);
285  if (is_dirichlet)
286  {
287  Real dapx = apx(i+1,j,0)-apx(i,j,0);
288  Real dapy = apy(i,j+1,0)-apy(i,j,0);
289  Real anorminv = Real(1.0)/std::sqrt(dapx*dapx+dapy*dapy);
290  Real anrmx = -dapx * anorminv;
291  Real anrmy = -dapy * anorminv;
292 
293  Real velb_0 = 0, velb_1 = 0;
294 
295  if (is_inhomog) {
296  velb_0 = velb(i,j,0,0);
297  velb_1 = velb(i,j,0,1);
298  }
299 
300  Real dx_eb = amrex::max(Real(0.3), (kappa*kappa-Real(0.25))/(Real(2.)*kappa));
301  Real dg = dx_eb / amrex::max(std::abs(anrmx),std::abs(anrmy));
302  Real dginv = Real(1.0)/dg;
303  Real gx = bc(i,j,0,0) - dg*anrmx;
304  Real gy = bc(i,j,0,1) - dg*anrmy;
305  int isx = (anrmx > Real(0.0)) ? 1 : -1;
306  int isy = (anrmy > Real(0.0)) ? 1 : -1;
307  int ii = i - isx;
308  int jj = j - isy;
309 
310  gx *= isx;
311  gy *= isy;
312  Real gxy = gx*gy;
313  Real oneggg = Real(1.0)+gx+gy+gxy;
314 
315  Real velg = oneggg * vel(i ,j ,0,0)
316  + (-gy - gxy) * vel(i ,jj,0,0)
317  + (-gx - gxy) * vel(ii,j ,0,0)
318  + gxy * vel(ii,jj,0,0);
319  Real dudn = (velb_0-velg) * dginv;
320 
321  velg = oneggg * vel(i ,j ,0,1)
322  + (-gy - gxy) * vel(i ,jj,0,1)
323  + (-gx - gxy) * vel(ii,j ,0,1)
324  + gxy * vel(ii,jj,0,1);
325  Real dvdn = (velb_1-velg) * dginv;
326 
327  // transverse derivatives are zero on the wall
328  Real dudx = dudn * anrmx;
329  Real dudy = dudn * anrmy;
330  Real dvdx = dvdn * anrmx;
331  Real dvdy = dvdn * anrmy;
332  Real divu = dudx + dvdy;
333  Real xi = kapb(i,j,0);
334  Real mu = etab(i,j,0);
335  Real tautmp = (xi-(2./3.)*mu)*divu;
336  // Note that mu*(grad v) has been included already in MLEBABecLap
337  Real tauxx = mu*dudx + tautmp;
338  Real tauyx = mu*dvdx;
339  Real tauyy = mu*dvdy + tautmp;
340  Real tauxy = mu*dudy;
341  // assuming dxi == dyi == dzi
342  feb_0 = dxi*(dapx*tauxx + dapy*tauyx);
343  feb_1 = dxi*(dapx*tauxy + dapy*tauyy);
344  }
345 
346  Real volinv = bscalar / kappa;
347  Ax(i,j,0,0) += volinv * (dxi*(apx(i,j,0)*fxm_0-apx(i+1,j,0)*fxp_0)
348  +dyi*(apy(i,j,0)*fym_0-apy(i,j+1,0)*fyp_0)
349  -dxi* feb_0);
350  Ax(i,j,0,1) += volinv * (dxi*(apx(i,j,0)*fxm_1-apx(i+1,j,0)*fxp_1)
351  +dyi*(apy(i,j,0)*fym_1-apy(i,j+1,0)*fyp_1)
352  -dxi* feb_1);
353  }
354  }
355  }
356 }
357 
358 
360 void mlebtensor_flux_0 (Box const& box,
361  Array4<Real> const& Ax,
362  Array4<Real const> const& fx,
363  Array4<Real const> const& ap,
364  Real bscalar) noexcept
365 {
366  const auto lo = amrex::lbound(box);
367  const auto hi = amrex::ubound(box);
368 
369  int k = 0;
370  for (int j = lo.y; j <= hi.y; ++j) {
372  for (int i = lo.x; i <= hi.x; ++i) {
373  if (ap(i,j,k) != Real(0.0)) {
374  for (int n=0; n<AMREX_SPACEDIM; n++) {
375  Ax(i,j,k,n) += bscalar*fx(i,j,k,n);
376  }
377  }
378  //else MLEBABec::compFlux already set Ax = 0
379  }
380  }
381 }
382 
383 
385 void mlebtensor_flux_x (Box const& box, Array4<Real> const& Ax,
386  Array4<Real const> const& fx, Array4<Real const> const& apx,
387  Array4<Real const> const& fcx,
388  Real const bscalar, Array4<int const> const& ccm,
389  int face_only, Box const& xbox) noexcept
390 {
391  int lof = xbox.smallEnd(0);
392  int hif = xbox.bigEnd(0);
393  amrex::LoopConcurrent(box, [=] (int i, int j, int k) noexcept
394  {
395  if (!face_only || lof == i || hif == i) {
396  if (apx(i,j,k) == Real(1.0)) {
397  for (int n=0; n<AMREX_SPACEDIM; n++) {
398  Ax(i,j,k,n) += bscalar*fx(i,j,k,n);
399  }
400  }
401  else if (apx(i,j,k) != 0.) {
402  Real fxm_0 = fx(i,j,0,0);
403  Real fxm_1 = fx(i,j,0,1);
404 
405  int jj = j + (fcx(i,j,0,0) >= Real(0.0) ? 1 : -1);
406  Real fracy = (ccm(i-1,jj,0) || ccm(i,jj,0)) ? std::abs(fcx(i,j,0,0)) : Real(0.0);
407  fxm_0 = (Real(1.0)-fracy)*fxm_0 + fracy*fx(i,jj,0,0);
408  fxm_1 = (Real(1.0)-fracy)*fxm_1 + fracy*fx(i,jj,0,1);
409 
410  Ax(i,j,k,0) += bscalar*fxm_0;
411  Ax(i,j,k,1) += bscalar*fxm_1;
412  }
413  //else MLEBABec::compFlux already set Ax = 0
414  }
415  });
416 }
417 
419 void mlebtensor_flux_y (Box const& box, Array4<Real> const& Ay,
420  Array4<Real const> const& fy, Array4<Real const> const& apy,
421  Array4<Real const> const& fcy,
422  Real const bscalar, Array4<int const> const& ccm,
423  int face_only, Box const& ybox) noexcept
424 {
425  int lof = ybox.smallEnd(1);
426  int hif = ybox.bigEnd(1);
427  amrex::LoopConcurrent(box, [=] (int i, int j, int k) noexcept
428  {
429  if (!face_only || lof == j || hif == j) {
430  if (apy(i,j,k) == Real(1.0)) {
431  for (int n=0; n<AMREX_SPACEDIM; n++) {
432  Ay(i,j,k,n) += bscalar*fy(i,j,k,n);
433  }
434  }
435  else if (apy(i,j,k) != 0.) {
436  Real fym_0 = fy(i,j,0,0);
437  Real fym_1 = fy(i,j,0,1);
438 
439  int ii = i + (fcy(i,j,0,0) >= Real(0.0) ? 1 : -1);
440  Real fracx = (ccm(ii,j-1,0) || ccm(ii,j,0)) ? std::abs(fcy(i,j,0,0)) : Real(0.0);
441  fym_0 = (Real(1.0)-fracx)*fym_0 + fracx*fy(ii,j,0,0);
442  fym_1 = (Real(1.0)-fracx)*fym_1 + fracx*fy(ii,j,0,1);
443 
444  Ay(i,j,k,0) += bscalar*fym_0;
445  Ay(i,j,k,1) += bscalar*fym_1;
446  }
447  //else MLEBABec::compFlux already set Ay = 0
448  }
449  });
450 }
451 
453 void mlebtensor_vel_grads_fx (Box const& box, Array4<Real> const& fx,
454  Array4<Real const> const& vel, Array4<Real const> const& apx,
455  Array4<EBCellFlag const> const& flag,
456  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
457 {
458  const Real dxi = dxinv[0];
459  const Real dyi = dxinv[1];
460  const auto lo = amrex::lbound(box);
461  const auto hi = amrex::ubound(box);
462 
463  int k = 0;
464  for (int j = lo.y; j <= hi.y; ++j) {
466  for (int i = lo.x; i <= hi.x; ++i) {
467  if (apx(i,j,0) == Real(0.0))
468  {
469  fx(i,j,0,0) = Real(0.0);
470  fx(i,j,0,1) = Real(0.0);
471  fx(i,j,0,2) = Real(0.0);
472  fx(i,j,0,3) = Real(0.0);
473  }
474  else
475  {
476  Real dudx = (vel(i,j,0,0) - vel(i-1,j,0,0))*dxi;
477  Real dvdx = (vel(i,j,0,1) - vel(i-1,j,0,1))*dxi;
478 
479  int jhip = j + flag(i ,j,0).isConnected(0, 1,0);
480  int jhim = j - flag(i ,j,0).isConnected(0,-1,0);
481  int jlop = j + flag(i-1,j,0).isConnected(0, 1,0);
482  int jlom = j - flag(i-1,j,0).isConnected(0,-1,0);
483  Real whi = mlebtensor_weight(jhip-jhim);
484  Real wlo = mlebtensor_weight(jlop-jlom);
485  Real dudy = mlebtensor_dy_on_xface(i,j,k,0,vel,dyi,
486  whi,wlo,jhip,jhim,jlop,jlom);
487  Real dvdy = mlebtensor_dy_on_xface(i,j,k,1,vel,dyi,
488  whi,wlo,jhip,jhim,jlop,jlom);
489  fx(i,j,0,0) = dudx;
490  fx(i,j,0,1) = dvdx;
491  fx(i,j,0,2) = dudy;
492  fx(i,j,0,3) = dvdy;
493  }
494  }
495  }
496 }
497 
499 void mlebtensor_vel_grads_fy (Box const& box, Array4<Real> const& fy,
500  Array4<Real const> const& vel, Array4<Real const> const& apy,
501  Array4<EBCellFlag const> const& flag,
502  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
503 {
504  const Real dxi = dxinv[0];
505  const Real dyi = dxinv[1];
506  const auto lo = amrex::lbound(box);
507  const auto hi = amrex::ubound(box);
508 
509  int k = 0;
510  for (int j = lo.y; j <= hi.y; ++j) {
512  for (int i = lo.x; i <= hi.x; ++i) {
513  if (apy(i,j,0) == 0.) {
514  fy(i,j,0,0) = Real(0.0);
515  fy(i,j,0,1) = Real(0.0);
516  fy(i,j,0,2) = Real(0.0);
517  fy(i,j,0,3) = Real(0.0);
518  }
519  else
520  {
521  int ihip = i + flag(i,j ,0).isConnected( 1,0,0);
522  int ihim = i - flag(i,j ,0).isConnected(-1,0,0);
523  int ilop = i + flag(i,j-1,0).isConnected( 1,0,0);
524  int ilom = i - flag(i,j-1,0).isConnected(-1,0,0);
525  Real whi = mlebtensor_weight(ihip-ihim);
526  Real wlo = mlebtensor_weight(ilop-ilom);
527  Real dudx = mlebtensor_dx_on_yface(i,j,k,0,vel,dxi,
528  whi,wlo,ihip,ihim,ilop,ilom);
529  Real dvdx = mlebtensor_dx_on_yface(i,j,k,1,vel,dxi,
530  whi,wlo,ihip,ihim,ilop,ilom);
531 
532  Real dudy = (vel(i,j,0,0) - vel(i,j-1,0,0))*dyi;
533  Real dvdy = (vel(i,j,0,1) - vel(i,j-1,0,1))*dyi;
534  fy(i,j,0,0) = dudx;
535  fy(i,j,0,1) = dvdx;
536  fy(i,j,0,2) = dudy;
537  fy(i,j,0,3) = dvdy;
538  }
539  }
540  }
541 }
542 
544 void mlebtensor_vel_grads_fx (Box const& box, Array4<Real> const& fx,
545  Array4<Real const> const& vel, Array4<Real const> const& apx,
546  Array4<EBCellFlag const> const& flag,
547  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
548  Array4<Real const> const& bvxlo,
549  Array4<Real const> const& bvxhi,
551  0,2*AMREX_SPACEDIM,
552  0,AMREX_SPACEDIM> const& bct,
553  Dim3 const& dlo, Dim3 const& dhi) noexcept
554 {
555  const Real dxi = dxinv[0];
556  const Real dyi = dxinv[1];
557  const auto lo = amrex::lbound(box);
558  const auto hi = amrex::ubound(box);
559 
560  int k = 0;
561  for (int j = lo.y; j <= hi.y; ++j) {
562  for (int i = lo.x; i <= hi.x; ++i) {
563  if (apx(i,j,0) == Real(0.0))
564  {
565  fx(i,j,0,0) = Real(0.0);
566  fx(i,j,0,1) = Real(0.0);
567  fx(i,j,0,2) = Real(0.0);
568  fx(i,j,0,3) = Real(0.0);
569  }
570  else
571  {
572  Real dudx = (vel(i,j,0,0) - vel(i-1,j,0,0))*dxi;
573  Real dvdx = (vel(i,j,0,1) - vel(i-1,j,0,1))*dxi;
574 
575  int jhip = j + flag(i ,j,0).isConnected(0, 1,0);
576  int jhim = j - flag(i ,j,0).isConnected(0,-1,0);
577  int jlop = j + flag(i-1,j,0).isConnected(0, 1,0);
578  int jlom = j - flag(i-1,j,0).isConnected(0,-1,0);
579  Real whi = mlebtensor_weight(jhip-jhim);
580  Real wlo = mlebtensor_weight(jlop-jlom);
581  Real dudy = mlebtensor_dy_on_xface(i,j,k,0,vel,dyi,
582  bvxlo,bvxhi,bct,dlo,dhi,
583  whi,wlo,jhip,jhim,jlop,jlom);
584  Real dvdy = mlebtensor_dy_on_xface(i,j,k,1,vel,dyi,
585  bvxlo,bvxhi,bct,dlo,dhi,
586  whi,wlo,jhip,jhim,jlop,jlom);
587  fx(i,j,0,0) = dudx;
588  fx(i,j,0,1) = dvdx;
589  fx(i,j,0,2) = dudy;
590  fx(i,j,0,3) = dvdy;
591  }
592  }
593  }
594 }
595 
597 void mlebtensor_vel_grads_fy (Box const& box, Array4<Real> const& fy,
598  Array4<Real const> const& vel, Array4<Real const> const& apy,
599  Array4<EBCellFlag const> const& flag,
600  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
601  Array4<Real const> const& bvylo,
602  Array4<Real const> const& bvyhi,
604  0,2*AMREX_SPACEDIM,
605  0,AMREX_SPACEDIM> const& bct,
606  Dim3 const& dlo, Dim3 const& dhi) noexcept
607 {
608  const Real dxi = dxinv[0];
609  const Real dyi = dxinv[1];
610  const auto lo = amrex::lbound(box);
611  const auto hi = amrex::ubound(box);
612 
613  int k = 0;
614  for (int j = lo.y; j <= hi.y; ++j) {
615  for (int i = lo.x; i <= hi.x; ++i) {
616  if (apy(i,j,0) == 0.) {
617  fy(i,j,0,0) = Real(0.0);
618  fy(i,j,0,1) = Real(0.0);
619  fy(i,j,0,2) = Real(0.0);
620  fy(i,j,0,3) = Real(0.0);
621  }
622  else
623  {
624  int ihip = i + flag(i,j ,0).isConnected( 1,0,0);
625  int ihim = i - flag(i,j ,0).isConnected(-1,0,0);
626  int ilop = i + flag(i,j-1,0).isConnected( 1,0,0);
627  int ilom = i - flag(i,j-1,0).isConnected(-1,0,0);
628  Real whi = mlebtensor_weight(ihip-ihim);
629  Real wlo = mlebtensor_weight(ilop-ilom);
630  Real dudx = mlebtensor_dx_on_yface(i,j,k,0,vel,dxi,
631  bvylo,bvyhi,bct,dlo,dhi,
632  whi,wlo,ihip,ihim,ilop,ilom);
633  Real dvdx = mlebtensor_dx_on_yface(i,j,k,1,vel,dxi,
634  bvylo,bvyhi,bct,dlo,dhi,
635  whi,wlo,ihip,ihim,ilop,ilom);
636 
637  Real dudy = (vel(i,j,0,0) - vel(i,j-1,0,0))*dyi;
638  Real dvdy = (vel(i,j,0,1) - vel(i,j-1,0,1))*dyi;
639  fy(i,j,0,0) = dudx;
640  fy(i,j,0,1) = dvdx;
641  fy(i,j,0,2) = dudy;
642  fy(i,j,0,3) = dvdy;
643  }
644  }
645  }
646 }
647 
649 void mlebtensor_vel_grads_fx (Box const& box, Array4<Real> const& fx,
650  Array4<Real const> const& vel, Array4<Real const> const& apx,
651  Array4<EBCellFlag const> const& flag,
652  Array4<int const> const& ccm,
653  Array4<Real const> const& fcx,
654  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
655 {
656  const Real dxi = dxinv[0];
657  const Real dyi = dxinv[1];
658  const auto lo = amrex::lbound(box);
659  const auto hi = amrex::ubound(box);
660 
661  int k = 0;
662  for (int j = lo.y; j <= hi.y; ++j) {
664  for (int i = lo.x; i <= hi.x; ++i) {
665  if (apx(i,j,0) == Real(0.0))
666  {
667  fx(i,j,0,0) = Real(0.0);
668  fx(i,j,0,1) = Real(0.0);
669  fx(i,j,0,2) = Real(0.0);
670  fx(i,j,0,3) = Real(0.0);
671  }
672  else
673  {
674  Real dudx = (vel(i,j,0,0) - vel(i-1,j,0,0))*dxi;
675  Real dvdx = (vel(i,j,0,1) - vel(i-1,j,0,1))*dxi;
676  if (apx(i,j,0) < Real(1.0)) {
677  int jj = j + (fcx(i,j,0,0) >= Real(0.0) ? 1 : -1);
678  Real fracy = (ccm(i-1,jj,0) || ccm(i,jj,0)) ? std::abs(fcx(i,j,0,0)) : Real(0.0);
679  dudx = (Real(1.0)-fracy)*dudx
680  + fracy *(vel(i,jj,0,0) - vel(i-1,jj,0,0))*dxi;
681  dvdx = (Real(1.0)-fracy)*dvdx
682  + fracy *(vel(i,jj,0,1) - vel(i-1,jj,0,1))*dxi;
683  }
684 
685  int jhip = j + flag(i ,j,0).isConnected(0, 1,0);
686  int jhim = j - flag(i ,j,0).isConnected(0,-1,0);
687  int jlop = j + flag(i-1,j,0).isConnected(0, 1,0);
688  int jlom = j - flag(i-1,j,0).isConnected(0,-1,0);
689  Real whi = mlebtensor_weight(jhip-jhim);
690  Real wlo = mlebtensor_weight(jlop-jlom);
691  Real dudy = mlebtensor_dy_on_xface(i,j,k,0,vel,dyi,
692  whi,wlo,jhip,jhim,jlop,jlom);
693  Real dvdy = mlebtensor_dy_on_xface(i,j,k,1,vel,dyi,
694  whi,wlo,jhip,jhim,jlop,jlom);
695  fx(i,j,0,0) = dudx;
696  fx(i,j,0,1) = dvdx;
697  fx(i,j,0,2) = dudy;
698  fx(i,j,0,3) = dvdy;
699  }
700  }
701  }
702 }
703 
705 void mlebtensor_vel_grads_fy (Box const& box, Array4<Real> const& fy,
706  Array4<Real const> const& vel, Array4<Real const> const& apy,
707  Array4<EBCellFlag const> const& flag,
708  Array4<int const> const& ccm,
709  Array4<Real const> const& fcy,
710  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
711 {
712  const Real dxi = dxinv[0];
713  const Real dyi = dxinv[1];
714  const auto lo = amrex::lbound(box);
715  const auto hi = amrex::ubound(box);
716 
717  int k = 0;
718  for (int j = lo.y; j <= hi.y; ++j) {
720  for (int i = lo.x; i <= hi.x; ++i) {
721  if (apy(i,j,0) == 0.) {
722  fy(i,j,0,0) = Real(0.0);
723  fy(i,j,0,1) = Real(0.0);
724  fy(i,j,0,2) = Real(0.0);
725  fy(i,j,0,3) = Real(0.0);
726  }
727  else
728  {
729  int ihip = i + flag(i,j ,0).isConnected( 1,0,0);
730  int ihim = i - flag(i,j ,0).isConnected(-1,0,0);
731  int ilop = i + flag(i,j-1,0).isConnected( 1,0,0);
732  int ilom = i - flag(i,j-1,0).isConnected(-1,0,0);
733  Real whi = mlebtensor_weight(ihip-ihim);
734  Real wlo = mlebtensor_weight(ilop-ilom);
735  Real dudx = mlebtensor_dx_on_yface(i,j,k,0,vel,dxi,
736  whi,wlo,ihip,ihim,ilop,ilom);
737  Real dvdx = mlebtensor_dx_on_yface(i,j,k,1,vel,dxi,
738  whi,wlo,ihip,ihim,ilop,ilom);
739 
740  Real dudy = (vel(i,j,0,0) - vel(i,j-1,0,0))*dyi;
741  Real dvdy = (vel(i,j,0,1) - vel(i,j-1,0,1))*dyi;
742  if (apy(i,j,0) < Real(1.0)) {
743  int ii = i + (fcy(i,j,0,0) >= Real(0.0) ? 1 : -1);
744  Real fracx = (ccm(ii,j-1,0) || ccm(ii,j,0)) ? std::abs(fcy(i,j,0,0)) : Real(0.0);
745  dudy = (Real(1.0)-fracx)*dudy
746  + fracx *(vel(ii,j,0,0) - vel(ii,j-1,0,0))*dyi;
747  dvdy = (Real(1.0)-fracx)*dvdy
748  + fracx *(vel(ii,j,0,1) - vel(ii,j-1,0,1))*dyi;
749  }
750 
751  fy(i,j,0,0) = dudx;
752  fy(i,j,0,1) = dvdx;
753  fy(i,j,0,2) = dudy;
754  fy(i,j,0,3) = dvdy;
755  }
756  }
757  }
758 }
759 
761 void mlebtensor_vel_grads_fx (Box const& box, Array4<Real> const& fx,
762  Array4<Real const> const& vel, Array4<Real const> const& apx,
763  Array4<EBCellFlag const> const& flag,
764  Array4<int const> const& ccm,
765  Array4<Real const> const& fcx,
766  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
767  Array4<Real const> const& bvxlo,
768  Array4<Real const> const& bvxhi,
770  0,2*AMREX_SPACEDIM,
771  0,AMREX_SPACEDIM> const& bct,
772  Dim3 const& dlo, Dim3 const& dhi) noexcept
773 {
774  const Real dxi = dxinv[0];
775  const Real dyi = dxinv[1];
776  const auto lo = amrex::lbound(box);
777  const auto hi = amrex::ubound(box);
778 
779  int k = 0;
780  for (int j = lo.y; j <= hi.y; ++j) {
781  for (int i = lo.x; i <= hi.x; ++i) {
782  if (apx(i,j,0) == Real(0.0))
783  {
784  fx(i,j,0,0) = Real(0.0);
785  fx(i,j,0,1) = Real(0.0);
786  fx(i,j,0,2) = Real(0.0);
787  fx(i,j,0,3) = Real(0.0);
788  }
789  else
790  {
791  Real dudx = (vel(i,j,0,0) - vel(i-1,j,0,0))*dxi;
792  Real dvdx = (vel(i,j,0,1) - vel(i-1,j,0,1))*dxi;
793  if (apx(i,j,0) < Real(1.0)) {
794  int jj = j + (fcx(i,j,0,0) >= Real(0.0) ? 1 : -1);
795  Real fracy = (ccm(i-1,jj,0) || ccm(i,jj,0)) ? std::abs(fcx(i,j,0,0)) : Real(0.0);
796  dudx = (Real(1.0)-fracy)*dudx
797  + fracy *(vel(i,jj,0,0) - vel(i-1,jj,0,0))*dxi;
798  dvdx = (Real(1.0)-fracy)*dvdx
799  + fracy *(vel(i,jj,0,1) - vel(i-1,jj,0,1))*dxi;
800  }
801 
802  int jhip = j + flag(i ,j,0).isConnected(0, 1,0);
803  int jhim = j - flag(i ,j,0).isConnected(0,-1,0);
804  int jlop = j + flag(i-1,j,0).isConnected(0, 1,0);
805  int jlom = j - flag(i-1,j,0).isConnected(0,-1,0);
806  Real whi = mlebtensor_weight(jhip-jhim);
807  Real wlo = mlebtensor_weight(jlop-jlom);
808  Real dudy = mlebtensor_dy_on_xface(i,j,k,0,vel,dyi,
809  bvxlo,bvxhi,bct,dlo,dhi,
810  whi,wlo,jhip,jhim,jlop,jlom);
811  Real dvdy = mlebtensor_dy_on_xface(i,j,k,1,vel,dyi,
812  bvxlo,bvxhi,bct,dlo,dhi,
813  whi,wlo,jhip,jhim,jlop,jlom);
814  fx(i,j,0,0) = dudx;
815  fx(i,j,0,1) = dvdx;
816  fx(i,j,0,2) = dudy;
817  fx(i,j,0,3) = dvdy;
818  }
819  }
820  }
821 }
822 
824 void mlebtensor_vel_grads_fy (Box const& box, Array4<Real> const& fy,
825  Array4<Real const> const& vel, Array4<Real const> const& apy,
826  Array4<EBCellFlag const> const& flag,
827  Array4<int const> const& ccm,
828  Array4<Real const> const& fcy,
829  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
830  Array4<Real const> const& bvylo,
831  Array4<Real const> const& bvyhi,
833  0,2*AMREX_SPACEDIM,
834  0,AMREX_SPACEDIM> const& bct,
835  Dim3 const& dlo, Dim3 const& dhi) noexcept
836 {
837  const Real dxi = dxinv[0];
838  const Real dyi = dxinv[1];
839  const auto lo = amrex::lbound(box);
840  const auto hi = amrex::ubound(box);
841 
842  int k = 0;
843  for (int j = lo.y; j <= hi.y; ++j) {
844  for (int i = lo.x; i <= hi.x; ++i) {
845  if (apy(i,j,0) == 0.) {
846  fy(i,j,0,0) = Real(0.0);
847  fy(i,j,0,1) = Real(0.0);
848  fy(i,j,0,2) = Real(0.0);
849  fy(i,j,0,3) = Real(0.0);
850  }
851  else
852  {
853  int ihip = i + flag(i,j ,0).isConnected( 1,0,0);
854  int ihim = i - flag(i,j ,0).isConnected(-1,0,0);
855  int ilop = i + flag(i,j-1,0).isConnected( 1,0,0);
856  int ilom = i - flag(i,j-1,0).isConnected(-1,0,0);
857  Real whi = mlebtensor_weight(ihip-ihim);
858  Real wlo = mlebtensor_weight(ilop-ilom);
859  Real dudx = mlebtensor_dx_on_yface(i,j,k,0,vel,dxi,
860  bvylo,bvyhi,bct,dlo,dhi,
861  whi,wlo,ihip,ihim,ilop,ilom);
862  Real dvdx = mlebtensor_dx_on_yface(i,j,k,1,vel,dxi,
863  bvylo,bvyhi,bct,dlo,dhi,
864  whi,wlo,ihip,ihim,ilop,ilom);
865 
866  Real dudy = (vel(i,j,0,0) - vel(i,j-1,0,0))*dyi;
867  Real dvdy = (vel(i,j,0,1) - vel(i,j-1,0,1))*dyi;
868  if (apy(i,j,0) < Real(1.0)) {
869  int ii = i + (fcy(i,j,0,0) >= Real(0.0) ? 1 : -1);
870  Real fracx = (ccm(ii,j-1,0) || ccm(ii,j,0)) ? std::abs(fcy(i,j,0,0)) : Real(0.0);
871  dudy = (Real(1.0)-fracx)*dudy
872  + fracx *(vel(ii,j,0,0) - vel(ii,j-1,0,0))*dyi;
873  dvdy = (Real(1.0)-fracx)*dvdy
874  + fracx *(vel(ii,j,0,1) - vel(ii,j-1,0,1))*dyi;
875  }
876 
877  fy(i,j,0,0) = dudx;
878  fy(i,j,0,1) = dvdx;
879  fy(i,j,0,2) = dudy;
880  fy(i,j,0,3) = dvdy;
881  }
882  }
883  }
884 }
885 
886 }
887 #endif
#define AMREX_PRAGMA_SIMD
Definition: AMReX_Extension.H:80
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
Maintain an identifier for boundary condition types.
Definition: AMReX_BoundCond.H:20
Definition: AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mlebtensor_weight(int d)
Definition: AMReX_MLEBTensor_K.H:10
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & max(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:35
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_vel_grads_fx(Box const &box, Array4< Real > const &fx, Array4< Real const > const &vel, Array4< Real const > const &apx, Array4< EBCellFlag const > const &flag, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:453
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE T abs(const GpuComplex< T > &a_z) noexcept
Return the absolute value of a complex number.
Definition: AMReX_GpuComplex.H:356
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_flux_y(Box const &box, Array4< Real > const &Ay, Array4< Real const > const &fy, Array4< Real const > const &apy, Array4< Real const > const &fcy, Real const bscalar, Array4< int const > const &ccm, int face_only, Box const &ybox) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:419
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_cross_terms_fy(Box const &box, Array4< Real > const &fy, Array4< Real const > const &vel, Array4< Real const > const &etay, Array4< Real const > const &kapy, Array4< Real const > const &apy, Array4< EBCellFlag const > const &flag, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:56
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_cross_terms(Box const &box, Array4< Real > const &Ax, Array4< Real const > const &fx, Array4< Real const > const &fy, Array4< Real const > const &vel, Array4< Real const > const &velb, Array4< Real const > const &etab, Array4< Real const > const &kapb, Array4< int const > const &ccm, Array4< EBCellFlag const > const &flag, Array4< Real const > const &vol, Array4< Real const > const &apx, Array4< Real const > const &apy, Array4< Real const > const &fcx, Array4< Real const > const &fcy, Array4< Real const > const &bc, bool is_dirichlet, bool is_inhomog, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Real bscalar) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:210
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_vel_grads_fy(Box const &box, Array4< Real > const &fy, Array4< Real const > const &vel, Array4< Real const > const &apy, Array4< EBCellFlag const > const &flag, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:499
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 Real mlebtensor_dx_on_yface(int, int j, int k, int n, Array4< Real const > const &vel, Real dxi, Real whi, Real wlo, int ihip, int ihim, int ilop, int ilom) noexcept
Definition: AMReX_MLEBTensor_K.H:25
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mlebtensor_dy_on_xface(int i, int, int k, int n, Array4< Real const > const &vel, Real dyi, Real whi, Real wlo, int jhip, int jhim, int jlop, int jlom) noexcept
Definition: AMReX_MLEBTensor_K.H:15
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_cross_terms_fx(Box const &box, Array4< Real > const &fx, Array4< Real const > const &vel, Array4< Real const > const &etax, Array4< Real const > const &kapx, Array4< Real const > const &apx, Array4< EBCellFlag const > const &flag, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:10
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_flux_x(Box const &box, Array4< Real > const &Ax, Array4< Real const > const &fx, Array4< Real const > const &apx, Array4< Real const > const &fcx, Real const bscalar, Array4< int const > const &ccm, int face_only, Box const &xbox) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:385
AMREX_GPU_HOST_DEVICE AMREX_ATTRIBUTE_FLATTEN_FOR void LoopConcurrent(Dim3 lo, Dim3 hi, F const &f) noexcept
Definition: AMReX_Loop.H:150
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mlebtensor_flux_0(Box const &box, Array4< Real > const &Ax, Array4< Real const > const &fx, Array4< Real const > const &ap, Real bscalar) noexcept
Definition: AMReX_MLEBTensor_2D_K.H:360
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
Definition: AMReX_Array.H:282
Definition: AMReX_Array4.H:61
Definition: AMReX_Dim3.H:12