ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_AdvectionSrcForMom_T.H
Go to the documentation of this file.
1 #include <ERF_IndexDefines.H>
2 #include <ERF_TerrainMetrics.H>
3 #include <ERF_Interpolation.H>
4 
5 /**
6  * Function for computing the advective tendency for the x-component of momentum
7  * with metric terms and for spatial order > 2
8  *
9  * @param[in] i,j,k indices of x-face at which to create tendency
10  * @param[in] rho_u x-component of momentum
11  * @param[in] rho_v y-component of momentum
12  * @param[in] Omega component of the momentum normal to the z-coordinate surface
13  * @param[in] z_nd height coordinate at nodes
14  * @param[in] ax Area fractions on x-faces
15  * @param[in] ay Area fractions on y-faces
16  * @param[in] az Area fractions on z-faces
17  * @param[in] detJ Jacobian of the metric transformation (= 1 if use_terrain is false)
18  * @param[in] cellSizeInv inverse of the mesh spacing
19  * @param[in] mf_u map factor on x-faces
20  * @param[in] mf_v map factor on y-faces
21  */
22 template<typename InterpType_H, typename InterpType_V>
23 AMREX_GPU_DEVICE
24 AMREX_FORCE_INLINE
25 amrex::Real
26 AdvectionSrcForXMom (int i, int j, int k,
27  const amrex::Array4<const amrex::Real>& rho_u,
28  const amrex::Array4<const amrex::Real>& rho_v,
29  const amrex::Array4<const amrex::Real>& Omega,
30  const amrex::Array4<const amrex::Real>& z_nd,
31  const amrex::Array4<const amrex::Real>& ax,
32  const amrex::Array4<const amrex::Real>& /*ay*/,
33  const amrex::Array4<const amrex::Real>& az,
34  const amrex::Array4<const amrex::Real>& detJ,
35  InterpType_H interp_u_h,
36  InterpType_V interp_u_v,
37  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
38  const amrex::Array4<const amrex::Real>& mf_u_inv,
39  const amrex::Array4<const amrex::Real>& mf_v_inv)
40 {
41  amrex::Real advectionSrc;
42  auto dxInv = cellSizeInv[0], dyInv = cellSizeInv[1], dzInv = cellSizeInv[2];
43 
44  amrex::Real rho_u_avg_lo, rho_u_avg_hi;
45  amrex::Real rho_v_avg_lo, rho_v_avg_hi;
46  amrex::Real Omega_avg_lo, Omega_avg_hi;
47 
48  amrex::Real interp_hi(0.), interp_lo(0.);
49 
50  // ****************************************************************************************
51  // X-fluxes (at cell centers)
52  // ****************************************************************************************
53 
54  rho_u_avg_hi = 0.5 * (rho_u(i+1, j, k) * mf_u_inv(i+1,j ,0) + rho_u(i, j, k) * mf_u_inv(i ,j ,0));
55  interp_u_h.InterpolateInX(i+1,j,k,0,interp_hi,rho_u_avg_hi);
56 
57  rho_u_avg_lo = 0.5 * (rho_u(i-1, j, k) * mf_u_inv(i-1,j ,0) + rho_u(i, j, k) * mf_u_inv(i ,j ,0));
58  interp_u_h.InterpolateInX(i,j,k,0,interp_lo,rho_u_avg_lo);
59 
60  amrex::Real centFluxXXNext = rho_u_avg_hi * interp_hi * 0.5 * (ax(i,j,k) + ax(i+1,j,k));
61  amrex::Real centFluxXXPrev = rho_u_avg_lo * interp_lo * 0.5 * (ax(i,j,k) + ax(i-1,j,k));
62 
63  // ****************************************************************************************
64  // Y-fluxes (at edges in k-direction)
65  // ****************************************************************************************
66  rho_v_avg_hi = 0.5 * (rho_v(i, j+1, k) * mf_v_inv(i ,j+1,0) + rho_v(i-1, j+1, k) * mf_v_inv(i-1,j+1,0));
67  interp_u_h.InterpolateInY(i,j+1,k,0,interp_hi,rho_v_avg_hi);
68 
69  rho_v_avg_lo = 0.5 * (rho_v(i, j , k) * mf_v_inv(i ,j ,0) + rho_v(i-1, j , k) * mf_v_inv(i-1,j ,0));
70  interp_u_h.InterpolateInY(i,j,k,0,interp_lo,rho_v_avg_lo);
71 
72  amrex::Real edgeFluxXYNext = rho_v_avg_hi * interp_hi * Compute_h_zeta_AtEdgeCenterK(i,j+1,k,cellSizeInv,z_nd);
73  amrex::Real edgeFluxXYPrev = rho_v_avg_lo * interp_lo * Compute_h_zeta_AtEdgeCenterK(i,j ,k,cellSizeInv,z_nd);
74 
75  // ****************************************************************************************
76  // Z-fluxes (at edges in j-direction)
77  // ****************************************************************************************
78  Omega_avg_hi = 0.5 * (Omega(i, j, k+1) + Omega(i-1, j, k+1)) * 0.5 * (az(i,j,k+1) + az(i-1,j,k+1));
79  Omega_avg_lo = 0.5 * (Omega(i, j, k ) + Omega(i-1, j, k )) * 0.5 * (az(i,j,k ) + az(i-1,j,k ));
80 
81  interp_u_v.InterpolateInZ(i,j,k+1,0,interp_hi,Omega_avg_hi);
82  interp_u_v.InterpolateInZ(i,j,k ,0,interp_lo,Omega_avg_lo);
83 
84  amrex::Real edgeFluxXZNext = Omega_avg_hi * interp_hi;
85  amrex::Real edgeFluxXZPrev = Omega_avg_lo * interp_lo;
86 
87  // ****************************************************************************************
88 
89  amrex::Real mfsq = 1 / (mf_u_inv(i,j,0) * mf_u_inv(i,j,0));
90 
91  advectionSrc = (centFluxXXNext - centFluxXXPrev) * dxInv * mfsq
92  + (edgeFluxXYNext - edgeFluxXYPrev) * dyInv * mfsq
93  + (edgeFluxXZNext - edgeFluxXZPrev) * dzInv;
94  advectionSrc /= 0.5*(detJ(i,j,k) + detJ(i-1,j,k));
95 
96  return advectionSrc;
97 }
98 
99 /**
100  * Function for computing the advective tendency for the y-component of momentum
101  * with metric terms and for spatial order > 2
102  *
103  * @param[in] i,j,k indices of y-face at which to create tendency
104  * @param[in] rho_u x-component of momentum
105  * @param[in] rho_v y-component of momentum
106  * @param[in] Omega component of the momentum normal to the z-coordinate surface
107  * @param[in] z_nd height coordinate at nodes
108  * @param[in] ax Area fractions on x-faces
109  * @param[in] ay Area fractions on y-faces
110  * @param[in] az Area fractions on z-faces
111  * @param[in] detJ Jacobian of the metric transformation (= 1 if use_terrain is false)
112  * @param[in] cellSizeInv inverse of the mesh spacing
113  * @param[in] mf_u map factor on x-faces
114  * @param[in] mf_v map factor on y-faces
115  */
116 template<typename InterpType_H, typename InterpType_V>
117 AMREX_GPU_DEVICE
118 AMREX_FORCE_INLINE
119 amrex::Real
120 AdvectionSrcForYMom (int i, int j, int k,
121  const amrex::Array4<const amrex::Real>& rho_u,
122  const amrex::Array4<const amrex::Real>& rho_v,
123  const amrex::Array4<const amrex::Real>& Omega,
124  const amrex::Array4<const amrex::Real>& z_nd,
125  const amrex::Array4<const amrex::Real>& /*ax*/,
126  const amrex::Array4<const amrex::Real>& ay,
127  const amrex::Array4<const amrex::Real>& az,
128  const amrex::Array4<const amrex::Real>& detJ,
129  InterpType_H interp_v_h,
130  InterpType_V interp_v_v,
131  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
132  const amrex::Array4<const amrex::Real>& mf_u_inv,
133  const amrex::Array4<const amrex::Real>& mf_v_inv)
134 {
135  amrex::Real advectionSrc;
136  auto dxInv = cellSizeInv[0], dyInv = cellSizeInv[1], dzInv = cellSizeInv[2];
137 
138  amrex::Real rho_u_avg_lo, rho_u_avg_hi;
139  amrex::Real rho_v_avg_lo, rho_v_avg_hi;
140  amrex::Real Omega_avg_lo, Omega_avg_hi;
141 
142  amrex::Real interp_hi(0.), interp_lo(0.);
143 
144  // ****************************************************************************************
145  // x-fluxes (at edges in k-direction)
146  // ****************************************************************************************
147  rho_u_avg_hi = 0.5 * (rho_u(i+1,j,k) * mf_u_inv(i+1,j,0) + rho_u(i+1,j-1,k) * mf_u_inv(i+1,j-1,0));
148  rho_u_avg_lo = 0.5 * (rho_u(i ,j,k) * mf_u_inv(i ,j,0) + rho_u(i ,j-1,k) * mf_u_inv(i ,j-1,0));
149 
150  interp_v_h.InterpolateInX(i+1,j,k,0,interp_hi,rho_u_avg_hi);
151  interp_v_h.InterpolateInX(i ,j,k,0,interp_lo,rho_u_avg_lo);
152 
153  amrex::Real edgeFluxYXNext = rho_u_avg_hi * interp_hi * Compute_h_zeta_AtEdgeCenterK(i+1,j,k,cellSizeInv,z_nd);
154  amrex::Real edgeFluxYXPrev = rho_u_avg_lo * interp_lo * Compute_h_zeta_AtEdgeCenterK(i ,j,k,cellSizeInv,z_nd);
155 
156  // ****************************************************************************************
157  // y-fluxes (at cell centers)
158  // ****************************************************************************************
159  rho_v_avg_hi = 0.5 * (rho_v(i,j,k) * mf_v_inv(i,j,0) + rho_v(i,j+1,k) * mf_v_inv(i,j+1,0));
160  rho_v_avg_lo = 0.5 * (rho_v(i,j,k) * mf_v_inv(i,j,0) + rho_v(i,j-1,k) * mf_v_inv(i,j-1,0));
161 
162  interp_v_h.InterpolateInY(i,j+1,k,0,interp_hi,rho_v_avg_hi);
163  interp_v_h.InterpolateInY(i,j ,k,0,interp_lo,rho_v_avg_lo);
164 
165  amrex::Real centFluxYYNext = rho_v_avg_hi * 0.5 * (ay(i,j,k) + ay(i,j+1,k)) * interp_hi;
166  amrex::Real centFluxYYPrev = rho_v_avg_lo * 0.5 * (ay(i,j,k) + ay(i,j-1,k)) * interp_lo;
167 
168  // ****************************************************************************************
169  // Z-fluxes (at edges in j-direction)
170  // ****************************************************************************************
171  Omega_avg_hi = 0.5 * (Omega(i, j, k+1) + Omega(i, j-1, k+1)) * 0.5 * (az(i,j,k+1) + az(i,j-1,k+1));
172  Omega_avg_lo = 0.5 * (Omega(i, j, k ) + Omega(i, j-1, k )) * 0.5 * (az(i,j,k ) + az(i,j-1,k ));
173 
174  interp_v_v.InterpolateInZ(i,j,k+1,0,interp_hi,Omega_avg_hi);
175  interp_v_v.InterpolateInZ(i,j,k ,0,interp_lo,Omega_avg_lo);
176 
177  amrex::Real edgeFluxYZNext = Omega_avg_hi * interp_hi;
178  amrex::Real edgeFluxYZPrev = Omega_avg_lo * interp_lo;
179 
180  // ****************************************************************************************
181 
182  amrex::Real mfsq = 1 / (mf_v_inv(i,j,0) * mf_v_inv(i,j,0));
183 
184  advectionSrc = (edgeFluxYXNext - edgeFluxYXPrev) * dxInv * mfsq
185  + (centFluxYYNext - centFluxYYPrev) * dyInv * mfsq
186  + (edgeFluxYZNext - edgeFluxYZPrev) * dzInv;
187  advectionSrc /= 0.5*(detJ(i,j,k) + detJ(i,j-1,k));
188 
189  return advectionSrc;
190 }
191 
192 /**
193  * Function for computing the advective tendency for the z-component of momentum
194  * with metric terms and for spatial order > 2
195  *
196  * @param[in] i,j,k indices of z-face at which to create tendency
197  * @param[in] rho_u x-component of momentum
198  * @param[in] rho_v y-component of momentum
199  * @param[in] Omega component of the momentum normal to the z-coordinate surface
200  * @param[in] w z-component of velocity
201  * @param[in] z_nd height coordinate at nodes
202  * @param[in] ax Area fractions on x-faces
203  * @param[in] ay Area fractions on y-faces
204  * @param[in] az Area fractions on z-faces
205  * @param[in] detJ Jacobian of the metric transformation (= 1 if use_terrain is false)
206  * @param[in] cellSizeInv inverse of the mesh spacing
207  * @param[in] mf_m map factor on cell centers
208  * @param[in] mf_u map factor on x-faces
209  * @param[in] mf_v map factor on y-faces
210  * @param[in] vert_adv_type int that defines advection stencil
211  * @param[in] lo_z_face minimum k value (z-face-centered)_in the domain at this level
212  * @param[in] hi_z_face maximum k value (z-face-centered) in the domain at this level
213  */
214 template<typename InterpType_H, typename InterpType_V, typename WallInterpType>
215 AMREX_GPU_DEVICE
216 AMREX_FORCE_INLINE
217 amrex::Real
218 AdvectionSrcForZMom (int i, int j, int k,
219  const amrex::Array4<const amrex::Real>& rho_u,
220  const amrex::Array4<const amrex::Real>& rho_v,
221  const amrex::Array4<const amrex::Real>& Omega,
222  const amrex::Array4<const amrex::Real>& w,
223  const amrex::Array4<const amrex::Real>& z_nd,
224  const amrex::Array4<const amrex::Real>& /*ax*/,
225  const amrex::Array4<const amrex::Real>& /*ay*/,
226  const amrex::Array4<const amrex::Real>& az,
227  const amrex::Array4<const amrex::Real>& detJ,
228  InterpType_H interp_omega_h,
229  InterpType_V interp_omega_v,
230  WallInterpType interp_omega_wall,
231  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
232  const amrex::Array4<const amrex::Real>& mf_m,
233  const amrex::Array4<const amrex::Real>& mf_u_inv,
234  const amrex::Array4<const amrex::Real>& mf_v_inv,
235  const AdvType vert_adv_type,
236  const int lo_z_face, const int hi_z_face)
237 {
238  amrex::Real advectionSrc;
239  auto dxInv = cellSizeInv[0], dyInv = cellSizeInv[1], dzInv = cellSizeInv[2];
240 
241  amrex::Real rho_u_avg_lo, rho_u_avg_hi;
242  amrex::Real rho_v_avg_lo, rho_v_avg_hi;
243  amrex::Real Omega_avg_lo, Omega_avg_hi;
244 
245  amrex::Real interp_hi(0.), interp_lo(0.);
246 
247  // ****************************************************************************************
248  // x-fluxes (at edges in j-direction)
249  // ****************************************************************************************
250  rho_u_avg_hi = 0.5 * (rho_u(i+1, j, k) + rho_u(i+1, j, k-1)) * mf_u_inv(i+1,j ,0);
251  rho_u_avg_lo = 0.5 * (rho_u(i , j, k) + rho_u(i , j, k-1)) * mf_u_inv(i ,j ,0);
252 
253  interp_omega_h.InterpolateInX(i+1,j,k,0,interp_hi,rho_u_avg_hi);
254  interp_omega_h.InterpolateInX(i,j,k,0,interp_lo,rho_u_avg_lo);
255 
256  amrex::Real edgeFluxZXNext = rho_u_avg_hi * interp_hi * Compute_h_zeta_AtEdgeCenterJ(i+1,j,k,cellSizeInv,z_nd);
257  amrex::Real edgeFluxZXPrev = rho_u_avg_lo * interp_lo * Compute_h_zeta_AtEdgeCenterJ(i ,j,k,cellSizeInv,z_nd);
258 
259  // ****************************************************************************************
260  // y-fluxes (at edges in i-direction)
261  // ****************************************************************************************
262  rho_v_avg_hi = 0.5 * (rho_v(i, j+1, k) + rho_v(i, j+1, k-1)) * mf_v_inv(i ,j+1,0);
263  interp_omega_h.InterpolateInY(i,j+1,k,0,interp_hi,rho_v_avg_hi);
264 
265  rho_v_avg_lo = 0.5 * (rho_v(i, j , k) + rho_v(i, j , k-1)) * mf_v_inv(i ,j ,0);
266  interp_omega_h.InterpolateInY(i,j,k,0,interp_lo,rho_v_avg_lo);
267 
268  amrex::Real edgeFluxZYNext = rho_v_avg_hi * interp_hi * Compute_h_zeta_AtEdgeCenterI(i,j+1,k,cellSizeInv,z_nd);
269  amrex::Real edgeFluxZYPrev = rho_v_avg_lo * interp_lo * Compute_h_zeta_AtEdgeCenterI(i,j ,k,cellSizeInv,z_nd);
270 
271  // ****************************************************************************************
272  // z-fluxes (at cell centers)
273  // ****************************************************************************************
274 
275  Omega_avg_hi = (k == hi_z_face) ? Omega(i,j,k) * az(i,j,k) :
276  0.25 * (Omega(i,j,k) + Omega(i,j,k+1)) * (az(i,j,k) + az(i,j,k+1));
277  amrex::Real centFluxZZNext = Omega_avg_hi;
278 
279  // int l_spatial_order_hi = std::min(std::min(vert_spatial_order, 2*(hi_z_face-k)), 2*(k+1));
280  // If k == hi_z_face-1, l_spatial_order_hi = 2
281  // If k == hi_z_face-2, l_spatial_order_hi = std::min(vert_spatial_order, 4);
282  // If k == lo_z_face+11 , l_spatial_order_hi = std::min(vert_spatial_order, 4);
283  if (k == hi_z_face) {
284  centFluxZZNext *= w(i,j,k);
285  } else {
286  if (k == hi_z_face-1)
287  {
288  interp_omega_wall.InterpolateInZ(i,j,k+1,0,interp_hi,Omega_avg_hi,AdvType::Centered_2nd);
289  } else if (k == hi_z_face-2 || k == lo_z_face+1) {
290  if (vert_adv_type != AdvType::Centered_2nd && vert_adv_type != AdvType::Upwind_3rd) {
291  interp_omega_wall.InterpolateInZ(i,j,k+1,0,interp_hi,Omega_avg_hi,AdvType::Centered_4th);
292  } else {
293  interp_omega_wall.InterpolateInZ(i,j,k+1,0,interp_hi,Omega_avg_hi,vert_adv_type);
294  }
295  } else {
296  interp_omega_v.InterpolateInZ(i,j,k+1,0,interp_hi,Omega_avg_hi);
297  }
298  centFluxZZNext *= interp_hi;
299  }
300 
301  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
302 
303  Omega_avg_lo = (k == 0) ? Omega(i,j,k) * az(i,j,k) :
304  0.25 * (Omega(i,j,k) + Omega(i,j,k-1)) * (az(i,j,k) + az(i,j,k-1));
305  amrex::Real centFluxZZPrev = Omega_avg_lo;
306 
307  // int l_spatial_order_lo = std::min(std::min(vert_spatial_order, 2*(hi_z_face+1-k)), 2*k);
308  // If k == hi_z_face-1, l_spatial_order_hi = 2
309  // If k == hi_z_face-2, l_spatial_order_hi = std::min(vert_spatial_order, 4);
310  // If k == lo_z_face+1, l_spatial_order_hi = std::min(vert_spatial_order, 4);
311  if (k == 0) {
312  centFluxZZPrev *= w(i,j,k);
313  } else {
314  if (k == lo_z_face+1) {
315  interp_omega_wall.InterpolateInZ(i,j,k,0,interp_lo,Omega_avg_lo,AdvType::Centered_2nd);
316  } else if (k == lo_z_face+2 || k == hi_z_face-1) {
317  if (vert_adv_type != AdvType::Centered_2nd && vert_adv_type != AdvType::Upwind_3rd) {
318  interp_omega_wall.InterpolateInZ(i,j,k,0,interp_lo,Omega_avg_lo,AdvType::Centered_4th);
319  } else {
320  interp_omega_wall.InterpolateInZ(i,j,k,0,interp_lo,Omega_avg_lo,vert_adv_type);
321  }
322  } else {
323  interp_omega_v.InterpolateInZ(i,j,k,0,interp_lo,Omega_avg_lo);
324  }
325  centFluxZZPrev *= interp_lo;
326  }
327 
328  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
329 
330  amrex::Real mfsq = mf_m(i,j,0) * mf_m(i,j,0);
331 
332  advectionSrc = (edgeFluxZXNext - edgeFluxZXPrev) * dxInv * mfsq
333  + (edgeFluxZYNext - edgeFluxZYPrev) * dyInv * mfsq
334  + (centFluxZZNext - centFluxZZPrev) * dzInv;
335 
336  amrex::Real denom = 0.5*(detJ(i,j,k) + detJ(i,j,k-1));
337  advectionSrc /= denom;
338 
339  return advectionSrc;
340 }
341 
342 /**
343  * Wrapper function for computing the advective tendency w/ spatial order > 2.
344  */
345 template<typename InterpType_H, typename InterpType_V, typename WallInterpType>
346 void
347 AdvectionSrcForMomWrapper (const amrex::Box& bxx, const amrex::Box& bxy, const amrex::Box& bxz,
348  const amrex::Array4<amrex::Real>& rho_u_rhs,
349  const amrex::Array4<amrex::Real>& rho_v_rhs,
350  const amrex::Array4<amrex::Real>& rho_w_rhs,
351  const amrex::Array4<const amrex::Real>& rho_u,
352  const amrex::Array4<const amrex::Real>& rho_v,
353  const amrex::Array4<const amrex::Real>& Omega,
354  const amrex::Array4<const amrex::Real>& u,
355  const amrex::Array4<const amrex::Real>& v,
356  const amrex::Array4<const amrex::Real>& w,
357  const amrex::Array4<const amrex::Real>& z_nd,
358  const amrex::Array4<const amrex::Real>& ax,
359  const amrex::Array4<const amrex::Real>& ay,
360  const amrex::Array4<const amrex::Real>& az,
361  const amrex::Array4<const amrex::Real>& detJ,
362  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
363  const amrex::Array4<const amrex::Real>& mf_m,
364  const amrex::Array4<const amrex::Real>& mf_u_inv,
365  const amrex::Array4<const amrex::Real>& mf_v_inv,
366  const amrex::Real upw_frac_h,
367  const amrex::Real upw_frac_v,
368  const AdvType vert_adv_type,
369  const int lo_z_face, const int hi_z_face)
370 {
371  // Instantiate the appropriate structs
372  InterpType_H interp_u_h(u, upw_frac_h); InterpType_V interp_u_v(u, upw_frac_v); // X-MOM
373  InterpType_H interp_v_h(v, upw_frac_h); InterpType_V interp_v_v(v, upw_frac_v); // Y-MOM
374  InterpType_H interp_w_h(w, upw_frac_h); InterpType_V interp_w_v(w, upw_frac_v); // Z-MOM
375  WallInterpType interp_w_wall(w, upw_frac_v); // Z-MOM @ wall
376 
377  amrex::ParallelFor(bxx,
378  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
379  {
380  rho_u_rhs(i, j, k) = -AdvectionSrcForXMom(i, j, k, rho_u, rho_v, Omega, z_nd, ax, ay, az, detJ,
381  interp_u_h, interp_u_v,
382  cellSizeInv, mf_u_inv, mf_v_inv);
383  });
384 
385  amrex::ParallelFor(bxy,
386  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
387  {
388  rho_v_rhs(i, j, k) = -AdvectionSrcForYMom(i, j, k, rho_u, rho_v, Omega, z_nd, ax, ay, az, detJ,
389  interp_v_h, interp_v_v,
390  cellSizeInv, mf_u_inv, mf_v_inv);
391  });
392 
393  amrex::ParallelFor(bxz,
394  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
395  {
396  rho_w_rhs(i, j, k) = -AdvectionSrcForZMom(i, j, k, rho_u, rho_v, Omega, w, z_nd, ax, ay, az, detJ,
397  interp_w_h, interp_w_v, interp_w_wall,
398  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
399  vert_adv_type, lo_z_face, hi_z_face);
400  });
401 }
402 
403 /**
404  * Wrapper function for computing the advective tendency w/ spatial order > 2.
405  */
406 template<typename InterpType_H>
407 void
408 AdvectionSrcForMomVert (const amrex::Box& bxx, const amrex::Box& bxy, const amrex::Box& bxz,
409  const amrex::Array4<amrex::Real>& rho_u_rhs,
410  const amrex::Array4<amrex::Real>& rho_v_rhs,
411  const amrex::Array4<amrex::Real>& rho_w_rhs,
412  const amrex::Array4<const amrex::Real>& rho_u,
413  const amrex::Array4<const amrex::Real>& rho_v,
414  const amrex::Array4<const amrex::Real>& Omega,
415  const amrex::Array4<const amrex::Real>& u,
416  const amrex::Array4<const amrex::Real>& v,
417  const amrex::Array4<const amrex::Real>& w,
418  const amrex::Array4<const amrex::Real>& z_nd,
419  const amrex::Array4<const amrex::Real>& ax,
420  const amrex::Array4<const amrex::Real>& ay,
421  const amrex::Array4<const amrex::Real>& az,
422  const amrex::Array4<const amrex::Real>& detJ,
423  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
424  const amrex::Array4<const amrex::Real>& mf_m,
425  const amrex::Array4<const amrex::Real>& mf_u_inv,
426  const amrex::Array4<const amrex::Real>& mf_v_inv,
427  const amrex::Real upw_frac_h,
428  const amrex::Real upw_frac_v,
429  const AdvType vert_adv_type,
430  const int lo_z_face, const int hi_z_face)
431 {
432  if (vert_adv_type == AdvType::Centered_2nd) {
433  AdvectionSrcForMomWrapper<InterpType_H,CENTERED2,UPWINDALL>(bxx, bxy, bxz,
434  rho_u_rhs, rho_v_rhs, rho_w_rhs,
435  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
436  cellSizeInv, mf_m,
437  mf_u_inv, mf_v_inv,
438  upw_frac_h, upw_frac_v,
439  vert_adv_type, lo_z_face, hi_z_face);
440  } else if (vert_adv_type == AdvType::Upwind_3rd) {
441  AdvectionSrcForMomWrapper<InterpType_H,UPWIND3,UPWINDALL>(bxx, bxy, bxz,
442  rho_u_rhs, rho_v_rhs, rho_w_rhs,
443  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
444  cellSizeInv, mf_m,
445  mf_u_inv, mf_v_inv,
446  upw_frac_h, upw_frac_v,
447  vert_adv_type, lo_z_face, hi_z_face);
448  } else if (vert_adv_type == AdvType::Centered_4th) {
449  AdvectionSrcForMomWrapper<InterpType_H,CENTERED4,UPWINDALL>(bxx, bxy, bxz,
450  rho_u_rhs, rho_v_rhs, rho_w_rhs,
451  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
452  cellSizeInv, mf_m,
453  mf_u_inv, mf_v_inv,
454  upw_frac_h, upw_frac_v,
455  vert_adv_type, lo_z_face, hi_z_face);
456  } else if (vert_adv_type == AdvType::Upwind_5th) {
457  AdvectionSrcForMomWrapper<InterpType_H,UPWIND5,UPWINDALL>(bxx, bxy, bxz,
458  rho_u_rhs, rho_v_rhs, rho_w_rhs,
459  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
460  cellSizeInv, mf_m,
461  mf_u_inv, mf_v_inv,
462  upw_frac_h, upw_frac_v,
463  vert_adv_type, lo_z_face, hi_z_face);
464  } else if (vert_adv_type == AdvType::Centered_6th) {
465  AdvectionSrcForMomWrapper<InterpType_H,CENTERED6,UPWINDALL>(bxx, bxy, bxz,
466  rho_u_rhs, rho_v_rhs, rho_w_rhs,
467  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
468  cellSizeInv, mf_m,
469  mf_u_inv, mf_v_inv,
470  upw_frac_h, upw_frac_v,
471  vert_adv_type, lo_z_face, hi_z_face);
472  } else {
473  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
474  }
475 }
AdvectionSrcForZMom
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real AdvectionSrcForZMom(int i, int j, int k, const amrex::Array4< const amrex::Real > &rho_u, const amrex::Array4< const amrex::Real > &rho_v, const amrex::Array4< const amrex::Real > &Omega, const amrex::Array4< const amrex::Real > &w, const amrex::Array4< const amrex::Real > &z_nd, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, InterpType_H interp_omega_h, InterpType_V interp_omega_v, WallInterpType interp_omega_wall, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &mf_m, const amrex::Array4< const amrex::Real > &mf_u_inv, const amrex::Array4< const amrex::Real > &mf_v_inv, const AdvType vert_adv_type, const int lo_z_face, const int hi_z_face)
Definition: ERF_AdvectionSrcForMom_T.H:218
AdvectionSrcForMomWrapper
void AdvectionSrcForMomWrapper(const amrex::Box &bxx, const amrex::Box &bxy, const amrex::Box &bxz, const amrex::Array4< amrex::Real > &rho_u_rhs, const amrex::Array4< amrex::Real > &rho_v_rhs, const amrex::Array4< amrex::Real > &rho_w_rhs, const amrex::Array4< const amrex::Real > &rho_u, const amrex::Array4< const amrex::Real > &rho_v, const amrex::Array4< const amrex::Real > &Omega, const amrex::Array4< const amrex::Real > &u, const amrex::Array4< const amrex::Real > &v, const amrex::Array4< const amrex::Real > &w, const amrex::Array4< const amrex::Real > &z_nd, const amrex::Array4< const amrex::Real > &ax, const amrex::Array4< const amrex::Real > &ay, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &mf_m, const amrex::Array4< const amrex::Real > &mf_u_inv, const amrex::Array4< const amrex::Real > &mf_v_inv, const amrex::Real upw_frac_h, const amrex::Real upw_frac_v, const AdvType vert_adv_type, const int lo_z_face, const int hi_z_face)
Definition: ERF_AdvectionSrcForMom_T.H:347
AdvectionSrcForYMom
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real AdvectionSrcForYMom(int i, int j, int k, const amrex::Array4< const amrex::Real > &rho_u, const amrex::Array4< const amrex::Real > &rho_v, const amrex::Array4< const amrex::Real > &Omega, const amrex::Array4< const amrex::Real > &z_nd, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &ay, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, InterpType_H interp_v_h, InterpType_V interp_v_v, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &mf_u_inv, const amrex::Array4< const amrex::Real > &mf_v_inv)
Definition: ERF_AdvectionSrcForMom_T.H:120
AdvectionSrcForXMom
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real AdvectionSrcForXMom(int i, int j, int k, const amrex::Array4< const amrex::Real > &rho_u, const amrex::Array4< const amrex::Real > &rho_v, const amrex::Array4< const amrex::Real > &Omega, const amrex::Array4< const amrex::Real > &z_nd, const amrex::Array4< const amrex::Real > &ax, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, InterpType_H interp_u_h, InterpType_V interp_u_v, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &mf_u_inv, const amrex::Array4< const amrex::Real > &mf_v_inv)
Definition: ERF_AdvectionSrcForMom_T.H:26
AdvectionSrcForMomVert
void AdvectionSrcForMomVert(const amrex::Box &bxx, const amrex::Box &bxy, const amrex::Box &bxz, const amrex::Array4< amrex::Real > &rho_u_rhs, const amrex::Array4< amrex::Real > &rho_v_rhs, const amrex::Array4< amrex::Real > &rho_w_rhs, const amrex::Array4< const amrex::Real > &rho_u, const amrex::Array4< const amrex::Real > &rho_v, const amrex::Array4< const amrex::Real > &Omega, const amrex::Array4< const amrex::Real > &u, const amrex::Array4< const amrex::Real > &v, const amrex::Array4< const amrex::Real > &w, const amrex::Array4< const amrex::Real > &z_nd, const amrex::Array4< const amrex::Real > &ax, const amrex::Array4< const amrex::Real > &ay, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Array4< const amrex::Real > &mf_m, const amrex::Array4< const amrex::Real > &mf_u_inv, const amrex::Array4< const amrex::Real > &mf_v_inv, const amrex::Real upw_frac_h, const amrex::Real upw_frac_v, const AdvType vert_adv_type, const int lo_z_face, const int hi_z_face)
Definition: ERF_AdvectionSrcForMom_T.H:408
ERF_IndexDefines.H
AdvType
AdvType
Definition: ERF_IndexDefines.H:203
AdvType::Upwind_3rd
@ Upwind_3rd
AdvType::Centered_4th
@ Centered_4th
AdvType::Upwind_5th
@ Upwind_5th
AdvType::Centered_6th
@ Centered_6th
AdvType::Centered_2nd
@ Centered_2nd
ERF_Interpolation.H
ERF_TerrainMetrics.H
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:273
Compute_h_zeta_AtEdgeCenterK
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Compute_h_zeta_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:228
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:317