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

Functions

void MomentumToVelocity (MultiFab &xvel, MultiFab &yvel, MultiFab &zvel, const MultiFab &density, const MultiFab &xmom_in, const MultiFab &ymom_in, const MultiFab &zmom_in, const Box &domain, const Vector< BCRec > &domain_bcs_type_h, const MultiFab *c_vfrac)
 

Function Documentation

◆ MomentumToVelocity()

void MomentumToVelocity ( MultiFab &  xvel,
MultiFab &  yvel,
MultiFab &  zvel,
const MultiFab &  density,
const MultiFab &  xmom_in,
const MultiFab &  ymom_in,
const MultiFab &  zmom_in,
const Box &  domain,
const Vector< BCRec > &  domain_bcs_type_h,
const MultiFab *  c_vfrac 
)

Convert momentum to velocity by dividing by density averaged onto faces

Parameters
[out]xvelx-component of velocity
[out]yvely-component of velocity
[out]zvelz-component of velocity
[in]densitydensity at cell centers
[in]xmom_inx-component of momentum
[in]ymom_iny-component of momentum
[in]zmom_inz-component of momentum
[in]domainDomain at this level
[in]domain_bcs_type_hhost vector for domain boundary conditions
31  {
32  BL_PROFILE_VAR("MomentumToVelocity()",MomentumToVelocity);
33 
34  const BCRec* bc_ptr_h = domain_bcs_type_h.data();
35 
36 #ifdef _OPENMP
37 #pragma omp parallel if (amrex::Gpu::notInLaunchRegion())
38 #endif
39  for ( MFIter mfi(density,TilingIfNotGPU()); mfi.isValid(); ++mfi)
40  {
41  // We need velocity in the interior ghost cells (init == real)
42  Box bx = mfi.tilebox();
43 
44  const Box& tbx = surroundingNodes(bx,0);
45  const Box& tby = surroundingNodes(bx,1);
46  const Box& tbz = surroundingNodes(bx,2);
47 
48  // Conserved variables on cell centers -- we use this for density
49  const Array4<const Real>& dens_arr = density.array(mfi);
50 
51  // Momentum on faces
52  Array4<Real const> const& momx = xmom_in.const_array(mfi);
53  Array4<Real const> const& momy = ymom_in.const_array(mfi);
54  Array4<Real const> const& momz = zmom_in.const_array(mfi);
55 
56  // Velocity on faces
57  const Array4<Real>& velx = xvel.array(mfi);
58  const Array4<Real>& vely = yvel.array(mfi);
59  const Array4<Real>& velz = zvel.array(mfi);
60 
61  if (c_vfrac==nullptr) {
62  ParallelFor(tbx, tby, tbz,
63  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
64  velx(i,j,k) = momx(i,j,k) * 2.0 / (dens_arr(i,j,k,Rho_comp) + dens_arr(i-1,j,k,Rho_comp));
65  },
66  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
67  vely(i,j,k) = momy(i,j,k) * 2.0 / (dens_arr(i,j,k,Rho_comp) + dens_arr(i,j-1,k,Rho_comp));
68  },
69  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
70  velz(i,j,k) = momz(i,j,k) * 2.0 / (dens_arr(i,j,k,Rho_comp) + dens_arr(i,j,k-1,Rho_comp));
71  });
72  } else {
73  // EB
74  const Array4<const Real>& c_vfrac_arr = c_vfrac->const_array(mfi);
75 
76  ParallelFor(tbx, tby, tbz,
77  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
78  if (c_vfrac_arr(i,j,k) > 0.0 || c_vfrac_arr(i-1,j,k) > 0.0) {
79  Real rho = (c_vfrac_arr(i,j,k) * dens_arr(i,j,k,Rho_comp)
80  + c_vfrac_arr(i-1,j,k) * dens_arr(i-1,j,k,Rho_comp))
81  / (c_vfrac_arr(i,j,k) + c_vfrac_arr(i-1,j,k));
82  velx(i,j,k) = momx(i,j,k) / rho;
83  }
84  },
85  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
86  if (c_vfrac_arr(i,j,k) > 0.0 || c_vfrac_arr(i,j-1,k) > 0.0) {
87  Real rho = (c_vfrac_arr(i,j,k) * dens_arr(i,j,k,Rho_comp)
88  + c_vfrac_arr(i,j-1,k) * dens_arr(i,j-1,k,Rho_comp))
89  / (c_vfrac_arr(i,j,k) + c_vfrac_arr(i,j-1,k));
90  vely(i,j,k) = momy(i,j,k) / rho;
91  }
92  },
93  [=] AMREX_GPU_DEVICE (int i, int j, int k) {
94  if (c_vfrac_arr(i,j,k) > 0.0 || c_vfrac_arr(i,j,k-1) > 0.0) {
95  Real rho = (c_vfrac_arr(i,j,k) * dens_arr(i,j,k,Rho_comp)
96  + c_vfrac_arr(i,j,k-1) * dens_arr(i,j,k-1,Rho_comp))
97  / (c_vfrac_arr(i,j,k) + c_vfrac_arr(i,j,k-1));
98  velz(i,j,k) = momz(i,j,k) / rho;
99  }
100  });
101  }
102 
103  if (bx.smallEnd(0) == domain.smallEnd(0)) {
104  if (bc_ptr_h[BCVars::cons_bc].lo(0) == ERFBCType::ext_dir)
105  {
106  ParallelFor(makeSlab(tbx,0,domain.smallEnd(0)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
107  velx(i,j,k) = momx(i,j,k) / dens_arr(i-1,j,k,Rho_comp);
108  });
109  }
110  else if (bc_ptr_h[BCVars::cons_bc].lo(0) == ERFBCType::ext_dir_upwind)
111  {
112  ParallelFor(makeSlab(tbx,0,domain.smallEnd(0)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
113  if (momx(i,j,k) >= 0.) {
114  velx(i,j,k) = momx(i,j,k) / dens_arr(i-1,j,k,Rho_comp);
115  }
116  });
117  }
118  }
119 
120  if (bx.bigEnd(0) == domain.bigEnd(0)) {
121  if (bc_ptr_h[BCVars::cons_bc].hi(0) == ERFBCType::ext_dir)
122  {
123  ParallelFor(makeSlab(tbx,0,domain.bigEnd(0)+1), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
124  velx(i,j,k) = momx(i,j,k) / dens_arr(i,j,k,Rho_comp);
125  });
126  }
127  else if (bc_ptr_h[BCVars::cons_bc].hi(0) == ERFBCType::ext_dir_upwind)
128  {
129  ParallelFor(makeSlab(tbx,0,domain.smallEnd(0)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
130  if (momx(i,j,k) <= 0.) {
131  velx(i,j,k) = momx(i,j,k) / dens_arr(i,j,k,Rho_comp);
132  }
133  });
134  }
135  }
136 
137  if (bx.smallEnd(1) == domain.smallEnd(1)) {
138  if (bc_ptr_h[BCVars::cons_bc].lo(1) == ERFBCType::ext_dir)
139  {
140  ParallelFor(makeSlab(tby,1,domain.smallEnd(1)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
141  vely(i,j,k) = momy(i,j,k) / dens_arr(i,j-1,k,Rho_comp);
142  });
143  }
144  else if (bc_ptr_h[BCVars::cons_bc].lo(1) == ERFBCType::ext_dir_upwind)
145  {
146  ParallelFor(makeSlab(tby,1,domain.smallEnd(1)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
147  if (momy(i,j,k) >= 0.) {
148  vely(i,j,k) = momy(i,j,k) / dens_arr(i,j-1,k,Rho_comp);
149  }
150  });
151  }
152  }
153 
154  if (bx.bigEnd(1) == domain.bigEnd(1)) {
155  if (bc_ptr_h[BCVars::cons_bc].hi(1) == ERFBCType::ext_dir)
156  {
157  ParallelFor(makeSlab(tby,1,domain.bigEnd(1)+1), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
158  vely(i,j,k) = momy(i,j,k) / dens_arr(i,j,k,Rho_comp);
159  });
160  }
161  else if (bc_ptr_h[BCVars::cons_bc].hi(1) == ERFBCType::ext_dir_upwind)
162  {
163  ParallelFor(makeSlab(tby,1,domain.smallEnd(1)), [=] AMREX_GPU_DEVICE (int i, int j, int k) {
164  if (momy(i,j,k) <= 0.) {
165  vely(i,j,k) = momy(i,j,k) / dens_arr(i,j,k,Rho_comp);
166  }
167  });
168  }
169  }
170  } // end MFIter
171 }
#define Rho_comp
Definition: ERF_IndexDefines.H:36
void MomentumToVelocity(MultiFab &xvel, MultiFab &yvel, MultiFab &zvel, const MultiFab &density, const MultiFab &xmom_in, const MultiFab &ymom_in, const MultiFab &zmom_in, const Box &domain, const Vector< BCRec > &domain_bcs_type_h, const MultiFab *c_vfrac)
Definition: ERF_MomentumToVelocity.cpp:25
amrex::Real Real
Definition: ERF_ShocInterface.H:19
@ cons_bc
Definition: ERF_IndexDefines.H:76
@ ext_dir
Definition: ERF_IndexDefines.H:209
@ ext_dir_upwind
Definition: ERF_IndexDefines.H:217
@ rho
Definition: ERF_Kessler.H:22
@ xvel
Definition: ERF_IndexDefines.H:141
@ zvel
Definition: ERF_IndexDefines.H:143
@ yvel
Definition: ERF_IndexDefines.H:142

Referenced by ERF::AverageDownTo(), ERF::FillCoarsePatch(), ERF::FillIntermediatePatch(), ERF::FillPatchFineLevel(), and ERF::project_velocity().

Here is the caller graph for this function: