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.
72 auto dxInv = cellSizeInv[0], dyInv = cellSizeInv[1], dzInv = cellSizeInv[2];
74 AMREX_ALWAYS_ASSERT(bxz.smallEnd(2) > 0);
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();
84 const bool use_terrain_fitted_coords = ( terrain_type == TerrainType::StaticFittedMesh ||
85 terrain_type == TerrainType::MovingFittedMesh);
87 ParallelFor(box2d_u, box2d_v,
88 [=] AMREX_GPU_DEVICE (
int i,
int j,
int) noexcept
90 mf_u_inv(i,j,0) = 1. / mf_u(i,j,0);
92 [=] AMREX_GPU_DEVICE (
int i,
int j,
int) noexcept
94 mf_v_inv(i,j,0) = 1. / mf_v(i,j,0);
97 if ( terrain_type == TerrainType::EB)
99 amrex::ignore_unused(use_terrain_fitted_coords);
100 ParallelFor(bxx, bxy, bxz,
101 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
103 rho_u_rhs(i, j, k) = 0.0;
105 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
107 rho_v_rhs(i, j, k) = 0.0;
109 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
111 rho_w_rhs(i, j, k) = 0.0;
114 if ( !use_terrain_fitted_coords) {
118 ParallelFor(bxx, bxy, bxz,
119 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
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));
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));
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));
130 Real mfsq = mf_u(i,j,0) * mf_u(i,j,0);
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;
137 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
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));
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));
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));
148 Real mfsq = mf_v(i,j,0) * mf_v(i,j,0);
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;
155 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
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));
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));
163 Real zflux_lo = 0.25 * (Omega(i,j,k) + Omega(i,j,k-1)) * (w(i,j,k) + w(i,j,k-1));
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));
168 Real mfsq = mf_m(i,j,0) * mf_m(i,j,0);
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;
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,
183 horiz_upw_frac, vert_upw_frac,
184 vert_adv_type, lo_z_face, hi_z_face);
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,
191 horiz_upw_frac, vert_upw_frac,
192 vert_adv_type, lo_z_face, hi_z_face);
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,
199 horiz_upw_frac, vert_upw_frac,
200 vert_adv_type, lo_z_face, hi_z_face);
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,
207 horiz_upw_frac, vert_upw_frac,
208 vert_adv_type, lo_z_face, hi_z_face);
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,
215 horiz_upw_frac, vert_upw_frac,
216 vert_adv_type, lo_z_face, hi_z_face);
218 AMREX_ASSERT_WITH_MESSAGE(
false,
"Unknown advection scheme!");
227 ParallelFor(bxx, bxy, bxz,
228 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
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));
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));
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;
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;
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 ));
249 Real mfsq = mf_u(i,j,0) * mf_u(i,j,0);
251 Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
252 + (yflux_hi - yflux_lo) * dyInv * mfsq
253 + (zflux_hi - zflux_lo) * dzInv;
255 rho_u_rhs(i, j, k) = -advectionSrc / (0.5 * (detJ(i,j,k) + detJ(i-1,j,k)));
257 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
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;
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;
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));
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));
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 ));
279 Real mfsq = mf_v(i,j,0) * mf_v(i,j,0);
281 Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
282 + (yflux_hi - yflux_lo) * dyInv * mfsq
283 + (zflux_hi - zflux_lo) * dzInv;
285 rho_v_rhs(i, j, k) = -advectionSrc / (0.5 * (detJ(i,j,k) + detJ(i,j-1,k)));
287 [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
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;
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;
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;
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;
305 Real zflux_lo = 0.25 * (Omega(i,j,k) + Omega(i,j,k-1)) * (w(i,j,k) + w(i,j,k-1));
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));
311 Real mfsq = mf_m(i,j,0) * mf_m(i,j,0);
313 Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
314 + (yflux_hi - yflux_lo) * dyInv * mfsq
315 + (zflux_hi - zflux_lo) * dzInv;
317 rho_w_rhs(i, j, k) = -advectionSrc / (0.5*(detJ(i,j,k) + detJ(i,j,k-1)));
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);
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);
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);
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);
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);
357 AMREX_ASSERT_WITH_MESSAGE(
false,
"Unknown advection scheme!");
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);
375 const int domhi_z = domain.bigEnd(2);
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));}
390 ay, az, detJ, cellSizeInv,
394 ax, ay, az, detJ, cellSizeInv,
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) );}
407 ay, az, detJ, cellSizeInv);
410 ax, ay, az, detJ, cellSizeInv,
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));}
423 ax, az, detJ, cellSizeInv,
428 ax, ay, az, detJ, cellSizeInv,
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) );}
440 ax, az, detJ, cellSizeInv);
444 ax, ay, az, detJ, cellSizeInv,
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)
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