ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Functions
ERF_AdvectionSrcForMom_StretchedDz.cpp File Reference
#include "AMReX_BCRec.H"
#include <ERF_Advection.H>
#include <ERF_AdvectionSrcForMom_N.H>
#include <ERF_AdvectionSrcForMom_T.H>
Include dependency graph for ERF_AdvectionSrcForMom_StretchedDz.cpp:

Functions

void AdvectionSrcForMom_StretchedDz (const Box &bxx, const Box &bxy, const Box &bxz, const Array4< Real > &rho_u_rhs, const Array4< Real > &rho_v_rhs, const Array4< Real > &rho_w_rhs, const Array4< const Real > &u, const Array4< const Real > &v, const Array4< const Real > &w, const Array4< const Real > &rho_u, const Array4< const Real > &rho_v, const Array4< const Real > &omega, const GpuArray< Real, AMREX_SPACEDIM > &cellSizeInv, const Gpu::DeviceVector< Real > &stretched_dz_d, 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 AdvType horiz_adv_type, const AdvType vert_adv_type, const Real horiz_upw_frac, const Real vert_upw_frac, const int lo_z_face, const int hi_z_face)
 

Function Documentation

◆ AdvectionSrcForMom_StretchedDz()

void AdvectionSrcForMom_StretchedDz ( const Box &  bxx,
const Box &  bxy,
const Box &  bxz,
const Array4< Real > &  rho_u_rhs,
const Array4< Real > &  rho_v_rhs,
const Array4< Real > &  rho_w_rhs,
const Array4< const Real > &  u,
const Array4< const Real > &  v,
const Array4< const Real > &  w,
const Array4< const Real > &  rho_u,
const Array4< const Real > &  rho_v,
const Array4< const Real > &  omega,
const GpuArray< Real, AMREX_SPACEDIM > &  cellSizeInv,
const Gpu::DeviceVector< Real > &  stretched_dz_d,
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 AdvType  horiz_adv_type,
const AdvType  vert_adv_type,
const Real  horiz_upw_frac,
const Real  vert_upw_frac,
const int  lo_z_face,
const int  hi_z_face 
)

Function for computing the advective tendency for the momentum equations when using stretched dz with no EB and no terrain-fitted coordinates.

Parameters
[in]bxxbox over which the x-momentum is updated
[in]bxybox over which the y-momentum is updated
[in]bxzbox over which the z-momentum is updated
[out]rho_u_rhstendency for the x-momentum equation
[out]rho_v_rhstendency for the y-momentum equation
[out]rho_w_rhstendency for the z-momentum equation
[in]ux-component of the velocity
[in]vy-component of the velocity
[in]wz-component of the velocity
[in]rho_ux-component of the momentum
[in]rho_vy-component of the momentum
[in]omegacomponent of the momentum normal to the z-coordinate surface
[in]cellSizeInvinverse of the grid spacing
[in]stretched_dz_ddevice vector of vertical cell spacing
[in]mf_mxx map factor at cell centers
[in]mf_uxx map factor at x-faces
[in]mf_vxx map factor at y-faces
[in]mf_myy map factor at cell centers
[in]mf_uyy map factor at x-faces
[in]mf_vyy map factor at y-faces
[in]horiz_adv_typesets the spatial order to be used for lateral derivatives
[in]vert_adv_typesets the spatial order to be used for vertical derivatives
[in]horiz_upw_frachorizontal upwind blending fraction
[in]vert_upw_fracvertical upwind blending fraction
[in]lo_z_faceminimum z-face k-index at this level
[in]hi_z_facemaximum z-face k-index at this level
64 {
65  BL_PROFILE_VAR("AdvectionSrcForMom_StretchedDz", AdvectionSrcForMom_StretchedDz);
66 
67  AMREX_ALWAYS_ASSERT(bxz.smallEnd(2) > 0);
68 
69  auto dxInv = cellSizeInv[0]; auto dyInv = cellSizeInv[1];
70 
71  // compute mapfactor inverses
72  Box box2d_u(bxx); box2d_u.setRange(2,0); box2d_u.grow({3,3,0});
73  Box box2d_v(bxy); box2d_v.setRange(2,0); box2d_v.grow({3,3,0});
74 
75  FArrayBox mf_ux_invFAB(box2d_u,1,The_Async_Arena());
76  FArrayBox mf_uy_invFAB(box2d_u,1,The_Async_Arena());
77  const Array4<Real>& mf_ux_inv = mf_ux_invFAB.array();
78  const Array4<Real>& mf_uy_inv = mf_uy_invFAB.array();
79 
80  FArrayBox mf_vx_invFAB(box2d_v,1,The_Async_Arena());
81  FArrayBox mf_vy_invFAB(box2d_v,1,The_Async_Arena());
82  const Array4<Real>& mf_vx_inv = mf_vx_invFAB.array();
83  const Array4<Real>& mf_vy_inv = mf_vy_invFAB.array();
84 
85  ParallelFor(box2d_u, box2d_v,
86  [=] AMREX_GPU_DEVICE (int i, int j, int) noexcept
87  {
88  mf_ux_inv(i,j,0) = one / mf_ux(i,j,0);
89  mf_uy_inv(i,j,0) = one / mf_uy(i,j,0);
90  },
91  [=] AMREX_GPU_DEVICE (int i, int j, int) noexcept
92  {
93  mf_vx_inv(i,j,0) = one / mf_vx(i,j,0);
94  mf_vy_inv(i,j,0) = one / mf_vy(i,j,0);
95  });
96 
97  auto dz_ptr = stretched_dz_d.data();
98 
99  // Inline with 2nd order for efficiency
100  if (horiz_adv_type == AdvType::Centered_2nd && vert_adv_type == AdvType::Centered_2nd)
101  {
102  ParallelFor(bxx, bxy, bxz,
103  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
104  {
105  Real xflux_hi = fourth * (rho_u(i, j , k) * mf_uy_inv(i,j,0) + rho_u(i+1, j , k) * mf_uy_inv(i+1,j,0)) * (u(i+1,j,k) + u(i,j,k));
106  Real xflux_lo = fourth * (rho_u(i, j , k) * mf_uy_inv(i,j,0) + rho_u(i-1, j , k) * mf_uy_inv(i-1,j,0)) * (u(i-1,j,k) + u(i,j,k));
107 
108  Real yflux_hi = fourth * (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)) * (u(i,j+1,k) + u(i,j,k));
109  Real yflux_lo = fourth * (rho_v(i, j , k) * mf_vx_inv(i,j ,0) + rho_v(i-1, j , k) * mf_vx_inv(i-1,j ,0)) * (u(i,j-1,k) + u(i,j,k));
110 
111  Real zflux_hi = fourth * (omega(i, j, k+1) + omega(i-1, j, k+1)) * (u(i,j,k+1) + u(i,j,k));
112  Real zflux_lo = fourth * (omega(i, j, k ) + omega(i-1, j, k )) * (u(i,j,k-1) + u(i,j,k));
113 
114  Real mfsq = mf_ux(i,j,0) * mf_uy(i,j,0);
115 
116  Real dzInv = one/dz_ptr[k];
117 
118  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
119  + (yflux_hi - yflux_lo) * dyInv * mfsq
120  + (zflux_hi - zflux_lo) * dzInv;
121  rho_u_rhs(i, j, k) = -advectionSrc;
122  },
123  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
124  {
125  Real xflux_hi = fourth * (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)) * (v(i+1,j,k) + v(i,j,k));
126  Real xflux_lo = fourth * (rho_u(i , j, k) * mf_uy_inv(i ,j,0) + rho_u(i , j-1, k) * mf_uy_inv(i ,j-1,0)) * (v(i-1,j,k) + v(i,j,k));
127 
128  Real yflux_hi = fourth * (rho_v(i ,j+1,k) * mf_vx_inv(i,j+1,0) + rho_v(i ,j ,k) * mf_vx_inv(i,j ,0)) * (v(i,j+1,k) + v(i,j,k));
129  Real yflux_lo = fourth * (rho_v(i ,j ,k) * mf_vx_inv(i,j ,0) + rho_v(i ,j-1,k) * mf_vx_inv(i,j-1,0)) * (v(i,j-1,k) + v(i,j,k));
130 
131  Real zflux_hi = fourth * (omega(i, j, k+1) + omega(i, j-1, k+1)) * (v(i,j,k+1) + v(i,j,k));
132  Real zflux_lo = fourth * (omega(i, j, k ) + omega(i, j-1, k )) * (v(i,j,k-1) + v(i,j,k));
133 
134  Real mfsq = mf_vx(i,j,0) * mf_vy(i,j,0);
135 
136  Real dzInv = one/dz_ptr[k];
137 
138  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
139  + (yflux_hi - yflux_lo) * dyInv * mfsq
140  + (zflux_hi - zflux_lo) * dzInv;
141  rho_v_rhs(i, j, k) = -advectionSrc;
142  },
143  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
144  {
145  Real xflux_hi = fourth*(rho_u(i+1,j ,k) + rho_u(i+1, j, k-1)) * mf_uy_inv(i+1,j ,0) * (w(i+1,j,k) + w(i,j,k));
146  Real xflux_lo = fourth*(rho_u(i ,j ,k) + rho_u(i , j, k-1)) * mf_uy_inv(i ,j ,0) * (w(i-1,j,k) + w(i,j,k));
147 
148  Real yflux_hi = fourth*(rho_v(i ,j+1,k) + rho_v(i, j+1, k-1)) * mf_vx_inv(i ,j+1,0) * (w(i,j+1,k) + w(i,j,k));
149  Real yflux_lo = fourth*(rho_v(i ,j ,k) + rho_v(i, j , k-1)) * mf_vx_inv(i ,j ,0) * (w(i,j-1,k) + w(i,j,k));
150 
151  Real zflux_lo = fourth * (omega(i,j,k) + omega(i,j,k-1)) * (w(i,j,k) + w(i,j,k-1));
152 
153  Real zflux_hi = (k == hi_z_face) ? omega(i,j,k) * w(i,j,k) :
154  fourth * (omega(i,j,k) + omega(i,j,k+1)) * (w(i,j,k) + w(i,j,k+1));
155 
156  Real mfsq = mf_mx(i,j,0) * mf_my(i,j,0);
157 
158  // dz_ptr holds cell widths, so the top face (k == hi_z_face) has only the cell below it
159  Real dzInv = (k == 0) ? one / dz_ptr[k] :
160  (k == hi_z_face) ? one / dz_ptr[k-1] : two/(dz_ptr[k] + dz_ptr[k-1]);
161 
162  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
163  + (yflux_hi - yflux_lo) * dyInv * mfsq
164  + (zflux_hi - zflux_lo) * dzInv;
165  rho_w_rhs(i, j, k) = -advectionSrc;
166  });
167 
168  // Template higher order methods
169  } else {
170  if (horiz_adv_type == AdvType::Centered_2nd) {
171  AdvectionSrcForMomVert_N<CENTERED2>(bxx, bxy, bxz,
172  rho_u_rhs, rho_v_rhs, rho_w_rhs,
173  rho_u, rho_v, omega, u, v, w,
174  cellSizeInv, stretched_dz_d,
175  mf_mx, mf_ux_inv, mf_vx_inv,
176  mf_my, mf_uy_inv, mf_vy_inv,
177  horiz_upw_frac, vert_upw_frac,
178  vert_adv_type, lo_z_face, hi_z_face);
179  } else if (horiz_adv_type == AdvType::Upwind_3rd) {
180  AdvectionSrcForMomVert_N<UPWIND3>(bxx, bxy, bxz,
181  rho_u_rhs, rho_v_rhs, rho_w_rhs,
182  rho_u, rho_v, omega, u, v, w,
183  cellSizeInv, stretched_dz_d,
184  mf_mx, mf_ux_inv, mf_vx_inv,
185  mf_my, mf_uy_inv, mf_vy_inv,
186  horiz_upw_frac, vert_upw_frac,
187  vert_adv_type, lo_z_face, hi_z_face);
188  } else if (horiz_adv_type == AdvType::Centered_4th) {
189  AdvectionSrcForMomVert_N<CENTERED4>(bxx, bxy, bxz,
190  rho_u_rhs, rho_v_rhs, rho_w_rhs,
191  rho_u, rho_v, omega, u, v, w,
192  cellSizeInv, stretched_dz_d,
193  mf_mx, mf_ux_inv, mf_vx_inv,
194  mf_my, mf_uy_inv, mf_vy_inv,
195  horiz_upw_frac, vert_upw_frac,
196  vert_adv_type, lo_z_face, hi_z_face);
197  } else if (horiz_adv_type == AdvType::Upwind_5th) {
198  AdvectionSrcForMomVert_N<UPWIND5>(bxx, bxy, bxz,
199  rho_u_rhs, rho_v_rhs, rho_w_rhs,
200  rho_u, rho_v, omega, u, v, w,
201  cellSizeInv, stretched_dz_d,
202  mf_mx, mf_ux_inv, mf_vx_inv,
203  mf_my, mf_uy_inv, mf_vy_inv,
204  horiz_upw_frac, vert_upw_frac,
205  vert_adv_type, lo_z_face, hi_z_face);
206  } else if (horiz_adv_type == AdvType::Centered_6th) {
207  AdvectionSrcForMomVert_N<CENTERED6>(bxx, bxy, bxz,
208  rho_u_rhs, rho_v_rhs, rho_w_rhs,
209  rho_u, rho_v, omega, u, v, w,
210  cellSizeInv, stretched_dz_d,
211  mf_mx, mf_ux_inv, mf_vx_inv,
212  mf_my, mf_uy_inv, mf_vy_inv,
213  horiz_upw_frac, vert_upw_frac,
214  vert_adv_type, lo_z_face, hi_z_face);
215  } else if (horiz_adv_type == AdvType::Weno_3) {
216  AdvectionSrcForMomVert_N<WENO3>(bxx, bxy, bxz,
217  rho_u_rhs, rho_v_rhs, rho_w_rhs,
218  rho_u, rho_v, omega, u, v, w,
219  cellSizeInv, stretched_dz_d,
220  mf_mx, mf_ux_inv, mf_vx_inv,
221  mf_my, mf_uy_inv, mf_vy_inv,
222  horiz_upw_frac, vert_upw_frac,
223  vert_adv_type, lo_z_face, hi_z_face);
224  } else if (horiz_adv_type == AdvType::Weno_3Z) {
225  AdvectionSrcForMomVert_N<WENO_Z3>(bxx, bxy, bxz,
226  rho_u_rhs, rho_v_rhs, rho_w_rhs,
227  rho_u, rho_v, omega, u, v, w,
228  cellSizeInv, stretched_dz_d,
229  mf_mx, mf_ux_inv, mf_vx_inv,
230  mf_my, mf_uy_inv, mf_vy_inv,
231  horiz_upw_frac, vert_upw_frac,
232  vert_adv_type, lo_z_face, hi_z_face);
233  } else if (horiz_adv_type == AdvType::Weno_3MZQ) {
234  AdvectionSrcForMomVert_N<WENO_MZQ3>(bxx, bxy, bxz,
235  rho_u_rhs, rho_v_rhs, rho_w_rhs,
236  rho_u, rho_v, omega, u, v, w,
237  cellSizeInv, stretched_dz_d,
238  mf_mx, mf_ux_inv, mf_vx_inv,
239  mf_my, mf_uy_inv, mf_vy_inv,
240  horiz_upw_frac, vert_upw_frac,
241  vert_adv_type, lo_z_face, hi_z_face);
242  } else if (horiz_adv_type == AdvType::Weno_5) {
243  AdvectionSrcForMomVert_N<WENO5>(bxx, bxy, bxz,
244  rho_u_rhs, rho_v_rhs, rho_w_rhs,
245  rho_u, rho_v, omega, u, v, w,
246  cellSizeInv, stretched_dz_d,
247  mf_mx, mf_ux_inv, mf_vx_inv,
248  mf_my, mf_uy_inv, mf_vy_inv,
249  horiz_upw_frac, vert_upw_frac,
250  vert_adv_type, lo_z_face, hi_z_face);
251  } else if (horiz_adv_type == AdvType::Weno_5Z) {
252  AdvectionSrcForMomVert_N<WENO_Z5>(bxx, bxy, bxz,
253  rho_u_rhs, rho_v_rhs, rho_w_rhs,
254  rho_u, rho_v, omega, u, v, w,
255  cellSizeInv, stretched_dz_d,
256  mf_mx, mf_ux_inv, mf_vx_inv,
257  mf_my, mf_uy_inv, mf_vy_inv,
258  horiz_upw_frac, vert_upw_frac,
259  vert_adv_type, lo_z_face, hi_z_face);
260  } else if (horiz_adv_type == AdvType::Weno_7) {
261  AdvectionSrcForMomVert_N<WENO7>(bxx, bxy, bxz,
262  rho_u_rhs, rho_v_rhs, rho_w_rhs,
263  rho_u, rho_v, omega, u, v, w,
264  cellSizeInv, stretched_dz_d,
265  mf_mx, mf_ux_inv, mf_vx_inv,
266  mf_my, mf_uy_inv, mf_vy_inv,
267  horiz_upw_frac, vert_upw_frac,
268  vert_adv_type, lo_z_face, hi_z_face);
269  } else if (horiz_adv_type == AdvType::Weno_7Z) {
270  AdvectionSrcForMomVert_N<WENO_Z7>(bxx, bxy, bxz,
271  rho_u_rhs, rho_v_rhs, rho_w_rhs,
272  rho_u, rho_v, omega, u, v, w,
273  cellSizeInv, stretched_dz_d,
274  mf_mx, mf_ux_inv, mf_vx_inv,
275  mf_my, mf_uy_inv, mf_vy_inv,
276  horiz_upw_frac, vert_upw_frac,
277  vert_adv_type, lo_z_face, hi_z_face);
278  } else {
279  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
280  }
281  }
282 }
AdvectionSrcForMom_StretchedDz
void AdvectionSrcForMom_StretchedDz(const Box &bxx, const Box &bxy, const Box &bxz, const Array4< Real > &rho_u_rhs, const Array4< Real > &rho_v_rhs, const Array4< Real > &rho_w_rhs, const Array4< const Real > &u, const Array4< const Real > &v, const Array4< const Real > &w, const Array4< const Real > &rho_u, const Array4< const Real > &rho_v, const Array4< const Real > &omega, const GpuArray< Real, AMREX_SPACEDIM > &cellSizeInv, const Gpu::DeviceVector< Real > &stretched_dz_d, 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 AdvType horiz_adv_type, const AdvType vert_adv_type, const Real horiz_upw_frac, const Real vert_upw_frac, const int lo_z_face, const int hi_z_face)
Definition: ERF_AdvectionSrcForMom_StretchedDz.cpp:41
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
AMREX_ALWAYS_ASSERT
AMREX_ALWAYS_ASSERT(bx.length()[2]==khi+1)
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);})
two
constexpr amrex::Real two
Definition: ERF_NumericalConstants.H:38
one
constexpr amrex::Real one
Definition: ERF_NumericalConstants.H:37
fourth
constexpr amrex::Real fourth
Definition: ERF_NumericalConstants.H:42
w
Real w
Definition: ERF_Plotfile2DInterpolator.cpp:23
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")
MicVar_Morr::omega
@ omega
Definition: ERF_Morrison.H:57
Here is the call graph for this function: