Skip to content

Commit

Permalink
Merge pull request #172 from AMReX-FHD/slice_SF_compressible_stag
Browse files Browse the repository at this point in the history 
Implement slice structure factor for compressible_stag
  • Loading branch information
@ajnonaka
ajnonaka authored Nov 26, 2024
2 parents 6041c07 + 81ae8ac commit 9cd3bd9
Show file tree
Hide file tree
Showing 2 changed files with 70 additions and 16 deletions.
41 changes: 39 additions & 2 deletions src_common/ComputeAverages.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -384,7 +384,44 @@ void ExtractSlice(const MultiFab& mf, MultiFab& mf_slice,

// create a new DistributionMapping and define the MultiFab
DistributionMapping dmap_slice(ba_slice);
mf_slice.define(ba_slice,dmap_slice,ncomp,0);
MultiFab mf_slice_tmp(ba_slice,dmap_slice,ncomp,0);

mf_slice.ParallelCopy(mf, incomp, 0, ncomp);
mf_slice_tmp.ParallelCopy(mf, incomp, 0, ncomp);

// now copy this into a multifab with index zero in the dir direction rather than slicepoint
// (structure factor code requires this)
dom_lo[dir] = 0;
dom_hi[dir] = 0;

Box domain_slice2(dom_lo,dom_hi);
BoxArray ba_slice2(domain_slice2);
ba_slice2.maxSize(IntVect(max_grid_slice));
mf_slice.define(ba_slice2,dmap_slice,ncomp,0);

for ( MFIter mfi(mf_slice_tmp,TilingIfNotGPU()); mfi.isValid(); ++mfi ) {

const Box& bx = mfi.tilebox();

const Array4<Real> & slice = mf_slice.array(mfi);
const Array4<Real> & slice_tmp = mf_slice_tmp.array(mfi);

if (dir == 0) {
amrex::ParallelFor(bx, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
{
slice(0,j,k,n) = slice_tmp(i,j,k,n);
});
}
if (dir == 1) {
amrex::ParallelFor(bx, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
{
slice(i,0,k,n) = slice_tmp(i,j,k,n);
});
}
if (dir == 2) {
amrex::ParallelFor(bx, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
{
slice(i,j,0,n) = slice_tmp(i,j,k,n);
});
}
}
}
45 changes: 31 additions & 14 deletions src_compressible_stag/main_driver.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -208,8 +208,13 @@ void main_driver(const char* argv)
if ((plot_cross) and ((cross_cell < 0) or (cross_cell > n_cells[0]-1))) {
Abort("Cross cell needs to be within the domain: 0 <= cross_cell <= n_cells[0] - 1");
}
if ((do_slab_sf) and ((membrane_cell <= 0) or (membrane_cell >= n_cells[0]-1))) {
Abort("Slab structure factor needs a membrane cell within the domain: 0 < cross_cell < n_cells[0] - 1");
if (project_dir >= 0) {
if (do_slab_sf and ((membrane_cell <= 0) or (membrane_cell >= n_cells[project_dir]-1))) {
Abort("Slab structure factor needs a membrane cell within the domain: 0 < cross_cell < n_cells[project_dir] - 1");
}
if (do_slab_sf and slicepoint >= 0) {
Abort("Cannot use do_slab_sf and slicepoint");
}
}
if ((project_dir >= 0) and ((do_1D) or (do_2D))) {
Abort("Projected structure factors (project_dir) works only for 3D case");
Expand Down Expand Up @@ -276,9 +281,9 @@ void main_driver(const char* argv)
MultiFab structFactPrimMF;
MultiFab structFactConsMF;

// Structure factor for 2D averaged data
StructFact structFactPrimVerticalAverage;
StructFact structFactConsVerticalAverage;
// Structure factor for vertically-averaged or sliced data
StructFact structFactPrimFlattened;
StructFact structFactConsFlattened;

// Structure factor for 2D averaged data (across a membrane)
StructFact structFactPrimVerticalAverage0;
Expand Down Expand Up @@ -781,7 +786,11 @@ void main_driver(const char* argv)

{
MultiFab X, XRot;
ComputeVerticalAverage(prim, X, geom, project_dir, 0, nprimvars);
if (slicepoint < 0) {
ComputeVerticalAverage(prim, X, geom, project_dir, 0, nprimvars);
} else {
ExtractSlice(prim, X, geom, project_dir, slicepoint, 0, 1);
}
XRot = RotateFlattenedMF(X);
ba_flat = XRot.boxArray();
dmap_flat = XRot.DistributionMap();
Expand Down Expand Up @@ -845,8 +854,8 @@ void main_driver(const char* argv)
}

if (do_slab_sf == 0) {
structFactPrimVerticalAverage.define(ba_flat,dmap_flat,prim_var_names,var_scaling_prim,2);
structFactConsVerticalAverage.define(ba_flat,dmap_flat,cons_var_names,var_scaling_cons,2);
structFactPrimFlattened.define(ba_flat,dmap_flat,prim_var_names,var_scaling_prim);
structFactConsFlattened.define(ba_flat,dmap_flat,cons_var_names,var_scaling_cons);
}
else {
structFactPrimVerticalAverage0.define(ba_flat,dmap_flat,prim_var_names,var_scaling_prim);
Expand Down Expand Up @@ -1468,19 +1477,27 @@ void main_driver(const char* argv)
{
MultiFab X, XRot;

ComputeVerticalAverage(structFactPrimMF, X, geom, project_dir, 0, structVarsPrim);
if (slicepoint < 0) {
ComputeVerticalAverage(structFactPrimMF, X, geom, project_dir, 0, structVarsPrim);
} else {
ExtractSlice(structFactPrimMF, X, geom, project_dir, slicepoint, 0, structVarsPrim);
}
XRot = RotateFlattenedMF(X);
master_project_rot_prim.ParallelCopy(XRot, 0, 0, structVarsPrim);
structFactPrimVerticalAverage.FortStructure(master_project_rot_prim,geom_flat);
structFactPrimFlattened.FortStructure(master_project_rot_prim,geom_flat);
}

{
MultiFab X, XRot;

ComputeVerticalAverage(structFactConsMF, X, geom, project_dir, 0, structVarsCons);
if (slicepoint < 0) {
ComputeVerticalAverage(structFactConsMF, X, geom, project_dir, 0, structVarsCons);
} else {
ExtractSlice(structFactConsMF, X, geom, project_dir, slicepoint, 0, structVarsCons);
}
XRot = RotateFlattenedMF(X);
master_project_rot_cons.ParallelCopy(XRot, 0, 0, structVarsCons);
structFactConsVerticalAverage.FortStructure(master_project_rot_cons,geom_flat);
structFactConsFlattened.FortStructure(master_project_rot_cons,geom_flat);
}

}
Expand Down Expand Up @@ -1576,8 +1593,8 @@ void main_driver(const char* argv)

if (project_dir >= 0) {
if (do_slab_sf == 0) {
structFactPrimVerticalAverage.WritePlotFile(step,time,geom_flat,"plt_SF_prim_VerticalAverage");
structFactConsVerticalAverage.WritePlotFile(step,time,geom_flat,"plt_SF_cons_VerticalAverage");
structFactPrimFlattened.WritePlotFile(step,time,geom_flat,"plt_SF_prim_Flattened");
structFactConsFlattened.WritePlotFile(step,time,geom_flat,"plt_SF_cons_Flattened");
}
else {
structFactPrimVerticalAverage0.WritePlotFile(step,time,geom_flat,"plt_SF_prim_VerticalAverageSlab0");
Expand Down

0 comments on commit 9cd3bd9

Please sign in to comment.