ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_ComputeStrain_EB.cpp File Reference
#include <ERF_Diffusion.H>
#include "ERF_EddyViscosity.H"
Include dependency graph for ERF_ComputeStrain_EB.cpp:

Functions

void ComputeStrain_EB (const MFIter &mfi, Box bxcc, Box tbxxy, Box tbxxz, Box tbxyz, Box domain, const Array4< const Real > &u, const Array4< const Real > &v, const Array4< const Real > &w, Array4< Real > &tau11, Array4< Real > &tau22, Array4< Real > &tau33, Array4< Real > &tau12, Array4< Real > &tau13, Array4< Real > &tau23, const GpuArray< Real, AMREX_SPACEDIM > &dxInv, const BCRec *bc_ptr, const eb_ &ebfact, Array4< Real > &tau13i, Array4< Real > &tau23i)
 

Function Documentation

◆ ComputeStrain_EB()

void ComputeStrain_EB ( const MFIter &  mfi,
Box  bxcc,
Box  tbxxy,
Box  tbxxz,
Box  tbxyz,
Box  domain,
const Array4< const Real > &  u,
const Array4< const Real > &  v,
const Array4< const Real > &  w,
Array4< Real > &  tau11,
Array4< Real > &  tau22,
Array4< Real > &  tau33,
Array4< Real > &  tau12,
Array4< Real > &  tau13,
Array4< Real > &  tau23,
const GpuArray< Real, AMREX_SPACEDIM > &  dxInv,
const BCRec *  bc_ptr,
const eb_ ebfact,
Array4< Real > &  tau13i,
Array4< Real > &  tau23i 
)

Function for computing the strain rates for EB.

