Block-Structured AMR Software Framework
AMReX_MLTensor_3D_K.H
Go to the documentation of this file.
1 #ifndef AMREX_MLTENSOR_3D_K_H_
2 #define AMREX_MLTENSOR_3D_K_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_MLLinOp_K.H>
6 
7 namespace amrex {
8 
9 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
10 void mltensor_fill_edges_xlo_ylo (int const i, int const j, int const k, Dim3 const& blen,
11  Array4<Real> const& vel,
12  Array4<int const> const& mxlo,
13  Array4<int const> const& mylo,
14  Array4<Real const> const& bcvalxlo,
15  Array4<Real const> const& bcvalylo,
16  Array2D<BoundCond,
17  0,2*AMREX_SPACEDIM,
18  0,AMREX_SPACEDIM> const& bct,
19  Array2D<Real,
20  0,2*AMREX_SPACEDIM,
21  0,AMREX_SPACEDIM> const& bcl,
22  int inhomog, int maxorder,
23  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
24  bool xlo_domain, bool ylo_domain) noexcept
25 {
26  if (mxlo(i,j,k) != BndryData::covered && (!xlo_domain || !ylo_domain)) {
27  bool x_interior = mylo(i+1,j ,k) == BndryData::covered;
28  bool x_exterior = mylo(i+1,j ,k) == BndryData::not_covered;
29  bool y_interior = mxlo(i ,j+1,k) == BndryData::covered;
30  bool y_exterior = mxlo(i ,j+1,k) == BndryData::not_covered;
31  if ((x_interior && y_interior) || (x_exterior && y_exterior)) {
32  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
33  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
34  bct(Orientation::xlo(), icomp),
35  bcl(Orientation::xlo(), icomp),
36  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
37  Real tmp = vel(i,j,k,icomp);
38  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
39  bct(Orientation::ylo(), icomp),
40  bcl(Orientation::ylo(), icomp),
41  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
42  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
43  }
44  } else if (x_interior || xlo_domain) {
45  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
46  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
47  bct(Orientation::xlo(), icomp),
48  bcl(Orientation::xlo(), icomp),
49  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
50  }
51  } else if (y_interior || ylo_domain) {
52  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
53  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
54  bct(Orientation::ylo(), icomp),
55  bcl(Orientation::ylo(), icomp),
56  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
57  }
58  }
59  }
60 }
61 
62 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
63 void mltensor_fill_edges_xhi_ylo (int const i, int const j, int const k, Dim3 const& blen,
64  Array4<Real> const& vel,
65  Array4<int const> const& mxhi,
66  Array4<int const> const& mylo,
67  Array4<Real const> const& bcvalxhi,
68  Array4<Real const> const& bcvalylo,
69  Array2D<BoundCond,
70  0,2*AMREX_SPACEDIM,
71  0,AMREX_SPACEDIM> const& bct,
72  Array2D<Real,
73  0,2*AMREX_SPACEDIM,
74  0,AMREX_SPACEDIM> const& bcl,
75  int inhomog, int maxorder,
76  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
77  bool xhi_domain, bool ylo_domain) noexcept
78 {
79  if (mxhi(i,j,k) != BndryData::covered && (!xhi_domain || !ylo_domain)) {
80  bool x_interior = mylo(i-1,j ,k) == BndryData::covered;
81  bool x_exterior = mylo(i-1,j ,k) == BndryData::not_covered;
82  bool y_interior = mxhi(i ,j+1,k) == BndryData::covered;
83  bool y_exterior = mxhi(i ,j+1,k) == BndryData::not_covered;
84  if ((x_interior && y_interior) || (x_exterior && y_exterior)) {
85  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
86  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
87  bct(Orientation::xhi(), icomp),
88  bcl(Orientation::xhi(), icomp),
89  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
90  Real tmp = vel(i,j,k,icomp);
91  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
92  bct(Orientation::ylo(), icomp),
93  bcl(Orientation::ylo(), icomp),
94  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
95  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
96  }
97  } else if (x_interior || xhi_domain) {
98  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
99  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
100  bct(Orientation::xhi(), icomp),
101  bcl(Orientation::xhi(), icomp),
102  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
103  }
104  } else if (y_interior || ylo_domain) {
105  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
106  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
107  bct(Orientation::ylo(), icomp),
108  bcl(Orientation::ylo(), icomp),
109  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
110  }
111  }
112  }
113 }
114 
115 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
116 void mltensor_fill_edges_xlo_yhi (int const i, int const j, int const k, Dim3 const& blen,
117  Array4<Real> const& vel,
118  Array4<int const> const& mxlo,
119  Array4<int const> const& myhi,
120  Array4<Real const> const& bcvalxlo,
121  Array4<Real const> const& bcvalyhi,
122  Array2D<BoundCond,
123  0,2*AMREX_SPACEDIM,
124  0,AMREX_SPACEDIM> const& bct,
125  Array2D<Real,
126  0,2*AMREX_SPACEDIM,
127  0,AMREX_SPACEDIM> const& bcl,
128  int inhomog, int maxorder,
129  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
130  bool xlo_domain, bool yhi_domain) noexcept
131 {
132  if (mxlo(i,j,k) != BndryData::covered && (!xlo_domain || !yhi_domain)) {
133  bool x_interior = myhi(i+1,j ,k) == BndryData::covered;
134  bool x_exterior = myhi(i+1,j ,k) == BndryData::not_covered;
135  bool y_interior = mxlo(i ,j-1,k) == BndryData::covered;
136  bool y_exterior = mxlo(i ,j-1,k) == BndryData::not_covered;
137  if ((x_interior && y_interior) || (x_exterior && y_exterior)) {
138  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
139  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
140  bct(Orientation::xlo(), icomp),
141  bcl(Orientation::xlo(), icomp),
142  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
143  Real tmp = vel(i,j,k,icomp);
144  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
145  bct(Orientation::yhi(), icomp),
146  bcl(Orientation::yhi(), icomp),
147  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
148  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
149  }
150  } else if (x_interior || xlo_domain) {
151  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
152  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
153  bct(Orientation::xlo(), icomp),
154  bcl(Orientation::xlo(), icomp),
155  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
156  }
157  } else if (y_interior || yhi_domain) {
158  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
159  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
160  bct(Orientation::yhi(), icomp),
161  bcl(Orientation::yhi(), icomp),
162  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
163  }
164  }
165  }
166 }
167 
168 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
169 void mltensor_fill_edges_xhi_yhi (int const i, int const j, int const k, Dim3 const& blen,
170  Array4<Real> const& vel,
171  Array4<int const> const& mxhi,
172  Array4<int const> const& myhi,
173  Array4<Real const> const& bcvalxhi,
174  Array4<Real const> const& bcvalyhi,
175  Array2D<BoundCond,
176  0,2*AMREX_SPACEDIM,
177  0,AMREX_SPACEDIM> const& bct,
178  Array2D<Real,
179  0,2*AMREX_SPACEDIM,
180  0,AMREX_SPACEDIM> const& bcl,
181  int inhomog, int maxorder,
182  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
183  bool xhi_domain, bool yhi_domain) noexcept
184 {
185  if (mxhi(i,j,k) != BndryData::covered && (!xhi_domain || !yhi_domain)) {
186  bool x_interior = myhi(i-1,j ,k) == BndryData::covered;
187  bool x_exterior = myhi(i-1,j ,k) == BndryData::not_covered;
188  bool y_interior = mxhi(i ,j-1,k) == BndryData::covered;
189  bool y_exterior = mxhi(i ,j-1,k) == BndryData::not_covered;
190  if ((x_interior && y_interior) || (x_exterior && y_exterior)) {
191  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
192  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
193  bct(Orientation::xhi(), icomp),
194  bcl(Orientation::xhi(), icomp),
195  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
196  Real tmp = vel(i,j,k,icomp);
197  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
198  bct(Orientation::yhi(), icomp),
199  bcl(Orientation::yhi(), icomp),
200  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
201  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
202  }
203  } else if (x_interior || xhi_domain) {
204  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
205  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
206  bct(Orientation::xhi(), icomp),
207  bcl(Orientation::xhi(), icomp),
208  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
209  }
210  } else if (y_interior || yhi_domain) {
211  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
212  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
213  bct(Orientation::yhi(), icomp),
214  bcl(Orientation::yhi(), icomp),
215  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
216  }
217  }
218  }
219 }
220 
221 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
222 void mltensor_fill_edges_xlo_zlo (int const i, int const j, int const k, Dim3 const& blen,
223  Array4<Real> const& vel,
224  Array4<int const> const& mxlo,
225  Array4<int const> const& mzlo,
226  Array4<Real const> const& bcvalxlo,
227  Array4<Real const> const& bcvalzlo,
228  Array2D<BoundCond,
229  0,2*AMREX_SPACEDIM,
230  0,AMREX_SPACEDIM> const& bct,
231  Array2D<Real,
232  0,2*AMREX_SPACEDIM,
233  0,AMREX_SPACEDIM> const& bcl,
234  int inhomog, int maxorder,
235  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
236  bool xlo_domain, bool zlo_domain) noexcept
237 {
238  if (mxlo(i,j,k) != BndryData::covered && (!xlo_domain || !zlo_domain)) {
239  bool x_interior = mzlo(i+1,j,k ) == BndryData::covered;
240  bool x_exterior = mzlo(i+1,j,k ) == BndryData::not_covered;
241  bool z_interior = mxlo(i ,j,k+1) == BndryData::covered;
242  bool z_exterior = mxlo(i ,j,k+1) == BndryData::not_covered;
243  if ((x_interior && z_interior) || (x_exterior && z_exterior)) {
244  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
245  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
246  bct(Orientation::xlo(), icomp),
247  bcl(Orientation::xlo(), icomp),
248  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
249  Real tmp = vel(i,j,k,icomp);
250  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
251  bct(Orientation::zlo(), icomp),
252  bcl(Orientation::zlo(), icomp),
253  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
254  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
255  }
256  } else if (x_interior || xlo_domain) {
257  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
258  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
259  bct(Orientation::xlo(), icomp),
260  bcl(Orientation::xlo(), icomp),
261  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
262  }
263  } else if (z_interior || zlo_domain) {
264  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
265  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
266  bct(Orientation::zlo(), icomp),
267  bcl(Orientation::zlo(), icomp),
268  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
269  }
270  }
271  }
272 }
273 
274 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
275 void mltensor_fill_edges_xhi_zlo (int const i, int const j, int const k, Dim3 const& blen,
276  Array4<Real> const& vel,
277  Array4<int const> const& mxhi,
278  Array4<int const> const& mzlo,
279  Array4<Real const> const& bcvalxhi,
280  Array4<Real const> const& bcvalzlo,
281  Array2D<BoundCond,
282  0,2*AMREX_SPACEDIM,
283  0,AMREX_SPACEDIM> const& bct,
284  Array2D<Real,
285  0,2*AMREX_SPACEDIM,
286  0,AMREX_SPACEDIM> const& bcl,
287  int inhomog, int maxorder,
288  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
289  bool xhi_domain, bool zlo_domain) noexcept
290 {
291  if (mxhi(i,j,k) != BndryData::covered && (!xhi_domain || !zlo_domain)) {
292  bool x_interior = mzlo(i-1,j,k ) == BndryData::covered;
293  bool x_exterior = mzlo(i-1,j,k ) == BndryData::not_covered;
294  bool z_interior = mxhi(i ,j,k+1) == BndryData::covered;
295  bool z_exterior = mxhi(i ,j,k+1) == BndryData::not_covered;
296  if ((x_interior && z_interior) || (x_exterior && z_exterior)) {
297  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
298  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
299  bct(Orientation::xhi(), icomp),
300  bcl(Orientation::xhi(), icomp),
301  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
302  Real tmp = vel(i,j,k,icomp);
303  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
304  bct(Orientation::zlo(), icomp),
305  bcl(Orientation::zlo(), icomp),
306  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
307  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
308  }
309  } else if (x_interior || xhi_domain) {
310  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
311  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
312  bct(Orientation::xhi(), icomp),
313  bcl(Orientation::xhi(), icomp),
314  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
315  }
316  } else if (z_interior || zlo_domain) {
317  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
318  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
319  bct(Orientation::zlo(), icomp),
320  bcl(Orientation::zlo(), icomp),
321  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
322  }
323  }
324  }
325 }
326 
327 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
328 void mltensor_fill_edges_xlo_zhi (int const i, int const j, int const k, Dim3 const& blen,
329  Array4<Real> const& vel,
330  Array4<int const> const& mxlo,
331  Array4<int const> const& mzhi,
332  Array4<Real const> const& bcvalxlo,
333  Array4<Real const> const& bcvalzhi,
334  Array2D<BoundCond,
335  0,2*AMREX_SPACEDIM,
336  0,AMREX_SPACEDIM> const& bct,
337  Array2D<Real,
338  0,2*AMREX_SPACEDIM,
339  0,AMREX_SPACEDIM> const& bcl,
340  int inhomog, int maxorder,
341  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
342  bool xlo_domain, bool zhi_domain) noexcept
343 {
344  if (mxlo(i,j,k) != BndryData::covered && (!xlo_domain || !zhi_domain)) {
345  bool x_interior = mzhi(i+1,j,k ) == BndryData::covered;
346  bool x_exterior = mzhi(i+1,j,k ) == BndryData::not_covered;
347  bool z_interior = mxlo(i ,j,k-1) == BndryData::covered;
348  bool z_exterior = mxlo(i ,j,k-1) == BndryData::not_covered;
349  if ((x_interior && z_interior) || (x_exterior && z_exterior)) {
350  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
351  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
352  bct(Orientation::xlo(), icomp),
353  bcl(Orientation::xlo(), icomp),
354  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
355  Real tmp = vel(i,j,k,icomp);
356  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
357  bct(Orientation::zhi(), icomp),
358  bcl(Orientation::zhi(), icomp),
359  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
360  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
361  }
362  } else if (x_interior || xlo_domain) {
363  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
364  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
365  bct(Orientation::xlo(), icomp),
366  bcl(Orientation::xlo(), icomp),
367  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
368  }
369  } else if (z_interior || zhi_domain) {
370  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
371  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
372  bct(Orientation::zhi(), icomp),
373  bcl(Orientation::zhi(), icomp),
374  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
375  }
376  }
377  }
378 }
379 
380 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
381 void mltensor_fill_edges_xhi_zhi (int const i, int const j, int const k, Dim3 const& blen,
382  Array4<Real> const& vel,
383  Array4<int const> const& mxhi,
384  Array4<int const> const& mzhi,
385  Array4<Real const> const& bcvalxhi,
386  Array4<Real const> const& bcvalzhi,
387  Array2D<BoundCond,
388  0,2*AMREX_SPACEDIM,
389  0,AMREX_SPACEDIM> const& bct,
390  Array2D<Real,
391  0,2*AMREX_SPACEDIM,
392  0,AMREX_SPACEDIM> const& bcl,
393  int inhomog, int maxorder,
394  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
395  bool xhi_domain, bool zhi_domain) noexcept
396 {
397  if (mxhi(i,j,k) != BndryData::covered && (!xhi_domain || !zhi_domain)) {
398  bool x_interior = mzhi(i-1,j,k ) == BndryData::covered;
399  bool x_exterior = mzhi(i-1,j,k ) == BndryData::not_covered;
400  bool z_interior = mxhi(i ,j,k-1) == BndryData::covered;
401  bool z_exterior = mxhi(i ,j,k-1) == BndryData::not_covered;
402  if ((x_interior && z_interior) || (x_exterior && z_exterior)) {
403  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
404  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
405  bct(Orientation::xhi(), icomp),
406  bcl(Orientation::xhi(), icomp),
407  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
408  Real tmp = vel(i,j,k,icomp);
409  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
410  bct(Orientation::zhi(), icomp),
411  bcl(Orientation::zhi(), icomp),
412  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
413  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
414  }
415  } else if (x_interior || xhi_domain) {
416  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
417  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
418  bct(Orientation::xhi(), icomp),
419  bcl(Orientation::xhi(), icomp),
420  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
421  }
422  } else if (z_interior || zhi_domain) {
423  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
424  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
425  bct(Orientation::zhi(), icomp),
426  bcl(Orientation::zhi(), icomp),
427  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
428  }
429  }
430  }
431 }
432 
433 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
434 void mltensor_fill_edges_ylo_zlo (int const i, int const j, int const k, Dim3 const& blen,
435  Array4<Real> const& vel,
436  Array4<int const> const& mylo,
437  Array4<int const> const& mzlo,
438  Array4<Real const> const& bcvalylo,
439  Array4<Real const> const& bcvalzlo,
440  Array2D<BoundCond,
441  0,2*AMREX_SPACEDIM,
442  0,AMREX_SPACEDIM> const& bct,
443  Array2D<Real,
444  0,2*AMREX_SPACEDIM,
445  0,AMREX_SPACEDIM> const& bcl,
446  int inhomog, int maxorder,
447  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
448  bool ylo_domain, bool zlo_domain) noexcept
449 {
450  if (mylo(i,j,k) != BndryData::covered && (!ylo_domain || !zlo_domain)) {
451  bool y_interior = mzlo(i,j+1,k ) == BndryData::covered;
452  bool y_exterior = mzlo(i,j+1,k ) == BndryData::not_covered;
453  bool z_interior = mylo(i,j ,k+1) == BndryData::covered;
454  bool z_exterior = mylo(i,j ,k+1) == BndryData::not_covered;
455  if ((y_interior && z_interior) || (y_exterior && z_exterior)) {
456  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
457  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
458  bct(Orientation::ylo(), icomp),
459  bcl(Orientation::ylo(), icomp),
460  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
461  Real tmp = vel(i,j,k,icomp);
462  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
463  bct(Orientation::zlo(), icomp),
464  bcl(Orientation::zlo(), icomp),
465  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
466  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
467  }
468  } else if (y_interior || ylo_domain) {
469  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
470  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
471  bct(Orientation::ylo(), icomp),
472  bcl(Orientation::ylo(), icomp),
473  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
474  }
475  } else if (z_interior || zlo_domain) {
476  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
477  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
478  bct(Orientation::zlo(), icomp),
479  bcl(Orientation::zlo(), icomp),
480  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
481  }
482  }
483  }
484 }
485 
486 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
487 void mltensor_fill_edges_yhi_zlo (int const i, int const j, int const k, Dim3 const& blen,
488  Array4<Real> const& vel,
489  Array4<int const> const& myhi,
490  Array4<int const> const& mzlo,
491  Array4<Real const> const& bcvalyhi,
492  Array4<Real const> const& bcvalzlo,
493  Array2D<BoundCond,
494  0,2*AMREX_SPACEDIM,
495  0,AMREX_SPACEDIM> const& bct,
496  Array2D<Real,
497  0,2*AMREX_SPACEDIM,
498  0,AMREX_SPACEDIM> const& bcl,
499  int inhomog, int maxorder,
500  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
501  bool yhi_domain, bool zlo_domain) noexcept
502 {
503  if (myhi(i,j,k) != BndryData::covered && (!yhi_domain || !zlo_domain)) {
504  bool y_interior = mzlo(i,j-1,k ) == BndryData::covered;
505  bool y_exterior = mzlo(i,j-1,k ) == BndryData::not_covered;
506  bool z_interior = myhi(i,j ,k+1) == BndryData::covered;
507  bool z_exterior = myhi(i,j ,k+1) == BndryData::not_covered;
508  if ((y_interior && z_interior) || (y_exterior && z_exterior)) {
509  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
510  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
511  bct(Orientation::yhi(), icomp),
512  bcl(Orientation::yhi(), icomp),
513  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
514  Real tmp = vel(i,j,k,icomp);
515  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
516  bct(Orientation::zlo(), icomp),
517  bcl(Orientation::zlo(), icomp),
518  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
519  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
520  }
521  } else if (y_interior || yhi_domain) {
522  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
523  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
524  bct(Orientation::yhi(), icomp),
525  bcl(Orientation::yhi(), icomp),
526  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
527  }
528  } else if (z_interior || zlo_domain) {
529  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
530  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
531  bct(Orientation::zlo(), icomp),
532  bcl(Orientation::zlo(), icomp),
533  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
534  }
535  }
536  }
537 }
538 
539 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
540 void mltensor_fill_edges_ylo_zhi (int const i, int const j, int const k, Dim3 const& blen,
541  Array4<Real> const& vel,
542  Array4<int const> const& mylo,
543  Array4<int const> const& mzhi,
544  Array4<Real const> const& bcvalylo,
545  Array4<Real const> const& bcvalzhi,
546  Array2D<BoundCond,
547  0,2*AMREX_SPACEDIM,
548  0,AMREX_SPACEDIM> const& bct,
549  Array2D<Real,
550  0,2*AMREX_SPACEDIM,
551  0,AMREX_SPACEDIM> const& bcl,
552  int inhomog, int maxorder,
553  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
554  bool ylo_domain, bool zhi_domain) noexcept
555 {
556  if (mylo(i,j,k) != BndryData::covered && (!ylo_domain || !zhi_domain)) {
557  bool y_interior = mzhi(i,j+1,k ) == BndryData::covered;
558  bool y_exterior = mzhi(i,j+1,k ) == BndryData::not_covered;
559  bool z_interior = mylo(i,j ,k-1) == BndryData::covered;
560  bool z_exterior = mylo(i,j ,k-1) == BndryData::not_covered;
561  if ((y_interior && z_interior) || (y_exterior && z_exterior)) {
562  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
563  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
564  bct(Orientation::ylo(), icomp),
565  bcl(Orientation::ylo(), icomp),
566  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
567  Real tmp = vel(i,j,k,icomp);
568  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
569  bct(Orientation::zhi(), icomp),
570  bcl(Orientation::zhi(), icomp),
571  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
572  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
573  }
574  } else if (y_interior || ylo_domain) {
575  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
576  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
577  bct(Orientation::ylo(), icomp),
578  bcl(Orientation::ylo(), icomp),
579  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
580  }
581  } else if (z_interior || zhi_domain) {
582  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
583  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
584  bct(Orientation::zhi(), icomp),
585  bcl(Orientation::zhi(), icomp),
586  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
587  }
588  }
589  }
590 }
591 
592 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
593 void mltensor_fill_edges_yhi_zhi (int const i, int const j, int const k, Dim3 const& blen,
594  Array4<Real> const& vel,
595  Array4<int const> const& myhi,
596  Array4<int const> const& mzhi,
597  Array4<Real const> const& bcvalyhi,
598  Array4<Real const> const& bcvalzhi,
599  Array2D<BoundCond,
600  0,2*AMREX_SPACEDIM,
601  0,AMREX_SPACEDIM> const& bct,
602  Array2D<Real,
603  0,2*AMREX_SPACEDIM,
604  0,AMREX_SPACEDIM> const& bcl,
605  int inhomog, int maxorder,
606  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
607  bool yhi_domain, bool zhi_domain) noexcept
608 {
609  if (myhi(i,j,k) != BndryData::covered && (!yhi_domain || !zhi_domain)) {
610  bool y_interior = mzhi(i,j-1,k ) == BndryData::covered;
611  bool y_exterior = mzhi(i,j-1,k ) == BndryData::not_covered;
612  bool z_interior = myhi(i,j ,k-1) == BndryData::covered;
613  bool z_exterior = myhi(i,j ,k-1) == BndryData::not_covered;
614  if ((y_interior && z_interior) || (y_exterior && z_exterior)) {
615  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
616  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
617  bct(Orientation::yhi(), icomp),
618  bcl(Orientation::yhi(), icomp),
619  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
620  Real tmp = vel(i,j,k,icomp);
621  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
622  bct(Orientation::zhi(), icomp),
623  bcl(Orientation::zhi(), icomp),
624  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
625  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
626  }
627  } else if (y_interior || yhi_domain) {
628  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
629  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
630  bct(Orientation::yhi(), icomp),
631  bcl(Orientation::yhi(), icomp),
632  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
633  }
634  } else if (z_interior || zhi_domain) {
635  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
636  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
637  bct(Orientation::zhi(), icomp),
638  bcl(Orientation::zhi(), icomp),
639  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
640  }
641  }
642  }
643 }
644 
645 #ifdef AMREX_USE_EB
646 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
647 void mltensor_fill_corners (int icorner, Box const& vbox, // vbox: the valid box
648  Array4<Real> const& vel,
649  Array4<int const> const& mxlo,
650  Array4<int const> const& mylo,
651  Array4<int const> const& mzlo,
652  Array4<int const> const& mxhi,
653  Array4<int const> const& myhi,
654  Array4<int const> const& mzhi,
655  Array4<Real const> const& bcvalxlo,
656  Array4<Real const> const& bcvalylo,
657  Array4<Real const> const& bcvalzlo,
658  Array4<Real const> const& bcvalxhi,
659  Array4<Real const> const& bcvalyhi,
660  Array4<Real const> const& bcvalzhi,
661  Array2D<BoundCond,
662  0,2*AMREX_SPACEDIM,
663  0,AMREX_SPACEDIM> const& bct,
664  Array2D<Real,
665  0,2*AMREX_SPACEDIM,
666  0,AMREX_SPACEDIM> const& bcl,
667  int inhomog, int maxorder,
668  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
669  Dim3 const& dlo, Dim3 const& dhi) noexcept
670 {
671  const auto blen = amrex::length(vbox);
672  const auto vlo = amrex::lbound(vbox);
673  const auto vhi = amrex::ubound(vbox);
674  bool xlo_domain = (vlo.x == dlo.x);
675  bool ylo_domain = (vlo.y == dlo.y);
676  bool zlo_domain = (vlo.z == dlo.z);
677  bool xhi_domain = (vhi.x == dhi.x);
678  bool yhi_domain = (vhi.y == dhi.y);
679  bool zhi_domain = (vhi.z == dhi.z);
680 
681  for (int icomp = 0; icomp < AMREX_SPACEDIM; ++icomp) {
682  switch (icorner) {
683  case 0: {
684  // xlo & ylo & zlo
685  int i = vlo.x-1;
686  int j = vlo.y-1;
687  int k = vlo.z-1;
688  if (mxlo(i,j,k) != BndryData::covered &&
689  (!xlo_domain || !ylo_domain || !zlo_domain)) {
690  bool x_interior = mylo(i+1,j ,k ) == BndryData::covered;
691  bool x_exterior = mylo(i+1,j ,k ) == BndryData::not_covered;
692  bool y_interior = mxlo(i ,j+1,k ) == BndryData::covered;
693  bool y_exterior = mxlo(i ,j+1,k ) == BndryData::not_covered;
694  bool z_interior = mxlo(i ,j ,k+1) == BndryData::covered;
695  bool z_exterior = mxlo(i ,j ,k+1) == BndryData::not_covered;
696  if ((x_interior && y_interior && z_interior) ||
697  (x_exterior && y_exterior && z_exterior)) {
698  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
699  bct(Orientation::xlo(), icomp),
700  bcl(Orientation::xlo(), icomp),
701  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
702  Real tmp = vel(i,j,k,icomp);
703  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
704  bct(Orientation::ylo(), icomp),
705  bcl(Orientation::ylo(), icomp),
706  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
707  tmp += vel(i,j,k,icomp);
708  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
709  bct(Orientation::zlo(), icomp),
710  bcl(Orientation::zlo(), icomp),
711  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
712  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
713  } else if (x_interior && y_interior) {
714  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
715  bct(Orientation::xlo(), icomp),
716  bcl(Orientation::xlo(), icomp),
717  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
718  Real tmp = vel(i,j,k,icomp);
719  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
720  bct(Orientation::ylo(), icomp),
721  bcl(Orientation::ylo(), icomp),
722  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
723  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
724  } else if (x_interior && z_interior) {
725  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
726  bct(Orientation::xlo(), icomp),
727  bcl(Orientation::xlo(), icomp),
728  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
729  Real tmp = vel(i,j,k,icomp);
730  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
731  bct(Orientation::zlo(), icomp),
732  bcl(Orientation::zlo(), icomp),
733  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
734  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
735  } else if (y_interior && z_interior) {
736  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
737  bct(Orientation::ylo(), icomp),
738  bcl(Orientation::ylo(), icomp),
739  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
740  Real tmp = vel(i,j,k,icomp);
741  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
742  bct(Orientation::zlo(), icomp),
743  bcl(Orientation::zlo(), icomp),
744  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
745  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
746  } else if (x_interior) {
747  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
748  bct(Orientation::xlo(), icomp),
749  bcl(Orientation::xlo(), icomp),
750  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
751  } else if (y_interior) {
752  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
753  bct(Orientation::ylo(), icomp),
754  bcl(Orientation::ylo(), icomp),
755  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
756  } else if (z_interior) {
757  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
758  bct(Orientation::zlo(), icomp),
759  bcl(Orientation::zlo(), icomp),
760  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
761  }
762  }
763  break;
764  }
765  case 1: {
766  // xhi & ylo & zlo
767  int i = vhi.x+1;
768  int j = vlo.y-1;
769  int k = vlo.z-1;
770  bool x_interior = mylo(i-1,j ,k ) == BndryData::covered;
771  bool x_exterior = mylo(i-1,j ,k ) == BndryData::not_covered;
772  bool y_interior = mxhi(i ,j+1,k ) == BndryData::covered;
773  bool y_exterior = mxhi(i ,j+1,k ) == BndryData::not_covered;
774  bool z_interior = mxhi(i ,j ,k+1) == BndryData::covered;
775  bool z_exterior = mxhi(i ,j ,k+1) == BndryData::not_covered;
776  if (mxhi(i,j,k) != BndryData::covered &&
777  (!xhi_domain || !ylo_domain || !zlo_domain)) {
778  if ((x_interior && y_interior && z_interior) ||
779  (x_exterior && y_exterior && z_exterior)) {
780  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
781  bct(Orientation::xhi(), icomp),
782  bcl(Orientation::xhi(), icomp),
783  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
784  Real tmp = vel(i,j,k,icomp);
785  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
786  bct(Orientation::ylo(), icomp),
787  bcl(Orientation::ylo(), icomp),
788  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
789  tmp += vel(i,j,k,icomp);
790  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
791  bct(Orientation::zlo(), icomp),
792  bcl(Orientation::zlo(), icomp),
793  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
794  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
795  } else if (x_interior && y_interior) {
796  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
797  bct(Orientation::xhi(), icomp),
798  bcl(Orientation::xhi(), icomp),
799  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
800  Real tmp = vel(i,j,k,icomp);
801  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
802  bct(Orientation::ylo(), icomp),
803  bcl(Orientation::ylo(), icomp),
804  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
805  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
806  } else if (x_interior && z_interior) {
807  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
808  bct(Orientation::xhi(), icomp),
809  bcl(Orientation::xhi(), icomp),
810  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
811  Real tmp = vel(i,j,k,icomp);
812  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
813  bct(Orientation::zlo(), icomp),
814  bcl(Orientation::zlo(), icomp),
815  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
816  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
817  } else if (y_interior && z_interior) {
818  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
819  bct(Orientation::ylo(), icomp),
820  bcl(Orientation::ylo(), icomp),
821  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
822  Real tmp = vel(i,j,k,icomp);
823  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
824  bct(Orientation::zlo(), icomp),
825  bcl(Orientation::zlo(), icomp),
826  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
827  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
828  } else if (x_interior) {
829  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
830  bct(Orientation::xhi(), icomp),
831  bcl(Orientation::xhi(), icomp),
832  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
833  } else if (y_interior) {
834  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
835  bct(Orientation::ylo(), icomp),
836  bcl(Orientation::ylo(), icomp),
837  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
838  } else if (z_interior) {
839  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
840  bct(Orientation::zlo(), icomp),
841  bcl(Orientation::zlo(), icomp),
842  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
843  }
844  }
845  break;
846  }
847  case 2: {
848  // xlo & yhi & zlo
849  int i = vlo.x-1;
850  int j = vhi.y+1;
851  int k = vlo.z-1;
852  bool x_interior = myhi(i+1,j ,k ) == BndryData::covered;
853  bool x_exterior = myhi(i+1,j ,k ) == BndryData::not_covered;
854  bool y_interior = mxlo(i ,j-1,k ) == BndryData::covered;
855  bool y_exterior = mxlo(i ,j-1,k ) == BndryData::not_covered;
856  bool z_interior = mxlo(i ,j ,k+1) == BndryData::covered;
857  bool z_exterior = mxlo(i ,j ,k+1) == BndryData::not_covered;
858  if (mxlo(i,j,k) != BndryData::covered &&
859  (!xlo_domain || !yhi_domain || !zlo_domain)) {
860  if ((x_interior && y_interior && z_interior) ||
861  (x_exterior && y_exterior && z_exterior)) {
862  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
863  bct(Orientation::xlo(), icomp),
864  bcl(Orientation::xlo(), icomp),
865  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
866  Real tmp = vel(i,j,k,icomp);
867  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
868  bct(Orientation::yhi(), icomp),
869  bcl(Orientation::yhi(), icomp),
870  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
871  tmp += vel(i,j,k,icomp);
872  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
873  bct(Orientation::zlo(), icomp),
874  bcl(Orientation::zlo(), icomp),
875  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
876  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
877  } else if (x_interior && y_interior) {
878  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
879  bct(Orientation::xlo(), icomp),
880  bcl(Orientation::xlo(), icomp),
881  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
882  Real tmp = vel(i,j,k,icomp);
883  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
884  bct(Orientation::yhi(), icomp),
885  bcl(Orientation::yhi(), icomp),
886  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
887  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
888  } else if (x_interior && z_interior) {
889  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
890  bct(Orientation::xlo(), icomp),
891  bcl(Orientation::xlo(), icomp),
892  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
893  Real tmp = vel(i,j,k,icomp);
894  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
895  bct(Orientation::zlo(), icomp),
896  bcl(Orientation::zlo(), icomp),
897  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
898  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
899  } else if (y_interior && z_interior) {
900  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
901  bct(Orientation::yhi(), icomp),
902  bcl(Orientation::yhi(), icomp),
903  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
904  Real tmp = vel(i,j,k,icomp);
905  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
906  bct(Orientation::zlo(), icomp),
907  bcl(Orientation::zlo(), icomp),
908  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
909  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
910  } else if (x_interior) {
911  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
912  bct(Orientation::xlo(), icomp),
913  bcl(Orientation::xlo(), icomp),
914  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
915  } else if (y_interior) {
916  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
917  bct(Orientation::yhi(), icomp),
918  bcl(Orientation::yhi(), icomp),
919  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
920  } else if (z_interior) {
921  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
922  bct(Orientation::zlo(), icomp),
923  bcl(Orientation::zlo(), icomp),
924  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
925  }
926  }
927  break;
928  }
929  case 3: {
930  // xhi & yhi & zlo
931  int i = vhi.x+1;
932  int j = vhi.y+1;
933  int k = vlo.z-1;
934  bool x_interior = myhi(i-1,j ,k ) == BndryData::covered;
935  bool x_exterior = myhi(i-1,j ,k ) == BndryData::not_covered;
936  bool y_interior = mxhi(i ,j-1,k ) == BndryData::covered;
937  bool y_exterior = mxhi(i ,j-1,k ) == BndryData::not_covered;
938  bool z_interior = mxhi(i ,j ,k+1) == BndryData::covered;
939  bool z_exterior = mxhi(i ,j ,k+1) == BndryData::not_covered;
940  if (mxhi(i,j,k) != BndryData::covered &&
941  (!xhi_domain || !yhi_domain || !zlo_domain)) {
942  if ((x_interior && y_interior && z_interior) ||
943  (x_exterior && y_exterior && z_exterior)) {
944  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
945  bct(Orientation::xhi(), icomp),
946  bcl(Orientation::xhi(), icomp),
947  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
948  Real tmp = vel(i,j,k,icomp);
949  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
950  bct(Orientation::yhi(), icomp),
951  bcl(Orientation::yhi(), icomp),
952  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
953  tmp += vel(i,j,k,icomp);
954  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
955  bct(Orientation::zlo(), icomp),
956  bcl(Orientation::zlo(), icomp),
957  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
958  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
959  } else if (x_interior && y_interior) {
960  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
961  bct(Orientation::xhi(), icomp),
962  bcl(Orientation::xhi(), icomp),
963  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
964  Real tmp = vel(i,j,k,icomp);
965  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
966  bct(Orientation::yhi(), icomp),
967  bcl(Orientation::yhi(), icomp),
968  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
969  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
970  } else if (x_interior && z_interior) {
971  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
972  bct(Orientation::xhi(), icomp),
973  bcl(Orientation::xhi(), icomp),
974  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
975  Real tmp = vel(i,j,k,icomp);
976  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
977  bct(Orientation::zlo(), icomp),
978  bcl(Orientation::zlo(), icomp),
979  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
980  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
981  } else if (y_interior && z_interior) {
982  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
983  bct(Orientation::yhi(), icomp),
984  bcl(Orientation::yhi(), icomp),
985  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
986  Real tmp = vel(i,j,k,icomp);
987  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
988  bct(Orientation::zlo(), icomp),
989  bcl(Orientation::zlo(), icomp),
990  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
991  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
992  } else if (x_interior) {
993  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
994  bct(Orientation::xhi(), icomp),
995  bcl(Orientation::xhi(), icomp),
996  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
997  } else if (y_interior) {
998  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
999  bct(Orientation::yhi(), icomp),
1000  bcl(Orientation::yhi(), icomp),
1001  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1002  } else if (z_interior) {
1003  mllinop_apply_bc_z(Orientation::low, i,j,k, blen.z, vel, mzlo,
1004  bct(Orientation::zlo(), icomp),
1005  bcl(Orientation::zlo(), icomp),
1006  bcvalzlo, maxorder, dxinv[2], inhomog, icomp);
1007  }
1008  }
1009  break;
1010  }
1011  case 4: {
1012  // xlo & ylo & zhi
1013  int i = vlo.x-1;
1014  int j = vlo.y-1;
1015  int k = vhi.z+1;
1016  bool x_interior = mylo(i+1,j ,k ) == BndryData::covered;
1017  bool x_exterior = mylo(i+1,j ,k ) == BndryData::not_covered;
1018  bool y_interior = mxlo(i ,j+1,k ) == BndryData::covered;
1019  bool y_exterior = mxlo(i ,j+1,k ) == BndryData::not_covered;
1020  bool z_interior = mxlo(i ,j ,k-1) == BndryData::covered;
1021  bool z_exterior = mxlo(i ,j ,k-1) == BndryData::not_covered;
1022  if (mxlo(i,j,k) != BndryData::covered &&
1023  (!xlo_domain || !ylo_domain || !zhi_domain)) {
1024  if ((x_interior && y_interior && z_interior) ||
1025  (x_exterior && y_exterior && z_exterior)) {
1026  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1027  bct(Orientation::xlo(), icomp),
1028  bcl(Orientation::xlo(), icomp),
1029  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1030  Real tmp = vel(i,j,k,icomp);
1031  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1032  bct(Orientation::ylo(), icomp),
1033  bcl(Orientation::ylo(), icomp),
1034  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1035  tmp += vel(i,j,k,icomp);
1036  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1037  bct(Orientation::zhi(), icomp),
1038  bcl(Orientation::zhi(), icomp),
1039  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1040  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
1041  } else if (x_interior && y_interior) {
1042  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1043  bct(Orientation::xlo(), icomp),
1044  bcl(Orientation::xlo(), icomp),
1045  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1046  Real tmp = vel(i,j,k,icomp);
1047  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1048  bct(Orientation::ylo(), icomp),
1049  bcl(Orientation::ylo(), icomp),
1050  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1051  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1052  } else if (x_interior && z_interior) {
1053  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1054  bct(Orientation::xlo(), icomp),
1055  bcl(Orientation::xlo(), icomp),
1056  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1057  Real tmp = vel(i,j,k,icomp);
1058  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1059  bct(Orientation::zhi(), icomp),
1060  bcl(Orientation::zhi(), icomp),
1061  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1062  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1063  } else if (y_interior && z_interior) {
1064  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1065  bct(Orientation::ylo(), icomp),
1066  bcl(Orientation::ylo(), icomp),
1067  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1068  Real tmp = vel(i,j,k,icomp);
1069  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1070  bct(Orientation::zhi(), icomp),
1071  bcl(Orientation::zhi(), icomp),
1072  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1073  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1074  } else if (x_interior) {
1075  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1076  bct(Orientation::xlo(), icomp),
1077  bcl(Orientation::xlo(), icomp),
1078  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1079  } else if (y_interior) {
1080  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1081  bct(Orientation::ylo(), icomp),
1082  bcl(Orientation::ylo(), icomp),
1083  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1084  } else if (z_interior) {
1085  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1086  bct(Orientation::zhi(), icomp),
1087  bcl(Orientation::zhi(), icomp),
1088  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1089  }
1090  }
1091  break;
1092  }
1093  case 5: {
1094  // xhi & ylo & zhi
1095  int i = vhi.x+1;
1096  int j = vlo.y-1;
1097  int k = vhi.z+1;
1098  bool x_interior = mylo(i-1,j ,k ) == BndryData::covered;
1099  bool x_exterior = mylo(i-1,j ,k ) == BndryData::not_covered;
1100  bool y_interior = mxhi(i ,j+1,k ) == BndryData::covered;
1101  bool y_exterior = mxhi(i ,j+1,k ) == BndryData::not_covered;
1102  bool z_interior = mxhi(i ,j ,k-1) == BndryData::covered;
1103  bool z_exterior = mxhi(i ,j ,k-1) == BndryData::not_covered;
1104  if (mxhi(i,j,k) != BndryData::covered &&
1105  (!xhi_domain || !ylo_domain || !zhi_domain)) {
1106  if ((x_interior && y_interior && z_interior) ||
1107  (x_exterior && y_exterior && z_exterior)) {
1108  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1109  bct(Orientation::xhi(), icomp),
1110  bcl(Orientation::xhi(), icomp),
1111  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1112  Real tmp = vel(i,j,k,icomp);
1113  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1114  bct(Orientation::ylo(), icomp),
1115  bcl(Orientation::ylo(), icomp),
1116  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1117  tmp += vel(i,j,k,icomp);
1118  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1119  bct(Orientation::zhi(), icomp),
1120  bcl(Orientation::zhi(), icomp),
1121  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1122  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
1123  } else if (x_interior && y_interior) {
1124  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1125  bct(Orientation::xhi(), icomp),
1126  bcl(Orientation::xhi(), icomp),
1127  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1128  Real tmp = vel(i,j,k,icomp);
1129  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1130  bct(Orientation::ylo(), icomp),
1131  bcl(Orientation::ylo(), icomp),
1132  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1133  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1134  } else if (x_interior && z_interior) {
1135  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1136  bct(Orientation::xhi(), icomp),
1137  bcl(Orientation::xhi(), icomp),
1138  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1139  Real tmp = vel(i,j,k,icomp);
1140  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1141  bct(Orientation::zhi(), icomp),
1142  bcl(Orientation::zhi(), icomp),
1143  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1144  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1145  } else if (y_interior && z_interior) {
1146  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1147  bct(Orientation::ylo(), icomp),
1148  bcl(Orientation::ylo(), icomp),
1149  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1150  Real tmp = vel(i,j,k,icomp);
1151  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1152  bct(Orientation::zhi(), icomp),
1153  bcl(Orientation::zhi(), icomp),
1154  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1155  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1156  } else if (x_interior) {
1157  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1158  bct(Orientation::xhi(), icomp),
1159  bcl(Orientation::xhi(), icomp),
1160  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1161  } else if (y_interior) {
1162  mllinop_apply_bc_y(Orientation::low, i,j,k, blen.y, vel, mylo,
1163  bct(Orientation::ylo(), icomp),
1164  bcl(Orientation::ylo(), icomp),
1165  bcvalylo, maxorder, dxinv[1], inhomog, icomp);
1166  } else if (z_interior) {
1167  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1168  bct(Orientation::zhi(), icomp),
1169  bcl(Orientation::zhi(), icomp),
1170  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1171  }
1172  }
1173  break;
1174  }
1175  case 6: {
1176  // xlo & yhi & zhi
1177  int i = vlo.x-1;
1178  int j = vhi.y+1;
1179  int k = vhi.z+1;
1180  bool x_interior = myhi(i+1,j ,k ) == BndryData::covered;
1181  bool x_exterior = myhi(i+1,j ,k ) == BndryData::not_covered;
1182  bool y_interior = mxlo(i ,j-1,k ) == BndryData::covered;
1183  bool y_exterior = mxlo(i ,j-1,k ) == BndryData::not_covered;
1184  bool z_interior = mxlo(i ,j ,k-1) == BndryData::covered;
1185  bool z_exterior = mxlo(i ,j ,k-1) == BndryData::not_covered;
1186  if (mxlo(i,j,k) != BndryData::covered &&
1187  (!xlo_domain || !yhi_domain || !zhi_domain)) {
1188  if ((x_interior && y_interior && z_interior) ||
1189  (x_exterior && y_exterior && z_exterior)) {
1190  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1191  bct(Orientation::xlo(), icomp),
1192  bcl(Orientation::xlo(), icomp),
1193  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1194  Real tmp = vel(i,j,k,icomp);
1195  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1196  bct(Orientation::yhi(), icomp),
1197  bcl(Orientation::yhi(), icomp),
1198  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1199  tmp += vel(i,j,k,icomp);
1200  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1201  bct(Orientation::zhi(), icomp),
1202  bcl(Orientation::zhi(), icomp),
1203  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1204  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
1205  } else if (x_interior && y_interior) {
1206  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1207  bct(Orientation::xlo(), icomp),
1208  bcl(Orientation::xlo(), icomp),
1209  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1210  Real tmp = vel(i,j,k,icomp);
1211  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1212  bct(Orientation::yhi(), icomp),
1213  bcl(Orientation::yhi(), icomp),
1214  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1215  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1216  } else if (x_interior && z_interior) {
1217  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1218  bct(Orientation::xlo(), icomp),
1219  bcl(Orientation::xlo(), icomp),
1220  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1221  Real tmp = vel(i,j,k,icomp);
1222  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1223  bct(Orientation::zhi(), icomp),
1224  bcl(Orientation::zhi(), icomp),
1225  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1226  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1227  } else if (y_interior && z_interior) {
1228  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1229  bct(Orientation::yhi(), icomp),
1230  bcl(Orientation::yhi(), icomp),
1231  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1232  Real tmp = vel(i,j,k,icomp);
1233  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1234  bct(Orientation::zhi(), icomp),
1235  bcl(Orientation::zhi(), icomp),
1236  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1237  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1238  } else if (x_interior) {
1239  mllinop_apply_bc_x(Orientation::low, i,j,k, blen.x, vel, mxlo,
1240  bct(Orientation::xlo(), icomp),
1241  bcl(Orientation::xlo(), icomp),
1242  bcvalxlo, maxorder, dxinv[0], inhomog, icomp);
1243  } else if (y_interior) {
1244  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1245  bct(Orientation::yhi(), icomp),
1246  bcl(Orientation::yhi(), icomp),
1247  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1248  } else if (z_interior) {
1249  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1250  bct(Orientation::zhi(), icomp),
1251  bcl(Orientation::zhi(), icomp),
1252  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1253  }
1254  }
1255  break;
1256  }
1257  case 7: {
1258  // xhi & yhi & zhi
1259  int i = vhi.x+1;
1260  int j = vhi.y+1;
1261  int k = vhi.z+1;
1262  bool x_interior = myhi(i-1,j ,k ) == BndryData::covered;
1263  bool x_exterior = myhi(i-1,j ,k ) == BndryData::not_covered;
1264  bool y_interior = mxhi(i ,j-1,k ) == BndryData::covered;
1265  bool y_exterior = mxhi(i ,j-1,k ) == BndryData::not_covered;
1266  bool z_interior = mxhi(i ,j ,k-1) == BndryData::covered;
1267  bool z_exterior = mxhi(i ,j ,k-1) == BndryData::not_covered;
1268  if (mxhi(i,j,k) != BndryData::covered &&
1269  (!xhi_domain || !yhi_domain || !zhi_domain)) {
1270  if ((x_interior && y_interior && z_interior) ||
1271  (x_exterior && y_exterior && z_exterior)) {
1272  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1273  bct(Orientation::xhi(), icomp),
1274  bcl(Orientation::xhi(), icomp),
1275  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1276  Real tmp = vel(i,j,k,icomp);
1277  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1278  bct(Orientation::yhi(), icomp),
1279  bcl(Orientation::yhi(), icomp),
1280  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1281  tmp += vel(i,j,k,icomp);
1282  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1283  bct(Orientation::zhi(), icomp),
1284  bcl(Orientation::zhi(), icomp),
1285  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1286  vel(i,j,k,icomp) = Real(1./3.)*(tmp+vel(i,j,k,icomp));
1287  } else if (x_interior && y_interior) {
1288  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1289  bct(Orientation::xhi(), icomp),
1290  bcl(Orientation::xhi(), icomp),
1291  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1292  Real tmp = vel(i,j,k,icomp);
1293  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1294  bct(Orientation::yhi(), icomp),
1295  bcl(Orientation::yhi(), icomp),
1296  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1297  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1298  } else if (x_interior && z_interior) {
1299  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1300  bct(Orientation::xhi(), icomp),
1301  bcl(Orientation::xhi(), icomp),
1302  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1303  Real tmp = vel(i,j,k,icomp);
1304  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1305  bct(Orientation::zhi(), icomp),
1306  bcl(Orientation::zhi(), icomp),
1307  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1308  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1309  } else if (y_interior && z_interior) {
1310  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1311  bct(Orientation::yhi(), icomp),
1312  bcl(Orientation::yhi(), icomp),
1313  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1314  Real tmp = vel(i,j,k,icomp);
1315  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1316  bct(Orientation::zhi(), icomp),
1317  bcl(Orientation::zhi(), icomp),
1318  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1319  vel(i,j,k,icomp) = 0.5_rt*(tmp+vel(i,j,k,icomp));
1320  } else if (x_interior) {
1321  mllinop_apply_bc_x(Orientation::high, i,j,k, blen.x, vel, mxhi,
1322  bct(Orientation::xhi(), icomp),
1323  bcl(Orientation::xhi(), icomp),
1324  bcvalxhi, maxorder, dxinv[0], inhomog, icomp);
1325  } else if (y_interior) {
1326  mllinop_apply_bc_y(Orientation::high, i,j,k, blen.y, vel, myhi,
1327  bct(Orientation::yhi(), icomp),
1328  bcl(Orientation::yhi(), icomp),
1329  bcvalyhi, maxorder, dxinv[1], inhomog, icomp);
1330  } else if (z_interior) {
1331  mllinop_apply_bc_z(Orientation::high, i,j,k, blen.z, vel, mzhi,
1332  bct(Orientation::zhi(), icomp),
1333  bcl(Orientation::zhi(), icomp),
1334  bcvalzhi, maxorder, dxinv[2], inhomog, icomp);
1335  }
1336  }
1337  break;
1338  }
1339  default: {}
1340  }
1341  }
1342 }
1343 #endif
1344 
1345 inline
1346 void mltensor_fill_edges (Box const& vbox, // vbox: the valid box
1347  Array4<Real> const& vel,
1348  Array4<int const> const& mxlo,
1349  Array4<int const> const& mylo,
1350  Array4<int const> const& mzlo,
1351  Array4<int const> const& mxhi,
1352  Array4<int const> const& myhi,
1353  Array4<int const> const& mzhi,
1354  Array4<Real const> const& bcvalxlo,
1355  Array4<Real const> const& bcvalylo,
1356  Array4<Real const> const& bcvalzlo,
1357  Array4<Real const> const& bcvalxhi,
1358  Array4<Real const> const& bcvalyhi,
1359  Array4<Real const> const& bcvalzhi,
1360  Array2D<BoundCond,
1361  0,2*AMREX_SPACEDIM,
1362  0,AMREX_SPACEDIM> const& bct,
1363  Array2D<Real,
1364  0,2*AMREX_SPACEDIM,
1365  0,AMREX_SPACEDIM> const& bcl,
1366  int inhomog, int maxorder,
1367  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1368  Dim3 const& dlo, Dim3 const& dhi) noexcept
1369 
1370 {
1371  const auto blen = amrex::length(vbox);
1372  const auto vlo = amrex::lbound(vbox);
1373  const auto vhi = amrex::ubound(vbox);
1374  bool xlo_domain = (vlo.x == dlo.x);
1375  bool ylo_domain = (vlo.y == dlo.y);
1376  bool zlo_domain = (vlo.z == dlo.z);
1377  bool xhi_domain = (vhi.x == dhi.x);
1378  bool yhi_domain = (vhi.y == dhi.y);
1379  bool zhi_domain = (vhi.z == dhi.z);
1380 
1381  for (int k = vlo.z; k <= vhi.z; ++k) {
1382  mltensor_fill_edges_xlo_ylo(vlo.x-1, vlo.y-1, k, blen, vel, mxlo, mylo, bcvalxlo, bcvalylo,
1383  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, ylo_domain);
1384  mltensor_fill_edges_xhi_ylo(vhi.x+1, vlo.y-1, k, blen, vel, mxhi, mylo, bcvalxhi, bcvalylo,
1385  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, ylo_domain);
1386  mltensor_fill_edges_xlo_yhi(vlo.x-1, vhi.y+1, k, blen, vel, mxlo, myhi, bcvalxlo, bcvalyhi,
1387  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, yhi_domain);
1388  mltensor_fill_edges_xhi_yhi(vhi.x+1, vhi.y+1, k, blen, vel, mxhi, myhi, bcvalxhi, bcvalyhi,
1389  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, yhi_domain);
1390  }
1391 
1392  for (int j = vlo.y; j <= vhi.y; ++j) {
1393  mltensor_fill_edges_xlo_zlo(vlo.x-1, j, vlo.z-1, blen, vel, mxlo, mzlo, bcvalxlo, bcvalzlo,
1394  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, zlo_domain);
1395  mltensor_fill_edges_xhi_zlo(vhi.x+1, j, vlo.z-1, blen, vel, mxhi, mzlo, bcvalxhi, bcvalzlo,
1396  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, zlo_domain);
1397  mltensor_fill_edges_xlo_zhi(vlo.x-1, j, vhi.z+1, blen, vel, mxlo, mzhi, bcvalxlo, bcvalzhi,
1398  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, zhi_domain);
1399  mltensor_fill_edges_xhi_zhi(vhi.x+1, j, vhi.z+1, blen, vel, mxhi, mzhi, bcvalxhi, bcvalzhi,
1400  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, zhi_domain);
1401  }
1402 
1403  for (int i = vlo.x; i <= vhi.x; ++i) {
1404  mltensor_fill_edges_ylo_zlo(i, vlo.y-1, vlo.z-1, blen, vel, mylo, mzlo, bcvalylo, bcvalzlo,
1405  bct, bcl, inhomog, maxorder, dxinv, ylo_domain, zlo_domain);
1406  mltensor_fill_edges_yhi_zlo(i, vhi.y+1, vlo.z-1, blen, vel, myhi, mzlo, bcvalyhi, bcvalzlo,
1407  bct, bcl, inhomog, maxorder, dxinv, yhi_domain, zlo_domain);
1408  mltensor_fill_edges_ylo_zhi(i, vlo.y-1, vhi.z+1, blen, vel, mylo, mzhi, bcvalylo, bcvalzhi,
1409  bct, bcl, inhomog, maxorder, dxinv, ylo_domain, zhi_domain);
1410  mltensor_fill_edges_yhi_zhi(i, vhi.y+1, vhi.z+1, blen, vel, myhi, mzhi, bcvalyhi, bcvalzhi,
1411  bct, bcl, inhomog, maxorder, dxinv, yhi_domain, zhi_domain);
1412  }
1413 }
1414 
1415 #ifdef AMREX_USE_GPU
1416 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
1417 void mltensor_fill_edges (int const bid, int const tid, int const bdim,
1418  Box const& vbox, // vbox: the valid box
1419  Array4<Real> const& vel,
1420  Array4<int const> const& mxlo,
1421  Array4<int const> const& mylo,
1422  Array4<int const> const& mzlo,
1423  Array4<int const> const& mxhi,
1424  Array4<int const> const& myhi,
1425  Array4<int const> const& mzhi,
1426  Array4<Real const> const& bcvalxlo,
1427  Array4<Real const> const& bcvalylo,
1428  Array4<Real const> const& bcvalzlo,
1429  Array4<Real const> const& bcvalxhi,
1430  Array4<Real const> const& bcvalyhi,
1431  Array4<Real const> const& bcvalzhi,
1432  Array2D<BoundCond,
1433  0,2*AMREX_SPACEDIM,
1434  0,AMREX_SPACEDIM> const& bct,
1435  Array2D<Real,
1436  0,2*AMREX_SPACEDIM,
1437  0,AMREX_SPACEDIM> const& bcl,
1438  int inhomog, int maxorder,
1439  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1440  Dim3 const& dlo, Dim3 const& dhi) noexcept
1441 {
1442  const auto blen = amrex::length(vbox);
1443  const auto vlo = amrex::lbound(vbox);
1444  const auto vhi = amrex::ubound(vbox);
1445  bool xlo_domain = (vlo.x == dlo.x);
1446  bool ylo_domain = (vlo.y == dlo.y);
1447  bool zlo_domain = (vlo.z == dlo.z);
1448  bool xhi_domain = (vhi.x == dhi.x);
1449  bool yhi_domain = (vhi.y == dhi.y);
1450  bool zhi_domain = (vhi.z == dhi.z);
1451  if (bid == 0) {
1452  for (int k = vlo.z + tid; k <= vhi.z; k += bdim) {
1453  mltensor_fill_edges_xlo_ylo(vlo.x-1, vlo.y-1, k, blen, vel, mxlo, mylo, bcvalxlo, bcvalylo,
1454  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, ylo_domain);
1455  }
1456  } else if (bid == 1) {
1457  for (int k = vlo.z + tid; k <= vhi.z; k += bdim) {
1458  mltensor_fill_edges_xhi_ylo(vhi.x+1, vlo.y-1, k, blen, vel, mxhi, mylo, bcvalxhi, bcvalylo,
1459  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, ylo_domain);
1460  }
1461  } else if (bid == 2) {
1462  for (int k = vlo.z + tid; k <= vhi.z; k += bdim) {
1463  mltensor_fill_edges_xlo_yhi(vlo.x-1, vhi.y+1, k, blen, vel, mxlo, myhi, bcvalxlo, bcvalyhi,
1464  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, yhi_domain);
1465  }
1466  } else if (bid == 3) {
1467  for (int k = vlo.z + tid; k <= vhi.z; k += bdim) {
1468  mltensor_fill_edges_xhi_yhi(vhi.x+1, vhi.y+1, k, blen, vel, mxhi, myhi, bcvalxhi, bcvalyhi,
1469  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, yhi_domain);
1470  }
1471  } else if (bid == 4) {
1472  for (int j = vlo.y + tid; j <= vhi.y; j += bdim) {
1473  mltensor_fill_edges_xlo_zlo(vlo.x-1, j, vlo.z-1, blen, vel, mxlo, mzlo, bcvalxlo, bcvalzlo,
1474  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, zlo_domain);
1475  }
1476  } else if (bid == 5) {
1477  for (int j = vlo.y + tid; j <= vhi.y; j += bdim) {
1478  mltensor_fill_edges_xhi_zlo(vhi.x+1, j, vlo.z-1, blen, vel, mxhi, mzlo, bcvalxhi, bcvalzlo,
1479  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, zlo_domain);
1480  }
1481  } else if (bid == 6) {
1482  for (int j = vlo.y + tid; j <= vhi.y; j += bdim) {
1483  mltensor_fill_edges_xlo_zhi(vlo.x-1, j, vhi.z+1, blen, vel, mxlo, mzhi, bcvalxlo, bcvalzhi,
1484  bct, bcl, inhomog, maxorder, dxinv, xlo_domain, zhi_domain);
1485  }
1486  } else if (bid == 7) {
1487  for (int j = vlo.y + tid; j <= vhi.y; j += bdim) {
1488  mltensor_fill_edges_xhi_zhi(vhi.x+1, j, vhi.z+1, blen, vel, mxhi, mzhi, bcvalxhi, bcvalzhi,
1489  bct, bcl, inhomog, maxorder, dxinv, xhi_domain, zhi_domain);
1490  }
1491  } else if (bid == 8) {
1492  for (int i = vlo.x + tid; i <= vhi.x; i += bdim) {
1493  mltensor_fill_edges_ylo_zlo(i, vlo.y-1, vlo.z-1, blen, vel, mylo, mzlo, bcvalylo, bcvalzlo,
1494  bct, bcl, inhomog, maxorder, dxinv, ylo_domain, zlo_domain);
1495  }
1496  } else if (bid == 9) {
1497  for (int i = vlo.x + tid; i <= vhi.x; i += bdim) {
1498  mltensor_fill_edges_yhi_zlo(i, vhi.y+1, vlo.z-1, blen, vel, myhi, mzlo, bcvalyhi, bcvalzlo,
1499  bct, bcl, inhomog, maxorder, dxinv, yhi_domain, zlo_domain);
1500  }
1501  } else if (bid == 10) {
1502  for (int i = vlo.x + tid; i <= vhi.x; i += bdim) {
1503  mltensor_fill_edges_ylo_zhi(i, vlo.y-1, vhi.z+1, blen, vel, mylo, mzhi, bcvalylo, bcvalzhi,
1504  bct, bcl, inhomog, maxorder, dxinv, ylo_domain, zhi_domain);
1505  }
1506  } else if (bid == 11) {
1507  for (int i = vlo.x + tid; i <= vhi.x; i += bdim) {
1508  mltensor_fill_edges_yhi_zhi(i, vhi.y+1, vhi.z+1, blen, vel, myhi, mzhi, bcvalyhi, bcvalzhi,
1509  bct, bcl, inhomog, maxorder, dxinv, yhi_domain, zhi_domain);
1510  }
1511  }
1512 }
1513 #endif
1514 
1515 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1516 Real mltensor_dz_on_xface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dzi) noexcept
1517 {
1518  return (vel(i,j,k+1,n)+vel(i-1,j,k+1,n)-vel(i,j,k-1,n)-vel(i-1,j,k-1,n))*(Real(0.25)*dzi);
1519 }
1520 
1521 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1522 Real mltensor_dz_on_yface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dzi) noexcept
1523 {
1524  return (vel(i,j,k+1,n)+vel(i,j-1,k+1,n)-vel(i,j,k-1,n)-vel(i,j-1,k-1,n))*(Real(0.25)*dzi);
1525 }
1526 
1527 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1528 Real mltensor_dx_on_zface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dxi) noexcept
1529 {
1530  return (vel(i+1,j,k,n)+vel(i+1,j,k-1,n)-vel(i-1,j,k,n)-vel(i-1,j,k-1,n))*(Real(0.25)*dxi);
1531 }
1532 
1533 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1534 Real mltensor_dy_on_zface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dyi) noexcept
1535 {
1536  return (vel(i,j+1,k,n)+vel(i,j+1,k-1,n)-vel(i,j-1,k,n)-vel(i,j-1,k-1,n))*(Real(0.25)*dyi);
1537 }
1538 
1539 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1540 void mltensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
1541  Array4<Real const> const& vel,
1542  Array4<Real const> const& etax,
1543  Array4<Real const> const& kapx,
1544  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
1545 {
1546  const Real dyi = dxinv[1];
1547  const Real dzi = dxinv[2];
1548  const auto lo = amrex::lbound(box);
1549  const auto hi = amrex::ubound(box);
1550  constexpr Real twoThirds = Real(2./3.);
1551 
1552  for (int k = lo.z; k <= hi.z; ++k) {
1553  for (int j = lo.y; j <= hi.y; ++j) {
1554  AMREX_PRAGMA_SIMD
1555  for (int i = lo.x; i <= hi.x; ++i) {
1556  Real dudy = mltensor_dy_on_xface(i,j,k,0,vel,dyi);
1557  Real dvdy = mltensor_dy_on_xface(i,j,k,1,vel,dyi);
1558  Real dudz = mltensor_dz_on_xface(i,j,k,0,vel,dzi);
1559  Real dwdz = mltensor_dz_on_xface(i,j,k,2,vel,dzi);
1560  Real divu = dvdy + dwdz;
1561  Real xif = kapx(i,j,k);
1562  Real mun = Real(0.75)*(etax(i,j,k,0)-xif); // restore the original eta
1563  Real mut = etax(i,j,k,1);
1564  fx(i,j,k,0) = -mun*(-twoThirds*divu) - xif*divu;
1565  fx(i,j,k,1) = -mut*(dudy);
1566  fx(i,j,k,2) = -mut*(dudz);
1567  }
1568  }
1569  }
1570 }
1571 
1572 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1573 void mltensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
1574  Array4<Real const> const& vel,
1575  Array4<Real const> const& etay,
1576  Array4<Real const> const& kapy,
1577  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
1578 {
1579  const Real dxi = dxinv[0];
1580  const Real dzi = dxinv[2];
1581  const auto lo = amrex::lbound(box);
1582  const auto hi = amrex::ubound(box);
1583  constexpr Real twoThirds = Real(2./3.);
1584 
1585  for (int k = lo.z; k <= hi.z; ++k) {
1586  for (int j = lo.y; j <= hi.y; ++j) {
1587  AMREX_PRAGMA_SIMD
1588  for (int i = lo.x; i <= hi.x; ++i) {
1589  Real dudx = mltensor_dx_on_yface(i,j,k,0,vel,dxi);
1590  Real dvdx = mltensor_dx_on_yface(i,j,k,1,vel,dxi);
1591  Real dvdz = mltensor_dz_on_yface(i,j,k,1,vel,dzi);
1592  Real dwdz = mltensor_dz_on_yface(i,j,k,2,vel,dzi);
1593  Real divu = dudx + dwdz;
1594  Real xif = kapy(i,j,k);
1595  Real mun = Real(0.75)*(etay(i,j,k,1)-xif); // restore the original eta
1596  Real mut = etay(i,j,k,0);
1597  fy(i,j,k,0) = -mut*(dvdx);
1598  fy(i,j,k,1) = -mun*(-twoThirds*divu) - xif*divu;
1599  fy(i,j,k,2) = -mut*(dvdz);
1600  }
1601  }
1602  }
1603 }
1604 
1605 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1606 void mltensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
1607  Array4<Real const> const& vel,
1608  Array4<Real const> const& etaz,
1609  Array4<Real const> const& kapz,
1610  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
1611 {
1612  const Real dxi = dxinv[0];
1613  const Real dyi = dxinv[1];
1614  const auto lo = amrex::lbound(box);
1615  const auto hi = amrex::ubound(box);
1616  constexpr Real twoThirds = Real(2./3.);
1617 
1618  for (int k = lo.z; k <= hi.z; ++k) {
1619  for (int j = lo.y; j <= hi.y; ++j) {
1620  AMREX_PRAGMA_SIMD
1621  for (int i = lo.x; i <= hi.x; ++i) {
1622  Real dudx = mltensor_dx_on_zface(i,j,k,0,vel,dxi);
1623  Real dwdx = mltensor_dx_on_zface(i,j,k,2,vel,dxi);
1624  Real dvdy = mltensor_dy_on_zface(i,j,k,1,vel,dyi);
1625  Real dwdy = mltensor_dy_on_zface(i,j,k,2,vel,dyi);
1626  Real divu = dudx + dvdy;
1627  Real xif = kapz(i,j,k);
1628  Real mun = Real(0.75)*(etaz(i,j,k,2)-xif); // restore the original eta
1629  Real mut = etaz(i,j,k,0);
1630  fz(i,j,k,0) = -mut*(dwdx);
1631  fz(i,j,k,1) = -mut*(dwdy);
1632  fz(i,j,k,2) = -mun*(-twoThirds*divu) - xif*divu;
1633  }
1634  }
1635  }
1636 }
1637 
1638 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1639 Real mltensor_dz_on_xface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dzi,
1640  Array4<Real const> const& bvxlo, Array4<Real const> const& bvxhi,
1641  Array2D<BoundCond,
1642  0,2*AMREX_SPACEDIM,
1643  0,AMREX_SPACEDIM> const& bct,
1644  Dim3 const& dlo, Dim3 const& dhi) noexcept
1645 {
1646  Real ddz;
1647  if (i == dlo.x) {
1648  if (bct(Orientation::xlo(),n) == AMREX_LO_DIRICHLET && bvxlo) {
1649  if (k == dlo.z) {
1650  ddz = (bvxlo(i-1,j,k ,n) * Real(-1.5) +
1651  bvxlo(i-1,j,k+1,n) * Real(2.) +
1652  bvxlo(i-1,j,k+2,n) * Real(-0.5)) * dzi;
1653  } else if (k == dhi.z) {
1654  ddz = -(bvxlo(i-1,j,k ,n) * Real(-1.5) +
1655  bvxlo(i-1,j,k-1,n) * Real(2.) +
1656  bvxlo(i-1,j,k-2,n) * Real(-0.5)) * dzi;
1657  } else {
1658  ddz = (bvxlo(i-1,j,k+1,n)-bvxlo(i-1,j,k-1,n))*(Real(0.5)*dzi);
1659  }
1660  } else if (bct(Orientation::xlo(),n) == AMREX_LO_NEUMANN) {
1661  ddz = (vel(i,j,k+1,n)-vel(i,j,k-1,n))*(Real(0.5)*dzi);
1662  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1663  ddz = Real(0.);
1664  }
1665  } else if (i == dhi.x+1) {
1666  if (bct(Orientation::xhi(),n) == AMREX_LO_DIRICHLET && bvxhi) {
1667  if (k == dlo.z) {
1668  ddz = (bvxhi(i,j,k ,n) * Real(-1.5) +
1669  bvxhi(i,j,k+1,n) * Real(2.) +
1670  bvxhi(i,j,k+2,n) * Real(-0.5)) * dzi;
1671  } else if (k == dhi.z) {
1672  ddz = -(bvxhi(i,j,k ,n) * Real(-1.5) +
1673  bvxhi(i,j,k-1,n) * Real(2.) +
1674  bvxhi(i,j,k-2,n) * Real(-0.5)) * dzi;
1675  } else {
1676  ddz = (bvxhi(i,j,k+1,n)-bvxhi(i,j,k-1,n))*(Real(0.5)*dzi);
1677  }
1678  } else if (bct(Orientation::xhi(),n) == AMREX_LO_NEUMANN) {
1679  ddz = (vel(i-1,j,k+1,n)-vel(i-1,j,k-1,n))*(Real(0.5)*dzi);
1680  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1681  ddz = Real(0.);
1682  }
1683  } else {
1684  ddz = mltensor_dz_on_xface(i,j,k,n,vel,dzi);
1685  }
1686  return ddz;
1687 }
1688 
1689 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1690 Real mltensor_dz_on_yface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dzi,
1691  Array4<Real const> const& bvylo, Array4<Real const> const& bvyhi,
1692  Array2D<BoundCond,
1693  0,2*AMREX_SPACEDIM,
1694  0,AMREX_SPACEDIM> const& bct,
1695  Dim3 const& dlo, Dim3 const& dhi) noexcept
1696 {
1697  Real ddz;
1698  if (j == dlo.y) {
1699  if (bct(Orientation::ylo(),n) == AMREX_LO_DIRICHLET && bvylo) {
1700  if (k == dlo.z) {
1701  ddz = (bvylo(i,j-1,k ,n) * Real(-1.5) +
1702  bvylo(i,j-1,k+1,n) * Real(2.) +
1703  bvylo(i,j-1,k+2,n) * Real(-0.5)) * dzi;
1704  } else if (k == dhi.z) {
1705  ddz = -(bvylo(i,j-1,k ,n) * Real(-1.5) +
1706  bvylo(i,j-1,k-1,n) * Real(2.) +
1707  bvylo(i,j-1,k-2,n) * Real(-0.5)) * dzi;
1708  } else {
1709  ddz = (bvylo(i,j-1,k+1,n)-bvylo(i,j-1,k-1,n))*(Real(0.5)*dzi);
1710  }
1711  } else if (bct(Orientation::ylo(),n) == AMREX_LO_NEUMANN) {
1712  ddz = (vel(i,j,k+1,n)-vel(i,j,k-1,n))*(Real(0.5)*dzi);
1713  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1714  ddz = Real(0.);
1715  }
1716  } else if (j == dhi.y+1) {
1717  if (bct(Orientation::yhi(),n) == AMREX_LO_DIRICHLET && bvyhi) {
1718  if (k == dlo.z) {
1719  ddz = (bvyhi(i,j,k ,n) * Real(-1.5) +
1720  bvyhi(i,j,k+1,n) * Real(2.) +
1721  bvyhi(i,j,k+2,n) * Real(-0.5)) * dzi;
1722  } else if (k == dhi.z) {
1723  ddz = -(bvyhi(i,j,k ,n) * Real(-1.5) +
1724  bvyhi(i,j,k-1,n) * Real(2.) +
1725  bvyhi(i,j,k-2,n) * Real(-0.5)) * dzi;
1726  } else {
1727  ddz = (bvyhi(i,j,k+1,n)-bvyhi(i,j,k-1,n))*(Real(0.5)*dzi);
1728  }
1729  } else if (bct(Orientation::yhi(),n) == AMREX_LO_NEUMANN) {
1730  ddz = (vel(i,j-1,k+1,n)-vel(i,j-1,k-1,n))*(Real(0.5)*dzi);
1731  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1732  ddz = Real(0.);
1733  }
1734  } else {
1735  ddz = mltensor_dz_on_yface(i,j,k,n,vel,dzi);
1736  }
1737  return ddz;
1738 }
1739 
1740 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1741 Real mltensor_dx_on_zface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dxi,
1742  Array4<Real const> const& bvzlo, Array4<Real const> const& bvzhi,
1743  Array2D<BoundCond,
1744  0,2*AMREX_SPACEDIM,
1745  0,AMREX_SPACEDIM> const& bct,
1746  Dim3 const& dlo, Dim3 const& dhi) noexcept
1747 {
1748  Real ddx;
1749  if (k == dlo.z) {
1750  if (bct(Orientation::zlo(),n) == AMREX_LO_DIRICHLET && bvzlo) {
1751  if (i == dlo.x) {
1752  ddx = (bvzlo(i ,j,k-1,n) * Real(-1.5) +
1753  bvzlo(i+1,j,k-1,n) * Real(2.) +
1754  bvzlo(i+2,j,k-1,n) * Real(-0.5)) * dxi;
1755  } else if (i == dhi.x) {
1756  ddx = -(bvzlo(i ,j,k-1,n) * Real(-1.5) +
1757  bvzlo(i-1,j,k-1,n) * Real(2.) +
1758  bvzlo(i-2,j,k-1,n) * Real(-0.5)) * dxi;
1759  } else {
1760  ddx = (bvzlo(i+1,j,k-1,n)-bvzlo(i-1,j,k-1,n))*(Real(0.5)*dxi);
1761  }
1762  } else if (bct(Orientation::zlo(),n) == AMREX_LO_NEUMANN) {
1763  ddx = (vel(i+1,j,k,n)-vel(i-1,j,k,n))*(Real(0.5)*dxi);
1764  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1765  ddx = Real(0.);
1766  }
1767  } else if (k == dhi.z+1) {
1768  if (bct(Orientation::zhi(),n) == AMREX_LO_DIRICHLET && bvzhi) {
1769  if (i == dlo.x) {
1770  ddx = (bvzhi(i ,j,k,n) * Real(-1.5) +
1771  bvzhi(i+1,j,k,n) * Real(2.) +
1772  bvzhi(i+2,j,k,n) * Real(-0.5)) * dxi;
1773  } else if (i == dhi.x) {
1774  ddx = -(bvzhi(i ,j,k,n) * Real(-1.5) +
1775  bvzhi(i-1,j,k,n) * Real(2.) +
1776  bvzhi(i-2,j,k,n) * Real(-0.5)) * dxi;
1777  } else {
1778  ddx = (bvzhi(i+1,j,k,n)-bvzhi(i-1,j,k,n))*(Real(0.5)*dxi);
1779  }
1780  } else if (bct(Orientation::zhi(),n) == AMREX_LO_NEUMANN) {
1781  ddx = (vel(i+1,j,k-1,n)-vel(i-1,j,k-1,n))*(Real(0.5)*dxi);
1782  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1783  ddx = Real(0.);
1784  }
1785  } else {
1786  ddx = mltensor_dx_on_zface(i,j,k,n,vel,dxi);
1787  }
1788  return ddx;
1789 }
1790 
1791 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1792 Real mltensor_dy_on_zface (int i, int j, int k, int n, Array4<Real const> const& vel, Real dyi,
1793  Array4<Real const> const& bvzlo, Array4<Real const> const& bvzhi,
1794  Array2D<BoundCond,
1795  0,2*AMREX_SPACEDIM,
1796  0,AMREX_SPACEDIM> const& bct,
1797  Dim3 const& dlo, Dim3 const& dhi) noexcept
1798 {
1799  Real ddy;
1800  if (k == dlo.z) {
1801  if (bct(Orientation::zlo(),n) == AMREX_LO_DIRICHLET && bvzlo) {
1802  if (j == dlo.y) {
1803  ddy = (bvzlo(i,j ,k-1,n) * Real(-1.5) +
1804  bvzlo(i,j+1,k-1,n) * Real(2.) +
1805  bvzlo(i,j+2,k-1,n) * Real(-0.5)) * dyi;
1806  } else if (j == dhi.y) {
1807  ddy = -(bvzlo(i,j ,k-1,n) * Real(-1.5) +
1808  bvzlo(i,j-1,k-1,n) * Real(2.) +
1809  bvzlo(i,j-2,k-1,n) * Real(-0.5)) * dyi;
1810  } else {
1811  ddy = (bvzlo(i,j+1,k-1,n)-bvzlo(i,j-1,k-1,n))*(Real(0.5)*dyi);
1812  }
1813  } else if (bct(Orientation::zlo(),n) == AMREX_LO_NEUMANN) {
1814  ddy = (vel(i,j+1,k,n)-vel(i,j-1,k,n))*(Real(0.5)*dyi);
1815  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1816  ddy = Real(0.);
1817  }
1818  } else if (k == dhi.z+1) {
1819  if (bct(Orientation::zhi(),n) == AMREX_LO_DIRICHLET && bvzhi) {
1820  if (j == dlo.y) {
1821  ddy = (bvzhi(i,j ,k,n) * Real(-1.5) +
1822  bvzhi(i,j+1,k,n) * Real(2.) +
1823  bvzhi(i,j+2,k,n) * Real(-0.5)) * dyi;
1824  } else if (j == dhi.y) {
1825  ddy = -(bvzhi(i,j ,k,n) * Real(-1.5) +
1826  bvzhi(i,j-1,k,n) * Real(2.) +
1827  bvzhi(i,j-2,k,n) * Real(-0.5)) * dyi;
1828  } else {
1829  ddy = (bvzhi(i,j+1,k,n)-bvzhi(i,j-1,k,n))*(Real(0.5)*dyi);
1830  }
1831  } else if (bct(Orientation::zhi(),n) == AMREX_LO_NEUMANN) {
1832  ddy = (vel(i,j+1,k-1,n)-vel(i,j-1,k-1,n))*(Real(0.5)*dyi);
1833  } else { // AMREX_LO_REFLECT_ODD or homogeneous Dirichlet
1834  ddy = Real(0.);
1835  }
1836  } else {
1837  ddy = mltensor_dy_on_zface(i,j,k,n,vel,dyi);
1838  }
1839  return ddy;
1840 }
1841 
1842 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1843 void mltensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
1844  Array4<Real const> const& vel,
1845  Array4<Real const> const& etax,
1846  Array4<Real const> const& kapx,
1847  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1848  Array4<Real const> const& bvxlo,
1849  Array4<Real const> const& bvxhi,
1850  Array2D<BoundCond,
1851  0,2*AMREX_SPACEDIM,
1852  0,AMREX_SPACEDIM> const& bct,
1853  Dim3 const& dlo, Dim3 const& dhi) noexcept
1854 {
1855  const Real dyi = dxinv[1];
1856  const Real dzi = dxinv[2];
1857  const auto lo = amrex::lbound(box);
1858  const auto hi = amrex::ubound(box);
1859  constexpr Real twoThirds = Real(2./3.);
1860 
1861  for (int k = lo.z; k <= hi.z; ++k) {
1862  for (int j = lo.y; j <= hi.y; ++j) {
1863  for (int i = lo.x; i <= hi.x; ++i) {
1864  Real dudy = mltensor_dy_on_xface(i,j,k,0,vel,dyi,bvxlo,bvxhi,bct,dlo,dhi);
1865  Real dvdy = mltensor_dy_on_xface(i,j,k,1,vel,dyi,bvxlo,bvxhi,bct,dlo,dhi);
1866  Real dudz = mltensor_dz_on_xface(i,j,k,0,vel,dzi,bvxlo,bvxhi,bct,dlo,dhi);
1867  Real dwdz = mltensor_dz_on_xface(i,j,k,2,vel,dzi,bvxlo,bvxhi,bct,dlo,dhi);
1868  Real divu = dvdy + dwdz;
1869  Real xif = kapx(i,j,k);
1870  Real mun = Real(0.75)*(etax(i,j,k,0)-xif); // restore the original eta
1871  Real mut = etax(i,j,k,1);
1872  fx(i,j,k,0) = -mun*(-twoThirds*divu) - xif*divu;
1873  fx(i,j,k,1) = -mut*(dudy);
1874  fx(i,j,k,2) = -mut*(dudz);
1875  }
1876  }
1877  }
1878 }
1879 
1880 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1881 void mltensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
1882  Array4<Real const> const& vel,
1883  Array4<Real const> const& etay,
1884  Array4<Real const> const& kapy,
1885  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1886  Array4<Real const> const& bvylo,
1887  Array4<Real const> const& bvyhi,
1888  Array2D<BoundCond,
1889  0,2*AMREX_SPACEDIM,
1890  0,AMREX_SPACEDIM> const& bct,
1891  Dim3 const& dlo, Dim3 const& dhi) noexcept
1892 {
1893  const Real dxi = dxinv[0];
1894  const Real dzi = dxinv[2];
1895  const auto lo = amrex::lbound(box);
1896  const auto hi = amrex::ubound(box);
1897  constexpr Real twoThirds = Real(2./3.);
1898 
1899  for (int k = lo.z; k <= hi.z; ++k) {
1900  for (int j = lo.y; j <= hi.y; ++j) {
1901  for (int i = lo.x; i <= hi.x; ++i) {
1902  Real dudx = mltensor_dx_on_yface(i,j,k,0,vel,dxi,bvylo,bvyhi,bct,dlo,dhi);
1903  Real dvdx = mltensor_dx_on_yface(i,j,k,1,vel,dxi,bvylo,bvyhi,bct,dlo,dhi);
1904  Real dvdz = mltensor_dz_on_yface(i,j,k,1,vel,dzi,bvylo,bvyhi,bct,dlo,dhi);
1905  Real dwdz = mltensor_dz_on_yface(i,j,k,2,vel,dzi,bvylo,bvyhi,bct,dlo,dhi);
1906  Real divu = dudx + dwdz;
1907  Real xif = kapy(i,j,k);
1908  Real mun = Real(0.75)*(etay(i,j,k,1)-xif); // restore the original eta
1909  Real mut = etay(i,j,k,0);
1910  fy(i,j,k,0) = -mut*(dvdx);
1911  fy(i,j,k,1) = -mun*(-twoThirds*divu) - xif*divu;
1912  fy(i,j,k,2) = -mut*(dvdz);
1913  }
1914  }
1915  }
1916 }
1917 
1918 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1919 void mltensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
1920  Array4<Real const> const& vel,
1921  Array4<Real const> const& etaz,
1922  Array4<Real const> const& kapz,
1923  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1924  Array4<Real const> const& bvzlo,
1925  Array4<Real const> const& bvzhi,
1926  Array2D<BoundCond,
1927  0,2*AMREX_SPACEDIM,
1928  0,AMREX_SPACEDIM> const& bct,
1929  Dim3 const& dlo, Dim3 const& dhi) noexcept
1930 {
1931  const Real dxi = dxinv[0];
1932  const Real dyi = dxinv[1];
1933  const auto lo = amrex::lbound(box);
1934  const auto hi = amrex::ubound(box);
1935  constexpr Real twoThirds = Real(2./3.);
1936 
1937  for (int k = lo.z; k <= hi.z; ++k) {
1938  for (int j = lo.y; j <= hi.y; ++j) {
1939  for (int i = lo.x; i <= hi.x; ++i) {
1940  Real dudx = mltensor_dx_on_zface(i,j,k,0,vel,dxi,bvzlo,bvzhi,bct,dlo,dhi);
1941  Real dwdx = mltensor_dx_on_zface(i,j,k,2,vel,dxi,bvzlo,bvzhi,bct,dlo,dhi);
1942  Real dvdy = mltensor_dy_on_zface(i,j,k,1,vel,dyi,bvzlo,bvzhi,bct,dlo,dhi);
1943  Real dwdy = mltensor_dy_on_zface(i,j,k,2,vel,dyi,bvzlo,bvzhi,bct,dlo,dhi);
1944  Real divu = dudx + dvdy;
1945  Real xif = kapz(i,j,k);
1946  Real mun = Real(0.75)*(etaz(i,j,k,2)-xif); // restore the original eta
1947  Real mut = etaz(i,j,k,0);
1948  fz(i,j,k,0) = -mut*(dwdx);
1949  fz(i,j,k,1) = -mut*(dwdy);
1950  fz(i,j,k,2) = -mun*(-twoThirds*divu) - xif*divu;
1951  }
1952  }
1953  }
1954 }
1955 
1956 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1957 void mltensor_cross_terms (Box const& box, Array4<Real> const& Ax,
1958  Array4<Real const> const& fx,
1959  Array4<Real const> const& fy,
1960  Array4<Real const> const& fz,
1961  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1962  Real bscalar) noexcept
1963 {
1964  const Real dxi = bscalar * dxinv[0];
1965  const Real dyi = bscalar * dxinv[1];
1966  const Real dzi = bscalar * dxinv[2];
1967  const auto lo = amrex::lbound(box);
1968  const auto hi = amrex::ubound(box);
1969 
1970  for (int k = lo.z; k <= hi.z; ++k) {
1971  for (int j = lo.y; j <= hi.y; ++j) {
1972  AMREX_PRAGMA_SIMD
1973  for (int i = lo.x; i <= hi.x; ++i) {
1974  Ax(i,j,k,0) += dxi*(fx(i+1,j ,k ,0) - fx(i,j,k,0))
1975  + dyi*(fy(i ,j+1,k ,0) - fy(i,j,k,0))
1976  + dzi*(fz(i ,j ,k+1,0) - fz(i,j,k,0));
1977  Ax(i,j,k,1) += dxi*(fx(i+1,j ,k ,1) - fx(i,j,k,1))
1978  + dyi*(fy(i ,j+1,k ,1) - fy(i,j,k,1))
1979  + dzi*(fz(i ,j ,k+1,1) - fz(i,j,k,1));
1980  Ax(i,j,k,2) += dxi*(fx(i+1,j ,k ,2) - fx(i,j,k,2))
1981  + dyi*(fy(i ,j+1,k ,2) - fy(i,j,k,2))
1982  + dzi*(fz(i ,j ,k+1,2) - fz(i,j,k,2));
1983  }
1984  }
1985  }
1986 }
1987 
1988 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1989 void mltensor_cross_terms_os (Box const& box, Array4<Real> const& Ax,
1990  Array4<Real const> const& fx,
1991  Array4<Real const> const& fy,
1992  Array4<Real const> const& fz,
1993  Array4<int const> const& osm,
1994  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
1995  Real bscalar) noexcept
1996 {
1997  const Real dxi = bscalar * dxinv[0];
1998  const Real dyi = bscalar * dxinv[1];
1999  const Real dzi = bscalar * dxinv[2];
2000  const auto lo = amrex::lbound(box);
2001  const auto hi = amrex::ubound(box);
2002 
2003  for (int k = lo.z; k <= hi.z; ++k) {
2004  for (int j = lo.y; j <= hi.y; ++j) {
2005  AMREX_PRAGMA_SIMD
2006  for (int i = lo.x; i <= hi.x; ++i) {
2007  if (osm(i,j,k) == 0) {
2008  Ax(i,j,k,0) = 0.0;
2009  Ax(i,j,k,1) = 0.0;
2010  Ax(i,j,k,2) = 0.0;
2011  } else {
2012  Ax(i,j,k,0) += dxi*(fx(i+1,j ,k ,0) - fx(i,j,k,0))
2013  + dyi*(fy(i ,j+1,k ,0) - fy(i,j,k,0))
2014  + dzi*(fz(i ,j ,k+1,0) - fz(i,j,k,0));
2015  Ax(i,j,k,1) += dxi*(fx(i+1,j ,k ,1) - fx(i,j,k,1))
2016  + dyi*(fy(i ,j+1,k ,1) - fy(i,j,k,1))
2017  + dzi*(fz(i ,j ,k+1,1) - fz(i,j,k,1));
2018  Ax(i,j,k,2) += dxi*(fx(i+1,j ,k ,2) - fx(i,j,k,2))
2019  + dyi*(fy(i ,j+1,k ,2) - fy(i,j,k,2))
2020  + dzi*(fz(i ,j ,k+1,2) - fz(i,j,k,2));
2021  }
2022  }
2023  }
2024  }
2025 }
2026 
2027 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
2028 void mltensor_vel_grads_fx (Box const& box, Array4<Real> const& fx,
2029  Array4<Real const> const& vel,
2030  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
2031 {
2032  const Real dxi = dxinv[0];
2033  const Real dyi = dxinv[1];
2034  const Real dzi = dxinv[2];
2035  const auto lo = amrex::lbound(box);
2036  const auto hi = amrex::ubound(box);
2037 
2038  for (int k = lo.z; k <= hi.z; ++k) {
2039  for (int j = lo.y; j <= hi.y; ++j) {
2040  AMREX_PRAGMA_SIMD
2041  for (int i = lo.x; i <= hi.x; ++i) {
2042 
2043  Real dudx = (vel(i,j,k,0) - vel(i-1,j,k,0))*dxi;
2044  Real dvdx = (vel(i,j,k,1) - vel(i-1,j,k,1))*dxi;
2045  Real dwdx = (vel(i,j,k,2) - vel(i-1,j,k,2))*dxi;
2046 
2047  Real dudy = mltensor_dy_on_xface(i,j,k,0,vel,dyi);
2048  Real dvdy = mltensor_dy_on_xface(i,j,k,1,vel,dyi);
2049  Real dwdy = mltensor_dy_on_xface(i,j,k,2,vel,dyi);
2050 
2051  Real dudz = mltensor_dz_on_xface(i,j,k,0,vel,dzi);
2052  Real dvdz = mltensor_dz_on_xface(i,j,k,1,vel,dzi);
2053  Real dwdz = mltensor_dz_on_xface(i,j,k,2,vel,dzi);
2054 
2055  fx(i,j,k,0) = dudx;
2056  fx(i,j,k,1) = dvdx;
2057  fx(i,j,k,2) = dwdx;
2058  fx(i,j,k,3) = dudy;
2059  fx(i,j,k,4) = dvdy;
2060  fx(i,j,k,5) = dwdy;
2061  fx(i,j,k,6) = dudz;
2062  fx(i,j,k,7) = dvdz;
2063  fx(i,j,k,8) = dwdz;
2064 
2065  }
2066  }
2067  }
2068 }
2069 
2070 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
2071 void mltensor_vel_grads_fy (Box const& box, Array4<Real> const& fy,
2072  Array4<Real const> const& vel,
2073  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
2074 {
2075  const Real dxi = dxinv[0];
2076  const Real dyi = dxinv[1];
2077  const Real dzi = dxinv[2];
2078  const auto lo = amrex::lbound(box);
2079  const auto hi = amrex::ubound(box);
2080 
2081  for (int k = lo.z; k <= hi.z; ++k) {
2082  for (int j = lo.y; j <= hi.y; ++j) {
2083  AMREX_PRAGMA_SIMD
2084  for (int i = lo.x; i <= hi.x; ++i) {
2085 
2086  Real dudx = mltensor_dx_on_yface(i,j,k,0,vel,dxi);
2087  Real dvdx = mltensor_dx_on_yface(i,j,k,1,vel,dxi);
2088  Real dwdx = mltensor_dx_on_yface(i,j,k,2,vel,dxi);
2089 
2090  Real dudy = (vel(i,j,k,0) - vel(i,j-1,k,0))*dyi;
2091  Real dvdy = (vel(i,j,k,1) - vel(i,j-1,k,1))*dyi;
2092  Real dwdy = (vel(i,j,k,2) - vel(i,j-1,k,2))*dyi;
2093 
2094  Real dudz = mltensor_dz_on_yface(i,j,k,0,vel,dzi);
2095  Real dvdz = mltensor_dz_on_yface(i,j,k,1,vel,dzi);
2096  Real dwdz = mltensor_dz_on_yface(i,j,k,2,vel,dzi);
2097 
2098  fy(i,j,k,0) = dudx;
2099  fy(i,j,k,1) = dvdx;
2100  fy(i,j,k,2) = dwdx;
2101  fy(i,j,k,3) = dudy;
2102  fy(i,j,k,4) = dvdy;
2103  fy(i,j,k,5) = dwdy;
2104  fy(i,j,k,6) = dudz;
2105  fy(i,j,k,7) = dvdz;
2106  fy(i,j,k,8) = dwdz;
2107 
2108  }
2109  }
2110  }
2111 }
2112 
2113 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
2114 void mltensor_vel_grads_fz (Box const& box, Array4<Real> const& fz,
2115  Array4<Real const> const& vel,
2116  GpuArray<Real,AMREX_SPACEDIM> const& dxinv) noexcept
2117 {
2118  const Real dxi = dxinv[0];
2119  const Real dyi = dxinv[1];
2120  const Real dzi = dxinv[2];
2121  const auto lo = amrex::lbound(box);
2122  const auto hi = amrex::ubound(box);
2123 
2124  for (int k = lo.z; k <= hi.z; ++k) {
2125  for (int j = lo.y; j <= hi.y; ++j) {
2126  AMREX_PRAGMA_SIMD
2127  for (int i = lo.x; i <= hi.x; ++i) {
2128 
2129  Real dudx = mltensor_dx_on_zface(i,j,k,0,vel,dxi);
2130  Real dvdx = mltensor_dx_on_zface(i,j,k,1,vel,dxi);
2131  Real dwdx = mltensor_dx_on_zface(i,j,k,2,vel,dxi);
2132 
2133  Real dudy = mltensor_dy_on_zface(i,j,k,0,vel,dyi);
2134  Real dvdy = mltensor_dy_on_zface(i,j,k,1,vel,dyi);
2135  Real dwdy = mltensor_dy_on_zface(i,j,k,2,vel,dyi);
2136 
2137  Real dudz = (vel(i,j,k,0) - vel(i,j,k-1,0))*dzi;
2138  Real dvdz = (vel(i,j,k,1) - vel(i,j,k-1,1))*dzi;
2139  Real dwdz = (vel(i,j,k,2) - vel(i,j,k-1,2))*dzi;
2140 
2141  fz(i,j,k,0) = dudx;
2142  fz(i,j,k,1) = dvdx;
2143  fz(i,j,k,2) = dwdx;
2144  fz(i,j,k,3) = dudy;
2145  fz(i,j,k,4) = dvdy;
2146  fz(i,j,k,5) = dwdy;
2147  fz(i,j,k,6) = dudz;
2148  fz(i,j,k,7) = dvdz;
2149  fz(i,j,k,8) = dwdz;
2150 
2151  }
2152  }
2153  }
2154 }
2155 
2156 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
2157 void mltensor_vel_grads_fx (Box const& box, Array4<Real> const& fx,
2158  Array4<Real const> const& vel,
2159  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
2160  Array4<Real const> const& bvxlo,
2161  Array4<Real const> const& bvxhi,
2162  Array2D<BoundCond,
2163  0,2*AMREX_SPACEDIM,
2164  0,AMREX_SPACEDIM> const& bct,
2165  Dim3 const& dlo, Dim3 const& dhi) noexcept
2166 {
2167  const Real dxi = dxinv[0];
2168  const Real dyi = dxinv[1];
2169  const Real dzi = dxinv[2];
2170  const auto lo = amrex::lbound(box);
2171  const auto hi = amrex::ubound(box);
2172 
2173  for (int k = lo.z; k <= hi.z; ++k) {
2174  for (int j = lo.y; j <= hi.y; ++j) {
2175  for (int i = lo.x; i <= hi.x; ++i) {
2176  Real dudx = (vel(i,j,k,0) - vel(i-1,j,k,0))*dxi;
2177  Real dvdx = (vel(i,j,k,1) - vel(i-1,j,k,1))*dxi;
2178  Real dwdx = (vel(i,j,k,2) - vel(i-1,j,k,2))*dxi;
2179  Real dudy = mltensor_dy_on_xface(i,j,k,0,vel,dyi,bvxlo,bvxhi,bct,dlo,dhi);
2180  Real dvdy = mltensor_dy_on_xface(i,j,k,1,vel,dyi,bvxlo,bvxhi,bct,dlo,dhi);
2181  Real dwdy = mltensor_dy_on_xface(i,j,k,2,vel,dyi,bvxlo,bvxhi,bct,dlo,dhi);
2182  Real dudz = mltensor_dz_on_xface(i,j,k,0,vel,dzi,bvxlo,bvxhi,bct,dlo,dhi);
2183  Real dvdz = mltensor_dz_on_xface(i,j,k,1,vel,dzi,bvxlo,bvxhi,bct,dlo,dhi);
2184  Real dwdz = mltensor_dz_on_xface(i,j,k,2,vel,dzi,bvxlo,bvxhi,bct,dlo,dhi);
2185  fx(i,j,k,0) = dudx;
2186  fx(i,j,k,1) = dvdx;
2187  fx(i,j,k,2) = dwdx;
2188  fx(i,j,k,3) = dudy;
2189  fx(i,j,k,4) = dvdy;
2190  fx(i,j,k,5) = dwdy;
2191  fx(i,j,k,6) = dudz;
2192  fx(i,j,k,7) = dvdz;
2193  fx(i,j,k,8) = dwdz;
2194 
2195  }
2196  }
2197  }
2198 }
2199 
2200 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
2201 void mltensor_vel_grads_fy (Box const& box, Array4<Real> const& fy,
2202  Array4<Real const> const& vel,
2203  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
2204  Array4<Real const> const& bvylo,
2205  Array4<Real const> const& bvyhi,
2206  Array2D<BoundCond,
2207  0,2*AMREX_SPACEDIM,
2208  0,AMREX_SPACEDIM> const& bct,
2209  Dim3 const& dlo, Dim3 const& dhi) noexcept
2210 {
2211  const Real dxi = dxinv[0];
2212  const Real dyi = dxinv[1];
2213  const Real dzi = dxinv[2];
2214  const auto lo = amrex::lbound(box);
2215  const auto hi = amrex::ubound(box);
2216 
2217  for (int k = lo.z; k <= hi.z; ++k) {
2218  for (int j = lo.y; j <= hi.y; ++j) {
2219  for (int i = lo.x; i <= hi.x; ++i) {
2220  Real dudx = mltensor_dx_on_yface(i,j,k,0,vel,dxi,bvylo,bvyhi,bct,dlo,dhi);
2221  Real dvdx = mltensor_dx_on_yface(i,j,k,1,vel,dxi,bvylo,bvyhi,bct,dlo,dhi);
2222  Real dwdx = mltensor_dx_on_yface(i,j,k,2,vel,dxi,bvylo,bvyhi,bct,dlo,dhi);
2223  Real dudy = (vel(i,j,k,0) - vel(i,j-1,k,0))*dyi;
2224  Real dvdy = (vel(i,j,k,1) - vel(i,j-1,k,1))*dyi;
2225  Real dwdy = (vel(i,j,k,2) - vel(i,j-1,k,2))*dyi;
2226  Real dudz = mltensor_dz_on_yface(i,j,k,0,vel,dzi,bvylo,bvyhi,bct,dlo,dhi);
2227  Real dvdz = mltensor_dz_on_yface(i,j,k,1,vel,dzi,bvylo,bvyhi,bct,dlo,dhi);
2228  Real dwdz = mltensor_dz_on_yface(i,j,k,2,vel,dzi,bvylo,bvyhi,bct,dlo,dhi);
2229  fy(i,j,k,0) = dudx;
2230  fy(i,j,k,1) = dvdx;
2231  fy(i,j,k,2) = dwdx;
2232  fy(i,j,k,3) = dudy;
2233  fy(i,j,k,4) = dvdy;
2234  fy(i,j,k,5) = dwdy;
2235  fy(i,j,k,6) = dudz;
2236  fy(i,j,k,7) = dvdz;
2237  fy(i,j,k,8) = dwdz;
2238 
2239  }
2240  }
2241  }
2242 }
2243 
2244 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
2245 void mltensor_vel_grads_fz (Box const& box, Array4<Real> const& fz,
2246  Array4<Real const> const& vel,
2247  GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
2248  Array4<Real const> const& bvzlo,
2249  Array4<Real const> const& bvzhi,
2250  Array2D<BoundCond,
2251  0,2*AMREX_SPACEDIM,
2252  0,AMREX_SPACEDIM> const& bct,
2253  Dim3 const& dlo, Dim3 const& dhi) noexcept
2254 {
2255  const Real dxi = dxinv[0];
2256  const Real dyi = dxinv[1];
2257  const Real dzi = dxinv[2];
2258  const auto lo = amrex::lbound(box);
2259  const auto hi = amrex::ubound(box);
2260 
2261  for (int k = lo.z; k <= hi.z; ++k) {
2262  for (int j = lo.y; j <= hi.y; ++j) {
2263  for (int i = lo.x; i <= hi.x; ++i) {
2264  Real dudx = mltensor_dx_on_zface(i,j,k,0,vel,dxi,bvzlo,bvzhi,bct,dlo,dhi);
2265  Real dvdx = mltensor_dx_on_zface(i,j,k,1,vel,dxi,bvzlo,bvzhi,bct,dlo,dhi);
2266  Real dwdx = mltensor_dx_on_zface(i,j,k,2,vel,dxi,bvzlo,bvzhi,bct,dlo,dhi);
2267  Real dudy = mltensor_dy_on_zface(i,j,k,0,vel,dyi,bvzlo,bvzhi,bct,dlo,dhi);
2268  Real dvdy = mltensor_dy_on_zface(i,j,k,1,vel,dyi,bvzlo,bvzhi,bct,dlo,dhi);
2269  Real dwdy = mltensor_dy_on_zface(i,j,k,2,vel,dyi,bvzlo,bvzhi,bct,dlo,dhi);
2270  Real dudz = (vel(i,j,k,0) - vel(i,j,k-1,0))*dzi;
2271  Real dvdz = (vel(i,j,k,1) - vel(i,j,k-1,1))*dzi;
2272  Real dwdz = (vel(i,j,k,2) - vel(i,j,k-1,2))*dzi;
2273  fz(i,j,k,0) = dudx;
2274  fz(i,j,k,1) = dvdx;
2275  fz(i,j,k,2) = dwdx;
2276  fz(i,j,k,3) = dudy;
2277  fz(i,j,k,4) = dvdy;
2278  fz(i,j,k,5) = dwdy;
2279  fz(i,j,k,6) = dudz;
2280  fz(i,j,k,7) = dvdz;
2281  fz(i,j,k,8) = dwdz;
2282 
2283  }
2284  }
2285  }
2286 }
2287 
2288 }
2289 
2290 #endif
AMREX_PRAGMA_SIMD
#define AMREX_PRAGMA_SIMD
Definition: AMReX_Extension.H:80
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
AMREX_GPU_DEVICE
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
AMREX_GPU_HOST_DEVICE
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
AMREX_LO_NEUMANN
#define AMREX_LO_NEUMANN
Definition: AMReX_LO_BCTYPES.H:6
AMREX_LO_DIRICHLET
#define AMREX_LO_DIRICHLET
Definition: AMReX_LO_BCTYPES.H:5
AMReX_MLLinOp_K.H
amrex::BndryDataT::not_covered
@ not_covered
Definition: AMReX_BndryData.H:44
amrex::BndryDataT::covered
@ covered
Definition: AMReX_BndryData.H:44
amrex::BoundCond
Maintain an identifier for boundary condition types.
Definition: AMReX_BoundCond.H:20
amrex::Orientation::zhi
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE int zhi() noexcept
Int value of the z-hi-face.
Definition: AMReX_Orientation.H:118
amrex::Orientation::yhi
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE int yhi() noexcept
Int value of the y-hi-face.
Definition: AMReX_Orientation.H:110
amrex::Orientation::ylo
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE int ylo() noexcept
Int value of the y-lo-face.
Definition: AMReX_Orientation.H:106
amrex::Orientation::xlo
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE int xlo() noexcept
Int value of the x-lo-face.
Definition: AMReX_Orientation.H:98
amrex::Orientation::low
@ low
Definition: AMReX_Orientation.H:34
amrex::Orientation::high
@ high
Definition: AMReX_Orientation.H:34
amrex::Orientation::xhi
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE int xhi() noexcept
Int value of the x-hi-face.
Definition: AMReX_Orientation.H:102
amrex::Orientation::zlo
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE int zlo() noexcept
Int value of the z-lo-face.
Definition: AMReX_Orientation.H:114
amrex
Definition: AMReX_Amr.cpp:49
amrex::mltensor_cross_terms_fz
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_cross_terms_fz(Box const &box, Array4< Real > const &fz, Array4< Real const > const &vel, Array4< Real const > const &etaz, Array4< Real const > const &kapz, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLTensor_3D_K.H:1606
amrex::mltensor_dx_on_yface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mltensor_dx_on_yface(int i, int j, int k, int n, Array4< Real const > const &vel, Real dxi) noexcept
Definition: AMReX_MLTensor_K.H:17
amrex::Box
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
amrex::mltensor_fill_edges
void mltensor_fill_edges(Box const &vbox, Array4< Real > const &vel, Array4< int const > const &mxlo, Array4< int const > const &mylo, Array4< int const > const &mzlo, Array4< int const > const &mxhi, Array4< int const > const &myhi, Array4< int const > const &mzhi, Array4< Real const > const &bcvalxlo, Array4< Real const > const &bcvalylo, Array4< Real const > const &bcvalzlo, Array4< Real const > const &bcvalxhi, Array4< Real const > const &bcvalyhi, Array4< Real const > const &bcvalzhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Dim3 const &dlo, Dim3 const &dhi) noexcept
Definition: AMReX_MLTensor_3D_K.H:1346
amrex::mltensor_fill_edges_xhi_ylo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xhi_ylo(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxhi, Array4< int const > const &mylo, Array4< Real const > const &bcvalxhi, Array4< Real const > const &bcvalylo, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xhi_domain, bool ylo_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:63
amrex::mltensor_vel_grads_fy
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_vel_grads_fy(Box const &box, Array4< Real > const &fy, Array4< Real const > const &vel, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLTensor_2D_K.H:389
amrex::mltensor_dx_on_zface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mltensor_dx_on_zface(int i, int j, int k, int n, Array4< Real const > const &vel, Real dxi) noexcept
Definition: AMReX_MLTensor_3D_K.H:1528
amrex::mltensor_fill_edges_xhi_yhi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xhi_yhi(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxhi, Array4< int const > const &myhi, Array4< Real const > const &bcvalxhi, Array4< Real const > const &bcvalyhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xhi_domain, bool yhi_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:169
amrex::mltensor_cross_terms_os
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_cross_terms_os(Box const &box, Array4< Real > const &Ax, Array4< Real const > const &fx, Array4< Real const > const &fy, Array4< int const > const &osm, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Real bscalar) noexcept
Definition: AMReX_MLTensor_2D_K.H:334
amrex::mltensor_fill_edges_xlo_yhi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xlo_yhi(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxlo, Array4< int const > const &myhi, Array4< Real const > const &bcvalxlo, Array4< Real const > const &bcvalyhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xlo_domain, bool yhi_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:116
amrex::mltensor_fill_edges_ylo_zhi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_ylo_zhi(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mylo, Array4< int const > const &mzhi, Array4< Real const > const &bcvalylo, Array4< Real const > const &bcvalzhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool ylo_domain, bool zhi_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:540
amrex::mltensor_fill_edges_xlo_ylo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xlo_ylo(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxlo, Array4< int const > const &mylo, Array4< Real const > const &bcvalxlo, Array4< Real const > const &bcvalylo, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xlo_domain, bool ylo_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:10
amrex::mltensor_fill_edges_xhi_zlo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xhi_zlo(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxhi, Array4< int const > const &mzlo, Array4< Real const > const &bcvalxhi, Array4< Real const > const &bcvalzlo, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xhi_domain, bool zlo_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:275
amrex::ubound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:315
amrex::mltensor_fill_edges_yhi_zlo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_yhi_zlo(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &myhi, Array4< int const > const &mzlo, Array4< Real const > const &bcvalyhi, Array4< Real const > const &bcvalzlo, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool yhi_domain, bool zlo_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:487
amrex::mltensor_fill_edges_xlo_zlo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xlo_zlo(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxlo, Array4< int const > const &mzlo, Array4< Real const > const &bcvalxlo, Array4< Real const > const &bcvalzlo, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xlo_domain, bool zlo_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:222
amrex::mltensor_cross_terms_fx
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_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, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLTensor_2D_K.H:187
amrex::lbound
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:308
amrex::mltensor_vel_grads_fx
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_vel_grads_fx(Box const &box, Array4< Real > const &fx, Array4< Real const > const &vel, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLTensor_2D_K.H:363
amrex::mltensor_dz_on_xface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mltensor_dz_on_xface(int i, int j, int k, int n, Array4< Real const > const &vel, Real dzi) noexcept
Definition: AMReX_MLTensor_3D_K.H:1516
amrex::mltensor_fill_edges_ylo_zlo
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_ylo_zlo(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mylo, Array4< int const > const &mzlo, Array4< Real const > const &bcvalylo, Array4< Real const > const &bcvalzlo, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool ylo_domain, bool zlo_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:434
amrex::mltensor_dy_on_zface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mltensor_dy_on_zface(int i, int j, int k, int n, Array4< Real const > const &vel, Real dyi) noexcept
Definition: AMReX_MLTensor_3D_K.H:1534
amrex::mltensor_cross_terms
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_cross_terms(Box const &box, Array4< Real > const &Ax, Array4< Real const > const &fx, Array4< Real const > const &fy, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Real bscalar) noexcept
Definition: AMReX_MLTensor_2D_K.H:311
amrex::mltensor_fill_edges_xlo_zhi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xlo_zhi(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxlo, Array4< int const > const &mzhi, Array4< Real const > const &bcvalxlo, Array4< Real const > const &bcvalzhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xlo_domain, bool zhi_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:328
amrex::mltensor_dz_on_yface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mltensor_dz_on_yface(int i, int j, int k, int n, Array4< Real const > const &vel, Real dzi) noexcept
Definition: AMReX_MLTensor_3D_K.H:1522
amrex::length
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 length(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:322
amrex::mllinop_apply_bc_x
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_apply_bc_x(int side, Box const &box, int blen, Array4< T > const &phi, Array4< int const > const &mask, BoundCond bct, T bcl, Array4< T const > const &bcval, int maxorder, T dxinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:14
amrex::mltensor_fill_corners
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_corners(int icorner, Box const &vbox, Array4< Real > const &vel, Array4< int const > const &mxlo, Array4< int const > const &mylo, Array4< int const > const &mxhi, Array4< int const > const &myhi, Array4< Real const > const &bcvalxlo, Array4< Real const > const &bcvalylo, Array4< Real const > const &bcvalxhi, Array4< Real const > const &bcvalyhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, Dim3 const &dlo, Dim3 const &dhi) noexcept
Definition: AMReX_MLTensor_2D_K.H:10
amrex::mltensor_fill_edges_xhi_zhi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_xhi_zhi(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &mxhi, Array4< int const > const &mzhi, Array4< Real const > const &bcvalxhi, Array4< Real const > const &bcvalzhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool xhi_domain, bool zhi_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:381
amrex::mllinop_apply_bc_z
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_apply_bc_z(int side, Box const &box, int blen, Array4< T > const &phi, Array4< int const > const &mask, BoundCond bct, T bcl, Array4< T const > const &bcval, int maxorder, T dzinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:224
amrex::mltensor_vel_grads_fz
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_vel_grads_fz(Box const &box, Array4< Real > const &fz, Array4< Real const > const &vel, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLTensor_3D_K.H:2114
amrex::mltensor_dy_on_xface
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real mltensor_dy_on_xface(int i, int j, int k, int n, Array4< Real const > const &vel, Real dyi) noexcept
Definition: AMReX_MLTensor_K.H:11
amrex::mllinop_apply_bc_y
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mllinop_apply_bc_y(int side, Box const &box, int blen, Array4< T > const &phi, Array4< int const > const &mask, BoundCond bct, T bcl, Array4< T const > const &bcval, int maxorder, T dyinv, int inhomog, int icomp) noexcept
Definition: AMReX_MLLinOp_K.H:119
amrex::mltensor_fill_edges_yhi_zhi
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_fill_edges_yhi_zhi(int const i, int const j, int const k, Dim3 const &blen, Array4< Real > const &vel, Array4< int const > const &myhi, Array4< int const > const &mzhi, Array4< Real const > const &bcvalyhi, Array4< Real const > const &bcvalzhi, Array2D< BoundCond, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bct, Array2D< Real, 0, 2 *AMREX_SPACEDIM, 0, AMREX_SPACEDIM > const &bcl, int inhomog, int maxorder, GpuArray< Real, AMREX_SPACEDIM > const &dxinv, bool yhi_domain, bool zhi_domain) noexcept
Definition: AMReX_MLTensor_3D_K.H:593
amrex::mltensor_cross_terms_fy
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void mltensor_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, GpuArray< Real, AMREX_SPACEDIM > const &dxinv) noexcept
Definition: AMReX_MLTensor_2D_K.H:215
amrex::Array2D
Definition: AMReX_Array.H:282
amrex::Array4
Definition: AMReX_Array4.H:61
amrex::Dim3
Definition: AMReX_Dim3.H:12