Skip to content

Commit

Permalink
MFsurfchem - k_beta (temperature dependency) and e_beta (energy update)
Browse files Browse the repository at this point in the history 
compressible_stag_mui - plot for surface covariance
  • Loading branch information
@HyuntaeJung
HyuntaeJung committed May 27, 2024
1 parent 9ef8836 commit a6ac210
Show file tree
Hide file tree
Showing 7 changed files with 78 additions and 23 deletions.
15 changes: 10 additions & 5 deletions src_MFsurfchem/MFsurfchem_functions.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,8 @@ AMREX_GPU_MANAGED GpuArray<amrex::Real, MAX_SPECIES> MFsurfchem::des_rate;
AMREX_GPU_MANAGED int MFsurfchem::stoch_surfcov0;
AMREX_GPU_MANAGED int MFsurfchem::stoch_MFsurfchem;

// temperature correction exponent
#define BETA 0.5
AMREX_GPU_MANAGED amrex::Real MFsurfchem::k_beta;
AMREX_GPU_MANAGED amrex::Real MFsurfchem::e_beta;

void InitializeMFSurfchemNamespace()
{
Expand Down Expand Up @@ -63,6 +63,12 @@ void InitializeMFSurfchemNamespace()

stoch_MFsurfchem = 1; // default value
pp.query("stoch_MFsurfchem",stoch_MFsurfchem);

k_beta = 0.5; // default value
pp.query("k_beta",k_beta);

e_beta = 0.5; // default value
pp.query("e_beta",e_beta);
return;
}

Expand Down Expand Up @@ -146,7 +152,7 @@ void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab

amrex::Real theta = surfcov_arr(i,j,k,m);

amrex::Real meanNads = ads_rate_const[m]*dens*(1-sumtheta)*Ntot*dt*pow(tempratio,BETA);
amrex::Real meanNads = ads_rate_const[m]*dens*(1-sumtheta)*Ntot*dt*pow(tempratio,k_beta);
amrex::Real meanNdes = des_rate[m]*theta*Ntot*dt;

amrex::Real Nads;
Expand Down Expand Up @@ -193,8 +199,7 @@ void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab
for (int m=0;m<n_ads_spec;m++) {
amrex::Real dN = dNadsdes_arr(i,j,k,m);
amrex::Real factor1 = molmass[m]/AVONUM/(dx[0]*dx[1]*dx[2]);
amrex::Real factor2 = (BETA*k_B*T_inst+(e0[m]+hcv[m]*T_inst)*molmass[m]/AVONUM)/(dx[0]*dx[1]*dx[2]);

amrex::Real factor2 = (e_beta*k_B*T_inst+(e0[m]+hcv[m]*T_inst)*molmass[m]/AVONUM)/(dx[0]*dx[1]*dx[2]);
surfcov_arr(i,j,k,m) += dN/Ntot;
cu_arr(i,j,k,0) -= factor1*dN;
cu_arr(i,j,k,5+m) -= factor1*dN;
Expand Down
3 changes: 3 additions & 0 deletions src_MFsurfchem/MFsurfchem_namespace.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,4 +11,7 @@ namespace MFsurfchem {
extern AMREX_GPU_MANAGED int stoch_surfcov0;
extern AMREX_GPU_MANAGED int stoch_MFsurfchem;

extern AMREX_GPU_MANAGED amrex::Real k_beta;
extern AMREX_GPU_MANAGED amrex::Real e_beta;

}
7 changes: 7 additions & 0 deletions src_compressible_stag/Checkpoint.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,7 @@ void WriteCheckPoint3D(int step,
const amrex::MultiFab& surfcov,
const amrex::MultiFab& surfcovMeans,
const amrex::MultiFab& surfcovVars,
const amrex::MultiFab& surfcovcoVars,
const Vector<Real>& spatialCross, int ncross,
TurbForcingComp& turbforce)

Expand Down Expand Up @@ -230,6 +231,8 @@ void WriteCheckPoint3D(int step,
amrex::MultiFabFileFullPrefix(0, checkpointname, "Level_", "surfcovMeans"));
VisMF::Write(surfcovVars,
amrex::MultiFabFileFullPrefix(0, checkpointname, "Level_", "surfcovVars"));
VisMF::Write(surfcovcoVars,
amrex::MultiFabFileFullPrefix(0, checkpointname, "Level_", "surfcovcoVars"));
}
}

