Skip to content

Commit

Permalink
fix some more benign bugs
Browse files Browse the repository at this point in the history
  • Loading branch information
isriva committed Sep 19, 2023
1 parent 9fdc958 commit 32cf520
Show file tree
Hide file tree
Showing 5 changed files with 88 additions and 95 deletions.
6 changes: 0 additions & 6 deletions src_compressible_stag/boundaryStag.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -484,7 +484,6 @@ void BCWallReservoirFluxStag(std::array< MultiFab, AMREX_SPACEDIM >& faceflux,
});
}
if ((bc_mass_lo[0] == 4) && (bx.smallEnd(0) < dom.smallEnd(0))) {
int lo = dom.smallEnd(0);
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (i < dom.smallEnd(0)) cenx(i,j,k) = cenx(i+1,j,k);
Expand All @@ -499,7 +498,6 @@ void BCWallReservoirFluxStag(std::array< MultiFab, AMREX_SPACEDIM >& faceflux,
});
}
if ((bc_mass_hi[0] == 4) && (bx.bigEnd(0) > dom.bigEnd(0))) {
int hi = dom.bigEnd(0);
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (i > dom.bigEnd(0)) cenx(i,j,k) = cenx(i-1,j,k);
Expand All @@ -514,7 +512,6 @@ void BCWallReservoirFluxStag(std::array< MultiFab, AMREX_SPACEDIM >& faceflux,
});
}
if ((bc_mass_lo[1] == 4) && (bx.smallEnd(1) < dom.smallEnd(1))) {
int lo = dom.smallEnd(1);
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (j < dom.smallEnd(1)) ceny(i,j,k) = ceny(i,j+1,k);
Expand All @@ -529,7 +526,6 @@ void BCWallReservoirFluxStag(std::array< MultiFab, AMREX_SPACEDIM >& faceflux,
});
}
if ((bc_mass_hi[1] == 4) && (bx.bigEnd(1) > dom.bigEnd(1))) {
int hi = dom.bigEnd(1);
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (j > dom.bigEnd(1)) ceny(i,j,k) = ceny(i,j-1,k);
Expand All @@ -544,7 +540,6 @@ void BCWallReservoirFluxStag(std::array< MultiFab, AMREX_SPACEDIM >& faceflux,
});
}
if ((bc_mass_lo[2] == 4) && (bx.smallEnd(2) < dom.smallEnd(2))) {
int lo = dom.smallEnd(2);
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (k < dom.smallEnd(2)) cenz(i,j,k) = cenz(i,j,k+1);;
Expand All @@ -559,7 +554,6 @@ void BCWallReservoirFluxStag(std::array< MultiFab, AMREX_SPACEDIM >& faceflux,
});
}
if ((bc_mass_hi[2] == 4) && (bx.bigEnd(2) > dom.bigEnd(2))) {
int hi = dom.bigEnd(2);
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (k > dom.bigEnd(2)) cenz(i,j,k) = cenz(i,j,k-1);;
Expand Down
12 changes: 6 additions & 6 deletions src_compressible_stag/compressible_functions_stag.H
Original file line number Diff line number Diff line change
Expand Up @@ -165,12 +165,12 @@ ResetReservoirMom(std::array<MultiFab, AMREX_SPACEDIM>& cumom,
const std::array<MultiFab, AMREX_SPACEDIM>& cumom_res,
const amrex::Geometry& geom);

void
ReFluxCons(MultiFab& cu, const MultiFab& cu0,
const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_res,
const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_cont,
const amrex::Geometry& geom,
const amrex::Real dt);
//void
//ReFluxCons(MultiFab& cu, const MultiFab& cu0,
// const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_res,
// const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_cont,
// const amrex::Geometry& geom,
// const amrex::Real dt);


void WriteCheckPoint3D(int step,
Expand Down
151 changes: 75 additions & 76 deletions src_compressible_stag/reservoirStag.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -22,9 +22,6 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
Box dom(geom.Domain());

Real N_A = Runiv/k_B; // Avagadro's number
Real vol = dx[0]*dx[1]*dx[2];
Real PI = 4.0*atan(1.0);
Real sqrtPI = sqrt(PI);

GpuArray<Real,MAX_SPECIES> mass;
for (int l=0;l<nspecies;++l) {
Expand All @@ -36,11 +33,13 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
faceflux_res[d].setVal(0.0);
}

Real area, vol;

// Reservoir in LO X
if (bc_mass_lo[0] == 4) {

Real area = dx[1]*dx[2];
Real vol = dx[0]*dx[1]*dx[2];
area = dx[1]*dx[2];
vol = dx[0]*dx[1]*dx[2];

// face-based flux (mass and energy) and normal momentum /////
//////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -177,8 +176,8 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
// Reservoir in HI X
if (bc_mass_hi[0] == 4) {

Real area = dx[1]*dx[2];
Real vol = dx[0]*dx[1]*dx[2];
area = dx[1]*dx[2];
vol = dx[0]*dx[1]*dx[2];

// face-based flux (mass and energy) and normal momentum /////
//////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -316,8 +315,8 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
// Reservoir in LO Y
if (bc_mass_lo[1] == 4) {

Real area = dx[0]*dx[2];
Real vol = dx[0]*dx[1]*dx[2];
area = dx[0]*dx[2];
vol = dx[0]*dx[1]*dx[2];

// face-based flux (mass and energy) and normal momentum /////
//////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -454,8 +453,8 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
// Reservoir in HI Y
if (bc_mass_hi[1] == 4) {

Real area = dx[0]*dx[2];
Real vol = dx[0]*dx[1]*dx[2];
area = dx[0]*dx[2];
vol = dx[0]*dx[1]*dx[2];

// face-based flux (mass and energy) and normal momentum /////
//////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -593,8 +592,8 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
// Reservoir in LO Z
if (bc_mass_lo[2] == 4) {

Real area = dx[0]*dx[1];
Real vol = dx[0]*dx[1]*dx[2];
area = dx[0]*dx[1];
vol = dx[0]*dx[1]*dx[2];

// face-based flux (mass and energy) and normal momentum /////
//////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -731,8 +730,8 @@ ComputeFluxMomReservoir(const MultiFab& cons0_in, const MultiFab& prim0_in,
// Reservoir in HI Z
if (bc_mass_hi[2] == 4) {

Real area = dx[0]*dx[1];
Real vol = dx[0]*dx[1]*dx[2];
area = dx[0]*dx[1];
vol = dx[0]*dx[1]*dx[2];

// face-based flux (mass and energy) and normal momentum /////
//////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -1565,65 +1564,65 @@ ResetReservoirMom(std::array<MultiFab, AMREX_SPACEDIM>& cumom,
// Reflux conserved qtys at the cell next to reservoir ////////////////////
// not used in the current implementation /////////////////////////////////
///////////////////////////////////////////////////////////////////////////
void
ReFluxCons(MultiFab& cu, const MultiFab& cu0,
const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_res,
const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_cont,
const amrex::Geometry& geom,
const amrex::Real dt)
{
BL_PROFILE_VAR("ReFluxCons()",ReFluxCons);

const GpuArray<Real, AMREX_SPACEDIM> dx = geom.CellSizeArray();
Box dom(geom.Domain());

for ( MFIter mfi(cu0); mfi.isValid(); ++mfi) {

const Box& bx = mfi.validbox();

AMREX_D_TERM(Array4<Real const> const& xflux_res = faceflux_res[0].array(mfi);,
Array4<Real const> const& yflux_res = faceflux_res[1].array(mfi);,
Array4<Real const> const& zflux_res = faceflux_res[2].array(mfi););
AMREX_D_TERM(Array4<Real const> const& xflux_cont = faceflux_cont[0].array(mfi);,
Array4<Real const> const& yflux_cont = faceflux_cont[1].array(mfi);,
Array4<Real const> const& zflux_cont = faceflux_cont[2].array(mfi););

const Array4<Real>& cons = cu.array(mfi);
const Array4<const Real>& cons0 = cu0.array(mfi);

// Reservoir in LO X
if ((bc_mass_lo[0] == 4) and (bx.smallEnd(0) <= dom.smallEnd(0))) {
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (i == dom.smallEnd(0)) {
cons(i,j,k,0) = cons0(i,j,k,0)
- (dt/dx[0])*(xflux_cont(i,j,k,0) - xflux_res(i,j,k,0)); // correct density
cons(i,j,k,4) = cons0(i,j,k,4)
- (dt/dx[0])*(xflux_cont(i,j,k,4) - xflux_res(i,j,k,4)); // correct en. density
for (int n=0;n<nspecies;++n) {
cons(i,j,k,n+5) = cons0(i,j,k,n+5)
- (dt/dx[0])*(xflux_cont(i,j,k,n+5) - xflux_res(i,j,k,n+5)); // correct species
}
}
});
}

// Reservoir in HI X
if ((bc_mass_hi[0] == 4) and (bx.bigEnd(0) >= dom.bigEnd(0))) {
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
if (i == dom.bigEnd(0)) {
cons(i,j,k,0) = cons0(i,j,k,0)
+ (dt/dx[0])*(xflux_cont(i+1,j,k,0) - xflux_res(i+1,j,k,0)); // correct density
cons(i,j,k,4) = cons0(i,j,k,4)
+ (dt/dx[0])*(xflux_cont(i+1,j,k,4) - xflux_res(i+1,j,k,4)); // correct en. density
for (int n=0;n<nspecies;++n) {
cons(i,j,k,n+5) = cons0(i,j,k,n+5)
+ (dt/dx[0])*(xflux_cont(i+1,j,k,n+5) - xflux_res(i+1,j,k,n+5)); // correct species
}
}
});
}
}
}
//void
//ReFluxCons(MultiFab& cu, const MultiFab& cu0,
// const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_res,
// const std::array<MultiFab, AMREX_SPACEDIM>& faceflux_cont,
// const amrex::Geometry& geom,
// const amrex::Real dt)
//{
// BL_PROFILE_VAR("ReFluxCons()",ReFluxCons);
//
// const GpuArray<Real, AMREX_SPACEDIM> dx = geom.CellSizeArray();
// Box dom(geom.Domain());
//
// for ( MFIter mfi(cu0); mfi.isValid(); ++mfi) {
//
// const Box& bx = mfi.validbox();
//
// AMREX_D_TERM(Array4<Real const> const& xflux_res = faceflux_res[0].array(mfi);,
// Array4<Real const> const& yflux_res = faceflux_res[1].array(mfi);,
// Array4<Real const> const& zflux_res = faceflux_res[2].array(mfi););
// AMREX_D_TERM(Array4<Real const> const& xflux_cont = faceflux_cont[0].array(mfi);,
// Array4<Real const> const& yflux_cont = faceflux_cont[1].array(mfi);,
// Array4<Real const> const& zflux_cont = faceflux_cont[2].array(mfi););
//
// const Array4<Real>& cons = cu.array(mfi);
// const Array4<const Real>& cons0 = cu0.array(mfi);
//
// // Reservoir in LO X
// if ((bc_mass_lo[0] == 4) and (bx.smallEnd(0) <= dom.smallEnd(0))) {
// amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
// {
// if (i == dom.smallEnd(0)) {
// cons(i,j,k,0) = cons0(i,j,k,0)
// - (dt/dx[0])*(xflux_cont(i,j,k,0) - xflux_res(i,j,k,0)); // correct density
// cons(i,j,k,4) = cons0(i,j,k,4)
// - (dt/dx[0])*(xflux_cont(i,j,k,4) - xflux_res(i,j,k,4)); // correct en. density
// for (int n=0;n<nspecies;++n) {
// cons(i,j,k,n+5) = cons0(i,j,k,n+5)
// - (dt/dx[0])*(xflux_cont(i,j,k,n+5) - xflux_res(i,j,k,n+5)); // correct species
// }
// }
// });
// }
//
// // Reservoir in HI X
// if ((bc_mass_hi[0] == 4) and (bx.bigEnd(0) >= dom.bigEnd(0))) {
// amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
// {
// if (i == dom.bigEnd(0)) {
// cons(i,j,k,0) = cons0(i,j,k,0)
// + (dt/dx[0])*(xflux_cont(i+1,j,k,0) - xflux_res(i+1,j,k,0)); // correct density
// cons(i,j,k,4) = cons0(i,j,k,4)
// + (dt/dx[0])*(xflux_cont(i+1,j,k,4) - xflux_res(i+1,j,k,4)); // correct en. density
// for (int n=0;n<nspecies;++n) {
// cons(i,j,k,n+5) = cons0(i,j,k,n+5)
// + (dt/dx[0])*(xflux_cont(i+1,j,k,n+5) - xflux_res(i+1,j,k,n+5)); // correct species
// }
// }
// });
// }
// }
//}

8 changes: 4 additions & 4 deletions src_compressible_stag/reservoirStag_K.H
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,6 @@ random_cross_vel(const amrex::Real& a, const amrex::Real& vT, bool tangential, a
}
}
else { // tangential velocity crossing
amrex::Real v;
v = Vtangential + vT*amrex::RandomNormal(0.0,sqrt(0.5),engine);
}
return v;
Expand Down Expand Up @@ -120,26 +119,27 @@ poisson_process_reservoir(const amrex::GpuArray<amrex::Real,MAX_SPECIES>& mass,
while(!stop) {

// add a particle for each species (if stop[n] = false)
amrex::Real vel_rand;
for (int n=0;n<nspec;++n) {
if (!stop_flag[n]) {

// normal component
mass_cross += mass[n]; // add a particle crossing for total mass
spec_mass_cross[n] += mass[n]; // add a particle crossing for species mass
amrex::Real vel_rand = random_cross_vel(speed_ratio[n],thermal_speed[n],false,0.0,engine); // get random velocity
vel_rand = random_cross_vel(speed_ratio[n],thermal_speed[n],false,0.0,engine); // get random velocity
mom_cross[0] += mass[n]*vel_rand; // add a particle crossing for total momentum
en_cross += 0.5*mass[n]*(vel_rand)*(vel_rand); // add a particle crossing for total energy

// tangential component 1
if (dim >= 2) {
amrex::Real vel_rand = random_cross_vel(speed_ratio[n],thermal_speed[n],true,VT1,engine); // get random velocity
vel_rand = random_cross_vel(speed_ratio[n],thermal_speed[n],true,VT1,engine); // get random velocity
mom_cross[1] += mass[n]*vel_rand; // add a particle crossing for total momentum
en_cross += 0.5*mass[n]*(vel_rand)*(vel_rand); // add a particle crossing for total energy
}

// tangential component 2
if (dim == 3) {
amrex::Real vel_rand = random_cross_vel(speed_ratio[n],thermal_speed[n],true,VT2,engine); // get random velocity
vel_rand = random_cross_vel(speed_ratio[n],thermal_speed[n],true,VT2,engine); // get random velocity
mom_cross[2] += mass[n]*vel_rand; // add a particle crossing for total momentum
en_cross += 0.5*mass[n]*(vel_rand)*(vel_rand); // add a particle crossing for total energy
}
Expand Down
6 changes: 3 additions & 3 deletions src_compressible_stag/timeStepStag.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -360,8 +360,8 @@ void RK3stepStag(MultiFab& cu,
compute_chemistry_source_CLE(dt, dx[0]*dx[1]*dx[2], prim, source, ranchem);
}

amrex::Real energy_in = amrex::Real(0.0);
#if defined(TURB)
amrex::Real energy_in = amrex::Real(0.0);
if (turbForcing == 2) { // random forcing tubulence : get average energy input
ReduceOps<ReduceOpSum> reduce_op;
ReduceData<Real> reduce_data(reduce_op);
Expand Down Expand Up @@ -687,8 +687,8 @@ void RK3stepStag(MultiFab& cu,
compute_chemistry_source_CLE(dt, dx[0]*dx[1]*dx[2], prim, source, ranchem);
}

amrex::Real energyp_in = amrex::Real(0.0);
#if defined(TURB)
amrex::Real energyp_in = amrex::Real(0.0);
if (turbForcing == 2) { // random forcing tubulence : get average energy input
ReduceOps<ReduceOpSum> reduce_op;
ReduceData<Real> reduce_data(reduce_op);
Expand Down Expand Up @@ -1021,8 +1021,8 @@ void RK3stepStag(MultiFab& cu,
compute_chemistry_source_CLE(dt, dx[0]*dx[1]*dx[2], prim, source, ranchem);
}

amrex::Real energyp2_in = amrex::Real(0.0);
#if defined(TURB)
amrex::Real energyp2_in = amrex::Real(0.0);
if (turbForcing == 2) { // random forcing tubulence : get average energy input
ReduceOps<ReduceOpSum> reduce_op;
ReduceData<Real> reduce_data(reduce_op);
Expand Down

0 comments on commit 32cf520

Please sign in to comment.