Skip to content

Commit

Permalink
Merge pull request #183 from AMReX-FHD/surfchemtest
Browse files Browse the repository at this point in the history
MFsurfchem - modification of desorption species
  • Loading branch information
ajnonaka authored Feb 11, 2025
2 parents a1a7552 + 11d1551 commit ecd7eb2
Show file tree
Hide file tree
Showing 4 changed files with 58 additions and 17 deletions.
4 changes: 2 additions & 2 deletions src_MFsurfchem/MFsurfchem_functions.H
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,8 @@ void InitializeMFSurfchemNamespace();

void init_surfcov(MultiFab& surfcov, const amrex::Geometry& geom);

void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes, const amrex::Geometry& geom, const amrex::Real dt);
void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes, MultiFab& dNads, MultiFab& dNdes, const amrex::Geometry& geom, const amrex::Real dt);

void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes, const amrex::Geometry& geom);
void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes, MultiFab& dNads, MultiFab& dNdes, const amrex::Geometry& geom);

#endif
54 changes: 44 additions & 10 deletions src_MFsurfchem/MFsurfchem_functions.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@ AMREX_GPU_MANAGED amrex::Real MFsurfchem::k_beta;
AMREX_GPU_MANAGED amrex::Real MFsurfchem::e_beta;

AMREX_GPU_MANAGED int MFsurfchem::splitting_MFsurfchem;
AMREX_GPU_MANAGED int MFsurfchem::conversion_MFsurfchem;

void InitializeMFSurfchemNamespace()
{
Expand Down Expand Up @@ -79,6 +80,13 @@ void InitializeMFSurfchemNamespace()
// get splitting type: first order (0), strang (1)
splitting_MFsurfchem = 0; // default value
pp.query("splitting_MFsurfchem",splitting_MFsurfchem);

// ads/des of species 1 -> ads of species 1 + desorption of species (1+n)
conversion_MFsurfchem = 0; // default value
pp.query("conversion_MFsurfchem",conversion_MFsurfchem);
if ( (n_ads_spec + conversion_MFsurfchem) > nspecies) {
Abort("ERROR: desorption species is not included in nspecies");
}
return;
}

Expand Down Expand Up @@ -130,7 +138,7 @@ void init_surfcov(MultiFab& surfcov, const amrex::Geometry& geom)
return;
}

void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes,
void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes, MultiFab& dNads, MultiFab& dNdes,
const amrex::Geometry& geom, const amrex::Real dt)
{

Expand All @@ -143,6 +151,8 @@ void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab
const Array4<Real> & prim_arr = prim.array(mfi);
const Array4<Real> & surfcov_arr = surfcov.array(mfi);
const Array4<Real> & dNadsdes_arr = dNadsdes.array(mfi);
const Array4<Real> & dNads_arr = dNads.array(mfi);
const Array4<Real> & dNdes_arr = dNdes.array(mfi);

amrex::Real Ntot = surf_site_num_dens*dx[0]*dx[1]; // total number of reactive sites

Expand Down Expand Up @@ -177,7 +187,13 @@ void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab
Ndes = RandomPoisson(meanNdes,engine);
}

dNadsdes_arr(i,j,k,m) = Nads-Ndes;
if (conversion_MFsurfchem > 0) {
dNads_arr(i,j,k,m) = Nads;
dNdes_arr(i,j,k,m) = Ndes;
}
else {
dNadsdes_arr(i,j,k,m) = Nads-Ndes;
}
}
}
});
Expand All @@ -186,7 +202,7 @@ void sample_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab
return;
}

