ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Functions
ERF_ComputeStrain_T.cpp File Reference
#include <ERF_Diffusion.H>
#include <ERF_TerrainMetrics.H>
Include dependency graph for ERF_ComputeStrain_T.cpp:

Functions

void ComputeStrain_T (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 > &tau21, Array4< Real > &tau13, Array4< Real > &tau31, Array4< Real > &tau23, Array4< Real > &tau32, const Array4< const Real > &z_nd, const Array4< const Real > &detJ, const GpuArray< Real, AMREX_SPACEDIM > &dxInv, const Array4< const Real > &mf_mx, const Array4< const Real > &mf_ux, const Array4< const Real > &mf_vx, const Array4< const Real > &mf_my, const Array4< const Real > &mf_uy, const Array4< const Real > &mf_vy, const BCRec *bc_ptr)
 

Function Documentation

◆ ComputeStrain_T()

void ComputeStrain_T ( 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 > &  tau21,
Array4< Real > &  tau13,
Array4< Real > &  tau31,
Array4< Real > &  tau23,
Array4< Real > &  tau32,
const Array4< const Real > &  z_nd,
const Array4< const Real > &  detJ,
const GpuArray< Real, AMREX_SPACEDIM > &  dxInv,
const Array4< const Real > &  mf_mx,
const Array4< const Real > &  mf_ux,
const Array4< const Real > &  mf_vx,
const Array4< const Real > &  mf_my,
const Array4< const Real > &  mf_uy,
const Array4< const Real > &  mf_vy,
const BCRec *  bc_ptr 
)

Function for computing the strain rates with terrain.

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]tau2121 strain
[out]tau2323 strain
[out]tau3131 strain
[out]tau3232 strain
[in]z_ndnodal array of physical z heights
[in]bc_ptrcontainer with boundary condition types
[in]dxInvinverse cell size array
[in]mf_mxmap factor at cell center
[in]mf_uxmap factor at x-face
[in]mf_vxmap factor at y-face
[in]mf_mymap factor at cell center
[in]mf_uymap factor at x-face
[in]mf_vymap factor at y-face
56 {
57  // Convert domain to each index type to test if we are on dirichlet boundary
58  Box domain_xy = convert(domain, tbxxy.ixType());
59  Box domain_xz = convert(domain, tbxxz.ixType());
60  Box domain_yz = convert(domain, tbxyz.ixType());
61 
62  const auto& dom_lo = lbound(domain);
63  const auto& dom_hi = ubound(domain);
64 
65  // Dirichlet on left or right plane
66  bool xl_v_dir = ( (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir) ||
67  (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ||
68  (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_ingested) );
69  xl_v_dir = ( xl_v_dir && (tbxxy.smallEnd(0) == domain_xy.smallEnd(0)) );
70 
71  bool xh_v_dir = ( (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir) ||
72  (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ||
73  (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_ingested) );
74  xh_v_dir = ( xh_v_dir && (tbxxy.bigEnd(0) == domain_xy.bigEnd(0)) );
75 
76  bool xl_w_dir = ( (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir) ||
77  (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ||
78  (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_ingested) );
79  xl_w_dir = ( xl_w_dir && (tbxxz.smallEnd(0) == domain_xz.smallEnd(0)) );
80 
81  bool xh_w_dir = ( (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir) ||
82  (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ||
83  (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_ingested) );
84  xh_w_dir = ( xh_w_dir && (tbxxz.bigEnd(0) == domain_xz.bigEnd(0)) );
85 
86  // Dirichlet on front or back plane
87  bool yl_u_dir = ( (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir) ||
88  (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ||
89  (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_ingested) );
90  yl_u_dir = ( yl_u_dir && (tbxxy.smallEnd(1) == domain_xy.smallEnd(1)) );
91 
92  bool yh_u_dir = ( (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir) ||
93  (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ||
94  (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_ingested) );
95  yh_u_dir = ( yh_u_dir && (tbxxy.bigEnd(1) == domain_xy.bigEnd(1)) );
96 
97  bool yl_w_dir = ( (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir) ||
98  (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ||
99  (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_ingested) );
100  yl_w_dir = ( yl_w_dir && (tbxyz.smallEnd(1) == domain_yz.smallEnd(1)) );
101 
102  bool yh_w_dir = ( (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir) ||
103  (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ||
104  (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_ingested) );
105  yh_w_dir = ( yh_w_dir && (tbxyz.bigEnd(1) == domain_yz.bigEnd(1)) );
106 
107  // Dirichlet on top or bottom plane
108  bool zl_u_dir = ( (bc_ptr[BCVars::xvel_bc].lo(2) == ERFBCType::ext_dir) ||
109  (bc_ptr[BCVars::xvel_bc].lo(2) == ERFBCType::ext_dir_ingested) );
110  zl_u_dir = ( zl_u_dir && (tbxxz.smallEnd(2) == domain_xz.smallEnd(2)) );
111 
112  bool zh_u_dir = ( (bc_ptr[BCVars::xvel_bc].hi(2) == ERFBCType::ext_dir) ||
113  (bc_ptr[BCVars::xvel_bc].hi(2) == ERFBCType::ext_dir_ingested) );
114  zh_u_dir = ( zh_u_dir && (tbxxz.bigEnd(2) == domain_xz.bigEnd(2)) );
115 
116  bool zl_v_dir = ( (bc_ptr[BCVars::yvel_bc].lo(2) == ERFBCType::ext_dir) ||
117  (bc_ptr[BCVars::yvel_bc].lo(2) == ERFBCType::ext_dir_ingested) );
118  zl_v_dir = ( zl_v_dir && (tbxyz.smallEnd(2) == domain_yz.smallEnd(2)) );
119 
120  bool zh_v_dir = ( (bc_ptr[BCVars::yvel_bc].hi(2) == ERFBCType::ext_dir) ||
121  (bc_ptr[BCVars::yvel_bc].hi(2) == ERFBCType::ext_dir_ingested) );
122  zh_v_dir = ( zh_v_dir && (tbxyz.bigEnd(2) == domain_yz.bigEnd(2)) );
123 
124  //***********************************************************************************
125  // X-Dirichlet
126  //***********************************************************************************
127  if (xl_v_dir)
128  {
129  Box planexy = tbxxy; planexy.setBig(0, planexy.smallEnd(0) );
130  tbxxy.growLo(0,-1);
131  bool need_to_test = (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ? true : false;
132 
133  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
134  Real dz0 = ( z_nd(i,j,k+1) - z_nd(i,j,k-1) );
135  Real idz0 = 1.0 / dz0;
136 
137  Real GradUz = 0.5 * idz0 * ( u(i ,j ,k+1) + u(i ,j-1,k+1)
138  - u(i ,j ,k-1) - u(i ,j-1,k-1) );
139  Real GradVz = 0.5 * idz0 * ( v(i ,j ,k+1) + v(i-1,j ,k+1)
140  - v(i ,j ,k-1) - v(i-1,j ,k-1) );
141 
142  Real mfy = 0.5 * (mf_uy(i,j,0) + mf_uy(i ,j-1,0));
143  Real mfx = 0.5 * (mf_vx(i,j,0) + mf_vx(i-1,j ,0));
144 
145  Real met_h_xi = Compute_h_xi_AtEdgeCenterK (i,j,k,dxInv,z_nd);
146  Real met_h_eta = Compute_h_eta_AtEdgeCenterK (i,j,k,dxInv,z_nd);
147 
148  if (!need_to_test || u(dom_lo.x,j,k) >= 0.) {
149  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k) )*dxInv[1]*mfy
150  + (-(8./3.) * v(i-1,j,k) + 3. * v(i,j,k) - (1./3.) * v(i+1,j,k))*dxInv[0]*mfx
151  - (met_h_eta)*GradUz*mfy
152  - (met_h_xi )*GradVz*mfx );
153  } else {
154  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j-1, k))*dxInv[1]*mfy
155  + (v(i, j, k) - v(i-1, j , k))*dxInv[0]*mfx
156  - (met_h_eta)*GradUz*mfy
157  - (met_h_xi )*GradVz*mfx );
158  }
159  tau21(i,j,k) = tau12(i,j,k);
160  });
161  }
162  if (xh_v_dir)
163  {
164  Box planexy = tbxxy; planexy.setSmall(0, planexy.bigEnd(0) );
165  tbxxy.growHi(0,-1);
166  bool need_to_test = (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ? true : false;
167 
168  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
169  Real dz0 = ( z_nd(i,j,k+1) - z_nd(i,j,k-1) );
170  Real idz0 = 1.0 / dz0;
171 
172  Real GradUz = 0.5 * idz0 * ( u(i ,j ,k+1) + u(i ,j-1,k+1)
173  - u(i ,j ,k-1) - u(i ,j-1,k-1) );
174  Real GradVz = 0.5 * idz0 * ( v(i ,j ,k+1) + v(i-1,j ,k+1)
175  - v(i ,j ,k-1) - v(i-1,j ,k-1) );
176 
177  Real mfy = 0.5 * (mf_uy(i,j,0) + mf_uy(i ,j-1,0));
178  Real mfx = 0.5 * (mf_vx(i,j,0) + mf_vx(i-1,j ,0));
179 
180  Real met_h_xi = Compute_h_xi_AtEdgeCenterK (i,j,k,dxInv,z_nd);
181  Real met_h_eta = Compute_h_eta_AtEdgeCenterK (i,j,k,dxInv,z_nd);
182 
183  if (!need_to_test || u(dom_hi.x+1,j,k) <= 0.) {
184  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j-1, k) )*dxInv[1]*mfy
185  - (-(8./3.) * v(i,j,k) + 3. * v(i-1,j,k) - (1./3.) * v(i-2,j,k))*dxInv[0]*mfx
186  - (met_h_eta)*GradUz*mfy
187  - (met_h_xi )*GradVz*mfx );
188  } else {
189  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j-1, k))*dxInv[1]*mfy
190  + (v(i, j, k) - v(i-1, j , k))*dxInv[0]*mfx
191  - (met_h_eta)*GradUz*mfy
192  - (met_h_xi )*GradVz*mfx );
193  }
194  tau21(i,j,k) = tau12(i,j,k);
195  });
196  }
197 
198  if (xl_w_dir)
199  {
200  Box planexz = tbxxz; planexz.setBig(0, planexz.smallEnd(0) );
201  planexz.setSmall(2, planexz.smallEnd(2)+1 ); planexz.setBig(2, planexz.bigEnd(2)-1 );
202  tbxxz.growLo(0,-1);
203  bool need_to_test = (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_upwind) ? true : false;
204 
205  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
206  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i,j+1,k+1)
207  - z_nd(i,j,k-1) - z_nd(i,j+1,k-1) );
208  Real idz0 = 1.0 / dz0;
209 
210  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i-1,j ,k+1)
211  - w(i ,j ,k-1) - w(i-1,j ,k-1) );
212  Real mfx = mf_ux(i,j,0);
213 
214  Real met_h_xi = Compute_h_xi_AtEdgeCenterJ (i,j,k,dxInv,z_nd);
215  Real met_h_zeta = Compute_h_zeta_AtEdgeCenterJ(i,j,k,dxInv,z_nd);
216 
217  if (!need_to_test || u(dom_lo.x,j,k) <= 0.) {
218  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j, k-1))*dxInv[2]/met_h_zeta
219  + ( (-(8./3.) * w(i-1,j,k) + 3. * w(i,j,k) - (1./3.) * w(i+1,j,k))*dxInv[0]
220  - (met_h_xi)*GradWz ) * mfx );
221  } else {
222  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j, k-1))*dxInv[2]/met_h_zeta
223  + ( (w(i, j, k) - w(i-1, j, k ))*dxInv[0]
224  - (met_h_xi)*GradWz ) * mfx );
225  }
226  tau31(i,j,k) = tau13(i,j,k);
227  });
228  }
229 
230  if (xh_w_dir)
231  {
232  Box planexz = tbxxz; planexz.setSmall(0, planexz.bigEnd(0) );
233  planexz.setSmall(2, planexz.smallEnd(2)+1 ); planexz.setBig(2, planexz.bigEnd(2)-1 );
234  tbxxz.growHi(0,-1);
235  bool need_to_test = (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_upwind) ? true : false;
236 
237  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
238  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i,j+1,k+1)
239  - z_nd(i,j,k-1) - z_nd(i,j+1,k-1) );
240  Real idz0 = 1.0 / dz0;
241 
242  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i-1,j ,k+1)
243  - w(i ,j ,k-1) - w(i-1,j ,k-1) );
244 
245  Real mfx = mf_ux(i,j,0);
246 
247  Real met_h_xi = Compute_h_xi_AtEdgeCenterJ (i,j,k,dxInv,z_nd);
248  Real met_h_zeta = Compute_h_zeta_AtEdgeCenterJ(i,j,k,dxInv,z_nd);
249 
250  if (!need_to_test || u(dom_hi.x+1,j,k) <= 0.) {
251  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j, k-1))*dxInv[2]/met_h_zeta
252  - ( (-(8./3.) * w(i,j,k) + 3. * w(i-1,j,k) - (1./3.) * w(i-2,j,k))*dxInv[0]
253  - (met_h_xi)*GradWz ) * mfx );
254  } else {
255  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j, k-1))*dxInv[2]/met_h_zeta
256  + ( (w(i, j, k) - w(i-1, j, k ))*dxInv[0]
257  - (met_h_xi)*GradWz ) * mfx );
258  }
259  tau31(i,j,k) = tau13(i,j,k);
260  });
261  }
262 
263  //***********************************************************************************
264  // Y-Dirichlet
265  //***********************************************************************************
266  if (yl_u_dir) {
267  Box planexy = tbxxy; planexy.setBig(1, planexy.smallEnd(1) );
268  tbxxy.growLo(1,-1);
269  bool need_to_test = (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ? true : false;
270 
271  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
272  Real dz0 = ( z_nd(i,j,k+1) - z_nd(i,j,k-1) );
273  Real idz0 = 1.0 / dz0;
274 
275  Real GradUz = 0.5 * idz0 * ( u(i ,j ,k+1) + u(i ,j-1,k+1)
276  - u(i ,j ,k-1) - u(i ,j-1,k-1) );
277  Real GradVz = 0.5 * idz0 * ( v(i ,j ,k+1) + v(i-1,j ,k+1)
278  - v(i ,j ,k-1) - v(i-1,j ,k-1) );
279 
280  Real mfy = 0.5 * (mf_uy(i,j,0) + mf_uy(i ,j-1,0));
281  Real mfx = 0.5 * (mf_vx(i,j,0) + mf_vx(i-1,j ,0));
282 
283  Real met_h_xi = Compute_h_xi_AtEdgeCenterK (i,j,k,dxInv,z_nd);
284  Real met_h_eta = Compute_h_eta_AtEdgeCenterK (i,j,k,dxInv,z_nd);
285 
286  if (!need_to_test || v(i,dom_lo.y,k) >= 0.) {
287  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]*mfy
288  + (v(i, j, k) - v(i-1, j, k))*dxInv[0]*mfx
289  - (met_h_eta)*GradUz*mfy
290  - (met_h_xi )*GradVz*mfx );
291  } else {
292  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j-1, k))*dxInv[1]*mfy
293  + (v(i, j, k) - v(i-1, j , k))*dxInv[0]*mfx
294  - (met_h_eta)*GradUz*mfy
295  - (met_h_xi )*GradVz*mfx );
296  }
297  tau21(i,j,k) = tau12(i,j,k);
298  });
299  }
300  if (yh_u_dir) {
301  Box planexy = tbxxy; planexy.setSmall(1, planexy.bigEnd(1) );
302  tbxxy.growHi(1,-1);
303  bool need_to_test = (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ? true : false;
304 
305  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
306  Real dz0 = ( z_nd(i,j,k+1) - z_nd(i,j,k-1) );
307  Real idz0 = 1.0 / dz0;
308 
309  Real GradUz = 0.5 * idz0 * ( u(i ,j ,k+1) + u(i ,j-1,k+1)
310  - u(i ,j ,k-1) - u(i ,j-1,k-1) );
311  Real GradVz = 0.5 * idz0 * ( v(i ,j ,k+1) + v(i-1,j ,k+1)
312  - v(i ,j ,k-1) - v(i-1,j ,k-1) );
313 
314  Real mfy = 0.5 * (mf_uy(i,j,0) + mf_uy(i ,j-1,0));
315  Real mfx = 0.5 * (mf_vx(i,j,0) + mf_vx(i-1,j ,0));
316 
317  Real met_h_xi = Compute_h_xi_AtEdgeCenterK (i,j,k,dxInv,z_nd);
318  Real met_h_eta = Compute_h_eta_AtEdgeCenterK (i,j,k,dxInv,z_nd);
319 
320  if (!need_to_test || v(i,dom_hi.y+1,k) >= 0.) {
321  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]*mfy +
322  + (v(i, j, k) - v(i-1, j, k))*dxInv[0]*mfx
323  - (met_h_eta)*GradUz*mfy
324  - (met_h_xi )*GradVz*mfx );
325  } else {
326  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j-1, k))*dxInv[1]*mfy
327  + (v(i, j, k) - v(i-1, j , k))*dxInv[0]*mfx
328  - (met_h_eta)*GradUz*mfy
329  - (met_h_xi )*GradVz*mfx );
330  }
331  tau21(i,j,k) = tau12(i,j,k);
332  });
333  }
334 
335  if (yl_w_dir) {
336  Box planeyz = tbxyz; planeyz.setBig(1, planeyz.smallEnd(1) );
337  planeyz.setSmall(2, planeyz.smallEnd(2)+1 ); planeyz.setBig(2, planeyz.bigEnd(2)-1 );
338  tbxyz.growLo(1,-1);
339  bool need_to_test = (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_upwind) ? true : false;
340 
341  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
342  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i+1,j,k+1)
343  - z_nd(i,j,k-1) - z_nd(i+1,j,k-1) );
344  Real idz0 = 1.0 / dz0;
345 
346  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i ,j-1,k+1)
347  - w(i ,j ,k-1) - w(i ,j-1,k-1) );
348 
349  Real mfy = mf_vy(i,j,0);
350 
351  Real met_h_eta = Compute_h_eta_AtEdgeCenterI (i,j,k,dxInv,z_nd);
352  Real met_h_zeta = Compute_h_zeta_AtEdgeCenterI(i,j,k,dxInv,z_nd);
353 
354  if (!need_to_test || v(i,dom_lo.y,k) >= 0.) {
355  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j, k-1))*dxInv[2]/met_h_zeta
356  + ( (-(8./3.) * w(i,j-1,k) + 3. * w(i,j ,k) - (1./3.) * w(i,j+1,k))*dxInv[1]
357  - (met_h_eta)*GradWz ) * mfy );
358  } else {
359  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j , k-1))*dxInv[2]/met_h_zeta
360  + ( (w(i, j, k) - w(i, j-1, k ))*dxInv[1]
361  - (met_h_eta)*GradWz ) * mfy );
362  }
363  tau32(i,j,k) = tau23(i,j,k);
364  });
365  }
366  if (yh_w_dir) {
367  Box planeyz = tbxyz; planeyz.setSmall(1, planeyz.bigEnd(1) );
368  planeyz.setSmall(2, planeyz.smallEnd(2)+1 ); planeyz.setBig(2, planeyz.bigEnd(2)-1 );
369  tbxyz.growHi(1,-1);
370  bool need_to_test = (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_upwind) ? true : false;
371 
372  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
373  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i+1,j,k+1)
374  - z_nd(i,j,k-1) - z_nd(i+1,j,k-1) );
375  Real idz0 = 1.0 / dz0;
376 
377  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i ,j-1,k+1)
378  - w(i ,j ,k-1) - w(i ,j-1,k-1) );
379 
380  Real mfy = mf_vy(i,j,0);
381 
382  Real met_h_eta = Compute_h_eta_AtEdgeCenterI (i,j,k,dxInv,z_nd);
383  Real met_h_zeta = Compute_h_zeta_AtEdgeCenterI(i,j,k,dxInv,z_nd);
384 
385  if (!need_to_test || v(i,dom_hi.y+1,k) >= 0.) {
386  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j, k-1))*dxInv[2]/met_h_zeta
387  - ( (-(8./3.) * w(i,j ,k) + 3. * w(i,j-1,k) - (1./3.) * w(i,j-2,k))*dxInv[1]
388  - (met_h_eta)*GradWz ) * mfy );
389  } else {
390  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j , k-1))*dxInv[2]/met_h_zeta
391  + ( (w(i, j, k) - w(i, j-1, k ))*dxInv[1]
392  - (met_h_eta)*GradWz ) * mfy );
393  }
394  tau32(i,j,k) = tau23(i,j,k);
395  });
396  }
397 
398  //***********************************************************************************
399  // Z-Dirichlet
400  //***********************************************************************************
401  if (zl_u_dir) {
402  Box planexz = tbxxz; planexz.setBig(2, planexz.smallEnd(2) );
403  tbxxz.growLo(2,-1);
404  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
405  // Third order stencil with variable dz
406  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i,j+1,k+1)
407  - z_nd(i,j,k ) - z_nd(i,j+1,k ) );
408  Real dz1 = 0.5 * ( z_nd(i,j,k+2) + z_nd(i,j+1,k+2)
409  - z_nd(i,j,k+1) - z_nd(i,j+1,k+1) );
410  Real idz0 = 1.0 / dz0;
411  Real f = (dz1 / dz0) + 2.0;
412  Real f2 = f*f;
413  Real c3 = 2.0 / (f - f2);
414  Real c2 = -f2*c3;
415  Real c1 = -(1.0-f2)*c3;
416 
417  Real GradWz = 0.5 * idz0 * ( w(i,j,k+1) + w(i-1,j,k+1)
418  - w(i,j,k ) - w(i-1,j,k ) );
419 
420  Real mfx = mf_ux(i,j,0);
421 
422  Real met_h_xi;
423  met_h_xi = Compute_h_xi_AtEdgeCenterJ (i,j,k,dxInv,z_nd);
424 
425  // Note: u(i,j,k) and u(i,j,k+1) are located at cell centers
426  // whereas u(i,j,k-1) is the value on the boundary
427  tau13(i,j,k) = 0.5 * ( (c1 * u(i,j,k-1) + c2 * u(i,j,k) + c3 * u(i,j,k+1))*idz0
428  + ( (w(i, j, k) - w(i-1, j, k))*dxInv[0]
429  - (met_h_xi)*GradWz ) * mfx );
430  tau31(i,j,k) = tau13(i,j,k);
431  });
432  }
433  if (zh_u_dir) {
434  // NOTE: h_xi = 0
435  Box planexz = tbxxz; planexz.setSmall(2, planexz.bigEnd(2) );
436  tbxxz.growHi(2,-1);
437  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
438  // Third order stencil with variable dz
439  Real dz0 = 0.5 * ( z_nd(i,j,k ) + z_nd(i,j+1,k )
440  - z_nd(i,j,k-1) - z_nd(i,j+1,k-1) );
441  Real dz1 = 0.5 * ( z_nd(i,j,k-1) + z_nd(i,j+1,k-1)
442  - z_nd(i,j,k-2) - z_nd(i,j+1,k-2) );
443  Real idz0 = 1.0 / dz0;
444  Real f = (dz1 / dz0) + 2.0;
445  Real f2 = f*f;
446  Real c3 = 2.0 / (f - f2);
447  Real c2 = -f2*c3;
448  Real c1 = -(1.0-f2)*c3;
449 
450  Real mfx = mf_ux(i,j,0);
451 
452  tau13(i,j,k) = 0.5 * ( -(c1 * u(i,j,k) + c2 * u(i,j,k-1) + c3 * u(i,j,k-2))*idz0
453  + (w(i, j, k) - w(i-1, j, k))*dxInv[0]*mfx );
454  tau31(i,j,k) = tau13(i,j,k);
455  });
456  }
457 
458  if (zl_v_dir) {
459  Box planeyz = tbxyz; planeyz.setBig(2, planeyz.smallEnd(2) );
460  tbxyz.growLo(2,-1);
461  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
462  // Third order stencil with variable dz
463  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i+1,j,k+1)
464  - z_nd(i,j,k ) - z_nd(i+1,j,k ) );
465  Real dz1 = 0.5 * ( z_nd(i,j,k+2) + z_nd(i+1,j,k+2)
466  - z_nd(i,j,k+1) - z_nd(i+1,j,k+1) );
467  Real idz0 = 1.0 / dz0;
468  Real f = (dz1 / dz0) + 2.0;
469  Real f2 = f*f;
470  Real c3 = 2.0 / (f - f2);
471  Real c2 = -f2*c3;
472  Real c1 = -(1.0-f2)*c3;
473 
474  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i ,j-1,k+1)
475  - w(i ,j ,k ) - w(i ,j-1,k ) );
476 
477  Real mfy = mf_vy(i,j,0);
478 
479  Real met_h_eta;
480  met_h_eta = Compute_h_eta_AtEdgeCenterI (i,j,k,dxInv,z_nd);
481 
482  // Note: v(i,j,k) and v(i,j,k+1) are located at cell centers
483  // whereas v(i,j,k-1) is the value on the boundary
484  tau23(i,j,k) = 0.5 * ( (c1 * v(i,j,k-1) + c2 * v(i,j,k ) + c3 * v(i,j,k+1))*idz0
485  + ( (w(i, j, k) - w(i, j-1, k))*dxInv[1]
486  - (met_h_eta)*GradWz ) * mfy );
487  tau32(i,j,k) = tau23(i,j,k);
488  });
489  }
490  if (zh_v_dir && (tbxyz.bigEnd(2) == domain_yz.bigEnd(2))) {
491  // NOTE: h_eta = 0
492  Box planeyz = tbxyz; planeyz.setSmall(2, planeyz.bigEnd(2) );
493  tbxyz.growHi(2,-1);
494  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
495  // Third order stencil with variable dz
496  Real dz0 = 0.5 * ( z_nd(i,j,k ) + z_nd(i+1,j,k )
497  - z_nd(i,j,k-1) - z_nd(i+1,j,k-1) );
498  Real dz1 = 0.5 * ( z_nd(i,j,k-1) + z_nd(i+1,j,k-1)
499  - z_nd(i,j,k-2) - z_nd(i+1,j,k-2) );
500  Real idz0 = 1.0 / dz0;
501  Real f = (dz1 / dz0) + 2.0;
502  Real f2 = f*f;
503  Real c3 = 2.0 / (f - f2);
504  Real c2 = -f2*c3;
505  Real c1 = -(1.0-f2)*c3;
506 
507  Real mfy = mf_vy(i,j,0);
508 
509  tau23(i,j,k) = 0.5 * ( -(c1 * v(i,j,k ) + c2 * v(i,j,k-1) + c3 * v(i,j,k-2))*idz0
510  + (w(i, j, k) - w(i, j-1, k))*dxInv[1]*mfy );
511  tau32(i,j,k) = tau23(i,j,k);
512  });
513  }
514 
515  // HO derivatives w/ Dirichlet BC (\partival <var> / \partial z from terrain transform)
516  if (zl_u_dir && zl_v_dir) {
517  Box planecc = bxcc; planecc.setBig(2, planecc.smallEnd(2) );
518  bxcc.growLo(2,-1);
519  ParallelFor(planecc, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
520  // Third order stencil with variable dz
521  Real dz0 = 0.25 * ( z_nd(i,j,k+1) + z_nd(i,j+1,k+1) + z_nd(i+1,j,k+1) + z_nd(i+1,j+1,k+1)
522  - z_nd(i,j,k ) - z_nd(i,j+1,k ) - z_nd(i+1,j,k ) - z_nd(i+1,j+1,k ) );
523  Real dz1 = 0.25 * ( z_nd(i,j,k+2) + z_nd(i,j+1,k+2) + z_nd(i+1,j,k+2) + z_nd(i+1,j+1,k+2)
524  - z_nd(i,j,k+1) - z_nd(i,j+1,k+1) - z_nd(i+1,j,k+1) - z_nd(i+1,j+1,k+1) );
525  Real idz0 = 1.0 / dz0;
526  Real f = (dz1 / dz0) + 2.0;
527  Real f2 = f*f;
528  Real c3 = 2.0 / (f - f2);
529  Real c2 = -f2*c3;
530  Real c1 = -(1.0-f2)*c3;
531 
532  Real GradUz = 0.5 * idz0 * ( (c1 * u(i ,j,k-1) + c2 * u(i ,j,k) + c3 * u(i ,j,k+1))
533  + (c1 * u(i-1,j,k-1) + c2 * u(i-1,j,k) + c3 * u(i-1,j,k+1)) );
534  Real GradVz = 0.5 * idz0 * ( (c1 * v(i,j ,k-1) + c2 * v(i,j ,k) + c3 * v(i,j ,k+1))
535  + (c1 * v(i,j-1,k-1) + c2 * v(i,j-1,k) + c3 * v(i,j-1,k+1)) );
536 
537  Real mfx = mf_mx(i,j,0);
538  Real mfy = mf_my(i,j,0);
539 
540  Real met_h_xi,met_h_eta;
541  met_h_xi = Compute_h_xi_AtCellCenter (i,j,k,dxInv,z_nd);
542  met_h_eta = Compute_h_eta_AtCellCenter (i,j,k,dxInv,z_nd);
543 
544  tau11(i,j,k) = ( (u(i+1, j, k) - u(i, j, k) )*dxInv[0]
545  - (met_h_xi)*GradUz ) * mfx;
546  tau22(i,j,k) = ( (v(i, j+1, k) - v(i, j, k) )*dxInv[1]
547  - (met_h_eta)*GradVz ) * mfy;
548  tau33(i,j,k) = ( w(i, j, k+1) - w(i, j, k) )*idz0;
549  });
550 
551  Box planexy = tbxxy; planexy.setBig(2, planexy.smallEnd(2) );
552  tbxxy.growLo(2,-1);
553  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
554  // Third order stencil with variable dz
555  Real dz0 = ( z_nd(i,j,k+1) - z_nd(i,j,k ) );
556  Real dz1 = ( z_nd(i,j,k+2) - z_nd(i,j,k+1) );
557  Real idz0 = 1.0 / dz0;
558  Real f = (dz1 / dz0) + 2.0;
559  Real f2 = f*f;
560  Real c3 = 2.0 / (f - f2);
561  Real c2 = -f2*c3;
562  Real c1 = -(1.0-f2)*c3;
563 
564  Real GradUz = 0.5 * idz0 * ( (c1 * u(i,j ,k-1) + c2 * u(i,j ,k) + c3 * u(i,j ,k+1))
565  + (c1 * u(i,j-1,k-1) + c2 * u(i,j-1,k) + c3 * u(i,j-1,k+1)) );
566  Real GradVz = 0.5 * idz0 * ( (c1 * v(i ,j,k-1) + c2 * v(i ,j,k) + c3 * v(i ,j,k+1))
567  + (c1 * v(i-1,j,k-1) + c2 * v(i-1,j,k) + c3 * v(i-1,j,k+1)) );
568 
569  Real mfy = 0.5 * (mf_uy(i,j,0) + mf_uy(i ,j-1,0));
570  Real mfx = 0.5 * (mf_vx(i,j,0) + mf_vx(i-1,j ,0));
571 
572  Real met_h_xi,met_h_eta;
573  met_h_xi = Compute_h_xi_AtEdgeCenterK (i,j,k,dxInv,z_nd);
574  met_h_eta = Compute_h_eta_AtEdgeCenterK (i,j,k,dxInv,z_nd);
575 
576  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j-1, k))*dxInv[1]*mfy
577  + (v(i, j, k) - v(i-1, j , k))*dxInv[0]*mfx
578  - (met_h_eta)*GradUz*mfy
579  - (met_h_xi )*GradVz*mfx );
580  tau21(i,j,k) = tau12(i,j,k);
581  });
582  }
583 
584  //***********************************************************************************
585  // Z-lo w/out Z-Dirichlet (GradWz extrapolation)
586  //***********************************************************************************
587  if (!zl_u_dir && (tbxxz.smallEnd(2) == domain_xz.smallEnd(2)) ) {
588  Box planexz = tbxxz; planexz.setBig(2, planexz.smallEnd(2));
589  tbxxz.growLo(2,-1);
590  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
591  Real GradWz = 0.5 * dxInv[2] * ( w(i ,j ,k+1) + w(i-1,j ,k+1)
592  - w(i ,j ,k ) - w(i-1,j ,k ) );
593 
594  Real mfx = mf_ux(i,j,0);
595 
596  Real met_h_xi,met_h_zeta;
597  met_h_xi = Compute_h_xi_AtEdgeCenterJ (i,j,k,dxInv,z_nd);
598  met_h_zeta = Compute_h_zeta_AtEdgeCenterJ(i,j,k,dxInv,z_nd);
599 
600  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j, k-1))*dxInv[2]/met_h_zeta
601  + ( (w(i, j, k) - w(i-1, j, k ))*dxInv[0]
602  - (met_h_xi/met_h_zeta)*GradWz ) * mfx);
603  tau31(i,j,k) = tau13(i,j,k);
604  });
605  }
606  if (!zl_v_dir && (tbxyz.smallEnd(2) == domain_yz.smallEnd(2))) {
607  Box planeyz = tbxyz; planeyz.setBig(2, planeyz.smallEnd(2) );
608  tbxyz.growLo(2,-1);
609  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
610  Real GradWz = 0.5 * dxInv[2] * ( w(i ,j ,k+1) + w(i ,j-1,k+1)
611  - w(i ,j ,k ) - w(i ,j-1,k ) );
612 
613  Real mfy = mf_vy(i,j,0);
614 
615  Real met_h_eta,met_h_zeta;
616  met_h_eta = Compute_h_eta_AtEdgeCenterI (i,j,k,dxInv,z_nd);
617  met_h_zeta = Compute_h_zeta_AtEdgeCenterI(i,j,k,dxInv,z_nd);
618 
619  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j , k-1))*dxInv[2]/met_h_zeta
620  + ( (w(i, j, k) - w(i, j-1, k ))*dxInv[1]
621  - (met_h_eta/met_h_zeta)*GradWz ) * mfy );
622  tau32(i,j,k) = tau23(i,j,k);
623  });
624  }
625 
626  //***********************************************************************************
627  // Z-hi w/out Z-Dirichlet (h_xi = h_eta = 0)
628  //***********************************************************************************
629  if (!zh_u_dir && (tbxxz.bigEnd(2) == domain_xz.bigEnd(2))) {
630  Box planexz = tbxxz; planexz.setSmall(2, planexz.bigEnd(2) );
631  tbxxz.growHi(2,-1);
632  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
633  Real mfx = mf_ux(i,j,0);
634 
635  Real met_h_zeta;
636  met_h_zeta = Compute_h_zeta_AtEdgeCenterJ(i,j,k,dxInv,z_nd);
637 
638  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j, k-1))*dxInv[2]/met_h_zeta
639  + (w(i, j, k) - w(i-1, j, k ))*dxInv[0]*mfx );
640  tau31(i,j,k) = tau13(i,j,k);
641  });
642  }
643  if (!zh_v_dir && (tbxyz.bigEnd(2) == domain_yz.bigEnd(2))) {
644  Box planeyz = tbxyz; planeyz.setSmall(2, planeyz.bigEnd(2) );
645  tbxyz.growHi(2,-1);
646  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
647  Real mfy = mf_vy(i,j,0);
648 
649  Real met_h_zeta;
650  met_h_zeta = Compute_h_zeta_AtEdgeCenterI(i,j,k,dxInv,z_nd);
651 
652  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j , k-1))*dxInv[2]/met_h_zeta
653  + (w(i, j, k) - w(i, j-1, k ))*dxInv[1]*mfy );
654  tau32(i,j,k) = tau23(i,j,k);
655  });
656  }
657 
658  //***********************************************************************************
659  // Fill the interior cells
660  //***********************************************************************************
661  // Cell centered strains
662  ParallelFor(bxcc, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
663  Real dz0 = 0.25 * ( z_nd(i,j,k+1) + z_nd(i,j+1,k+1) + z_nd(i+1,j,k+1) + z_nd(i+1,j+1,k+1)
664  - z_nd(i,j,k-1) - z_nd(i,j+1,k-1) - z_nd(i+1,j,k-1) - z_nd(i+1,j+1,k-1) );
665  Real idz0 = 1.0 / dz0;
666 
667  Real GradUz = 0.5 * idz0 * ( u(i ,j ,k+1) + u(i-1,j ,k+1)
668  - u(i ,j ,k-1) - u(i-1,j ,k-1) );
669  Real GradVz = 0.5 * idz0 * ( v(i ,j ,k+1) + v(i ,j-1,k+1)
670  - v(i ,j ,k-1) - v(i ,j-1,k-1) );
671 
672  Real mfx = mf_mx(i,j,0);
673  Real mfy = mf_my(i,j,0);
674 
675  Real met_h_xi,met_h_eta,met_h_zeta;
676  met_h_xi = Compute_h_xi_AtCellCenter (i,j,k,dxInv,z_nd);
677  met_h_eta = Compute_h_eta_AtCellCenter (i,j,k,dxInv,z_nd);
678  met_h_zeta = detJ(i,j,k);
679 
680  tau11(i,j,k) = ( (u(i+1, j, k) - u(i, j, k))*dxInv[0] - met_h_xi*GradUz ) * mfx;
681  tau22(i,j,k) = ( (v(i, j+1, k) - v(i, j, k))*dxInv[1] - met_h_eta*GradVz ) * mfy;
682  tau33(i,j,k) = ( w(i, j, k+1) - w(i, j, k) )*dxInv[2]/met_h_zeta;
683  });
684 
685  // Off-diagonal strains
686  ParallelFor(tbxxy,tbxxz,tbxyz,
687  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
688  Real dz0 = ( z_nd(i,j,k+1) - z_nd(i,j,k-1) );
689  Real idz0 = 1.0 / dz0;
690 
691  Real GradUz = 0.5 * idz0 * ( u(i ,j ,k+1) + u(i ,j-1,k+1)
692  - u(i ,j ,k-1) - u(i ,j-1,k-1) );
693  Real GradVz = 0.5 * idz0 * ( v(i ,j ,k+1) + v(i-1,j ,k+1)
694  - v(i ,j ,k-1) - v(i-1,j ,k-1) );
695 
696  Real mfy = 0.5 * (mf_uy(i,j,0) + mf_uy(i ,j-1,0));
697  Real mfx = 0.5 * (mf_vx(i,j,0) + mf_vx(i-1,j ,0));
698 
699  Real met_h_xi,met_h_eta;
700  met_h_xi = Compute_h_xi_AtEdgeCenterK (i,j,k,dxInv,z_nd);
701  met_h_eta = Compute_h_eta_AtEdgeCenterK (i,j,k,dxInv,z_nd);
702 
703  tau12(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j-1, k))*dxInv[1]*mfy
704  + (v(i, j, k) - v(i-1, j , k))*dxInv[0]*mfx
705  - (met_h_eta)*GradUz*mfy
706  - (met_h_xi)*GradVz*mfx );
707  tau21(i,j,k) = tau12(i,j,k);
708  },
709  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
710  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i,j+1,k+1)
711  - z_nd(i,j,k-1) - z_nd(i,j+1,k-1) );
712  Real idz0 = 1.0 / dz0;
713 
714  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i-1,j ,k+1)
715  - w(i ,j ,k-1) - w(i-1,j ,k-1) );
716 
717  Real mfx = mf_ux(i,j,0);
718 
719  Real met_h_xi,met_h_zeta;
720  met_h_xi = Compute_h_xi_AtEdgeCenterJ (i,j,k,dxInv,z_nd);
721  met_h_zeta = Compute_h_zeta_AtEdgeCenterJ(i,j,k,dxInv,z_nd);
722 
723  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i , j, k-1))*dxInv[2]/met_h_zeta
724  + ( (w(i, j, k) - w(i-1, j, k ))*dxInv[0]
725  - (met_h_xi)*GradWz ) * mfx );
726  tau31(i,j,k) = tau13(i,j,k);
727  },
728  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
729  Real dz0 = 0.5 * ( z_nd(i,j,k+1) + z_nd(i+1,j,k+1)
730  - z_nd(i,j,k-1) - z_nd(i+1,j,k-1) );
731  Real idz0 = 1.0 / dz0;
732 
733  Real GradWz = 0.5 * idz0 * ( w(i ,j ,k+1) + w(i ,j-1,k+1)
734  - w(i ,j ,k-1) - w(i ,j-1,k-1) );
735 
736  Real mfy = mf_vy(i,j,0);
737 
738  Real met_h_eta,met_h_zeta;
739  met_h_eta = Compute_h_eta_AtEdgeCenterI (i,j,k,dxInv,z_nd);
740  met_h_zeta = Compute_h_zeta_AtEdgeCenterI(i,j,k,dxInv,z_nd);
741 
742  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j , k-1))*dxInv[2]/met_h_zeta
743  + ( (w(i, j, k) - w(i, j-1, k ))*dxInv[1]
744  - (met_h_eta)*GradWz ) * mfy );
745  tau32(i,j,k) = tau23(i,j,k);
746  });
747 }
tau12
@ tau12
Definition: ERF_DataStruct.H:30
tau23
@ tau23
Definition: ERF_DataStruct.H:30
tau33
@ tau33
Definition: ERF_DataStruct.H:30
tau22
@ tau22
Definition: ERF_DataStruct.H:30
tau11
@ tau11
Definition: ERF_DataStruct.H:30
tau32
@ tau32
Definition: ERF_DataStruct.H:30
tau31
@ tau31
Definition: ERF_DataStruct.H:30
tau21
@ tau21
Definition: ERF_DataStruct.H:30
tau13
@ tau13
Definition: ERF_DataStruct.H:30
dom_hi
const auto & dom_hi
Definition: ERF_DiffSetup.H:10
dom_lo
const auto & dom_lo
Definition: ERF_DiffSetup.H:9
Compute_h_xi_AtEdgeCenterK
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_xi_AtEdgeCenterK(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:244
Compute_h_xi_AtEdgeCenterJ
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_xi_AtEdgeCenterJ(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:290
Compute_h_zeta_AtEdgeCenterJ
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_zeta_AtEdgeCenterJ(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:276
Compute_h_zeta_AtEdgeCenterI
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_zeta_AtEdgeCenterI(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:320
Compute_h_eta_AtCellCenter
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_eta_AtCellCenter(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:77
Compute_h_eta_AtEdgeCenterK
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_eta_AtEdgeCenterK(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:259
Compute_h_xi_AtCellCenter
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_xi_AtCellCenter(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:62
Compute_h_eta_AtEdgeCenterI
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_eta_AtEdgeCenterI(const int &i, const int &j, const int &k, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &z_nd)
Definition: ERF_TerrainMetrics.H:347
BCVars::zvel_bc
@ zvel_bc
Definition: ERF_IndexDefines.H:89
BCVars::yvel_bc
@ yvel_bc
Definition: ERF_IndexDefines.H:88
BCVars::xvel_bc
@ xvel_bc
Definition: ERF_IndexDefines.H:87
ERFBCType::ext_dir_ingested
@ ext_dir_ingested
Definition: ERF_IndexDefines.H:212
ERFBCType::ext_dir
@ ext_dir
Definition: ERF_IndexDefines.H:209
ERFBCType::ext_dir_upwind
@ ext_dir_upwind
Definition: ERF_IndexDefines.H:217
module_model_constants::c2
real(c_double), parameter c2
Definition: ERF_module_model_constants.F90:35
module_mp_morr_two_moment::c1
real(c_double), private c1
Definition: ERF_module_mp_morr_two_moment.F90:211

Referenced by ERF::advance_dycore(), and erf_make_tau_terms().

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