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

Functions

void ComputeStrain_N (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 BCRec *bc_ptr, const GpuArray< Real, AMREX_SPACEDIM > &dxInv, const Array4< const Real > &mf_m, const Array4< const Real > &mf_u, const Array4< const Real > &mf_v)
 

Function Documentation

◆ ComputeStrain_N()

void ComputeStrain_N ( 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 BCRec *  bc_ptr,
const GpuArray< Real, AMREX_SPACEDIM > &  dxInv,
const Array4< const Real > &  mf_m,
const Array4< const Real > &  mf_u,
const Array4< const Real > &  mf_v 
)

Function for computing the strain rates without 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]tau2323 strain
[in]bc_ptrcontainer with boundary condition types
[in]dxInvinverse cell size array
[in]mf_mmap factor at cell center
[in]mf_umap factor at x-face
[in]mf_vmap factor at y-face
34 {
35  // Convert domain to each index type to test if we are on dirichlet boundary
36  Box domain_xy = convert(domain, tbxxy.ixType());
37  Box domain_xz = convert(domain, tbxxz.ixType());
38  Box domain_yz = convert(domain, tbxyz.ixType());
39 
40  // Dirichlet on left or right plane
41  bool xl_v_dir = ( (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir) ||
42  (bc_ptr[BCVars::yvel_bc].lo(0) == ERFBCType::ext_dir_ingested) );
43  xl_v_dir = ( xl_v_dir && (tbxxy.smallEnd(0) == domain_xy.smallEnd(0)) );
44 
45  bool xh_v_dir = ( (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir) ||
46  (bc_ptr[BCVars::yvel_bc].hi(0) == ERFBCType::ext_dir_ingested) );
47  xh_v_dir = ( xh_v_dir && (tbxxy.bigEnd(0) == domain_xy.bigEnd(0)) );
48 
49  bool xl_w_dir = ( (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir) ||
50  (bc_ptr[BCVars::zvel_bc].lo(0) == ERFBCType::ext_dir_ingested) );
51  xl_w_dir = ( xl_w_dir && (tbxxz.smallEnd(0) == domain_xz.smallEnd(0)) );
52 
53  bool xh_w_dir = ( (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir) ||
54  (bc_ptr[BCVars::zvel_bc].hi(0) == ERFBCType::ext_dir_ingested) );
55  xh_w_dir = ( xh_w_dir && (tbxxz.bigEnd(0) == domain_xz.bigEnd(0)) );
56 
57  // Dirichlet on front or back plane
58  bool yl_u_dir = ( (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir) ||
59  (bc_ptr[BCVars::xvel_bc].lo(1) == ERFBCType::ext_dir_ingested) );
60  yl_u_dir = ( yl_u_dir && (tbxxy.smallEnd(1) == domain_xy.smallEnd(1)) );
61 
62  bool yh_u_dir = ( (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir) ||
63  (bc_ptr[BCVars::xvel_bc].hi(1) == ERFBCType::ext_dir_ingested) );
64  yh_u_dir = ( yh_u_dir && (tbxxy.bigEnd(1) == domain_xy.bigEnd(1)) );
65 
66  bool yl_w_dir = ( (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir) ||
67  (bc_ptr[BCVars::zvel_bc].lo(1) == ERFBCType::ext_dir_ingested) );
68  yl_w_dir = ( yl_w_dir && (tbxyz.smallEnd(1) == domain_yz.smallEnd(1)) );
69 
70  bool yh_w_dir = ( (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir) ||
71  (bc_ptr[BCVars::zvel_bc].hi(1) == ERFBCType::ext_dir_ingested) );
72  yh_w_dir = ( yh_w_dir && (tbxyz.bigEnd(1) == domain_yz.bigEnd(1)) );
73 
74  // Dirichlet on top or bottom plane
75  bool zl_u_dir = ( (bc_ptr[BCVars::xvel_bc].lo(2) == ERFBCType::ext_dir) ||
76  (bc_ptr[BCVars::xvel_bc].lo(2) == ERFBCType::ext_dir_ingested) );
77  zl_u_dir = ( zl_u_dir && (tbxxz.smallEnd(2) == domain_xz.smallEnd(2)) );
78 
79  bool zh_u_dir = ( (bc_ptr[BCVars::xvel_bc].hi(2) == ERFBCType::ext_dir) ||
80  (bc_ptr[BCVars::xvel_bc].hi(2) == ERFBCType::ext_dir_ingested) );
81  zh_u_dir = ( zh_u_dir && (tbxxz.bigEnd(2) == domain_xz.bigEnd(2)) );
82 
83  bool zl_v_dir = ( (bc_ptr[BCVars::yvel_bc].lo(2) == ERFBCType::ext_dir) ||
84  (bc_ptr[BCVars::yvel_bc].lo(2) == ERFBCType::ext_dir_ingested) );
85  zl_v_dir = ( zl_v_dir && (tbxyz.smallEnd(2) == domain_yz.smallEnd(2)) );
86 
87  bool zh_v_dir = ( (bc_ptr[BCVars::yvel_bc].hi(2) == ERFBCType::ext_dir) ||
88  (bc_ptr[BCVars::yvel_bc].hi(2) == ERFBCType::ext_dir_ingested) );
89  zh_v_dir = ( zh_v_dir && (tbxyz.bigEnd(2) == domain_yz.bigEnd(2)) );
90 
91  // X-Dirichlet
92  //***********************************************************************************
93  if (xl_v_dir) {
94  Box planexy = tbxxy; planexy.setBig(0, planexy.smallEnd(0) );
95  tbxxy.growLo(0,-1);
96  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
97  tau12(i,j,k) = 0.5 * ( (u(i, j, k)/mf_u(i,j,0) - u(i, j-1, k)/mf_u(i,j-1,0))*dxInv[1] +
98  (-(8./3.) * v(i-1,j,k)/mf_v(i-1,j,0) + 3. * v(i,j,k)/mf_v(i,j,0) -
99  (1./3.) * v(i+1,j,k)/mf_v(i+1,j,0))*dxInv[0] ) * mf_u(i,j,0)*mf_u(i,j,0);
100  });
101  }
102  if (xh_v_dir) {
103  // 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
104  Box planexy = tbxxy; planexy.setSmall(0, planexy.bigEnd(0) );
105  tbxxy.growHi(0,-1);
106  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
107  tau12(i,j,k) = 0.5 * ( (u(i, j, k)/mf_u(i,j,0) - u(i, j-1, k)/mf_u(i,j-1,0))*dxInv[1] +
108  -(-(8./3.) * v(i,j,k)/mf_v(i,j,0) + 3. * v(i-1,j,k)/mf_v(i-1,j,0) -
109  (1./3.) * v(i-2,j,k)/mf_v(i-2,j,0))*dxInv[0] ) * mf_u(i,j,0)*mf_u(i,j,0);
110  });
111  }
112 
113  if (xl_w_dir) {
114  Box planexz = tbxxz; planexz.setBig(0, planexz.smallEnd(0) );
115  tbxxz.growLo(0,-1);
116  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
117  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j, k-1))*dxInv[2] +
118  (-(8./3.) * w(i-1,j,k) + 3. * w(i,j,k) - (1./3.) * w(i+1,j,k))*dxInv[0]*mf_u(i,j,0) );
119  });
120  }
121  if (xh_w_dir) {
122  // 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
123  Box planexz = tbxxz; planexz.setSmall(0, planexz.bigEnd(0) );
124  tbxxz.growHi(0,-1);
125  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
126  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j, k-1))*dxInv[2] +
127  -(-(8./3.) * w(i,j,k) + 3. * w(i-1,j,k) - (1./3.) * w(i-2,j,k))*dxInv[0]*mf_u(i,j,0) );
128  });
129  }
130 
131  // Y-Dirichlet
132  //***********************************************************************************
133  if (yl_u_dir) {
134  Box planexy = tbxxy; planexy.setBig(1, planexy.smallEnd(1) );
135  tbxxy.growLo(1,-1);
136  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
137  tau12(i,j,k) = 0.5 * ( (-(8./3.) * u(i,j-1,k)/mf_u(i,j-1,0) + 3. * u(i,j,k)/mf_u(i,j,0) - (1./3.) * u(i,j+1,k)/mf_m(i,j+1,0))*dxInv[1] +
138  (v(i, j, k)/mf_v(i,j,0) - v(i-1, j, k)/mf_v(i-1,j,0))*dxInv[0] ) * mf_u(i,j,0)*mf_u(i,j,0);
139  });
140  }
141  if (yh_u_dir) {
142  // 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
143  Box planexy = tbxxy; planexy.setSmall(1, planexy.bigEnd(1) );
144  tbxxy.growHi(1,-1);
145  ParallelFor(planexy,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
146  tau12(i,j,k) = 0.5 * ( -(-(8./3.) * u(i,j,k)/mf_u(i,j,0) + 3. * u(i,j-1,k)/mf_u(i,j-1,0) - (1./3.) * u(i,j-2,k)/mf_u(i,j-2,0))*dxInv[1] +
147  (v(i, j, k)/mf_v(i,j,0) - v(i-1, j, k)/mf_v(i-1,j,0))*dxInv[0] ) * mf_u(i,j,0)*mf_u(i,j,0);
148  });
149  }
150 
151  if (yl_w_dir) {
152  Box planeyz = tbxyz; planeyz.setBig(1, planeyz.smallEnd(1) );
153  tbxyz.growLo(1,-1);
154  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
155  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j, k-1))*dxInv[2] +
156  (-(8./3.) * w(i,j-1,k) + 3. * w(i,j ,k) - (1./3.) * w(i,j+1,k))*dxInv[1]*mf_v(i,j,0)*mf_v(i,j,0) );
157  });
158  }
159  if (yh_w_dir) {
160  // 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
161  Box planeyz = tbxyz; planeyz.setSmall(1, planeyz.bigEnd(1) );
162  tbxyz.growHi(1,-1);
163  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
164  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j, k-1))*dxInv[2] +
165  -(-(8./3.) * w(i,j ,k) + 3. * w(i,j-1,k) - (1./3.) * w(i,j-2,k))*dxInv[1]*mf_v(i,j,0)*mf_v(i,j,0) );
166  });
167  }
168 
169  // Z-Dirichlet
170  //***********************************************************************************
171  if (zl_u_dir) {
172  Box planexz = tbxxz; planexz.setBig(2, planexz.smallEnd(2) );
173  tbxxz.growLo(2,-1);
174  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
175  tau13(i,j,k) = 0.5 * ( (-(8./3.) * u(i,j,k-1) + 3. * u(i,j,k) - (1./3.) * u(i,j,k+1))*dxInv[2] +
176  (w(i, j, k) - w(i-1, j, k))*dxInv[0]*mf_u(i,j,0) );
177  });
178  }
179  if (zh_u_dir) {
180  // 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
181  Box planexz = tbxxz; planexz.setSmall(2, planexz.bigEnd(2) );
182  tbxxz.growHi(2,-1);
183  ParallelFor(planexz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
184  tau13(i,j,k) = 0.5 * ( -(-(8./3.) * u(i,j,k) + 3. * u(i,j,k-1) - (1./3.) * u(i,j,k-2))*dxInv[2] +
185  (w(i, j, k) - w(i-1, j, k))*dxInv[0]*mf_u(i,j,0) );
186  });
187  }
188 
189  if (zl_v_dir) {
190  Box planeyz = tbxyz; planeyz.setBig(2, planeyz.smallEnd(2) );
191  tbxyz.growLo(2,-1);
192  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
193  tau23(i,j,k) = 0.5 * ( (-(8./3.) * v(i,j,k-1) + 3. * v(i,j,k ) - (1./3.) * v(i,j,k+1))*dxInv[2] +
194  (w(i, j, k) - w(i, j-1, k))*dxInv[1]*mf_v(i,j,0) );
195  });
196  }
197  if (zh_v_dir) {
198  // 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
199  Box planeyz = tbxyz; planeyz.setSmall(2, planeyz.bigEnd(2) );
200  tbxyz.growHi(2,-1);
201  ParallelFor(planeyz,[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
202  tau23(i,j,k) = 0.5 * ( -(-(8./3.) * v(i,j,k ) + 3. * v(i,j,k-1) - (1./3.) * v(i,j,k-2))*dxInv[2] +
203  (w(i, j, k) - w(i, j-1, k))*dxInv[1]*mf_v(i,j,0) );
204  });
205  }
206 
207  // Fill the remaining cells
208  //***********************************************************************************
209  // Cell centered strains
210  ParallelFor(bxcc, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
211  tau11(i,j,k) = (u(i+1, j , k )/mf_u(i+1,j,0) - u(i, j, k)/mf_u(i,j,0))*dxInv[0]*mf_m(i,j,0)*mf_m(i,j,0);
212  tau22(i,j,k) = (v(i , j+1, k )/mf_v(i,j+1,0) - v(i, j, k)/mf_v(i,j,0))*dxInv[1]*mf_m(i,j,0)*mf_m(i,j,0);
213  tau33(i,j,k) = (w(i , j , k+1) - w(i, j, k))*dxInv[2];
214  });
215 
216  // Off-diagonal strains
217  ParallelFor(tbxxy,tbxxz,tbxyz,
218  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
219  tau12(i,j,k) = 0.5 * ( (u(i, j, k)/mf_u(i,j,0) - u(i, j-1, k)/mf_u(i,j-1,0))*dxInv[1] +
220  (v(i, j, k)/mf_v(i,j,0) - v(i-1, j, k)/mf_v(i-1,j,0))*dxInv[0] ) * mf_u(i,j,0)*mf_u(i,j,0);
221  },
222  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
223  tau13(i,j,k) = 0.5 * ( (u(i, j, k) - u(i, j, k-1))*dxInv[2] + (w(i, j, k) - w(i-1, j, k))*dxInv[0]*mf_u(i,j,0) );
224  },
225  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
226  tau23(i,j,k) = 0.5 * ( (v(i, j, k) - v(i, j, k-1))*dxInv[2] + (w(i, j, k) - w(i, j-1, k))*dxInv[1]*mf_v(i,j,0) );
227  });
228 
229 }
BCVars::zvel_bc
@ zvel_bc
Definition: ERF_IndexDefines.H:90
BCVars::yvel_bc
@ yvel_bc
Definition: ERF_IndexDefines.H:89
BCVars::xvel_bc
@ xvel_bc
Definition: ERF_IndexDefines.H:88
ERFBCType::ext_dir_ingested
@ ext_dir_ingested
Definition: ERF_IndexDefines.H:183
ERFBCType::ext_dir
@ ext_dir
Definition: ERF_IndexDefines.H:180

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

Here is the caller graph for this function: