ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_AdvectionSrcForMom_N.H
Go to the documentation of this file.
1 #include "ERF_IndexDefines.H"
2 #include "ERF_Interpolation.H"
3 
4 /**
5  * Function for computing the advective tendency for the x-component of momentum
6  * without metric terms (constant OR stretched dz) and for higher-order stencils
7  *
8  * @tparam InterpType_H horizontal interpolation type
9  * @tparam InterpType_V vertical interpolation type
10  * @param[in] i i-index of x-face at which to create tendency
11  * @param[in] j j-index of x-face at which to create tendency
12  * @param[in] k k-index of x-face at which to create tendency
13  * @param[in] rho_u x-component of momentum
14  * @param[in] rho_v y-component of momentum
15  * @param[in] rho_w z-component of momentum
16  * @param[in] interp_u_h horizontal interpolator for x-component velocity
17  * @param[in] interp_u_v vertical interpolator for x-component velocity
18  * @param[in] mf_ux_inv inverse x map factor on x-faces
19  * @param[in] mf_uy_inv inverse y map factor on x-faces
20  * @param[in] mf_vx_inv inverse x map factor on y-faces
21  * @param[in] dxInv inverse grid spacing in x
22  * @param[in] dyInv inverse grid spacing in y
23  * @param[in] dzInv inverse grid spacing in z
24  * @return advective tendency for x-momentum
25  */
26 template<typename InterpType_H, typename InterpType_V>
27 AMREX_GPU_DEVICE
28 AMREX_FORCE_INLINE
29 amrex::Real
30 AdvectionSrcForXMom_N (int i, int j, int k,
31  const amrex::Array4<const amrex::Real>& rho_u,
32  const amrex::Array4<const amrex::Real>& rho_v,
33  const amrex::Array4<const amrex::Real>& rho_w,
34  InterpType_H interp_u_h,
35  InterpType_V interp_u_v,
36  const amrex::Array4<const amrex::Real>& mf_ux_inv,
37  const amrex::Array4<const amrex::Real>& mf_uy_inv,
38  const amrex::Array4<const amrex::Real>& mf_vx_inv,
39  const amrex::Real dxInv, const amrex::Real dyInv, const amrex::Real dzInv)
40 {
41  amrex::Real advectionSrc;
42 
43  amrex::Real rho_u_avg_lo, rho_u_avg_hi;
44  amrex::Real rho_v_avg_lo, rho_v_avg_hi;
45  amrex::Real rho_w_avg_lo, rho_w_avg_hi;
46 
47  amrex::Real xflux_hi; amrex::Real xflux_lo;
48  amrex::Real yflux_hi; amrex::Real yflux_lo;
49  amrex::Real zflux_hi; amrex::Real zflux_lo;
50 
51  amrex::Real interp_hi(zero), interp_lo(zero);
52 
53  rho_u_avg_hi = myhalf * (rho_u(i+1, j, k) * mf_uy_inv(i+1,j,0) + rho_u(i, j, k) * mf_uy_inv(i,j,0));
54  interp_u_h.InterpolateInX(i+1,j,k,0,interp_hi,rho_u_avg_hi);
55  xflux_hi = rho_u_avg_hi * interp_hi;
56 
57  rho_u_avg_lo = myhalf * (rho_u(i-1, j, k) * mf_uy_inv(i-1,j,0) + rho_u(i, j, k) * mf_uy_inv(i,j,0));
58  interp_u_h.InterpolateInX(i,j,k,0,interp_lo,rho_u_avg_lo);
59  xflux_lo = rho_u_avg_lo * interp_lo;
60 
61  rho_v_avg_hi = myhalf * (rho_v(i, j+1, k) * mf_vx_inv(i,j+1,0) + rho_v(i-1, j+1, k) * mf_vx_inv(i-1,j+1,0));
62  interp_u_h.InterpolateInY(i,j+1,k,0,interp_hi,rho_v_avg_hi);
63  yflux_hi = rho_v_avg_hi * interp_hi;
64 
65  rho_v_avg_lo = myhalf * (rho_v(i, j , k) * mf_vx_inv(i,j ,0) + rho_v(i-1, j , k) * mf_vx_inv(i-1,j ,0));
66  interp_u_h.InterpolateInY(i,j,k,0,interp_lo,rho_v_avg_lo);
67  yflux_lo = rho_v_avg_lo * interp_lo;
68 
69  rho_w_avg_hi = myhalf * (rho_w(i, j, k+1) + rho_w(i-1, j, k+1));
70  interp_u_v.InterpolateInZ(i,j,k+1,0,interp_hi,rho_w_avg_hi);
71  zflux_hi = rho_w_avg_hi * interp_hi;
72 
73  rho_w_avg_lo = myhalf * (rho_w(i, j, k ) + rho_w(i-1, j, k ));
74  interp_u_v.InterpolateInZ(i,j,k,0,interp_lo,rho_w_avg_lo);
75  zflux_lo = rho_w_avg_lo * interp_lo;
76 
77  amrex::Real mfsq = 1 / (mf_ux_inv(i,j,0) * mf_uy_inv(i,j,0));
78 
79  advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
80  + (yflux_hi - yflux_lo) * dyInv * mfsq
81  + (zflux_hi - zflux_lo) * dzInv;
82 
83  return advectionSrc;
84 }
85 
86 /**
87  * Function for computing the advective tendency for the y-component of momentum
88  * without metric terms and for higher-order stencils
89  *
90  * @tparam InterpType_H horizontal interpolation type
91  * @tparam InterpType_V vertical interpolation type
92  * @param[in] i i-index of y-face at which to create tendency
93  * @param[in] j j-index of y-face at which to create tendency
94  * @param[in] k k-index of y-face at which to create tendency
95  * @param[in] rho_u x-component of momentum
96  * @param[in] rho_v y-component of momentum
97  * @param[in] rho_w z-component of momentum
98  * @param[in] interp_v_h horizontal interpolator for y-component velocity
99  * @param[in] interp_v_v vertical interpolator for y-component velocity
100  * @param[in] mf_uy_inv inverse y map factor on x-faces
101  * @param[in] mf_vx_inv inverse x map factor on y-faces
102  * @param[in] mf_vy_inv inverse y map factor on y-faces
103  * @param[in] dxInv inverse grid spacing in x
104  * @param[in] dyInv inverse grid spacing in y
105  * @param[in] dzInv inverse grid spacing in z
106  * @return advective tendency for y-momentum
107  */
108 template<typename InterpType_H, typename InterpType_V>
109 AMREX_GPU_DEVICE
110 AMREX_FORCE_INLINE
111 amrex::Real
112 AdvectionSrcForYMom_N (int i, int j, int k,
113  const amrex::Array4<const amrex::Real>& rho_u,
114  const amrex::Array4<const amrex::Real>& rho_v,
115  const amrex::Array4<const amrex::Real>& rho_w,
116  InterpType_H interp_v_h,
117  InterpType_V interp_v_v,
118  const amrex::Array4<const amrex::Real>& mf_uy_inv,
119  const amrex::Array4<const amrex::Real>& mf_vx_inv,
120  const amrex::Array4<const amrex::Real>& mf_vy_inv,
121  const amrex::Real dxInv, const amrex::Real dyInv, const amrex::Real dzInv)
122 {
123  amrex::Real advectionSrc;
124 
125  amrex::Real rho_u_avg_lo, rho_u_avg_hi;
126  amrex::Real rho_v_avg_lo, rho_v_avg_hi;
127  amrex::Real rho_w_avg_lo, rho_w_avg_hi;
128 
129  amrex::Real xflux_hi; amrex::Real xflux_lo;
130  amrex::Real yflux_hi; amrex::Real yflux_lo;
131  amrex::Real zflux_hi; amrex::Real zflux_lo;
132 
133  amrex::Real interp_hi(zero), interp_lo(zero);
134 
135  rho_u_avg_hi = myhalf * (rho_u(i+1, j, k) * mf_uy_inv(i+1,j ,0) + rho_u(i+1, j-1, k) * mf_uy_inv(i+1,j-1,0));
136  interp_v_h.InterpolateInX(i+1,j,k,0,interp_hi,rho_u_avg_hi);
137  xflux_hi = rho_u_avg_hi * interp_hi;
138 
139  rho_u_avg_lo = myhalf * (rho_u(i , j, k) * mf_uy_inv(i ,j ,0) + rho_u(i , j-1, k) * mf_uy_inv(i ,j-1,0));
140  interp_v_h.InterpolateInX(i,j,k,0,interp_lo,rho_u_avg_lo);
141  xflux_lo = rho_u_avg_lo * interp_lo;
142 
143  rho_v_avg_hi = myhalf * (rho_v(i, j, k) * mf_vx_inv(i ,j ,0) + rho_v(i, j+1, k) * mf_vx_inv(i ,j+1,0));
144  interp_v_h.InterpolateInY(i,j+1,k,0,interp_hi,rho_v_avg_hi);
145  yflux_hi = rho_v_avg_hi * interp_hi;
146 
147  rho_v_avg_lo = myhalf * (rho_v(i, j, k) * mf_vx_inv(i ,j ,0) + rho_v(i, j-1, k) * mf_vx_inv(i ,j-1,0));
148  interp_v_h.InterpolateInY(i,j,k,0,interp_lo,rho_v_avg_lo);
149  yflux_lo = rho_v_avg_lo * interp_lo;
150 
151  rho_w_avg_hi = myhalf * (rho_w(i, j, k+1) + rho_w(i, j-1, k+1));
152  interp_v_v.InterpolateInZ(i,j,k+1,0,interp_hi,rho_w_avg_hi);
153  zflux_hi = rho_w_avg_hi * interp_hi;
154 
155  rho_w_avg_lo = myhalf * (rho_w(i, j, k ) + rho_w(i, j-1, k ));
156  interp_v_v.InterpolateInZ(i,j,k ,0,interp_lo,rho_w_avg_lo);
157  zflux_lo = rho_w_avg_lo * interp_lo;
158 
159  amrex::Real mfsq = 1 / (mf_vx_inv(i,j,0) * mf_vy_inv(i,j,0));
160 
161  advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
162  + (yflux_hi - yflux_lo) * dyInv * mfsq
163  + (zflux_hi - zflux_lo) * dzInv;
164 
165  return advectionSrc;
166 }
167 
168 /**
169  * Function for computing the advective tendency for the z-component of momentum
170  * without metric terms and for higher-order stencils
171  *
172  * @tparam InterpType_H horizontal interpolation type
173  * @tparam InterpType_V vertical interpolation type
174  * @tparam WallInterpType wall-aware vertical interpolation type
175  * @param[in] i i-index of z-face at which to create tendency
176  * @param[in] j j-index of z-face at which to create tendency
177  * @param[in] k k-index of z-face at which to create tendency
178  * @param[in] rho_u x-component of momentum
179  * @param[in] rho_v y-component of momentum
180  * @param[in] rho_w z-component of momentum
181  * @param[in] w z-component of velocity
182  * @param[in] interp_w_h horizontal interpolator for z-component velocity
183  * @param[in] interp_w_v vertical interpolator for z-component velocity
184  * @param[in] interp_w_wall wall-aware vertical interpolator for z-component velocity
185  * @param[in] mf_mx x map factor on cell centers
186  * @param[in] mf_my y map factor on cell centers
187  * @param[in] mf_uy_inv inverse y map factor on x-faces
188  * @param[in] mf_vx_inv inverse x map factor on y-faces
189  * @param[in] vert_adv_type vertical advection stencil
190  * @param[in] lo_z_face minimum z-face k-index at this level
191  * @param[in] hi_z_face maximum z-face k-index at this level
192  * @param[in] dxInv inverse grid spacing in x
193  * @param[in] dyInv inverse grid spacing in y
194  * @param[in] dzInv inverse grid spacing in z
195  * @return advective tendency for z-momentum
196  */
197 template<typename InterpType_H, typename InterpType_V, typename WallInterpType>
198 AMREX_GPU_DEVICE
199 AMREX_FORCE_INLINE
200 amrex::Real
201 AdvectionSrcForZMom_N (int i, int j, int k,
202  const amrex::Array4<const amrex::Real>& rho_u,
203  const amrex::Array4<const amrex::Real>& rho_v,
204  const amrex::Array4<const amrex::Real>& rho_w,
205  const amrex::Array4<const amrex::Real>& w,
206  InterpType_H interp_w_h,
207  InterpType_V interp_w_v,
208  WallInterpType interp_w_wall,
209  const amrex::Array4<const amrex::Real>& mf_mx,
210  const amrex::Array4<const amrex::Real>& mf_my,
211  const amrex::Array4<const amrex::Real>& mf_uy_inv,
212  const amrex::Array4<const amrex::Real>& mf_vx_inv,
213  const AdvType vert_adv_type,
214  const int lo_z_face, const int hi_z_face,
215  const amrex::Real dxInv, const amrex::Real dyInv, const amrex::Real dzInv)
216 {
217 
218  amrex::Real advectionSrc;
219 
220  amrex::Real rho_u_avg_lo, rho_u_avg_hi;
221  amrex::Real rho_v_avg_lo, rho_v_avg_hi;
222  amrex::Real rho_w_avg_lo, rho_w_avg_hi;
223 
224  amrex::Real xflux_hi; amrex::Real xflux_lo;
225  amrex::Real yflux_hi; amrex::Real yflux_lo;
226  amrex::Real zflux_hi; amrex::Real zflux_lo;
227 
228  amrex::Real interp_hi(zero), interp_lo(zero);
229 
230  rho_u_avg_hi = myhalf * (rho_u(i+1, j, k) + rho_u(i+1, j, k-1)) * mf_uy_inv(i+1,j ,0);
231  interp_w_h.InterpolateInX(i+1,j,k,0,interp_hi,rho_u_avg_hi);
232  xflux_hi = rho_u_avg_hi * interp_hi;
233 
234  rho_u_avg_lo = myhalf * (rho_u(i , j, k) + rho_u(i , j, k-1)) * mf_uy_inv(i ,j ,0);
235  interp_w_h.InterpolateInX(i,j,k,0,interp_lo,rho_u_avg_lo);
236  xflux_lo = rho_u_avg_lo * interp_lo;
237 
238  rho_v_avg_hi = myhalf * (rho_v(i, j+1, k) + rho_v(i, j+1, k-1)) * mf_vx_inv(i ,j+1,0);
239  interp_w_h.InterpolateInY(i,j+1,k,0,interp_hi,rho_v_avg_hi);
240  yflux_hi = rho_v_avg_hi * interp_hi;
241 
242  rho_v_avg_lo = myhalf * (rho_v(i, j , k) + rho_v(i, j , k-1)) * mf_vx_inv(i ,j ,0);
243  interp_w_h.InterpolateInY(i,j,k,0,interp_lo,rho_v_avg_lo);
244  yflux_lo = rho_v_avg_lo * interp_lo;
245 
246  // int l_spatial_order_hi = std::min(std::min(vert_spatial_order, 2*(hi_z_face-k)), 2*(k+1));
247  // If k == hi_z_face-1, l_spatial_order_hi = 2
248  // If k == hi_z_face-2, l_spatial_order_hi = std::min(vert_spatial_order, 4);
249  // If k == lo_z+face+1, l_spatial_order_hi = std::min(vert_spatial_order, 4);
250 
251  if (k == hi_z_face) {
252  zflux_hi = rho_w(i,j,k) * w(i,j,k);
253  } else {
254  rho_w_avg_hi = myhalf * (rho_w(i,j,k) + rho_w(i,j,k+1));
255  if (k == hi_z_face-1) {
256  interp_w_wall.InterpolateInZ(i,j,k+1,0,interp_hi,rho_w_avg_hi,AdvType::Centered_2nd);
257  } else if (k == hi_z_face-2 || k == lo_z_face+1) {
258  if (vert_adv_type != AdvType::Centered_2nd && vert_adv_type != AdvType::Upwind_3rd) {
259  interp_w_wall.InterpolateInZ(i,j,k+1,0,interp_hi,rho_w_avg_hi,AdvType::Centered_4th);
260  } else {
261  interp_w_wall.InterpolateInZ(i,j,k+1,0,interp_hi,rho_w_avg_hi,vert_adv_type);
262  }
263  } else {
264  interp_w_v.InterpolateInZ(i,j,k+1,0,interp_hi,rho_w_avg_hi);
265  }
266  zflux_hi = rho_w_avg_hi * interp_hi;
267  }
268 
269  // int l_spatial_order_lo = std::min(std::min(vert_spatial_order, 2*(hi_z_face+1-k)), 2*k);
270  // If k == lo_z_face+1, l_spatial_order_hi = 2
271  // If k == lo_z_face+2, l_spatial_order_hi = std::min(vert_spatial_order, 4);
272  // If k == hi_z_face-1, l_spatial_order_hi = std::min(vert_spatial_order, 4);
273 
274  if (k == lo_z_face) {
275  zflux_lo = rho_w(i,j,k) * w(i,j,k);
276  } else {
277  rho_w_avg_lo = myhalf * (rho_w(i,j,k) + rho_w(i,j,k-1));
278  if (k == lo_z_face+1) {
279  interp_w_wall.InterpolateInZ(i,j,k,0,interp_lo,rho_w_avg_lo,AdvType::Centered_2nd);
280  } else if (k == lo_z_face+2 || k == hi_z_face-1) {
281  if (vert_adv_type != AdvType::Centered_2nd && vert_adv_type != AdvType::Upwind_3rd) {
282  interp_w_wall.InterpolateInZ(i,j,k,0,interp_lo,rho_w_avg_lo,AdvType::Centered_4th);
283  } else {
284  interp_w_wall.InterpolateInZ(i,j,k,0,interp_lo,rho_w_avg_lo,vert_adv_type);
285  }
286  } else {
287  interp_w_v.InterpolateInZ(i,j,k,0,interp_lo,rho_w_avg_lo);
288  }
289  zflux_lo = rho_w_avg_lo * interp_lo;
290  }
291 
292  amrex::Real mfsq = mf_mx(i,j,0) * mf_my(i,j,0);
293 
294  advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
295  + (yflux_hi - yflux_lo) * dyInv * mfsq
296  + (zflux_hi - zflux_lo) * dzInv;
297 
298  return advectionSrc;
299 }
300 
301 /**
302  * Wrapper function for computing the advective tendency w/ spatial order > two
303  *
304  * @tparam InterpType_H horizontal interpolation type
305  * @tparam InterpType_V vertical interpolation type
306  * @tparam WallInterpType wall-aware vertical interpolation type
307  * @param[in] bxx box over which x-momentum is updated
308  * @param[in] bxy box over which y-momentum is updated
309  * @param[in] bxz box over which z-momentum is updated
310  * @param[out] rho_u_rhs tendency for the x-momentum equation
311  * @param[out] rho_v_rhs tendency for the y-momentum equation
312  * @param[out] rho_w_rhs tendency for the z-momentum equation
313  * @param[in] rho_u x-component of momentum
314  * @param[in] rho_v y-component of momentum
315  * @param[in] rho_w z-component of momentum
316  * @param[in] u x-component of velocity
317  * @param[in] v y-component of velocity
318  * @param[in] w z-component of velocity
319  * @param[in] cellSizeInv inverse grid spacing
320  * @param[in] stretched_dz_d device vector of vertical cell spacing
321  * @param[in] mf_mx x map factor on cell centers
322  * @param[in] mf_ux_inv inverse x map factor on x-faces
323  * @param[in] mf_vx_inv inverse x map factor on y-faces
324  * @param[in] mf_my y map factor on cell centers
325  * @param[in] mf_uy_inv inverse y map factor on x-faces
326  * @param[in] mf_vy_inv inverse y map factor on y-faces
327  * @param[in] upw_frac_h horizontal upwind blending fraction
328  * @param[in] upw_frac_v vertical upwind blending fraction
329  * @param[in] vert_adv_type vertical advection stencil
330  * @param[in] lo_z_face minimum z-face k-index at this level
331  * @param[in] hi_z_face maximum z-face k-index at this level
332  */
333 template<typename InterpType_H, typename InterpType_V, typename WallInterpType>
334 void
335 AdvectionSrcForMomWrapper_N (const amrex::Box& bxx, const amrex::Box& bxy, const amrex::Box& bxz,
336  const amrex::Array4<amrex::Real>& rho_u_rhs,
337  const amrex::Array4<amrex::Real>& rho_v_rhs,
338  const amrex::Array4<amrex::Real>& rho_w_rhs,
339  const amrex::Array4<const amrex::Real>& rho_u,
340  const amrex::Array4<const amrex::Real>& rho_v,
341  const amrex::Array4<const amrex::Real>& rho_w,
342  const amrex::Array4<const amrex::Real>& u,
343  const amrex::Array4<const amrex::Real>& v,
344  const amrex::Array4<const amrex::Real>& w,
345  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
346  const amrex::Gpu::DeviceVector<amrex::Real>& stretched_dz_d,
347  const amrex::Array4<const amrex::Real>& mf_mx,
348  const amrex::Array4<const amrex::Real>& mf_ux_inv,
349  const amrex::Array4<const amrex::Real>& mf_vx_inv,
350  const amrex::Array4<const amrex::Real>& mf_my,
351  const amrex::Array4<const amrex::Real>& mf_uy_inv,
352  const amrex::Array4<const amrex::Real>& mf_vy_inv,
353  const amrex::Real upw_frac_h,
354  const amrex::Real upw_frac_v,
355  const AdvType vert_adv_type,
356  const int lo_z_face, const int hi_z_face)
357 {
358  // Instantiate the appropriate structs
359  InterpType_H interp_u_h(u, upw_frac_h); InterpType_V interp_u_v(u, upw_frac_v); // X-MOM
360  InterpType_H interp_v_h(v, upw_frac_h); InterpType_V interp_v_v(v, upw_frac_v); // Y-MOM
361  InterpType_H interp_w_h(w, upw_frac_h); InterpType_V interp_w_v(w, upw_frac_v); // Z-MOM
362  WallInterpType interp_w_wall(w, upw_frac_v); // Z-MOM @ wall
363 
364  auto dxInv = cellSizeInv[0];
365  auto dyInv = cellSizeInv[1];
366 
367  auto dz_ptr = stretched_dz_d.data();
368 
369  amrex::ParallelFor(bxx,
370  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
371  {
372  auto dzInv = one/dz_ptr[k];
373  rho_u_rhs(i, j, k) = -AdvectionSrcForXMom_N(i, j, k, rho_u, rho_v, rho_w,
374  interp_u_h, interp_u_v,
375  mf_ux_inv, mf_uy_inv, mf_vx_inv,
376  dxInv, dyInv, dzInv);
377  });
378 
379  amrex::ParallelFor(bxy,
380  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
381  {
382  auto dzInv = one/dz_ptr[k];
383  rho_v_rhs(i, j, k) = -AdvectionSrcForYMom_N(i, j, k, rho_u, rho_v, rho_w,
384  interp_v_h, interp_v_v,
385  mf_uy_inv, mf_vx_inv, mf_vy_inv,
386  dxInv, dyInv, dzInv);
387  });
388 
389  amrex::ParallelFor(bxz,
390  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
391  {
392  auto dzInv = (k == 0) ? one / dz_ptr[k] : two/(dz_ptr[k] + dz_ptr[k-1]);
393  rho_w_rhs(i, j, k) = -AdvectionSrcForZMom_N(i, j, k, rho_u, rho_v, rho_w, w,
394  interp_w_h, interp_w_v, interp_w_wall,
395  mf_mx, mf_my, mf_uy_inv, mf_vx_inv,
396  vert_adv_type, lo_z_face, hi_z_face,
397  dxInv, dyInv, dzInv);
398  });
399 }
400 
401 /**
402  * Wrapper function for computing the advective tendency w/ spatial order > two
403  *
404  * @tparam InterpType_H horizontal interpolation type
405  * @param[in] bxx box over which x-momentum is updated
406  * @param[in] bxy box over which y-momentum is updated
407  * @param[in] bxz box over which z-momentum is updated
408  * @param[out] rho_u_rhs tendency for the x-momentum equation
409  * @param[out] rho_v_rhs tendency for the y-momentum equation
410  * @param[out] rho_w_rhs tendency for the z-momentum equation
411  * @param[in] rho_u x-component of momentum
412  * @param[in] rho_v y-component of momentum
413  * @param[in] rho_w z-component of momentum
414  * @param[in] u x-component of velocity
415  * @param[in] v y-component of velocity
416  * @param[in] w z-component of velocity
417  * @param[in] cellSizeInv inverse grid spacing
418  * @param[in] stretched_dz_d device vector of vertical cell spacing
419  * @param[in] mf_mx x map factor on cell centers
420  * @param[in] mf_ux_inv inverse x map factor on x-faces
421  * @param[in] mf_vx_inv inverse x map factor on y-faces
422  * @param[in] mf_my y map factor on cell centers
423  * @param[in] mf_uy_inv inverse y map factor on x-faces
424  * @param[in] mf_vy_inv inverse y map factor on y-faces
425  * @param[in] upw_frac_h horizontal upwind blending fraction
426  * @param[in] upw_frac_v vertical upwind blending fraction
427  * @param[in] vert_adv_type vertical advection stencil
428  * @param[in] lo_z_face minimum z-face k-index at this level
429  * @param[in] hi_z_face maximum z-face k-index at this level
430  */
431 template<typename InterpType_H>
432 void
433 AdvectionSrcForMomVert_N (const amrex::Box& bxx, const amrex::Box& bxy, const amrex::Box& bxz,
434  const amrex::Array4<amrex::Real>& rho_u_rhs,
435  const amrex::Array4<amrex::Real>& rho_v_rhs,
436  const amrex::Array4<amrex::Real>& rho_w_rhs,
437  const amrex::Array4<const amrex::Real>& rho_u,
438  const amrex::Array4<const amrex::Real>& rho_v,
439  const amrex::Array4<const amrex::Real>& rho_w,
440  const amrex::Array4<const amrex::Real>& u,
441  const amrex::Array4<const amrex::Real>& v,
442  const amrex::Array4<const amrex::Real>& w,
443  const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>& cellSizeInv,
444  const amrex::Gpu::DeviceVector<amrex::Real>& stretched_dz_d,
445  const amrex::Array4<const amrex::Real>& mf_mx,
446  const amrex::Array4<const amrex::Real>& mf_ux_inv,
447  const amrex::Array4<const amrex::Real>& mf_vx_inv,
448  const amrex::Array4<const amrex::Real>& mf_my,
449  const amrex::Array4<const amrex::Real>& mf_uy_inv,
450  const amrex::Array4<const amrex::Real>& mf_vy_inv,
451  const amrex::Real upw_frac_h,
452  const amrex::Real upw_frac_v,
453  const AdvType vert_adv_type,
454  const int lo_z_face, const int hi_z_face)
455 {
456  if (vert_adv_type == AdvType::Centered_2nd) {
457  AdvectionSrcForMomWrapper_N<InterpType_H,CENTERED2,UPWINDALL>(bxx, bxy, bxz,
458  rho_u_rhs, rho_v_rhs, rho_w_rhs,
459  rho_u, rho_v, rho_w, u, v, w,
460  cellSizeInv, stretched_dz_d,
461  mf_mx, mf_ux_inv, mf_vx_inv,
462  mf_my, mf_uy_inv, mf_vy_inv,
463  upw_frac_h, upw_frac_v,
464  vert_adv_type,
465  lo_z_face, hi_z_face);
466  } else if (vert_adv_type == AdvType::Upwind_3rd) {
467  AdvectionSrcForMomWrapper_N<InterpType_H,UPWIND3,UPWINDALL>(bxx, bxy, bxz,
468  rho_u_rhs, rho_v_rhs, rho_w_rhs,
469  rho_u, rho_v, rho_w, u, v, w,
470  cellSizeInv, stretched_dz_d,
471  mf_mx, mf_ux_inv, mf_vx_inv,
472  mf_my, mf_uy_inv, mf_vy_inv,
473  upw_frac_h, upw_frac_v,
474  vert_adv_type,
475  lo_z_face, hi_z_face);
476  } else if (vert_adv_type == AdvType::Centered_4th) {
477  AdvectionSrcForMomWrapper_N<InterpType_H,CENTERED4,UPWINDALL>(bxx, bxy, bxz,
478  rho_u_rhs, rho_v_rhs, rho_w_rhs,
479  rho_u, rho_v, rho_w, u, v, w,
480  cellSizeInv, stretched_dz_d,
481  mf_mx, mf_ux_inv, mf_vx_inv,
482  mf_my, mf_uy_inv, mf_vy_inv,
483  upw_frac_h, upw_frac_v,
484  vert_adv_type,
485  lo_z_face, hi_z_face);
486  } else if (vert_adv_type == AdvType::Upwind_5th) {
487  AdvectionSrcForMomWrapper_N<InterpType_H,UPWIND5,UPWINDALL>(bxx, bxy, bxz,
488  rho_u_rhs, rho_v_rhs, rho_w_rhs,
489  rho_u, rho_v, rho_w, u, v, w,
490  cellSizeInv, stretched_dz_d,
491  mf_mx, mf_ux_inv, mf_vx_inv,
492  mf_my, mf_uy_inv, mf_vy_inv,
493  upw_frac_h, upw_frac_v,
494  vert_adv_type,
495  lo_z_face, hi_z_face);
496  } else if (vert_adv_type == AdvType::Centered_6th) {
497  AdvectionSrcForMomWrapper_N<InterpType_H,CENTERED6,UPWINDALL>(bxx, bxy, bxz,
498  rho_u_rhs, rho_v_rhs, rho_w_rhs,
499  rho_u, rho_v, rho_w, u, v, w,
500  cellSizeInv, stretched_dz_d,
501  mf_mx, mf_ux_inv, mf_vx_inv,
502  mf_my, mf_uy_inv, mf_vy_inv,
503  upw_frac_h, upw_frac_v,
504  vert_adv_type,
505  lo_z_face, hi_z_face);
506  } else if (vert_adv_type == AdvType::Weno_3) {
507  AdvectionSrcForMomWrapper_N<InterpType_H,WENO3,UPWINDALL>(bxx, bxy, bxz,
508  rho_u_rhs, rho_v_rhs, rho_w_rhs,
509  rho_u, rho_v, rho_w, u, v, w,
510  cellSizeInv, stretched_dz_d,
511  mf_mx, mf_ux_inv, mf_vx_inv,
512  mf_my, mf_uy_inv, mf_vy_inv,
513  upw_frac_h, upw_frac_v,
514  vert_adv_type,
515  lo_z_face, hi_z_face);
516  } else if (vert_adv_type == AdvType::Weno_3Z) {
517  AdvectionSrcForMomWrapper_N<InterpType_H,WENO_Z3,UPWINDALL>(bxx, bxy, bxz,
518  rho_u_rhs, rho_v_rhs, rho_w_rhs,
519  rho_u, rho_v, rho_w, u, v, w,
520  cellSizeInv, stretched_dz_d,
521  mf_mx, mf_ux_inv, mf_vx_inv,
522  mf_my, mf_uy_inv, mf_vy_inv,
523  upw_frac_h, upw_frac_v,
524  vert_adv_type,
525  lo_z_face, hi_z_face);
526  } else if (vert_adv_type == AdvType::Weno_3MZQ) {
527  AdvectionSrcForMomWrapper_N<InterpType_H,WENO_MZQ3,UPWINDALL>(bxx, bxy, bxz,
528  rho_u_rhs, rho_v_rhs, rho_w_rhs,
529  rho_u, rho_v, rho_w, u, v, w,
530  cellSizeInv, stretched_dz_d,
531  mf_mx, mf_ux_inv, mf_vx_inv,
532  mf_my, mf_uy_inv, mf_vy_inv,
533  upw_frac_h, upw_frac_v,
534  vert_adv_type,
535  lo_z_face, hi_z_face);
536  } else if (vert_adv_type == AdvType::Weno_5) {
537  AdvectionSrcForMomWrapper_N<InterpType_H,WENO5,UPWINDALL>(bxx, bxy, bxz,
538  rho_u_rhs, rho_v_rhs, rho_w_rhs,
539  rho_u, rho_v, rho_w, u, v, w,
540  cellSizeInv, stretched_dz_d,
541  mf_mx, mf_ux_inv, mf_vx_inv,
542  mf_my, mf_uy_inv, mf_vy_inv,
543  upw_frac_h, upw_frac_v,
544  vert_adv_type,
545  lo_z_face, hi_z_face);
546  } else if (vert_adv_type == AdvType::Weno_5Z) {
547  AdvectionSrcForMomWrapper_N<InterpType_H,WENO_Z5,UPWINDALL>(bxx, bxy, bxz,
548  rho_u_rhs, rho_v_rhs, rho_w_rhs,
549  rho_u, rho_v, rho_w, u, v, w,
550  cellSizeInv, stretched_dz_d,
551  mf_mx, mf_ux_inv, mf_vx_inv,
552  mf_my, mf_uy_inv, mf_vy_inv,
553  upw_frac_h, upw_frac_v,
554  vert_adv_type,
555  lo_z_face, hi_z_face);
556  } else if (vert_adv_type == AdvType::Weno_7) {
557  AdvectionSrcForMomWrapper_N<InterpType_H,WENO7,UPWINDALL>(bxx, bxy, bxz,
558  rho_u_rhs, rho_v_rhs, rho_w_rhs,
559  rho_u, rho_v, rho_w, u, v, w,
560  cellSizeInv, stretched_dz_d,
561  mf_mx, mf_ux_inv, mf_vx_inv,
562  mf_my, mf_uy_inv, mf_vy_inv,
563  upw_frac_h, upw_frac_v,
564  vert_adv_type,
565  lo_z_face, hi_z_face);
566  } else if (vert_adv_type == AdvType::Weno_7Z) {
567  AdvectionSrcForMomWrapper_N<InterpType_H,WENO_Z7,UPWINDALL>(bxx, bxy, bxz,
568  rho_u_rhs, rho_v_rhs, rho_w_rhs,
569  rho_u, rho_v, rho_w, u, v, w,
570  cellSizeInv, stretched_dz_d,
571  mf_mx, mf_ux_inv, mf_vx_inv,
572  mf_my, mf_uy_inv, mf_vy_inv,
573  upw_frac_h, upw_frac_v,
574  vert_adv_type,
575  lo_z_face, hi_z_face);
576  } else {
577  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
578  }
579 }
AdvectionSrcForMomVert_N
void AdvectionSrcForMomVert_N(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 > &rho_w, const amrex::Array4< const amrex::Real > &u, const amrex::Array4< const amrex::Real > &v, const amrex::Array4< const amrex::Real > &w, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Gpu::DeviceVector< amrex::Real > &stretched_dz_d, const amrex::Array4< const amrex::Real > &mf_mx, const amrex::Array4< const amrex::Real > &mf_ux_inv, const amrex::Array4< const amrex::Real > &mf_vx_inv, const amrex::Array4< const amrex::Real > &mf_my, const amrex::Array4< const amrex::Real > &mf_uy_inv, const amrex::Array4< const amrex::Real > &mf_vy_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_N.H:433
AdvectionSrcForYMom_N
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real AdvectionSrcForYMom_N(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 > &rho_w, InterpType_H interp_v_h, InterpType_V interp_v_v, const amrex::Array4< const amrex::Real > &mf_uy_inv, const amrex::Array4< const amrex::Real > &mf_vx_inv, const amrex::Array4< const amrex::Real > &mf_vy_inv, const amrex::Real dxInv, const amrex::Real dyInv, const amrex::Real dzInv)
Definition: ERF_AdvectionSrcForMom_N.H:112
AdvectionSrcForMomWrapper_N
void AdvectionSrcForMomWrapper_N(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 > &rho_w, const amrex::Array4< const amrex::Real > &u, const amrex::Array4< const amrex::Real > &v, const amrex::Array4< const amrex::Real > &w, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const amrex::Gpu::DeviceVector< amrex::Real > &stretched_dz_d, const amrex::Array4< const amrex::Real > &mf_mx, const amrex::Array4< const amrex::Real > &mf_ux_inv, const amrex::Array4< const amrex::Real > &mf_vx_inv, const amrex::Array4< const amrex::Real > &mf_my, const amrex::Array4< const amrex::Real > &mf_uy_inv, const amrex::Array4< const amrex::Real > &mf_vy_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_N.H:335
AdvectionSrcForZMom_N
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real AdvectionSrcForZMom_N(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 > &rho_w, const amrex::Array4< const amrex::Real > &w, InterpType_H interp_w_h, InterpType_V interp_w_v, WallInterpType interp_w_wall, const amrex::Array4< const amrex::Real > &mf_mx, const amrex::Array4< const amrex::Real > &mf_my, const amrex::Array4< const amrex::Real > &mf_uy_inv, const amrex::Array4< const amrex::Real > &mf_vx_inv, const AdvType vert_adv_type, const int lo_z_face, const int hi_z_face, const amrex::Real dxInv, const amrex::Real dyInv, const amrex::Real dzInv)
Definition: ERF_AdvectionSrcForMom_N.H:201
AdvectionSrcForXMom_N
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real AdvectionSrcForXMom_N(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 > &rho_w, InterpType_H interp_u_h, InterpType_V interp_u_v, const amrex::Array4< const amrex::Real > &mf_ux_inv, const amrex::Array4< const amrex::Real > &mf_uy_inv, const amrex::Array4< const amrex::Real > &mf_vx_inv, const amrex::Real dxInv, const amrex::Real dyInv, const amrex::Real dzInv)
Definition: ERF_AdvectionSrcForMom_N.H:30
two
constexpr amrex::Real two
Definition: ERF_Constants.H:10
one
constexpr amrex::Real one
Definition: ERF_Constants.H:9
zero
constexpr amrex::Real zero
Definition: ERF_Constants.H:8
myhalf
constexpr amrex::Real myhalf
Definition: ERF_Constants.H:13
ERF_IndexDefines.H
AdvType
AdvType
Definition: ERF_IndexDefines.H:261
AdvType::Weno_7
@ Weno_7
AdvType::Upwind_3rd
@ Upwind_3rd
AdvType::Weno_7Z
@ Weno_7Z
AdvType::Centered_4th
@ Centered_4th
AdvType::Upwind_5th
@ Upwind_5th
AdvType::Centered_6th
@ Centered_6th
AdvType::Weno_3Z
@ Weno_3Z
AdvType::Weno_3
@ Weno_3
AdvType::Centered_2nd
@ Centered_2nd
AdvType::Weno_5Z
@ Weno_5Z
AdvType::Weno_5
@ Weno_5
AdvType::Weno_3MZQ
@ Weno_3MZQ
dxInv
amrex::GpuArray< Real, AMREX_SPACEDIM > dxInv
Definition: ERF_InitCustomPertVels_ParticleTests.H:17
ERF_Interpolation.H
ParallelFor
ParallelFor(fab_box, [=] AMREX_GPU_DEVICE(int i, int j, int k) { qrcuten_arr(i, j, k)=Real(0);qscuten_arr(i, j, k)=Real(0);qicuten_arr(i, j, k)=Real(0);})
w
Real w
Definition: ERF_Plotfile2DInterpolator.cpp:22
Real
amrex::Real Real
Definition: ERF_ShocInterface.H:19
AMREX_ASSERT_WITH_MESSAGE
AMREX_ASSERT_WITH_MESSAGE(wbar_cutoff_min > wbar_cutoff_max, "ERROR: wbar_cutoff_min < wbar_cutoff_max")