Expand Down Expand Up @@ -612,6 +615,7 @@ void ReadCheckPoint3D(int& step,
amrex::MultiFab& surfcov,
amrex::MultiFab& surfcovMeans,
amrex::MultiFab& surfcovVars,
amrex::MultiFab& surfcovcoVars,
Vector<Real>& spatialCross,
int ncross,
TurbForcingComp& turbforce,
Expand Down Expand Up @@ -733,6 +737,7 @@ void ReadCheckPoint3D(int& step,
surfcov.define(ba,dmap,n_ads_spec,0);
surfcovMeans.define(ba,dmap,n_ads_spec,0);
surfcovVars.define(ba,dmap,n_ads_spec,0);
surfcovcoVars.define(ba,dmap,n_ads_spec*6,0);
}

}
Expand Down Expand Up @@ -830,6 +835,7 @@ void ReadCheckPoint3D(int& step,
for (int m=0;m<n_ads_spec;m++) {
surfcovMeans.setVal(0.0);
surfcovVars.setVal(0.0);
surfcovcoVars.setVal(0.0);
}
}
}
Expand Down Expand Up @@ -859,6 +865,7 @@ void ReadCheckPoint3D(int& step,
if (n_ads_spec>0) {
Read_Copy_MF_Checkpoint(surfcovMeans,"surfcovMeans",checkpointname,ba_old,dmap_old,n_ads_spec,0);
Read_Copy_MF_Checkpoint(surfcovVars,"surfcovVars",checkpointname,ba_old,dmap_old,n_ads_spec,0);
Read_Copy_MF_Checkpoint(surfcovcoVars,"surfcovcoVars",checkpointname,ba_old,dmap_old,n_ads_spec*6,0);
}
}

Expand Down
11 changes: 7 additions & 4 deletions src_compressible_stag/compressible_functions_stag.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,7 @@ void WritePlotFileStag(int step,
const MultiFab& surfcov,
const MultiFab& surfcovMeans,
const MultiFab& surfcovVars,
const MultiFab& surfcovcoVars,
const MultiFab& eta,
const MultiFab& kappa,
const MultiFab& zeta);
Expand Down Expand Up @@ -194,6 +195,7 @@ void WriteCheckPoint3D(int step,
const amrex::MultiFab& surfcov,
const amrex::MultiFab& surfcovMeans,
const amrex::MultiFab& surfcovVars,
const amrex::MultiFab& surfcovcoVars,
const Vector<Real>& spatialCross, int ncross,
TurbForcingComp& turbforce);

Expand All @@ -218,6 +220,7 @@ void ReadCheckPoint3D(int& step,
amrex::MultiFab& surfcov,
amrex::MultiFab& surfcovMeans,
amrex::MultiFab& surfcovVars,
amrex::MultiFab& surfcovcoVars,
Vector<Real>& spatialCross, int ncross,
TurbForcingComp& turbforce,
BoxArray& ba, DistributionMapping& dmap);
Expand Down Expand Up @@ -322,7 +325,7 @@ void evaluateStatsStag3D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
Vector<Real>& dataSliceMeans_xcross,
Vector<Real>& spatialCross3D, const int ncross,
const amrex::Box& domain,
Expand All @@ -338,7 +341,7 @@ void evaluateStatsStag2D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
MultiFab& spatialCross2D, const int ncross,
const int steps,
const Geometry& geom);
Expand All @@ -352,7 +355,7 @@ void evaluateStatsStag1D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
MultiFab& spatialCross1D, const int ncross,
const int steps,
const Geometry& geom);
Expand All @@ -374,7 +377,7 @@ void EvaluateVarsCoVars(const MultiFab& cons, const MultiFab& consMean, MultiFab
const std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
const MultiFab& theta, const MultiFab& thetaMean, MultiFab& thetaVar,
const MultiFab& theta, const MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
const int steps);