void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes,
void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab& dNadsdes, MultiFab& dNads, MultiFab& dNdes,
const amrex::Geometry& geom)
{
const GpuArray<Real, AMREX_SPACEDIM> dx = geom.CellSizeArray();
Expand All @@ -198,6 +214,8 @@ void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab
const Array4<Real> & prim_arr = prim.array(mfi);
const Array4<Real> & surfcov_arr = surfcov.array(mfi);
const Array4<Real> & dNadsdes_arr = dNadsdes.array(mfi);
const Array4<Real> & dNads_arr = dNads.array(mfi);
const Array4<Real> & dNdes_arr = dNdes.array(mfi);

amrex::Real Ntot = surf_site_num_dens*dx[0]*dx[1]; // total number of reactive sites
amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
Expand All @@ -207,13 +225,29 @@ void update_MFsurfchem(MultiFab& cu, MultiFab& prim, MultiFab& surfcov, MultiFab
amrex::Real T_inst = prim_arr(i,j,k,4);

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 = (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;
cu_arr(i,j,k,4) -= factor2*dN;
if (conversion_MFsurfchem > 0) {
amrex::Real mconv = m + conversion_MFsurfchem;
amrex::Real dNads = dNads_arr(i,j,k,m);
amrex::Real dNdes = dNdes_arr(i,j,k,m);
amrex::Real factor1 = molmass[m]/AVONUM/(dx[0]*dx[1]*dx[2]);
amrex::Real factor1conv = molmass[mconv]/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]);
amrex::Real factor2conv = (e_beta*k_B*T_inst+(e0[mconv]+hcv[mconv]*T_inst)*molmass[mconv]/AVONUM)/(dx[0]*dx[1]*dx[2]);
surfcov_arr(i,j,k,m) += (dNads-dNdes)/Ntot;
cu_arr(i,j,k,0) -= factor1*dNads - factor1conv*dNdes;
cu_arr(i,j,k,5+m) -= factor1*dNads;
cu_arr(i,j,k,5+mconv) += factor1conv*dNdes;
cu_arr(i,j,k,4) -= factor2*dNads - factor2conv*dNdes;
}
else {
amrex::Real dN = dNadsdes_arr(i,j,k,m);
amrex::Real factor1 = 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;
cu_arr(i,j,k,4) -= factor2*dN;
}
}
}
});
Expand Down
1 change: 1 addition & 0 deletions src_MFsurfchem/MFsurfchem_namespace.H
Original file line number Diff line number Diff line change
Expand Up @@ -16,5 +16,6 @@ namespace MFsurfchem {
extern AMREX_GPU_MANAGED amrex::Real e_beta;

extern AMREX_GPU_MANAGED int splitting_MFsurfchem;
extern AMREX_GPU_MANAGED int conversion_MFsurfchem;

}
16 changes: 11 additions & 5 deletions src_compressible_stag/main_driver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -214,6 +214,8 @@ void main_driver(const char* argv)
// MFsurfchem
MultiFab surfcov; // also used in surfchem_mui for stats and plotfiles
MultiFab dNadsdes;
MultiFab dNads;
MultiFab dNdes;

#if defined(MUI) || defined(USE_AMREX_MPMD)
MultiFab Ntot; // saves total number of sites
Expand Down Expand Up @@ -536,6 +538,8 @@ void main_driver(const char* argv)

if (n_ads_spec>0) {
dNadsdes.define(ba,dmap,n_ads_spec,0);
dNads.define(ba,dmap,n_ads_spec,0);
dNdes.define(ba,dmap,n_ads_spec,0);
nspec_surfcov = n_ads_spec;
}

Expand Down Expand Up @@ -606,6 +610,8 @@ void main_driver(const char* argv)
if (n_ads_spec>0) {
surfcov.define(ba,dmap,n_ads_spec,0);
dNadsdes.define(ba,dmap,n_ads_spec,0);
dNads.define(ba,dmap,n_ads_spec,0);
dNdes.define(ba,dmap,n_ads_spec,0);
nspec_surfcov = n_ads_spec;
}

Expand Down Expand Up @@ -997,10 +1003,10 @@ void main_driver(const char* argv)
#endif
if (n_ads_spec>0) {
if (splitting_MFsurfchem == 0) {
sample_MFsurfchem(cu, prim, surfcov, dNadsdes, geom, dt);
sample_MFsurfchem(cu, prim, surfcov, dNadsdes, dNads, dNdes, geom, dt);
} else if (splitting_MFsurfchem == 1) {
sample_MFsurfchem(cu, prim, surfcov, dNadsdes, geom, dt/2.0);
update_MFsurfchem(cu, prim, surfcov, dNadsdes, geom);
sample_MFsurfchem(cu, prim, surfcov, dNadsdes, dNads, dNdes, geom, dt/2.0);
update_MFsurfchem(cu, prim, surfcov, dNadsdes, dNads, dNdes, geom);

for (int d=0; d<AMREX_SPACEDIM; d++) {
cumom[d].FillBoundary(geom.periodicity());
Expand Down Expand Up @@ -1031,7 +1037,7 @@ void main_driver(const char* argv)
}

if (n_ads_spec>0 && splitting_MFsurfchem == 1) {
sample_MFsurfchem(cu, prim, surfcov, dNadsdes, geom, dt/2.0);
sample_MFsurfchem(cu, prim, surfcov, dNadsdes, dNads, dNdes, geom, dt/2.0);
}

// update surface chemistry (via either surfchem_mui or MFsurfchem)
Expand Down Expand Up @@ -1066,7 +1072,7 @@ void main_driver(const char* argv)

if (n_ads_spec>0) {

update_MFsurfchem(cu, prim, surfcov, dNadsdes, geom);
update_MFsurfchem(cu, prim, surfcov, dNadsdes, dNads, dNdes, geom);

for (int d=0; d<AMREX_SPACEDIM; d++) {
cumom[d].FillBoundary(geom.periodicity());
Expand Down

0 comments on commit ecd7eb2

Please sign in to comment.