Parameters
[in]bxcccell center box for tau_ii
[in]tbxxynodal xy box for tau_12
[in]tbxxznodal xz box for tau_13
[in]tbxyznodal yz box for tau_23
[in]ux-direction velocity
[in]vy-direction velocity
[in]wz-direction velocity
[out]tau1111 strain
[out]tau2222 strain
[out]tau3333 strain
[out]tau1212 strain
[out]tau1313 strain
[out]tau2323 strain
[in]bc_ptrcontainer with boundary condition types
[in]dxInvinverse cell size array
[in]tau13icontribution to strain from du/dz
[in]tau23icontribution to strain from dv/dz
43 {
44  // Convert domain to each index type to test if we are on Dirichlet boundary
45  Box domain_xy = convert(domain, tbxxy.ixType());
46  Box domain_xz = convert(domain, tbxxz.ixType());
47  Box domain_yz = convert(domain, tbxyz.ixType());
48 
49  const auto& dom_lo = lbound(domain);
50  const auto& dom_hi = ubound(domain);
51 
52  // EB
53  // Array4<const EBCellFlag> cflag = (ebfact.get_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
54  Array4<const EBCellFlag> u_cflag = (ebfact.get_u_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
55  Array4<const EBCellFlag> v_cflag = (ebfact.get_v_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
56  Array4<const EBCellFlag> w_cflag = (ebfact.get_w_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
57 
58  // Dirichlet on left or right plane
59  bool xl_v_dir = ( (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir) ||
60  (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ||
61  (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_ingested) );
62  xl_v_dir = ( xl_v_dir && (tbxxy.smallEnd(0) == domain_xy.smallEnd(0)) );
63 
64  bool xh_v_dir = ( (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir) ||
65  (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ||
66  (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_ingested) );
67  xh_v_dir = ( xh_v_dir && (tbxxy.bigEnd(0) == domain_xy.bigEnd(0)) );
68 
69  bool xl_w_dir = ( (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir) ||
70  (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ||
71  (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_ingested) );
72  xl_w_dir = ( xl_w_dir && (tbxxz.smallEnd(0) == domain_xz.smallEnd(0)) );
73 
74  bool xh_w_dir = ( (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir) ||
75  (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ||
76  (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_ingested) );
77  xh_w_dir = ( xh_w_dir && (tbxxz.bigEnd(0) == domain_xz.bigEnd(0)) );
78 
79  // Dirichlet on front or back plane
80  bool yl_u_dir = ( (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir) ||
81  (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ||
82  (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_ingested) );
83  yl_u_dir = ( yl_u_dir && (tbxxy.smallEnd(1) == domain_xy.smallEnd(1)) );
84 
85  bool yh_u_dir = ( (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir) ||
86  (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ||
87  (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_ingested) );
88  yh_u_dir = ( yh_u_dir && (tbxxy.bigEnd(1) == domain_xy.bigEnd(1)) );
89 
90  bool yl_w_dir = ( (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir) ||
91  (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ||
92  (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_ingested) );
93  yl_w_dir = ( yl_w_dir && (tbxyz.smallEnd(1) == domain_yz.smallEnd(1)) );
94 
95  bool yh_w_dir = ( (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir) ||
96  (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ||
97  (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_ingested) );
98  yh_w_dir = ( yh_w_dir && (tbxyz.bigEnd(1) == domain_yz.bigEnd(1)) );
99 
100  // Dirichlet on top or bottom plane
101  bool zl_u_dir = ( (bc_ptr[BCVars::xvel_bc].lo(2) == ERFBCType::ext_dir) ||
102  (bc_ptr[BCVars::xvel_bc].lo(2) == ERFBCType::ext_dir_ingested) );
103  zl_u_dir = ( zl_u_dir && (tbxxz.smallEnd(2) == domain_xz.smallEnd(2)) );
104 
105  bool zh_u_dir = ( (bc_ptr[BCVars::xvel_bc].hi(2) == ERFBCType::ext_dir) ||
106  (bc_ptr[BCVars::xvel_bc].hi(2) == ERFBCType::ext_dir_ingested) );
107  zh_u_dir = ( zh_u_dir && (tbxxz.bigEnd(2) == domain_xz.bigEnd(2)) );
108 
109  bool zl_v_dir = ( (bc_ptr[BCVars::yvel_bc].lo(2) == ERFBCType::ext_dir) ||
110  (bc_ptr[BCVars::yvel_bc].lo(2) == ERFBCType::ext_dir_ingested) );
111  zl_v_dir = ( zl_v_dir && (tbxyz.smallEnd(2) == domain_yz.smallEnd(2)) );
112 
113  bool zh_v_dir = ( (bc_ptr[BCVars::yvel_bc].hi(2) == ERFBCType::ext_dir) ||
114  (bc_ptr[BCVars::yvel_bc].hi(2) == ERFBCType::ext_dir_ingested) );
115  zh_v_dir = ( zh_v_dir && (tbxyz.bigEnd(2) == domain_yz.bigEnd(2)) );
116 
117  //***********************************************************************************
118  // X-Dirichlet
119  //***********************************************************************************
120  if (xl_v_dir) {
121  Box planexy = tbxxy; planexy.setBig(0, planexy.smallEnd(0) );
122  tbxxy.growLo(0,-1);
123  bool need_to_test = (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ? true : false;
124 
125  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
126  if (!need_to_test || u(dom_lo.x,j,k) >= 0.) {
127  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k))*dxInv[1]
128  + (-(8./3.) * v(i-1,j,k) + 3. * v(i,j,k) - (1./3.) * v(i+1,j,k))*dxInv[0] );
129  } else {
130  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k))*dxInv[1] +
131  (v(i, j, k) - v(i-1, j, k))*dxInv[0] );
132  }
133  });
134  }
135  if (xh_v_dir) {
136  // note: tilebox xy should be nodal, so i|i-1|i-2 at the bigEnd is analogous to i-1|i|i+1 at the smallEnd
137  Box planexy = tbxxy; planexy.setSmall(0, planexy.bigEnd(0) );
138  tbxxy.growHi(0,-1);
139  bool need_to_test = (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ? true : false;
140 
141  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
142  if (!need_to_test || u(dom_hi.x+1,j,k) <= 0.) {
143  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k))*dxInv[1]
144  - (-(8./3.) * v(i,j,k) + 3. * v(i-1,j,k) - (1./3.) * v(i-2,j,k))*dxInv[0] );
145  } else {
146  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k))*dxInv[1] +
147  (v(i, j, k) - v(i-1, j, k))*dxInv[0] );
148  }
149  });
150  }
151 
152  if (xl_w_dir) {
153  Box planexz = tbxxz; planexz.setBig(0, planexz.smallEnd(0) );
154  tbxxz.growLo(0,-1);
155  bool need_to_test = (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ? true : false;
156 
157  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
158  Real du_dz = (u(i, j, k) - u(i, j, k-1))*dxInv[2];
159  if (!need_to_test || u(dom_lo.x,j,k) >= 0.) {
160  tau13(i,j,k) = 0.5 * ( du_dz
161  + (-(8./3.) * w(i-1,j,k) + 3. * w(i,j,k) - (1./3.) * w(i+1,j,k))*dxInv[0] );
162  } else {
163  tau13(i,j,k) = 0.5 * ( du_dz
164  + (w(i, j, k) - w(i-1, j, k))*dxInv[0] );
165  }
166 
167  if (tau13i) tau13i(i,j,k) = 0.5 * du_dz;
168  });
169  }
170  if (xh_w_dir) {
171  // note: tilebox xz should be nodal, so i|i-1|i-2 at the bigEnd is analogous to i-1|i|i+1 at the smallEnd
172  Box planexz = tbxxz; planexz.setSmall(0, planexz.bigEnd(0) );
173  tbxxz.growHi(0,-1);
174  bool need_to_test = (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ? true : false;
175 
176  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
177  Real du_dz = (u(i, j, k) - u(i, j, k-1))*dxInv[2];
178  if (!need_to_test || u(dom_hi.x+1,j,k) <= 0.) {
179  tau13(i,j,k) = 0.5 * ( du_dz
180  - (-(8./3.) * w(i,j,k) + 3. * w(i-1,j,k) - (1./3.) * w(i-2,j,k))*dxInv[0] );
181  } else {
182  tau13(i,j,k) = 0.5 * ( du_dz
183  + (w(i, j, k) - w(i-1, j, k))*dxInv[0] );
184  }
185 
186  if (tau13i) tau13i(i,j,k) = 0.5 * du_dz;
187  });
188  }
189 
190  //***********************************************************************************
191  // Y-Dirichlet
192  //***********************************************************************************
193  if (yl_u_dir) {
194  Box planexy = tbxxy; planexy.setBig(1, planexy.smallEnd(1) );
195  tbxxy.growLo(1,-1);
196  bool need_to_test = (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ? true : false;
197 
198  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
199  if (!need_to_test || v(i,dom_lo.y,k) >= 0.) {
200  tau12(i,j,k) = 0.5 * ( (-(8./3.) * u(i,j-1,k) + 3. * u(i,j,k) - (1./3.) * u(i,j+1,k))*dxInv[1]
201  + (v(i, j, k) - v(i-1, j, k))*dxInv[0] );
202  } else {
203  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k))*dxInv[1]
204  + (v(i, j, k) - v(i-1, j, k))*dxInv[0] );
205  }
206  });
207  }
208  if (yh_u_dir) {
209  // note: tilebox xy should be nodal, so j|j-1|j-2 at the bigEnd is analogous to j-1|j|j+1 at the smallEnd
210  Box planexy = tbxxy; planexy.setSmall(1, planexy.bigEnd(1) );
211  tbxxy.growHi(1,-1);
212  bool need_to_test = (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ? true : false;
213 
214  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
215  if (!need_to_test || v(i,dom_hi.y+1,k) <= 0.) {
216  tau12(i,j,k) = 0.5 * ( -(-(8./3.) * u(i,j,k) + 3. * u(i,j-1,k) - (1./3.) * u(i,j-2,k))*dxInv[1]
217  + (v(i, j, k) - v(i-1, j, k))*dxInv[0] );
218  } else {
219  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k))*dxInv[1]
220  + (v(i, j, k) - v(i-1, j, k))*dxInv[0] );
221  }
222  });
223  }
224 
225  if (yl_w_dir) {
226  Box planeyz = tbxyz; planeyz.setBig(1, planeyz.smallEnd(1) );
227  tbxyz.growLo(1,-1);
228  bool need_to_test = (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ? true : false;
229 
230  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
231  Real dv_dz = (v(i, j, k) - v(i, j, k-1))*dxInv[2];
232  if (!need_to_test || v(i,dom_lo.y,k) >= 0.) {
233  tau23(i,j,k) = 0.5 * ( dv_dz
234  + (-(8./3.) * w(i,j-1,k) + 3. * w(i,j ,k) - (1./3.) * w(i,j+1,k))*dxInv[1] );
235  } else {
236  tau23(i,j,k) = 0.5 * ( dv_dz
237  + (w(i, j, k) - w(i, j-1, k))*dxInv[1] );
238  }
239 
240  if (tau23i) tau23i(i,j,k) = 0.5 * dv_dz;
241  });
242  }
243  if (yh_w_dir) {
244  // note: tilebox yz should be nodal, so j|j-1|j-2 at the bigEnd is analogous to j-1|j|j+1 at the smallEnd
245  Box planeyz = tbxyz; planeyz.setSmall(1, planeyz.bigEnd(1) );
246  tbxyz.growHi(1,-1);
247  bool need_to_test = (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ? true : false;
248 
249  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
250  Real dv_dz = (v(i, j, k) - v(i, j, k-1))*dxInv[2];
251  if (!need_to_test || v(i,dom_hi.y+1,k) <= 0.) {
252  tau23(i,j,k) = 0.5 * ( dv_dz
253  - (-(8./3.) * w(i,j ,k) + 3. * w(i,j-1,k) - (1./3.) * w(i,j-2,k))*dxInv[1] );
254  } else {
255  tau23(i,j,k) = 0.5 * ( dv_dz
256  + (w(i, j, k) - w(i, j-1, k))*dxInv[1] );
257  }
258 
259  if (tau23i) tau23i(i,j,k) = 0.5 * dv_dz;
260  });
261  }
262 
263  //***********************************************************************************
264  // Z-Dirichlet
265  //***********************************************************************************
266  if (zl_u_dir) {
267  Box planexz = tbxxz; planexz.setBig(2, planexz.smallEnd(2) );
268  tbxxz.growLo(2,-1);
269 
270  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
271  Real du_dz = (-(8./3.) * u(i,j,k-1) + 3. * u(i,j,k) - (1./3.) * u(i,j,k+1))*dxInv[2];
272  tau13(i,j,k) = 0.5 * ( du_dz
273  + (w(i, j, k) - w(i-1, j, k))*dxInv[0] );
274 
275  if (tau13i) tau13i(i,j,k) = 0.5 * du_dz;
276  });
277  }
278  if (zh_u_dir) {
279  // note: tilebox xz should be nodal, so k|k-1|k-2 at the bigEnd is analogous to k-1|k|k+1 at the smallEnd
280  Box planexz = tbxxz; planexz.setSmall(2, planexz.bigEnd(2) );
281  tbxxz.growHi(2,-1);
282 
283  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
284  Real du_dz = -(-(8./3.) * u(i,j,k) + 3. * u(i,j,k-1) - (1./3.) * u(i,j,k-2))*dxInv[2];
285  tau13(i,j,k) = 0.5 * ( du_dz
286  + (w(i, j, k) - w(i-1, j, k))*dxInv[0] );
287 
288  if (tau13i) tau13i(i,j,k) = 0.5 * du_dz;
289  });
290  }
291 
292  if (zl_v_dir) {
293  Box planeyz = tbxyz; planeyz.setBig(2, planeyz.smallEnd(2) );
294  tbxyz.growLo(2,-1);
295 
296  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
297  Real dv_dz = (-(8./3.) * v(i,j,k-1) + 3. * v(i,j,k ) - (1./3.) * v(i,j,k+1))*dxInv[2];
298  tau23(i,j,k) = 0.5 * ( dv_dz
299  + (w(i, j, k) - w(i, j-1, k))*dxInv[1] );
300 
301  if (tau23i) tau23i(i,j,k) = 0.5 * dv_dz;
302  });
303  }
304  if (zh_v_dir) {
305  // note: tilebox yz should be nodal, so k|k-1|k-2 at the bigEnd is analogous to k-1|k|k+1 at the smallEnd
306  Box planeyz = tbxyz; planeyz.setSmall(2, planeyz.bigEnd(2) );
307  tbxyz.growHi(2,-1);
308 
309  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
310  Real dv_dz = -(-(8./3.) * v(i,j,k ) + 3. * v(i,j,k-1) - (1./3.) * v(i,j,k-2))*dxInv[2];
311  tau23(i,j,k) = 0.5 * ( dv_dz
312  + (w(i, j, k) - w(i, j-1, k))*dxInv[1] );
313 
314  if (tau23i) tau23i(i,j,k) = 0.5 * dv_dz;
315  });
316  }
317 
318  // Fill the remaining cells
319  //***********************************************************************************
320  // Cell centered strains
321  ParallelFor(bxcc, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
322 
323  Real du_dx{0.0};
324  if (u_cflag(i+1,j,k).isCovered() && !u_cflag(i,j,k).isCovered()) {
325  du_dx = ( 2.0*u(i, j, k) - 3.0*u(i-1, j, k) + u(i-2, j, k))*dxInv[0];
326  } else if (u_cflag(i,j,k).isCovered() && !u_cflag(i+1,j,k).isCovered()) {
327  du_dx = (- 2.0*u(i+1, j, k) + 3.0*u(i+2, j, k) - u(i+3, j, k))*dxInv[0];
328  } else if (!u_cflag(i+1,j,k).isCovered() && !u_cflag(i,j,k).isCovered()) {
329  du_dx = (u(i+1, j, k) - u(i, j, k))*dxInv[0];
330  }
331 
332  Real dv_dy{0.0};
333  if (v_cflag(i,j+1,k).isCovered() && !v_cflag(i,j,k).isCovered()) {
334  dv_dy = ( 2.0*v(i, j, k) - 3.0*v(i, j-1, k) + v(i, j-2, k))*dxInv[1];
335  } else if (v_cflag(i,j,k).isCovered() && !v_cflag(i,j+1,k).isCovered()) {
336  dv_dy = (- 2.0*v(i, j+1, k) + 3.0*v(i, j+2, k) - v(i, j+3, k))*dxInv[1];
337  } else if (!v_cflag(i,j+1,k).isCovered() && !v_cflag(i,j,k).isCovered()) {
338  dv_dy = (v(i, j+1, k) - v(i, j, k))*dxInv[1];
339  }
340 
341  Real dw_dz{0.0};
342  if (w_cflag(i,j,k+1).isCovered() && !w_cflag(i,j,k).isCovered()) {
343  dw_dz = ( 2.0*w(i, j, k) - 3.0*w(i, j, k-1) + w(i, j, k-2))*dxInv[2];
344  } else if (w_cflag(i,j,k).isCovered() && !w_cflag(i,j,k+1).isCovered()) {
345  dw_dz = (- 2.0*w(i, j, k+1) + 3.0*w(i, j, k+2) - w(i, j, k+3))*dxInv[2];
346  } else if (!w_cflag(i,j,k+1).isCovered() && !w_cflag(i,j,k).isCovered()) {
347  dw_dz = (w(i, j, k+1) - w(i, j, k))*dxInv[2];
348  }
349 
350  tau11(i,j,k) = du_dx;
351  tau22(i,j,k) = dv_dy;
352  tau33(i,j,k) = dw_dz;
353  });
354 
355  // Off-diagonal strains
356  ParallelFor(tbxxy,tbxxz,tbxyz,
357  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
358 
359  Real du_dy{0.0};
360  if (u_cflag(i,j,k).isCovered() && !u_cflag(i,j-1,k).isCovered()) {
361  du_dy = ( 2.0*u(i, j-1, k) - 3.0*u(i, j-2, k) + u(i, j-3, k))*dxInv[1];
362  } else if (u_cflag(i,j-1,k).isCovered() && !u_cflag(i,j,k).isCovered()) {
363  du_dy = (- 2.0*u(i, j, k) + 3.0*u(i, j+1, k) - u(i, j+2, k))*dxInv[1];
364  } else if (!u_cflag(i,j,k).isCovered() && !u_cflag(i,j-1,k).isCovered()) {
365  du_dy = (u(i, j, k) - u(i, j-1, k))*dxInv[1];
366  }
367 
368  Real dv_dx{0.0};
369  if (v_cflag(i,j,k).isCovered() && !v_cflag(i-1,j,k).isCovered()) {
370  dv_dx = ( 2.0*v(i-1, j, k) - 3.0*v(i-2, j, k) + v(i-3, j, k))*dxInv[0];
371  } else if (v_cflag(i-1,j,k).isCovered() && !v_cflag(i,j,k).isCovered()) {
372  dv_dx = (- 2.0*v(i, j, k) + 3.0*v(i+1, j, k) - v(i+2, j, k))*dxInv[0];
373  } else if (!v_cflag(i,j,k).isCovered() && !v_cflag(i-1,j,k).isCovered()) {
374  dv_dx = (v(i, j, k) - v(i-1, j, k))*dxInv[0];
375  }
376 
377  tau12(i,j,k) = 0.5 * ( du_dy + dv_dx );
378  },
379  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
380 
381  Real du_dz{0.0};
382  if (u_cflag(i,j,k).isCovered() && !u_cflag(i,j,k-1).isCovered()) {
383  du_dz = ( 2.0*u(i, j, k-1) - 3.0*u(i, j, k-2) + u(i, j, k-3))*dxInv[2];
384  } else if (u_cflag(i,j,k-1).isCovered() && !u_cflag(i,j,k).isCovered()) {
385  du_dz = (- 2.0*u(i, j, k) + 3.0*u(i, j, k+1) - u(i, j, k+2))*dxInv[2];
386  } else if (!u_cflag(i,j,k).isCovered() && !u_cflag(i,j,k-1).isCovered()) {
387  du_dz = (u(i, j, k) - u(i, j, k-1))*dxInv[2];
388  }
389 
390  Real dw_dx{0.0};
391  if (w_cflag(i,j,k).isCovered() && !w_cflag(i-1,j,k).isCovered()) {
392  dw_dx = ( 2.0*w(i-1, j, k) - 3.0*w(i-2, j, k) + w(i-3, j, k))*dxInv[0];
393  } else if (w_cflag(i-1,j,k).isCovered() && !w_cflag(i,j,k).isCovered()) {
394  dw_dx = (- 2.0*w(i, j, k) + 3.0*w(i+1, j, k) - w(i+2, j, k))*dxInv[0];
395  } else if (!w_cflag(i,j,k).isCovered() && !w_cflag(i-1,j,k).isCovered()) {
396  dw_dx = (w(i, j, k) - w(i-1, j, k))*dxInv[0];
397  }
398 
399  tau13(i,j,k) = 0.5 * ( du_dz + dw_dx );
400 
401  if (tau13i) tau13i(i,j,k) = 0.5 * du_dz;
402  },
403  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
404 
405  Real dv_dz{0.0};
406  if (v_cflag(i,j,k).isCovered() && !v_cflag(i,j,k-1).isCovered()) {
407  dv_dz = ( 2.0*v(i, j, k-1) - 3.0*v(i, j, k-2) + v(i, j, k-3))*dxInv[2];
408  } else if (v_cflag(i,j,k-1).isCovered() && !v_cflag(i,j,k).isCovered()) {
409  dv_dz = (- 2.0*v(i, j, k) + 3.0*v(i, j, k+1) - v(i, j, k+2))*dxInv[2];
410  } else if (!v_cflag(i,j,k).isCovered() && !v_cflag(i,j,k-1).isCovered()) {
411  dv_dz = (v(i, j, k) - v(i, j, k-1))*dxInv[2];
412  }
413 
414  Real dw_dy{0.0};
415  if (w_cflag(i,j,k).isCovered() && !w_cflag(i,j-1,k).isCovered()) {
416  dw_dy = ( 2.0*w(i, j-1, k) - 3.0*w(i, j-2, k) + w(i, j-3, k))*dxInv[1];
417  } else if (w_cflag(i,j-1,k).isCovered() && !w_cflag(i,j,k).isCovered()) {
418  dw_dy = (- 2.0*w(i, j, k) + 3.0*w(i, j+1, k) - w(i, j+2, k))*dxInv[1];
419  } else if (!w_cflag(i,j,k).isCovered() && !w_cflag(i,j-1,k).isCovered()) {
420  dw_dy = (w(i, j, k) - w(i, j-1, k))*dxInv[1];
421  }
422 
423  tau23(i,j,k) = 0.5 * ( dv_dz + dw_dy );
424 
425  if (tau23i) tau23i(i,j,k) = 0.5 * dv_dz;
426  });
427 }
@ tau12
Definition: ERF_DataStruct.H:32
@ tau23
Definition: ERF_DataStruct.H:32
@ tau33
Definition: ERF_DataStruct.H:32
@ tau22
Definition: ERF_DataStruct.H:32
@ tau11
Definition: ERF_DataStruct.H:32
@ tau13
Definition: ERF_DataStruct.H:32
const auto & dom_hi
Definition: ERF_SetupVertDiff.H:2
const auto & dom_lo
Definition: ERF_SetupVertDiff.H:1
amrex::Real Real
Definition: ERF_ShocInterface.H:19
eb_aux_ const * get_w_const_factory() const noexcept
Definition: ERF_EB.H:52
eb_aux_ const * get_v_const_factory() const noexcept
Definition: ERF_EB.H:51
eb_aux_ const * get_u_const_factory() const noexcept
Definition: ERF_EB.H:50
@ zvel_bc
Definition: ERF_IndexDefines.H:89
@ yvel_bc
Definition: ERF_IndexDefines.H:88
@ xvel_bc
Definition: ERF_IndexDefines.H:87
@ ext_dir_ingested
Definition: ERF_IndexDefines.H:212
@ ext_dir
Definition: ERF_IndexDefines.H:209
@ ext_dir_upwind
Definition: ERF_IndexDefines.H:217

Referenced by erf_make_tau_terms().

Here is the call graph for this function:
Here is the caller graph for this function: