ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_AdvectionSrcForMom.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.cpp:

Functions

void AdvectionSrcForMom (const Box &bx, 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 > &cell_data, 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 Array4< const Real > &z_nd, const Array4< const Real > &ax, const Array4< const Real > &ay, const Array4< const Real > &az, const Array4< const Real > &detJ, const GpuArray< Real, AMREX_SPACEDIM > &cellSizeInv, const Array4< const Real > &mf_m, const Array4< const Real > &mf_u, const Array4< const Real > &mf_v, const AdvType horiz_adv_type, const AdvType vert_adv_type, const Real horiz_upw_frac, const Real vert_upw_frac, TerrainType &terrain_type, const int lo_z_face, const int hi_z_face, const Box &domain, const BCRec *bc_ptr_h)
 

Function Documentation

◆ AdvectionSrcForMom()

void AdvectionSrcForMom ( const Box &  bx,
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 > &  cell_data,
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 Array4< const Real > &  z_nd,
const Array4< const Real > &  ax,
const Array4< const Real > &  ay,
const Array4< const Real > &  az,
const Array4< const Real > &  detJ,
const GpuArray< Real, AMREX_SPACEDIM > &  cellSizeInv,
const Array4< const Real > &  mf_m,
const Array4< const Real > &  mf_u,
const Array4< const Real > &  mf_v,
const AdvType  horiz_adv_type,
const AdvType  vert_adv_type,
const Real  horiz_upw_frac,
const Real  vert_upw_frac,
TerrainType &  terrain_type,
const int  lo_z_face,
const int  hi_z_face,
const Box &  domain,
const BCRec *  bc_ptr_h 
)

