ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Functions
ERF_VelocityToMomentum.cpp File Reference
#include <AMReX.H>
#include <AMReX_MultiFab.H>
#include <ERF_Utils.H>
#include <AMReX_MultiFabUtil.H>
Include dependency graph for ERF_VelocityToMomentum.cpp:

Functions

void VelocityToMomentum (const MultiFab &xvel_in, const IntVect &xvel_ngrow, const MultiFab &yvel_in, const IntVect &yvel_ngrow, const MultiFab &zvel_in, const IntVect &zvel_ngrow, const MultiFab &density, MultiFab &xmom, MultiFab &ymom, MultiFab &zmom, const Box &domain, const Vector< BCRec > &domain_bcs_type_h)
 

Function Documentation

◆ VelocityToMomentum()

void VelocityToMomentum ( const MultiFab &  xvel_in,
const IntVect &  xvel_ngrow,
const MultiFab &  yvel_in,
const IntVect &  yvel_ngrow,
const MultiFab &  zvel_in,
const IntVect &  zvel_ngrow,
const MultiFab &  density,
MultiFab &  xmom,
MultiFab &  ymom,
MultiFab &  zmom,
const Box &  domain,
const Vector< BCRec > &  domain_bcs_type_h 
)

Convert velocity to momentum by multiplying by density averaged onto faces.