void GetSliceAverageCross(Vector<Real>& dataAvMeans_x,
Expand Down
18 changes: 11 additions & 7 deletions src_compressible_stag/main_driver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -198,6 +198,7 @@ void main_driver(const char* argv)

MultiFab surfcovMeans; // used in either MFsurfchem or surfchem_mui
MultiFab surfcovVars; // used in either MFsurfchem or surfchem_mui
MultiFab surfcovcoVars; // used in either MFsurfchem or surfchem_mui

std::array< MultiFab, AMREX_SPACEDIM > velMeans;
std::array< MultiFab, AMREX_SPACEDIM > velVars;
Expand Down Expand Up @@ -482,7 +483,7 @@ void main_driver(const char* argv)
else {
ReadCheckPoint3D(step_start, time, statsCount, geom, domain, cu, cuMeans, cuVars, prim,
primMeans, primVars, cumom, cumomMeans, cumomVars,
vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, spatialCross3D, ncross, turbforce, ba, dmap);
vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, surfcovcoVars, spatialCross3D, ncross, turbforce, ba, dmap);
}

if (reset_stats == 1) statsCount = 1;
Expand Down Expand Up @@ -626,8 +627,10 @@ void main_driver(const char* argv)
if (nspec_surfcov>0) {
surfcovMeans.define(ba,dmap,nspec_surfcov,0);
surfcovVars.define(ba,dmap,nspec_surfcov,0);
surfcovcoVars.define(ba,dmap,nspec_surfcov*6,0);
surfcovMeans.setVal(0.0);
surfcovVars.setVal(0.0);
surfcovcoVars.setVal(0.0);
}