Function for computing the advective tendency for the momentum equations This routine has explicit expressions for all cases (terrain or not) when the horizontal and vertical spatial orders are <= 2, and calls more specialized functions when either (or both) spatial order(s) is greater than 2.

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]z_ndheight coordinate at nodes
[in]axArea fraction of x-faces
[in]ayArea fraction of y-faces
[in]azArea fraction of z-faces
[in]detJJacobian of the metric transformation (= 1 if use_terrain_fitted_coords is false)
[in]cellSizeInvinverse of the mesh spacing
[in]mf_mmap factor at cell centers
[in]mf_umap factor at x-faces
[in]mf_vmap 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
69 {
70  BL_PROFILE_VAR("AdvectionSrcForMom", AdvectionSrcForMom);
71 
72  auto dxInv = cellSizeInv[0], dyInv = cellSizeInv[1], dzInv = cellSizeInv[2];
73 
74  AMREX_ALWAYS_ASSERT(bxz.smallEnd(2) > 0);
75 
76  // compute mapfactor inverses
77  Box box2d_u(bxx); box2d_u.setRange(2,0); box2d_u.grow({3,3,0});
78  Box box2d_v(bxy); box2d_v.setRange(2,0); box2d_v.grow({3,3,0});
79  FArrayBox mf_u_invFAB(box2d_u,1,The_Async_Arena());
80  FArrayBox mf_v_invFAB(box2d_v,1,The_Async_Arena());
81  const Array4<Real>& mf_u_inv = mf_u_invFAB.array();
82  const Array4<Real>& mf_v_inv = mf_v_invFAB.array();
83 
84  const bool use_terrain_fitted_coords = ( terrain_type == TerrainType::StaticFittedMesh ||
85  terrain_type == TerrainType::MovingFittedMesh);
86 
87  ParallelFor(box2d_u, box2d_v,
88  [=] AMREX_GPU_DEVICE (int i, int j, int) noexcept
89  {
90  mf_u_inv(i,j,0) = 1. / mf_u(i,j,0);
91  },
92  [=] AMREX_GPU_DEVICE (int i, int j, int) noexcept
93  {
94  mf_v_inv(i,j,0) = 1. / mf_v(i,j,0);
95  });
96 
97  if ( terrain_type == TerrainType::EB)
98  {
99  amrex::ignore_unused(use_terrain_fitted_coords);
100  ParallelFor(bxx, bxy, bxz,
101  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
102  {
103  rho_u_rhs(i, j, k) = 0.0;
104  },
105  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
106  {
107  rho_v_rhs(i, j, k) = 0.0;
108  },
109  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
110  {
111  rho_w_rhs(i, j, k) = 0.0;
112  });
113  } else {
114  if ( !use_terrain_fitted_coords) {
115  // Inline with 2nd order for efficiency
116  if (horiz_adv_type == AdvType::Centered_2nd && vert_adv_type == AdvType::Centered_2nd)
117  {
118  ParallelFor(bxx, bxy, bxz,
119  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
120  {
121  Real xflux_hi = 0.25 * (rho_u(i, j , k) * mf_u_inv(i,j,0) + rho_u(i+1, j , k) * mf_u_inv(i+1,j,0)) * (u(i+1,j,k) + u(i,j,k));
122  Real xflux_lo = 0.25 * (rho_u(i, j , k) * mf_u_inv(i,j,0) + rho_u(i-1, j , k) * mf_u_inv(i-1,j,0)) * (u(i-1,j,k) + u(i,j,k));
123 
124  Real yflux_hi = 0.25 * (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)) * (u(i,j+1,k) + u(i,j,k));
125  Real yflux_lo = 0.25 * (rho_v(i, j , k) * mf_v_inv(i,j ,0) + rho_v(i-1, j , k) * mf_v_inv(i-1,j ,0)) * (u(i,j-1,k) + u(i,j,k));
126 
127  Real zflux_hi = 0.25 * (Omega(i, j, k+1) + Omega(i-1, j, k+1)) * (u(i,j,k+1) + u(i,j,k));
128  Real zflux_lo = 0.25 * (Omega(i, j, k ) + Omega(i-1, j, k )) * (u(i,j,k-1) + u(i,j,k));
129 
130  Real mfsq = mf_u(i,j,0) * mf_u(i,j,0);
131 
132  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
133  + (yflux_hi - yflux_lo) * dyInv * mfsq
134  + (zflux_hi - zflux_lo) * dzInv;
135  rho_u_rhs(i, j, k) = -advectionSrc;
136  },
137  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
138  {
139  Real xflux_hi = 0.25 * (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)) * (v(i+1,j,k) + v(i,j,k));
140  Real xflux_lo = 0.25 * (rho_u(i , j, k) * mf_u_inv(i ,j,0) + rho_u(i , j-1, k) * mf_u_inv(i ,j-1,0)) * (v(i-1,j,k) + v(i,j,k));
141 
142  Real yflux_hi = 0.25 * (rho_v(i ,j+1,k) * mf_v_inv(i,j+1,0) + rho_v(i ,j ,k) * mf_v_inv(i,j ,0)) * (v(i,j+1,k) + v(i,j,k));
143  Real yflux_lo = 0.25 * (rho_v(i ,j ,k) * mf_v_inv(i,j ,0) + rho_v(i ,j-1,k) * mf_v_inv(i,j-1,0) ) * (v(i,j-1,k) + v(i,j,k));
144 
145  Real zflux_hi = 0.25 * (Omega(i, j, k+1) + Omega(i, j-1, k+1)) * (v(i,j,k+1) + v(i,j,k));
146  Real zflux_lo = 0.25 * (Omega(i, j, k ) + Omega(i, j-1, k )) * (v(i,j,k-1) + v(i,j,k));
147 
148  Real mfsq = mf_v(i,j,0) * mf_v(i,j,0);
149 
150  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
151  + (yflux_hi - yflux_lo) * dyInv * mfsq
152  + (zflux_hi - zflux_lo) * dzInv;
153  rho_v_rhs(i, j, k) = -advectionSrc;
154  },
155  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
156  {
157  Real xflux_hi = 0.25*(rho_u(i+1,j ,k) + rho_u(i+1, j, k-1)) * mf_u_inv(i+1,j ,0) * (w(i+1,j,k) + w(i,j,k));
158  Real xflux_lo = 0.25*(rho_u(i ,j ,k) + rho_u(i , j, k-1)) * mf_u_inv(i ,j ,0) * (w(i-1,j,k) + w(i,j,k));
159 
160  Real yflux_hi = 0.25*(rho_v(i ,j+1,k) + rho_v(i, j+1, k-1)) * mf_v_inv(i ,j+1,0) * (w(i,j+1,k) + w(i,j,k));
161  Real yflux_lo = 0.25*(rho_v(i ,j ,k) + rho_v(i, j , k-1)) * mf_v_inv(i ,j ,0) * (w(i,j-1,k) + w(i,j,k));
162 
163  Real zflux_lo = 0.25 * (Omega(i,j,k) + Omega(i,j,k-1)) * (w(i,j,k) + w(i,j,k-1));
164 
165  Real zflux_hi = (k == hi_z_face) ? Omega(i,j,k) * w(i,j,k) :
166  0.25 * (Omega(i,j,k) + Omega(i,j,k+1)) * (w(i,j,k) + w(i,j,k+1));
167 
168  Real mfsq = mf_m(i,j,0) * mf_m(i,j,0);
169 
170  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
171  + (yflux_hi - yflux_lo) * dyInv * mfsq
172  + (zflux_hi - zflux_lo) * dzInv;
173  rho_w_rhs(i, j, k) = -advectionSrc;
174  });
175  // Template higher order methods
176  } else {
177  if (horiz_adv_type == AdvType::Centered_2nd) {
178  AdvectionSrcForMomVert_N<CENTERED2>(bxx, bxy, bxz,
179  rho_u_rhs, rho_v_rhs, rho_w_rhs,
180  rho_u, rho_v, Omega, u, v, w,
181  cellSizeInv, mf_m,
182  mf_u_inv, mf_v_inv,
183  horiz_upw_frac, vert_upw_frac,
184  vert_adv_type, lo_z_face, hi_z_face);
185  } else if (horiz_adv_type == AdvType::Upwind_3rd) {
186  AdvectionSrcForMomVert_N<UPWIND3>(bxx, bxy, bxz,
187  rho_u_rhs, rho_v_rhs, rho_w_rhs,
188  rho_u, rho_v, Omega, u, v, w,
189  cellSizeInv, mf_m,
190  mf_u_inv, mf_v_inv,
191  horiz_upw_frac, vert_upw_frac,
192  vert_adv_type, lo_z_face, hi_z_face);
193  } else if (horiz_adv_type == AdvType::Centered_4th) {
194  AdvectionSrcForMomVert_N<CENTERED4>(bxx, bxy, bxz,
195  rho_u_rhs, rho_v_rhs, rho_w_rhs,
196  rho_u, rho_v, Omega, u, v, w,
197  cellSizeInv, mf_m,
198  mf_u_inv, mf_v_inv,
199  horiz_upw_frac, vert_upw_frac,
200  vert_adv_type, lo_z_face, hi_z_face);
201  } else if (horiz_adv_type == AdvType::Upwind_5th) {
202  AdvectionSrcForMomVert_N<UPWIND5>(bxx, bxy, bxz,
203  rho_u_rhs, rho_v_rhs, rho_w_rhs,
204  rho_u, rho_v, Omega, u, v, w,
205  cellSizeInv, mf_m,
206  mf_u_inv, mf_v_inv,
207  horiz_upw_frac, vert_upw_frac,
208  vert_adv_type, lo_z_face, hi_z_face);
209  } else if (horiz_adv_type == AdvType::Centered_6th) {
210  AdvectionSrcForMomVert_N<CENTERED6>(bxx, bxy, bxz,
211  rho_u_rhs, rho_v_rhs, rho_w_rhs,
212  rho_u, rho_v, Omega, u, v, w,
213  cellSizeInv, mf_m,
214  mf_u_inv, mf_v_inv,
215  horiz_upw_frac, vert_upw_frac,
216  vert_adv_type, lo_z_face, hi_z_face);
217  } else {
218  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
219  }
220  }
221  } // end of use_terrain_fitted_coords == false
222  else
223  { // now do use_terrain_fitted_coords = true
224  // Inline with 2nd order for efficiency
225  if (horiz_adv_type == AdvType::Centered_2nd && vert_adv_type == AdvType::Centered_2nd)
226  {
227  ParallelFor(bxx, bxy, bxz,
228  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
229  {
230  Real xflux_hi = 0.25 * (rho_u(i,j,k) * mf_u_inv(i,j,0) + rho_u(i+1,j,k) * mf_u_inv(i+1,j,0)) *
231  (u(i+1,j,k) + u(i,j,k)) * 0.5 * (ax(i,j,k) + ax(i+1,j,k));
232 
233  Real xflux_lo = 0.25 * (rho_u(i,j,k) * mf_u_inv(i,j,0) + rho_u(i-1,j,k) * mf_u_inv(i-1,j,0)) *
234  (u(i-1,j,k) + u(i,j,k)) * 0.5 * (ax(i,j,k) + ax(i-1,j,k));
235 
236  Real met_h_zeta_yhi = Compute_h_zeta_AtEdgeCenterK(i,j+1,k,cellSizeInv,z_nd);
237  Real yflux_hi = 0.25 * (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)) *
238  (u(i,j+1,k) + u(i,j,k)) * met_h_zeta_yhi;
239 
240  Real met_h_zeta_ylo = Compute_h_zeta_AtEdgeCenterK(i,j ,k,cellSizeInv,z_nd);
241  Real yflux_lo = 0.25 * (rho_v(i,j ,k)*mf_v_inv(i,j ,0) + rho_v(i-1,j ,k)*mf_v_inv(i-1,j ,0)) *
242  (u(i,j-1,k) + u(i,j,k)) * met_h_zeta_ylo;
243 
244  Real zflux_hi = 0.25 * (Omega(i,j,k+1) + Omega(i-1,j,k+1)) * (u(i,j,k+1) + u(i,j,k)) *
245  0.5 * (az(i,j,k+1) + az(i-1,j,k+1));
246  Real zflux_lo = 0.25 * (Omega(i,j,k ) + Omega(i-1,j,k )) * (u(i,j,k-1) + u(i,j,k)) *
247  0.5 * (az(i,j,k ) + az(i-1,j,k ));
248 
249  Real mfsq = mf_u(i,j,0) * mf_u(i,j,0);
250 
251  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
252  + (yflux_hi - yflux_lo) * dyInv * mfsq
253  + (zflux_hi - zflux_lo) * dzInv;
254 
255  rho_u_rhs(i, j, k) = -advectionSrc / (0.5 * (detJ(i,j,k) + detJ(i-1,j,k)));
256  },
257  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
258  {
259 
260  Real met_h_zeta_xhi = Compute_h_zeta_AtEdgeCenterK(i+1,j,k,cellSizeInv,z_nd);
261  Real xflux_hi = 0.25 * (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)) *
262  (v(i+1,j,k) + v(i,j,k)) * met_h_zeta_xhi;
263 
264  Real met_h_zeta_xlo = Compute_h_zeta_AtEdgeCenterK(i ,j,k,cellSizeInv,z_nd);
265  Real xflux_lo = 0.25 * (rho_u(i, j, k)*mf_u_inv(i ,j,0) + rho_u(i ,j-1,k)*mf_u_inv(i-1,j ,0)) *
266  (v(i-1,j,k) + v(i,j,k)) * met_h_zeta_xlo;
267 
268  Real yflux_hi = 0.25 * (rho_v(i,j+1,k)*mf_v_inv(i,j+1,0) + rho_v(i,j ,k) * mf_v_inv(i,j ,0)) *
269  (v(i,j+1,k) + v(i,j,k)) * 0.5 * (ay(i,j,k) + ay(i,j+1,k));
270 
271  Real yflux_lo = 0.25 * (rho_v(i,j ,k)*mf_v_inv(i,j ,0) + rho_v(i,j-1,k) * mf_v_inv(i,j-1,0)) *
272  (v(i,j-1,k) + v(i,j,k)) * 0.5 * (ay(i,j,k) + ay(i,j-1,k));
273 
274  Real zflux_hi = 0.25 * (Omega(i,j,k+1) + Omega(i, j-1, k+1)) * (v(i,j,k+1) + v(i,j,k)) *
275  0.5 * (az(i,j,k+1) + az(i,j-1,k+1));
276  Real zflux_lo = 0.25 * (Omega(i,j,k ) + Omega(i, j-1, k )) * (v(i,j,k-1) + v(i,j,k)) *
277  0.5 * (az(i,j,k ) + az(i,j-1,k ));
278 
279  Real mfsq = mf_v(i,j,0) * mf_v(i,j,0);
280 
281  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
282  + (yflux_hi - yflux_lo) * dyInv * mfsq
283  + (zflux_hi - zflux_lo) * dzInv;
284 
285  rho_v_rhs(i, j, k) = -advectionSrc / (0.5 * (detJ(i,j,k) + detJ(i,j-1,k)));
286  },
287  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
288  {
289  Real met_h_zeta_xhi = Compute_h_zeta_AtEdgeCenterJ(i+1,j ,k ,cellSizeInv,z_nd);
290  Real xflux_hi = 0.25*(rho_u(i+1,j ,k) + rho_u(i+1,j,k-1)) * mf_u_inv(i+1,j,0) *
291  (w(i+1,j,k) + w(i,j,k)) * met_h_zeta_xhi;
292 
293  Real met_h_zeta_xlo = Compute_h_zeta_AtEdgeCenterJ(i ,j ,k ,cellSizeInv,z_nd);
294  Real xflux_lo = 0.25*(rho_u(i ,j ,k) + rho_u(i ,j,k-1)) * mf_u_inv(i ,j,0) *
295  (w(i-1,j,k) + w(i,j,k)) * met_h_zeta_xlo;
296 
297  Real met_h_zeta_yhi = Compute_h_zeta_AtEdgeCenterI(i ,j+1,k ,cellSizeInv,z_nd);
298  Real yflux_hi = 0.25*(rho_v(i,j+1,k) + rho_v(i,j+1,k-1)) * mf_v_inv(i,j+1,0) *
299  (w(i,j+1,k) + w(i,j,k)) * met_h_zeta_yhi;
300 
301  Real met_h_zeta_ylo = Compute_h_zeta_AtEdgeCenterI(i ,j ,k ,cellSizeInv,z_nd);
302  Real yflux_lo = 0.25*(rho_v(i,j ,k) + rho_v(i,j ,k-1)) * mf_v_inv(i,j ,0) *
303  (w(i,j-1,k) + w(i,j,k)) * met_h_zeta_ylo;
304 
305  Real zflux_lo = 0.25 * (Omega(i,j,k) + Omega(i,j,k-1)) * (w(i,j,k) + w(i,j,k-1));
306 
307  Real zflux_hi = (k == hi_z_face) ? Omega(i,j,k) * w(i,j,k) * az(i,j,k):
308  0.25 * (Omega(i,j,k) + Omega(i,j,k+1)) * (w(i,j,k) + w(i,j,k+1)) *
309  0.5 * (az(i,j,k) + az(i,j,k+1));
310 
311  Real mfsq = mf_m(i,j,0) * mf_m(i,j,0);
312 
313  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
314  + (yflux_hi - yflux_lo) * dyInv * mfsq
315  + (zflux_hi - zflux_lo) * dzInv;
316 
317  rho_w_rhs(i, j, k) = -advectionSrc / (0.5*(detJ(i,j,k) + detJ(i,j,k-1)));
318  });
319  // Template higher order methods
320  } else {
321  if (horiz_adv_type == AdvType::Centered_2nd) {
322  AdvectionSrcForMomVert<CENTERED2>(bxx, bxy, bxz,
323  rho_u_rhs, rho_v_rhs, rho_w_rhs,
324  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
325  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
326  horiz_upw_frac, vert_upw_frac,
327  vert_adv_type, lo_z_face, hi_z_face);
328  } else if (horiz_adv_type == AdvType::Upwind_3rd) {
329  AdvectionSrcForMomVert<UPWIND3>(bxx, bxy, bxz,
330  rho_u_rhs, rho_v_rhs, rho_w_rhs,
331  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
332  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
333  horiz_upw_frac, vert_upw_frac,
334  vert_adv_type, lo_z_face, hi_z_face);
335  } else if (horiz_adv_type == AdvType::Centered_4th) {
336  AdvectionSrcForMomVert<CENTERED4>(bxx, bxy, bxz,
337  rho_u_rhs, rho_v_rhs, rho_w_rhs,
338  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
339  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
340  horiz_upw_frac, vert_upw_frac,
341  vert_adv_type, lo_z_face, hi_z_face);
342  } else if (horiz_adv_type == AdvType::Upwind_5th) {
343  AdvectionSrcForMomVert<UPWIND5>(bxx, bxy, bxz,
344  rho_u_rhs, rho_v_rhs, rho_w_rhs,
345  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
346  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
347  horiz_upw_frac, vert_upw_frac,
348  vert_adv_type, lo_z_face, hi_z_face);
349  } else if (horiz_adv_type == AdvType::Centered_6th) {
350  AdvectionSrcForMomVert<CENTERED6>(bxx, bxy, bxz,
351  rho_u_rhs, rho_v_rhs, rho_w_rhs,
352  rho_u, rho_v, Omega, u, v, w, z_nd, ax, ay, az, detJ,
353  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
354  horiz_upw_frac, vert_upw_frac,
355  vert_adv_type, lo_z_face, hi_z_face);
356  } else {
357  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
358  }
359  } // higher order
360  } // terrain
361  } // Not EB
362 
363  // Open bc will be imposed upon all vars (we only access cons here for simplicity)
364  const bool xlo_open = (bc_ptr_h[BCVars::cons_bc].lo(0) == ERFBCType::open);
365  const bool xhi_open = (bc_ptr_h[BCVars::cons_bc].hi(0) == ERFBCType::open);
366  const bool ylo_open = (bc_ptr_h[BCVars::cons_bc].lo(1) == ERFBCType::open);
367  const bool yhi_open = (bc_ptr_h[BCVars::cons_bc].hi(1) == ERFBCType::open);
368 
369  // We recreate tbx, tbz, tbz here rather than using bxx, bxy, bxz because those
370  // have already been shrunk by one in the case of open BCs.
371  Box tbx(surroundingNodes(bx,0));
372  Box tby(surroundingNodes(bx,1));
373  Box tbz(surroundingNodes(bx,2)); tbz.growLo(2,-1); tbz.growHi(2,-1);
374 
375  const int domhi_z = domain.bigEnd(2);
376 
377  // Special advection operator for open BC (bndry normal/tangent operations)
378  if (xlo_open)
379  {
380  Box tbx_xlo, tby_xlo, tbz_xlo;
381  if (tbx.smallEnd(0) == domain.smallEnd(0)) { tbx_xlo = makeSlab(tbx,0,domain.smallEnd(0));}
382  if (tby.smallEnd(0) == domain.smallEnd(0)) { tby_xlo = makeSlab(tby,0,domain.smallEnd(0));}
383  if (tbz.smallEnd(0) == domain.smallEnd(0)) { tbz_xlo = makeSlab(tbz,0,domain.smallEnd(0));}
384 
385  bool do_lo = true;
386 
387  AdvectionSrcForOpenBC_Normal(tbx_xlo, 0, rho_u_rhs, u, cell_data, cellSizeInv, do_lo);
388  AdvectionSrcForOpenBC_Tangent_Ymom(tby_xlo, 0, rho_v_rhs, v,
389  rho_u, rho_v, Omega,
390  ay, az, detJ, cellSizeInv,
391  do_lo);
392  AdvectionSrcForOpenBC_Tangent_Zmom(tbz_xlo, 0, rho_w_rhs, w,
393  rho_u, rho_v, Omega,
394  ax, ay, az, detJ, cellSizeInv,
395  domhi_z, do_lo);
396  }
397  if (xhi_open)
398  {
399  Box tbx_xhi, tby_xhi, tbz_xhi;
400  if (tbx.bigEnd(0) == domain.bigEnd(0)+1) { tbx_xhi = makeSlab(tbx,0,domain.bigEnd(0)+1);}
401  if (tby.bigEnd(0) == domain.bigEnd(0)) { tby_xhi = makeSlab(tby,0,domain.bigEnd(0) );}
402  if (tbz.bigEnd(0) == domain.bigEnd(0)) { tbz_xhi = makeSlab(tbz,0,domain.bigEnd(0) );}
403 
404  AdvectionSrcForOpenBC_Normal(tbx_xhi, 0, rho_u_rhs, u, cell_data, cellSizeInv);
405  AdvectionSrcForOpenBC_Tangent_Ymom(tby_xhi, 0, rho_v_rhs, v,
406  rho_u, rho_v, Omega,
407  ay, az, detJ, cellSizeInv);
408  AdvectionSrcForOpenBC_Tangent_Zmom(tbz_xhi, 0, rho_w_rhs, w,
409  rho_u, rho_v, Omega,
410  ax, ay, az, detJ, cellSizeInv,
411  domhi_z);
412  }
413  if (ylo_open)
414  {
415  Box tbx_ylo, tby_ylo, tbz_ylo;
416  if (tbx.smallEnd(1) == domain.smallEnd(1)) { tbx_ylo = makeSlab(tbx,1,domain.smallEnd(1));}
417  if (tby.smallEnd(1) == domain.smallEnd(1)) { tby_ylo = makeSlab(tby,1,domain.smallEnd(1));}
418  if (tbz.smallEnd(1) == domain.smallEnd(1)) { tbz_ylo = makeSlab(tbz,1,domain.smallEnd(1));}
419 
420  bool do_lo = true;
421  AdvectionSrcForOpenBC_Tangent_Xmom(tbx_ylo, 1, rho_u_rhs, u,
422  rho_u, rho_v, Omega,
423  ax, az, detJ, cellSizeInv,
424  do_lo);
425  AdvectionSrcForOpenBC_Normal(tby_ylo, 1, rho_v_rhs, v, cell_data, cellSizeInv, do_lo);
426  AdvectionSrcForOpenBC_Tangent_Zmom(tbz_ylo, 1, rho_w_rhs, w,
427  rho_u, rho_v, Omega,
428  ax, ay, az, detJ, cellSizeInv,
429  domhi_z, do_lo);
430  }
431  if (yhi_open)
432  {
433  Box tbx_yhi, tby_yhi, tbz_yhi;
434  if (tbx.bigEnd(1) == domain.bigEnd(1)) { tbx_yhi = makeSlab(tbx,1,domain.bigEnd(1) );}
435  if (tby.bigEnd(1) == domain.bigEnd(1)+1) { tby_yhi = makeSlab(tby,1,domain.bigEnd(1)+1);}
436  if (tbz.bigEnd(1) == domain.bigEnd(1)) { tbz_yhi = makeSlab(tbz,1,domain.bigEnd(1) );}
437 
438  AdvectionSrcForOpenBC_Tangent_Xmom(tbx_yhi, 1, rho_u_rhs, u,
439  rho_u, rho_v, Omega,
440  ax, az, detJ, cellSizeInv);
441  AdvectionSrcForOpenBC_Normal(tby_yhi, 1, rho_v_rhs, v, cell_data, cellSizeInv);
442  AdvectionSrcForOpenBC_Tangent_Zmom(tbz_yhi, 1, rho_w_rhs, w,
443  rho_u, rho_v, Omega,
444  ax, ay, az, detJ, cellSizeInv,
445  domhi_z);
446  }
447 }
void AdvectionSrcForMom(const Box &bx, 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 > &cell_data, 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 Array4< const Real > &z_nd, const Array4< const Real > &ax, const Array4< const Real > &ay, const Array4< const Real > &az, const Array4< const Real > &detJ, const GpuArray< Real, AMREX_SPACEDIM > &cellSizeInv, const Array4< const Real > &mf_m, const Array4< const Real > &mf_u, const Array4< const Real > &mf_v, const AdvType horiz_adv_type, const AdvType vert_adv_type, const Real horiz_upw_frac, const Real vert_upw_frac, TerrainType &terrain_type, const int lo_z_face, const int hi_z_face, const Box &domain, const BCRec *bc_ptr_h)
Definition: ERF_AdvectionSrcForMom.cpp:40
void AdvectionSrcForOpenBC_Tangent_Ymom(const amrex::Box &bxy, const int &dir, const amrex::Array4< amrex::Real > &rho_v_rhs, const amrex::Array4< const amrex::Real > &v, 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 > &ay, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const bool do_lo=false)
void AdvectionSrcForOpenBC_Tangent_Xmom(const amrex::Box &bxx, const int &dir, const amrex::Array4< amrex::Real > &rho_u_rhs, const amrex::Array4< const amrex::Real > &u, 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 > &ax, const amrex::Array4< const amrex::Real > &az, const amrex::Array4< const amrex::Real > &detJ, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &cellSizeInv, const bool do_lo=false)
void AdvectionSrcForOpenBC_Normal(const amrex::Box &bx, const int &dir, const amrex::Array4< amrex::Real > &rhs_arr, const amrex::Array4< const amrex::Real > &vel_norm_arr, const amrex::Array4< const amrex::Real > &cell_data_arr, const amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > &dxInv, const bool do_lo=false)
void AdvectionSrcForOpenBC_Tangent_Zmom(const amrex::Box &bxz, const int &dir, const amrex::Array4< amrex::Real > &rho_w_rhs, const amrex::Array4< const amrex::Real > &w, 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 > &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 int domhi_z, const bool do_lo=false)
@ Centered_4th
@ Centered_6th
@ Centered_2nd
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:266
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:221
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:310
@ cons_bc
Definition: ERF_IndexDefines.H:76
@ open
Definition: ERF_IndexDefines.H:186
Here is the call graph for this function: