34 const BCRec* bc_ptr_h = domain_bcs_type_h.data();
37 #pragma omp parallel if (amrex::Gpu::notInLaunchRegion())
39 for ( MFIter mfi(density,TilingIfNotGPU()); mfi.isValid(); ++mfi)
42 Box bx = mfi.tilebox();
44 const Box& tbx = surroundingNodes(bx,0);
45 const Box& tby = surroundingNodes(bx,1);
46 const Box& tbz = surroundingNodes(bx,2);
49 const Array4<const Real>& dens_arr = density.array(mfi);
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);
57 const Array4<Real>& velx =
xvel.array(mfi);
58 const Array4<Real>& vely =
yvel.array(mfi);
59 const Array4<Real>& velz =
zvel.array(mfi);
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));
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));
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));
74 const Array4<const Real>& c_vfrac_arr = c_vfrac->const_array(mfi);
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) {
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;
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) {
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;
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) {
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;
103 if (bx.smallEnd(0) == domain.smallEnd(0)) {
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);
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);
120 if (bx.bigEnd(0) == domain.bigEnd(0)) {
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);
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);
137 if (bx.smallEnd(1) == domain.smallEnd(1)) {
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);
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);
154 if (bx.bigEnd(1) == domain.bigEnd(1)) {
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);
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);
#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