Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Implement slice structure factor for compressible_stag #172

Merged
merged 4 commits into from
Nov 26, 2024
Merged

Implement slice structure factor for compressible_stag #172

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
13f9d43
Implement slice structure factor for compressible_stag
ajnonaka Nov 5, 2024
0ff568f
Merge branch 'main' into slice_SF_compressible_stag
ajnonaka Nov 20, 2024
0cac562
introduce code that forces the slice multifab to have index 0:0 in th…
ajnonaka Nov 20, 2024
81ae8ac
fix copy of slice MF components
ajnonaka Nov 21, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
41 changes: 39 additions & 2 deletions src_common/ComputeAverages.cpp
View file
Open in desktop
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
Open in desktop
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
Loading