for (int d=0; d<AMREX_SPACEDIM; d++) {
Expand Down Expand Up @@ -717,7 +720,7 @@ void main_driver(const char* argv)

if (plot_int > 0) {
WritePlotFileStag(0, 0.0, geom, cu, cuMeans, cuVars, cumom, cumomMeans, cumomVars,
prim, primMeans, primVars, vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, eta, kappa, zeta);
prim, primMeans, primVars, vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, surfcovcoVars, eta, kappa, zeta);
#if defined(TURB)
if (turbForcing > 0) {
EvaluateWritePlotFileVelGrad(0, 0.0, geom, vel, vel_decomp);
Expand Down Expand Up @@ -1079,6 +1082,7 @@ void main_driver(const char* argv)
if (nspec_surfcov>0) {
surfcovMeans.setVal(0.0);
surfcovVars.setVal(0.0);
surfcovcoVars.setVal(0.0);
}

if (do_1D) {
Expand All @@ -1101,19 +1105,19 @@ void main_driver(const char* argv)
if (do_1D) {
evaluateStatsStag1D(cu, cuMeans, cuVars, prim, primMeans, primVars, vel,
velMeans, velVars, cumom, cumomMeans, cumomVars, coVars,
surfcov, surfcovMeans, surfcovVars,
surfcov, surfcovMeans, surfcovVars, surfcovcoVars,
spatialCross1D, ncross, statsCount, geom);
}
else if (do_2D) {
evaluateStatsStag2D(cu, cuMeans, cuVars, prim, primMeans, primVars, vel,
velMeans, velVars, cumom, cumomMeans, cumomVars, coVars,
surfcov, surfcovMeans, surfcovVars,
surfcov, surfcovMeans, surfcovVars, surfcovcoVars,
spatialCross2D, ncross, statsCount, geom);
}
else {
evaluateStatsStag3D(cu, cuMeans, cuVars, prim, primMeans, primVars, vel,
velMeans, velVars, cumom, cumomMeans, cumomVars, coVars,
surfcov, surfcovMeans, surfcovVars,
surfcov, surfcovMeans, surfcovVars, surfcovcoVars,
dataSliceMeans_xcross, spatialCross3D, ncross, domain,
statsCount, geom);
}
Expand Down Expand Up @@ -1145,7 +1149,7 @@ void main_driver(const char* argv)
//yzAverage(cuMeans, cuVars, primMeans, primVars, spatialCross,
// cuMeansAv, cuVarsAv, primMeansAv, primVarsAv, spatialCrossAv);
WritePlotFileStag(step, time, geom, cu, cuMeans, cuVars, cumom, cumomMeans, cumomVars,
prim, primMeans, primVars, vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, eta, kappa, zeta);
prim, primMeans, primVars, vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, surfcovcoVars, eta, kappa, zeta);

#if defined(TURB)
if (turbForcing > 0) {
Expand Down Expand Up @@ -1485,7 +1489,7 @@ void main_driver(const char* argv)
else {
WriteCheckPoint3D(step, time, statsCount, geom, cu, cuMeans, cuVars, prim,
primMeans, primVars, cumom, cumomMeans, cumomVars,
vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, spatialCross3D, ncross, turbforce);
vel, velMeans, velVars, coVars, surfcov, surfcovMeans, surfcovVars, surfcovcoVars, spatialCross3D, ncross, turbforce);
}
}

Expand Down
26 changes: 19 additions & 7 deletions src_compressible_stag/statsStag.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ void evaluateStatsStag3D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
Vector<Real>& dataSliceMeans_xcross,
Vector<Real>& spatialCross3D, const int ncross,
const amrex::Box& domain,
Expand All @@ -35,7 +35,7 @@ void evaluateStatsStag3D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
//// Evaluate Variances and Covariances
if ((plot_vars) or (plot_covars)) {
EvaluateVarsCoVars(cons,consMean,consVar,prim_in,primMean,primVar,velMean,velVar,
cumom,cumomMean,cumomVar,coVar,theta,thetaMean,thetaVar,steps);
cumom,cumomMean,cumomVar,coVar,theta,thetaMean,thetaVar,thetacoVar,steps);
consVar.FillBoundary(geom.periodicity());
primVar.FillBoundary(geom.periodicity());
}
Expand Down Expand Up @@ -81,7 +81,7 @@ void evaluateStatsStag2D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
MultiFab& /*spatialCross2D*/, const int /*ncross*/,
const int steps,
const Geometry& geom)
Expand All @@ -98,7 +98,7 @@ void evaluateStatsStag2D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
//// Evaluate Variances and Covariances
if ((plot_vars) or (plot_covars)) {
EvaluateVarsCoVars(cons,consMean,consVar,prim_in,primMean,primVar,velMean,velVar,
cumom,cumomMean,cumomVar,coVar,theta,thetaMean,thetaVar,steps);
cumom,cumomMean,cumomVar,coVar,theta,thetaMean,thetaVar,thetacoVar,steps);
consVar.FillBoundary(geom.periodicity());
primVar.FillBoundary(geom.periodicity());
}
Expand Down Expand Up @@ -130,7 +130,7 @@ void evaluateStatsStag1D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar,
MultiFab& theta, MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
MultiFab& spatialCross1D, const int ncross,
const int steps,
const Geometry& geom)
Expand All @@ -147,7 +147,7 @@ void evaluateStatsStag1D(MultiFab& cons, MultiFab& consMean, MultiFab& consVar,
//// Evaluate Variances and Covariances
if ((plot_vars) or (plot_covars)) {
EvaluateVarsCoVars(cons,consMean,consVar,prim_in,primMean,primVar,velMean,velVar,
cumom,cumomMean,cumomVar,coVar,theta,thetaMean,thetaVar,steps);
cumom,cumomMean,cumomVar,coVar,theta,thetaMean,thetaVar,thetacoVar,steps);
consVar.FillBoundary(geom.periodicity());
primVar.FillBoundary(geom.periodicity());
}
Expand Down Expand Up @@ -362,7 +362,7 @@ void EvaluateVarsCoVars(const MultiFab& cons, const MultiFab& consMean, MultiFab
const std::array<MultiFab, AMREX_SPACEDIM>& cumomMean,
std::array<MultiFab, AMREX_SPACEDIM>& cumomVar,
MultiFab& coVar,
const MultiFab& theta, const MultiFab& thetaMean, MultiFab& thetaVar,
const MultiFab& theta, const MultiFab& thetaMean, MultiFab& thetaVar, MultiFab& thetacoVar,
const int steps)
{
BL_PROFILE_VAR("EvaluateVarsCoVars()",EvaluateVarsCoVars);
Expand Down Expand Up @@ -397,6 +397,11 @@ void EvaluateVarsCoVars(const MultiFab& cons, const MultiFab& consMean, MultiFab
const Array4< Real> momyvars = cumomVar[1].array(mfi);
const Array4< Real> momzvars = cumomVar[2].array(mfi);

const Array4<const Real> surfcov = (nspec_surfcov>0) ? theta.array(mfi) : cons.array(mfi);
const Array4<const Real> surfcovmeans = (nspec_surfcov>0) ? thetaMean.array(mfi) : consMean.array(mfi);
const Array4< Real> surfcovvars = (nspec_surfcov>0) ? thetaVar.array(mfi) : consVar.array(mfi);
const Array4< Real> surfcovcoVars= (nspec_surfcov>0) ? thetacoVar.array(mfi) : consVar.array(mfi);

// update momentum and velocity variances
amrex::ParallelFor(tbx, tby, tbz,
[=] AMREX_GPU_DEVICE (int i, int j, int k)
Expand Down Expand Up @@ -460,6 +465,7 @@ void EvaluateVarsCoVars(const MultiFab& cons, const MultiFab& consMean, MultiFab
const Array4<const Real> surfcov = (nspec_surfcov>0) ? theta.array(mfi) : cons.array(mfi);
const Array4<const Real> surfcovmeans = (nspec_surfcov>0) ? thetaMean.array(mfi) : consMean.array(mfi);
const Array4< Real> surfcovvars = (nspec_surfcov>0) ? thetaVar.array(mfi) : consVar.array(mfi);
const Array4< Real> surfcovcoVars= (nspec_surfcov>0) ? thetacoVar.array(mfi) : consVar.array(mfi);

// other cell-centered variances and covariances (temperature fluctuation and variance)
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
Expand Down Expand Up @@ -570,6 +576,12 @@ void EvaluateVarsCoVars(const MultiFab& cons, const MultiFab& consMean, MultiFab
for (int m=0; m<nspec_surfcov; ++m) {
Real delsurfcov = surfcov(i,j,k,m) - surfcovmeans(i,j,k,m);
surfcovvars(i,j,k,m) = (surfcovvars(i,j,k,m)*stepsminusone + delsurfcov*delsurfcov)*stepsinv;
surfcovcoVars(i,j,k,6*m) = (surfcovcoVars(i,j,k,6*m)*stepsminusone + delrhoYk[m]*deltemp)*stepsinv; // <rhoYk T>
surfcovcoVars(i,j,k,6*m+1) = (surfcovcoVars(i,j,k,6*m+1)*stepsminusone + delsurfcov*delrhoYk[m])*stepsinv; // <theta rhoYk>
surfcovcoVars(i,j,k,6*m+2) = (surfcovcoVars(i,j,k,6*m+2)*stepsminusone + delsurfcov*delvelx)*stepsinv; // <theta vx>
surfcovcoVars(i,j,k,6*m+3) = (surfcovcoVars(i,j,k,6*m+3)*stepsminusone + delsurfcov*delvely)*stepsinv; // <theta vy>
surfcovcoVars(i,j,k,6*m+4) = (surfcovcoVars(i,j,k,6*m+4)*stepsminusone + delsurfcov*delvelz)*stepsinv; // <theta vz>
surfcovcoVars(i,j,k,6*m+5) = (surfcovcoVars(i,j,k,6*m+5)*stepsminusone + delsurfcov*deltemp)*stepsinv; // <theta T>
}
}
});
Expand Down
Loading

0 comments on commit a6ac210

Please sign in to comment.