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

Functions

void AdvectionSrcForMom_TF (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 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, const int lo_z_face, const int hi_z_face)
 

Function Documentation

◆ AdvectionSrcForMom_TF()

void AdvectionSrcForMom_TF ( 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 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,
const int  lo_z_face,
const int  hi_z_face 
)

Function for computing the advective tendency for the momentum equations ONLY when using 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]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
[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
62 {
63  BL_PROFILE_VAR("AdvectionSrcForMom_TF", AdvectionSrcForMom_TF);
64 
65  AMREX_ALWAYS_ASSERT(bxz.smallEnd(2) > 0);
66 
67  auto dxInv = cellSizeInv[0], dyInv = cellSizeInv[1], dzInv = cellSizeInv[2];
68 
69  // compute mapfactor inverses
70  Box box2d_u(bxx); box2d_u.setRange(2,0); box2d_u.grow({3,3,0});
71  Box box2d_v(bxy); box2d_v.setRange(2,0); box2d_v.grow({3,3,0});
72  FArrayBox mf_u_invFAB(box2d_u,1,The_Async_Arena());
73  FArrayBox mf_v_invFAB(box2d_v,1,The_Async_Arena());
74  const Array4<Real>& mf_u_inv = mf_u_invFAB.array();
75  const Array4<Real>& mf_v_inv = mf_v_invFAB.array();
76 
77  ParallelFor(box2d_u, box2d_v,
78  [=] AMREX_GPU_DEVICE (int i, int j, int) noexcept
79  {
80  mf_u_inv(i,j,0) = 1. / mf_u(i,j,0);
81  },
82  [=] AMREX_GPU_DEVICE (int i, int j, int) noexcept
83  {
84  mf_v_inv(i,j,0) = 1. / mf_v(i,j,0);
85  });
86 
87  // Inline with 2nd order for efficiency
88  if (horiz_adv_type == AdvType::Centered_2nd && vert_adv_type == AdvType::Centered_2nd)
89  {
90  ParallelFor(bxx, bxy, bxz,
91  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
92  {
93  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)) *
94  (u(i+1,j,k) + u(i,j,k)) * 0.5 * (ax(i,j,k) + ax(i+1,j,k));
95 
96  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)) *
97  (u(i-1,j,k) + u(i,j,k)) * 0.5 * (ax(i,j,k) + ax(i-1,j,k));
98 
99  Real met_h_zeta_yhi = Compute_h_zeta_AtEdgeCenterK(i,j+1,k,cellSizeInv,z_nd);
100  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)) *
101  (u(i,j+1,k) + u(i,j,k)) * met_h_zeta_yhi;
102 
103  Real met_h_zeta_ylo = Compute_h_zeta_AtEdgeCenterK(i,j ,k,cellSizeInv,z_nd);
104  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)) *
105  (u(i,j-1,k) + u(i,j,k)) * met_h_zeta_ylo;
106 
107  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)) *
108  0.5 * (az(i,j,k+1) + az(i-1,j,k+1));
109  Real zflux_lo = 0.25 * (Omega(i,j,k ) + Omega(i-1,j,k )) * (u(i,j,k-1) + u(i,j,k)) *
110  0.5 * (az(i,j,k ) + az(i-1,j,k ));
111 
112  Real mfsq = mf_u(i,j,0) * mf_u(i,j,0);
113 
114  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
115  + (yflux_hi - yflux_lo) * dyInv * mfsq
116  + (zflux_hi - zflux_lo) * dzInv;
117  rho_u_rhs(i, j, k) = -advectionSrc / (0.5 * (detJ(i,j,k) + detJ(i-1,j,k)));
118  },
119  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
120  {
121 
122  Real met_h_zeta_xhi = Compute_h_zeta_AtEdgeCenterK(i+1,j,k,cellSizeInv,z_nd);
123  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)) *
124  (v(i+1,j,k) + v(i,j,k)) * met_h_zeta_xhi;
125 
126  Real met_h_zeta_xlo = Compute_h_zeta_AtEdgeCenterK(i ,j,k,cellSizeInv,z_nd);
127  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)) *
128  (v(i-1,j,k) + v(i,j,k)) * met_h_zeta_xlo;
129 
130  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)) *
131  (v(i,j+1,k) + v(i,j,k)) * 0.5 * (ay(i,j,k) + ay(i,j+1,k));
132 
133  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)) *
134  (v(i,j-1,k) + v(i,j,k)) * 0.5 * (ay(i,j,k) + ay(i,j-1,k));
135 
136  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)) *
137  0.5 * (az(i,j,k+1) + az(i,j-1,k+1));
138  Real zflux_lo = 0.25 * (Omega(i,j,k ) + Omega(i, j-1, k )) * (v(i,j,k-1) + v(i,j,k)) *
139  0.5 * (az(i,j,k ) + az(i,j-1,k ));
140 
141  Real mfsq = mf_v(i,j,0) * mf_v(i,j,0);
142 
143  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
144  + (yflux_hi - yflux_lo) * dyInv * mfsq
145  + (zflux_hi - zflux_lo) * dzInv;
146  rho_v_rhs(i, j, k) = -advectionSrc / (0.5 * (detJ(i,j,k) + detJ(i,j-1,k)));
147  },
148  [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
149  {
150  Real met_h_zeta_xhi = Compute_h_zeta_AtEdgeCenterJ(i+1,j ,k ,cellSizeInv,z_nd);
151  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) *
152  (w(i+1,j,k) + w(i,j,k)) * met_h_zeta_xhi;
153 
154  Real met_h_zeta_xlo = Compute_h_zeta_AtEdgeCenterJ(i ,j ,k ,cellSizeInv,z_nd);
155  Real xflux_lo = 0.25*(rho_u(i ,j ,k) + rho_u(i ,j,k-1)) * mf_u_inv(i ,j,0) *
156  (w(i-1,j,k) + w(i,j,k)) * met_h_zeta_xlo;
157 
158  Real met_h_zeta_yhi = Compute_h_zeta_AtEdgeCenterI(i ,j+1,k ,cellSizeInv,z_nd);
159  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) *
160  (w(i,j+1,k) + w(i,j,k)) * met_h_zeta_yhi;
161 
162  Real met_h_zeta_ylo = Compute_h_zeta_AtEdgeCenterI(i ,j ,k ,cellSizeInv,z_nd);
163  Real yflux_lo = 0.25*(rho_v(i,j ,k) + rho_v(i,j ,k-1)) * mf_v_inv(i,j ,0) *
164  (w(i,j-1,k) + w(i,j,k)) * met_h_zeta_ylo;
165 
166  Real zflux_lo = 0.25 * (Omega(i,j,k) + Omega(i,j,k-1)) * (w(i,j,k) + w(i,j,k-1));
167 
168  Real zflux_hi = (k == hi_z_face) ? Omega(i,j,k) * w(i,j,k) * az(i,j,k):
169  0.25 * (Omega(i,j,k) + Omega(i,j,k+1)) * (w(i,j,k) + w(i,j,k+1)) *
170  0.5 * (az(i,j,k) + az(i,j,k+1));
171 
172  Real mfsq = mf_m(i,j,0) * mf_m(i,j,0);
173 
174  Real advectionSrc = (xflux_hi - xflux_lo) * dxInv * mfsq
175  + (yflux_hi - yflux_lo) * dyInv * mfsq
176  + (zflux_hi - zflux_lo) * dzInv;
177  rho_w_rhs(i, j, k) = -advectionSrc / (0.5*(detJ(i,j,k) + detJ(i,j,k-1)));
178  });
179 
180  // Template higher order methods
181  } else {
182  if (horiz_adv_type == AdvType::Centered_2nd) {
183  AdvectionSrcForMomVert<CENTERED2>(bxx, bxy, bxz,
184  rho_u_rhs, rho_v_rhs, rho_w_rhs,
185  rho_u, rho_v, Omega, u, v, w,
186  z_nd, ax, ay, az, detJ,
187  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
188  horiz_upw_frac, vert_upw_frac,
189  vert_adv_type, lo_z_face, hi_z_face);
190  } else if (horiz_adv_type == AdvType::Upwind_3rd) {
191  AdvectionSrcForMomVert<UPWIND3>(bxx, bxy, bxz,
192  rho_u_rhs, rho_v_rhs, rho_w_rhs,
193  rho_u, rho_v, Omega, u, v, w,
194  z_nd, ax, ay, az, detJ,
195  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
196  horiz_upw_frac, vert_upw_frac,
197  vert_adv_type, lo_z_face, hi_z_face);
198  } else if (horiz_adv_type == AdvType::Centered_4th) {
199  AdvectionSrcForMomVert<CENTERED4>(bxx, bxy, bxz,
200  rho_u_rhs, rho_v_rhs, rho_w_rhs,
201  rho_u, rho_v, Omega, u, v, w,
202  z_nd, ax, ay, az, detJ,
203  cellSizeInv, mf_m, mf_u_inv, mf_v_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::Upwind_5th) {
207  AdvectionSrcForMomVert<UPWIND5>(bxx, bxy, bxz,
208  rho_u_rhs, rho_v_rhs, rho_w_rhs,
209  rho_u, rho_v, Omega, u, v, w,
210  z_nd, ax, ay, az, detJ,
211  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
212  horiz_upw_frac, vert_upw_frac,
213  vert_adv_type, lo_z_face, hi_z_face);
214  } else if (horiz_adv_type == AdvType::Centered_6th) {
215  AdvectionSrcForMomVert<CENTERED6>(bxx, bxy, bxz,
216  rho_u_rhs, rho_v_rhs, rho_w_rhs,
217  rho_u, rho_v, Omega, u, v, w,
218  z_nd, ax, ay, az, detJ,
219  cellSizeInv, mf_m, mf_u_inv, mf_v_inv,
220  horiz_upw_frac, vert_upw_frac,
221  vert_adv_type, lo_z_face, hi_z_face);
222  } else {
223  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown advection scheme!");
224  }
225  }
226 }
void AdvectionSrcForMom_TF(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 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, const int lo_z_face, const int hi_z_face)
Definition: ERF_AdvectionSrcForMom_TF.cpp:38
@ 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
Here is the call graph for this function: