From 45b79e362b8d3707b35b7422d364380e3536a69e Mon Sep 17 00:00:00 2001 From: Andy Nonaka Date: Fri, 3 Jan 2025 12:59:16 -0800 Subject: [PATCH] cleanup flattened code, simplify a bit --- exec/hydro/main_driver.cpp | 93 +++------- src_common/Make.package | 1 - src_common/RotateFlattenedMF.cpp | 51 ------ src_common/common_functions.H | 5 - src_compressible/main_driver.cpp | 112 ++++-------- src_compressible_stag/main_driver.cpp | 249 +++++++++----------------- 6 files changed, 133 insertions(+), 378 deletions(-) delete mode 100644 src_common/RotateFlattenedMF.cpp diff --git a/exec/hydro/main_driver.cpp b/exec/hydro/main_driver.cpp index 2fdf97be8..d3edef023 100644 --- a/exec/hydro/main_driver.cpp +++ b/exec/hydro/main_driver.cpp @@ -350,7 +350,7 @@ void main_driver(const char* argv) /////////////////////////////////////////// StructFact structFactFlattened; - MultiFab FlattenedRotMaster; + MultiFab FlattenedMaster; Geometry geom_flat; @@ -365,76 +365,30 @@ void main_driver(const char* argv) } else { ExtractSlice(structFactMF, Flattened, geom, project_dir, slicepoint, 0, 1); } - // we rotate this flattened MultiFab to have normal in the z-direction since - // our structure factor class assumes this for flattened - MultiFab FlattenedRot = RotateFlattenedMF(Flattened); - BoxArray ba_flat = FlattenedRot.boxArray(); - const DistributionMapping& dmap_flat = FlattenedRot.DistributionMap(); - FlattenedRotMaster.define(ba_flat,dmap_flat,structVars,0); + BoxArray ba_flat = Flattened.boxArray(); + const DistributionMapping& dmap_flat = Flattened.DistributionMap(); + FlattenedMaster.define(ba_flat,dmap_flat,structVars,0); { - IntVect dom_lo(AMREX_D_DECL(0,0,0)); - IntVect dom_hi; - - // yes you could simplify this code but for now - // these are written out fully to better understand what is happening - // we wanted dom_hi[AMREX_SPACEDIM-1] to be equal to 0 - // and need to transmute the other indices depending on project_dir -#if (AMREX_SPACEDIM == 2) - if (project_dir == 0) { - dom_hi[0] = n_cells[1]-1; - } - else if (project_dir == 1) { - dom_hi[0] = n_cells[0]-1; - } - dom_hi[1] = 0; -#elif (AMREX_SPACEDIM == 3) - if (project_dir == 0) { - dom_hi[0] = n_cells[1]-1; - dom_hi[1] = n_cells[2]-1; - } else if (project_dir == 1) { - dom_hi[0] = n_cells[0]-1; - dom_hi[1] = n_cells[2]-1; - } else if (project_dir == 2) { - dom_hi[0] = n_cells[0]-1; - dom_hi[1] = n_cells[1]-1; - } - dom_hi[2] = 0; -#endif - Box domain(dom_lo, dom_hi); - + Box domain_flat = FlattenedMaster.boxArray().minimalBox(); + // This defines the physical box + // we retain prob_lo and prob_hi in all directions except project_dir, + // where the physical size is 0 to dx[project_dir] + Vector projected_lo(AMREX_SPACEDIM); Vector projected_hi(AMREX_SPACEDIM); - // yes you could simplify this code but for now - // these are written out fully to better understand what is happening - // we wanted projected_hi[AMREX_SPACEDIM-1] to be equal to dx[projected_dir] - // and need to transmute the other indices depending on project_dir -#if (AMREX_SPACEDIM == 2) - if (project_dir == 0) { - projected_hi[0] = prob_hi[1]; - } else if (project_dir == 1) { - projected_hi[0] = prob_hi[0]; + for (int d=0; d1 cells in the other directions) -// returns a flattened multifab that is now flattened in the AMREX_SPACEDIM-1 direction -// (z direction in 3D, y direction in 2D) -MultiFab RotateFlattenedMF(MultiFab const& mf) -{ - BoxArray const& old_ba = mf.boxArray(); - DistributionMapping const& dm = mf.DistributionMap(); - Box const& domain_box = old_ba.minimalBox(); - int short_direction; - int short_size = domain_box.shortside(short_direction); - if (short_size != 1) { - Print() << "RotateFlattenedMF needs a MF with short_size==1; returning the original input MultiFab\n"; - return MultiFab(mf, amrex::make_alias, 0, mf.nComp()); - } else if (short_direction == AMREX_SPACEDIM-1) { - return MultiFab(mf, amrex::make_alias, 0, mf.nComp()); - } else { - IntVect old_ng = mf.nGrowVect(); - AMREX_ALWAYS_ASSERT_WITH_MESSAGE(old_ng[short_direction] == 0, - "Not supposed to have ghost cells in the shortest direction"); - IntVect ng; - if (short_direction == 0) { - ng = IntVect(AMREX_D_DECL(old_ng[1],old_ng[2],0)); - } else { - ng = IntVect(AMREX_D_DECL(old_ng[0],old_ng[2],0)); - } - BoxList bl = old_ba.boxList(); - for (auto& b : bl) { - const auto lo = b.smallEnd(); - const auto hi = b.bigEnd(); - if (short_direction == 0) { - b = Box(IntVect(AMREX_D_DECL(lo[1],lo[2],0)), - IntVect(AMREX_D_DECL(hi[1],hi[2],0)), - b.ixType()); - } else { - b = Box(IntVect(AMREX_D_DECL(lo[0],lo[2],0)), - IntVect(AMREX_D_DECL(hi[0],hi[2],0)), - b.ixType()); - } - } - BoxArray new_ba(std::move(bl)); - const int ncomp = mf.nComp(); - MultiFab new_mf(new_ba, dm, ncomp, ng, MFInfo().SetAlloc(false)); - for (MFIter mfi(new_mf); mfi.isValid(); ++mfi) { - new_mf.setFab(mfi, FArrayBox(mfi.fabbox(), ncomp, mf[mfi.index()].dataPtr())); - } - return new_mf; - } -} diff --git a/src_common/common_functions.H b/src_common/common_functions.H index 7cfeba665..2e3a6d0c9 100644 --- a/src_common/common_functions.H +++ b/src_common/common_functions.H @@ -243,11 +243,6 @@ void CCL2Norm(const MultiFab & m1, amrex::MultiFab& mscr, Real & norm_l2); -/////////////////////////// -// in RotateFlattenedMF.cpp - -MultiFab RotateFlattenedMF(MultiFab const& mf); - /////////////////////////// // in InterpCoarsen.cpp void FaceFillCoarse(Vector>& mf, int map); diff --git a/src_compressible/main_driver.cpp b/src_compressible/main_driver.cpp index 73c8eec97..bf97b00e9 100644 --- a/src_compressible/main_driver.cpp +++ b/src_compressible/main_driver.cpp @@ -406,7 +406,7 @@ void main_driver(const char* argv) // structure factor class for flattened dataset StructFact structFactPrimFlattened; - MultiFab primFlattenedRotMaster; + MultiFab primFlattenedMaster; ////////////////////////////////////////////// @@ -465,7 +465,7 @@ void main_driver(const char* argv) // structure factor class for flattened dataset StructFact structFactConsFlattened; - MultiFab consFlattenedRotMaster; + MultiFab consFlattenedMaster; ////////////////////////////////////////////// @@ -482,77 +482,31 @@ void main_driver(const char* argv) } else { ExtractSlice(structFactPrimMF, Flattened, geom, project_dir, slicepoint, 0, 1); } - // we rotate this flattened MultiFab to have normal in the z-direction since - // our structure factor class assumes this for flattened - MultiFab FlattenedRot = RotateFlattenedMF(Flattened); - BoxArray ba_flat = FlattenedRot.boxArray(); - const DistributionMapping& dmap_flat = FlattenedRot.DistributionMap(); - primFlattenedRotMaster.define(ba_flat,dmap_flat,structVarsPrim,0); - consFlattenedRotMaster.define(ba_flat,dmap_flat,structVarsCons,0); + BoxArray ba_flat = Flattened.boxArray(); + const DistributionMapping& dmap_flat = Flattened.DistributionMap(); + primFlattenedMaster.define(ba_flat,dmap_flat,structVarsPrim,0); + consFlattenedMaster.define(ba_flat,dmap_flat,structVarsCons,0); { - IntVect dom_lo(AMREX_D_DECL(0,0,0)); - IntVect dom_hi; - - // yes you could simplify this code but for now - // these are written out fully to better understand what is happening - // we wanted dom_hi[AMREX_SPACEDIM-1] to be equal to 0 - // and need to transmute the other indices depending on project_dir -#if (AMREX_SPACEDIM == 2) - if (project_dir == 0) { - dom_hi[0] = n_cells[1]-1; - } - else if (project_dir == 1) { - dom_hi[0] = n_cells[0]-1; - } - dom_hi[1] = 0; -#elif (AMREX_SPACEDIM == 3) - if (project_dir == 0) { - dom_hi[0] = n_cells[1]-1; - dom_hi[1] = n_cells[2]-1; - } else if (project_dir == 1) { - dom_hi[0] = n_cells[0]-1; - dom_hi[1] = n_cells[2]-1; - } else if (project_dir == 2) { - dom_hi[0] = n_cells[0]-1; - dom_hi[1] = n_cells[1]-1; - } - dom_hi[2] = 0; -#endif - Box domain(dom_lo, dom_hi); - - // This defines the physical box - Vector projected_hi(AMREX_SPACEDIM); - - // yes you could simplify this code but for now - // these are written out fully to better understand what is happening - // we wanted projected_hi[AMREX_SPACEDIM-1] to be equal to dx[projected_dir] - // and need to transmute the other indices depending on project_dir -#if (AMREX_SPACEDIM == 2) - if (project_dir == 0) { - projected_hi[0] = prob_hi[1]; - } else if (project_dir == 1) { - projected_hi[0] = prob_hi[0]; - } - projected_hi[1] = prob_hi[project_dir] / n_cells[project_dir]; -#elif (AMREX_SPACEDIM == 3) - if (project_dir == 0) { - projected_hi[0] = prob_hi[1]; - projected_hi[1] = prob_hi[2]; - } else if (project_dir == 1) { - projected_hi[0] = prob_hi[0]; - projected_hi[1] = prob_hi[2]; - } else if (project_dir == 2) { - projected_hi[0] = prob_hi[0]; - projected_hi[1] = prob_hi[1]; - } - projected_hi[2] = prob_hi[project_dir] / n_cells[project_dir]; -#endif - - RealBox real_box({AMREX_D_DECL( prob_lo[0], prob_lo[1], prob_lo[2])}, - {AMREX_D_DECL(projected_hi[0],projected_hi[1],projected_hi[2])}); - - // This defines a Geometry object - geom_flat.define(domain,&real_box,CoordSys::cartesian,is_periodic.data()); + Box domain_flat = primFlattenedMaster.boxArray().minimalBox(); + + // This defines the physical box + // we retain prob_lo and prob_hi in all directions except project_dir, + // where the physical size is 0 to dx[project_dir] + Vector projected_lo(AMREX_SPACEDIM); + Vector projected_hi(AMREX_SPACEDIM); + + for (int d=0; d structFactPrimArray; Vector < StructFact > structFactConsArray; - MultiFab master_2D_rot_prim; - MultiFab master_2D_rot_cons; + MultiFab prim2DFlattenedMaster; + MultiFab cons2DFlattenedMaster; #if defined(TURB) // Structure factor for compressible turbulence @@ -315,7 +315,7 @@ void main_driver(const char* argv) // surface coverage structure factor StructFact surfcovStructFact; - MultiFab surfcovFlattenedRotMaster; + MultiFab surfcovFlattenedMaster; Geometry surfcov_geom_flat; BoxArray surfcov_ba_flat; DistributionMapping surfcov_dmap_flat; @@ -766,75 +766,43 @@ void main_driver(const char* argv) // structure factor class for vertically-averaged dataset if (project_dir >= 0) { - + MultiFab Flattened; // flattened multifab define below + + // we are only calling ComputeVerticalAverage or ExtractSlice here to obtain + // a built version of Flattened so can obtain what we need to build the + // structure factor and geometry objects for flattened data + if (slicepoint < 0) { + ComputeVerticalAverage(prim, Flattened, geom, project_dir, 0, nprimvars); + } else { + ExtractSlice(prim, Flattened, geom, project_dir, slicepoint, 0, 1); + } + ba_flat = Flattened.boxArray(); + dmap_flat = Flattened.DistributionMap(); + primFlattenedMaster.define(ba_flat,dmap_flat,structVarsPrim,0); + consFlattenedMaster.define(ba_flat,dmap_flat,structVarsCons,0); { - MultiFab X, XRot; - 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(); - master_project_rot_prim.define(ba_flat,dmap_flat,structVarsPrim,0); - master_project_rot_cons.define(ba_flat,dmap_flat,structVarsCons,0); - - IntVect dom_lo_flat(AMREX_D_DECL(0,0,0)); - IntVect dom_hi_flat; -#if (AMREX_SPACEDIM == 2) - if (project_dir == 0) { - dom_hi_flat[0] = n_cells[1]-1; - dom_hi_flat[1] = 0; - } - else if (project_dir == 1) { - dom_hi_flat[0] = n_cells[0]-1; - dom_hi_flat[1] = 0; - } -#elif (AMREX_SPACEDIM == 3) - if (project_dir == 0) { - dom_hi_flat[0] = n_cells[1]-1; - dom_hi_flat[1] = n_cells[2]-1; - dom_hi_flat[2] = 0; - } else if (project_dir == 1) { - dom_hi_flat[0] = n_cells[0]-1; - dom_hi_flat[1] = n_cells[2]-1; - dom_hi_flat[2] = 0; - } else if (project_dir == 2) { - dom_hi_flat[0] = n_cells[0]-1; - dom_hi_flat[1] = n_cells[1]-1; - dom_hi_flat[2] = 0; - } -#endif - Box domain_flat(dom_lo_flat, dom_hi_flat); + Box domain_flat = primFlattenedMaster.boxArray().minimalBox(); // This defines the physical box + // we retain prob_lo and prob_hi in all directions except project_dir, + // where the physical size is 0 to dx[project_dir] + Vector projected_lo(AMREX_SPACEDIM); Vector projected_hi(AMREX_SPACEDIM); - for (int d=0; d projected_lo(AMREX_SPACEDIM); Vector projected_hi(AMREX_SPACEDIM); - // yes you could simplify this code but for now - // these are written out fully to better understand what is happening - // we wanted projected_hi[AMREX_SPACEDIM-1] to be equal to dx[projected_dir] - // and need to transmute the other indices depending on surfcov_dir -#if (AMREX_SPACEDIM == 2) - if (surfcov_dir == 0) { - projected_hi[0] = prob_hi[1]; - } else if (surfcov_dir == 1) { - projected_hi[0] = prob_hi[0]; - } - projected_hi[1] = prob_hi[surfcov_dir] / n_cells[surfcov_dir]; -#elif (AMREX_SPACEDIM == 3) - if (surfcov_dir == 0) { - projected_hi[0] = prob_hi[1]; - projected_hi[1] = prob_hi[2]; - } else if (surfcov_dir == 1) { - projected_hi[0] = prob_hi[0]; - projected_hi[1] = prob_hi[2]; - } else if (surfcov_dir == 2) { - projected_hi[0] = prob_hi[0]; - projected_hi[1] = prob_hi[1]; + for (int d=0; d