Parameters
[in]xvel_inx-component of velocity
[in]xvel_ngrowhow many cells to grow the tilebox for the x-momentum
[in]yvel_iny-component of velocity
[in]yvel_ngrowhow many cells to grow the tilebox for the y-momentum
[in]zvel_inz-component of velocity
[in]zvel_ngrowhow many cells to grow the tilebox for the z-momentum
[in]densitydensity at cell centers
[out]xmomx-component of momentum
[out]ymomy-component of momentum
[out]zmomz-component of momentum
[in]domainDomain at this level
[in]domain_bcs_type_hhost vector for domain boundary conditions
[in]l_use_ndiffflag describing whether we will later add explicit numerical diffusion
38 {
39  BL_PROFILE_VAR("VelocityToMomentum()",VelocityToMomentum);
40 
41  const BCRec* bc_ptr_h = domain_bcs_type_h.data();
42 
43  // Loop over boxes = valid boxes grown by ngrow
44 #ifdef _OPENMP
45 #pragma omp parallel if (amrex::Gpu::notInLaunchRegion())
46 #endif
47  for ( MFIter mfi(density,TilingIfNotGPU()); mfi.isValid(); ++mfi)
48  {
49  // We need momentum in the interior ghost cells (init == real)
50  const Box& bx = mfi.tilebox();
51  Box tbx, tby, tbz;
52 
53  tbx = mfi.tilebox(IntVect(1,0,0),xvel_ngrow);
54  tby = mfi.tilebox(IntVect(0,1,0),yvel_ngrow);
55  tbz = mfi.tilebox(IntVect(0,0,1),zvel_ngrow);
56 
57  // Don't actually try to fill w above or below the domain
58  if (tbz.smallEnd(2) < domain.smallEnd(2)) tbz.setSmall(2,domain.smallEnd(2));
59  if (tbz.bigEnd(2) > domain.bigEnd(2)+1) tbz.setBig(2,domain.bigEnd(2)+1);
60 
61  // Conserved/state variables on cell centers -- we use this for density
62  const Array4<const Real>& dens_arr = density.array(mfi);
63 
64  // Momentum on faces, to be computed
65  Array4<Real> const& momx = xmom.array(mfi);
66  Array4<Real> const& momy = ymom.array(mfi);
67  Array4<Real> const& momz = zmom.array(mfi);
68 
69  // Velocity on faces, used in computation
70  const Array4<Real const>& velx = xvel_in.const_array(mfi);
71  const Array4<Real const>& vely = yvel_in.const_array(mfi);
72  const Array4<Real const>& velz = zvel_in.const_array(mfi);
73 
74  // ********************************************************************************************
75 
76  ParallelFor(tbx, tby, tbz,
77  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
78  momx(i,j,k) = velx(i,j,k) * 0.5 * (dens_arr(i,j,k,Rho_comp) + dens_arr(i-1,j,k,Rho_comp));
79  },
80  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
81  momy(i,j,k) = vely(i,j,k) * 0.5 * (dens_arr(i,j,k,Rho_comp) + dens_arr(i,j-1,k,Rho_comp));
82  },
83  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
84  momz(i,j,k) = velz(i,j,k) * 0.5 * (dens_arr(i,j,k,Rho_comp) + dens_arr(i,j,k-1,Rho_comp));
85  });
86 
87  // ********************************************************************************************
88 
89  if (bx.smallEnd(0) == domain.smallEnd(0)) {
90  if (bc_ptr_h[BCVars::cons_bc].lo(0) == ERFBCType::ext_dir)
91  {
92  ParallelFor(makeSlab(tbx,0,domain.smallEnd(0)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
93  momx(i,j,k) = velx(i,j,k) * dens_arr(i-1,j,k,Rho_comp);
94  });
95  }
96  else if (bc_ptr_h[BCVars::cons_bc].lo(0) == ERFBCType::ext_dir_upwind)
97  {
98  ParallelFor(makeSlab(tbx,0,domain.smallEnd(0)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
99  if (velx(i,j,k) >= 0.) {
100  momx(i,j,k) = velx(i,j,k) * dens_arr(i-1,j,k,Rho_comp);
101  }
102  });
103  }
104  }
105 
106  if (bx.bigEnd(0) == domain.bigEnd(0)) {
107  if (bc_ptr_h[BCVars::cons_bc].hi(0) == ERFBCType::ext_dir)
108  {
109  ParallelFor(makeSlab(tbx,0,domain.bigEnd(0)+1), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
110  momx(i,j,k) = velx(i,j,k) * dens_arr(i,j,k,Rho_comp);
111  });
112  }
113  else if (bc_ptr_h[BCVars::cons_bc].hi(0) == ERFBCType::ext_dir_upwind)
114  {
115  ParallelFor(makeSlab(tbx,0,domain.smallEnd(0)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
116  if (velx(i,j,k) <= 0.) {
117  momx(i,j,k) = velx(i,j,k) * dens_arr(i,j,k,Rho_comp);
118  }
119  });
120  }
121  }
122 
123  if (bx.smallEnd(1) == domain.smallEnd(1)) {
124  if (bc_ptr_h[BCVars::cons_bc].lo(1) == ERFBCType::ext_dir)
125  {
126  ParallelFor(makeSlab(tby,1,domain.smallEnd(1)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
127  momy(i,j,k) = vely(i,j,k) / dens_arr(i,j-1,k,Rho_comp);
128  });
129  }
130  else if (bc_ptr_h[BCVars::cons_bc].lo(1) == ERFBCType::ext_dir_upwind)
131  {
132  ParallelFor(makeSlab(tby,1,domain.smallEnd(1)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
133  if (vely(i,j,k) >= 0.) {
134  momy(i,j,k) = vely(i,j,k) / dens_arr(i,j-1,k,Rho_comp);
135  }
136  });
137  }
138  }
139 
140  if (bx.bigEnd(1) == domain.bigEnd(1)) {
141  if (bc_ptr_h[BCVars::cons_bc].hi(1) == ERFBCType::ext_dir)
142  {
143  ParallelFor(makeSlab(tby,1,domain.bigEnd(1)+1), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
144  momy(i,j,k) = vely(i,j,k) / dens_arr(i,j,k,Rho_comp);
145  });
146  }
147  else if (bc_ptr_h[BCVars::cons_bc].hi(1) == ERFBCType::ext_dir_upwind)
148  {
149  ParallelFor(makeSlab(tby,1,domain.smallEnd(1)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
150  if (vely(i,j,k) <= 0.) {
151  momy(i,j,k) = vely(i,j,k) / dens_arr(i,j,k,Rho_comp);
152  }
153  });
154  }
155  }
156 
157  // Extrapolate momenta at top of domain (in case they need to be used in OmegaFromW)
158  if (bx.bigEnd(2) == domain.bigEnd(2))
159  {
160  ParallelFor(makeSlab(tbx,2,domain.bigEnd(2)+1), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
161  momx(i,j,k) = momx(i,j,k-1);
162  });
163  ParallelFor(makeSlab(tby,2,domain.bigEnd(2)+1), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
164  momy(i,j,k) = momy(i,j,k-1);
165  });
166  }
167  } // end MFIter
168 }
Rho_comp
#define Rho_comp
Definition: ERF_IndexDefines.H:36
VelocityToMomentum
void VelocityToMomentum(const MultiFab &xvel_in, const IntVect &xvel_ngrow, const MultiFab &yvel_in, const IntVect &yvel_ngrow, const MultiFab &zvel_in, const IntVect &zvel_ngrow, const MultiFab &density, MultiFab &xmom, MultiFab &ymom, MultiFab &zmom, const Box &domain, const Vector< BCRec > &domain_bcs_type_h)
Definition: ERF_VelocityToMomentum.cpp:28
BCVars::cons_bc
@ cons_bc
Definition: ERF_IndexDefines.H:76
ERFBCType::ext_dir
@ ext_dir
Definition: ERF_IndexDefines.H:209
ERFBCType::ext_dir_upwind
@ ext_dir_upwind
Definition: ERF_IndexDefines.H:217
IntVars::ymom
@ ymom
Definition: ERF_IndexDefines.H:160
IntVars::zmom
@ zmom
Definition: ERF_IndexDefines.H:161
IntVars::xmom
@ xmom
Definition: ERF_IndexDefines.H:159