diff options
Diffstat (limited to 'Source')
| -rw-r--r-- | Source/Diagnostics/SliceDiagnostic.H | 9 | ||||
| -rw-r--r-- | Source/Diagnostics/SliceDiagnostic.cpp | 166 | ||||
| -rw-r--r-- | Source/Utils/WarpXMovingWindow.cpp | 20 |
3 files changed, 104 insertions, 91 deletions
diff --git a/Source/Diagnostics/SliceDiagnostic.H b/Source/Diagnostics/SliceDiagnostic.H index 86c1794e6..31eea83be 100644 --- a/Source/Diagnostics/SliceDiagnostic.H +++ b/Source/Diagnostics/SliceDiagnostic.H @@ -1,8 +1,6 @@ #ifndef WARPX_SliceDiagnostic_H_ #define WARPX_SliceDiagnostic_H_ -#include <vector> - #include <AMReX_VisMF.H> #include <AMReX_PlotFileUtil.H> #include <AMReX_ParallelDescriptor.H> @@ -25,10 +23,9 @@ std::unique_ptr<amrex::MultiFab> CreateSlice( const amrex::MultiFab& mf, void CheckSliceInput( const amrex::RealBox real_box, amrex::RealBox &slice_cc_nd_box, amrex::RealBox &slice_realbox, - amrex::IntVect &cr_ratio, amrex::IntVect slice_cr_ratio, - amrex::Vector<amrex::Geometry> dom_geom, amrex::IntVect const SliceType, - amrex::IntVect &slice_lo, amrex::IntVect &slice_hi, - amrex::IntVect &interp_lo); + amrex::IntVect &slice_cr_ratio, amrex::Vector<amrex::Geometry> dom_geom, + amrex::IntVect const SliceType, amrex::IntVect &slice_lo, + amrex::IntVect &slice_hi, amrex::IntVect &interp_lo); void InterpolateSliceValues( amrex::MultiFab& smf, amrex::IntVect interp_lo, amrex::RealBox slice_realbox, diff --git a/Source/Diagnostics/SliceDiagnostic.cpp b/Source/Diagnostics/SliceDiagnostic.cpp index 623d94f90..ce4427fb0 100644 --- a/Source/Diagnostics/SliceDiagnostic.cpp +++ b/Source/Diagnostics/SliceDiagnostic.cpp @@ -5,7 +5,7 @@ #include <WarpX.H> - +using namespace amrex; /* \brief @@ -18,17 +18,16 @@ * If interpolation is required, then on the smf, using data points stored in the ghost cells, the data in interpolated. * If coarsening is required, then a coarse slice multifab is generated (cs_mf) and the * values of the refined slice (smf) is averaged down to obtain the coarse slice. + * \param mf is the source multifab containing the field data + * \param dom_geom is the geometry of the domain and used in the function to obtain the + * CellSize of the underlying grid. * \param slice_realbox defines the extent of the slice * \param slice_cr_ratio provides the coarsening ratio for diagnostics - * \param SliceType is the index type of the slice. By default is slice is cell-centered. - * If the source multifab (mf) has a different index type then SliceType is converted. */ -using namespace amrex; -std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, - const amrex::Vector<Geometry> &dom_geom, - amrex::RealBox &slice_realbox, - amrex::IntVect &slice_cr_ratio ) +std::unique_ptr<MultiFab> +CreateSlice( const MultiFab& mf, const Vector<Geometry> &dom_geom, + RealBox &slice_realbox, IntVect &slice_cr_ratio ) { std::unique_ptr<MultiFab> smf; std::unique_ptr<MultiFab> cs_mf; @@ -38,16 +37,11 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, int nlevels = dom_geom.size(); int ncomp = (mf).nComp(); - IntVect cr_ratio(AMREX_D_DECL(0,0,0)); const auto conversionType = (mf).ixType(); IntVect SliceType(AMREX_D_DECL(0,0,0)); - for (int idim = 0; idim < AMREX_SPACEDIM; ++idim ) { - if ( conversionType.nodeCentered(idim)) - { - SliceType[idim] = 1; - } + SliceType[idim] = conversionType.nodeCentered(idim); } const RealBox& real_box = dom_geom[0].ProbDomain(); @@ -63,12 +57,12 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, IntVect slice_lo2(AMREX_D_DECL(0,0,0)); IntVect interp_lo(AMREX_D_DECL(0,0,0)); - CheckSliceInput(real_box, slice_cc_nd_box, slice_realbox, cr_ratio, - slice_cr_ratio, dom_geom, SliceType, slice_lo, + CheckSliceInput(real_box, slice_cc_nd_box, slice_realbox, slice_cr_ratio, + dom_geom, SliceType, slice_lo, slice_hi, interp_lo); - // Determine if interpolation is required and number of cells in slice // for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) { + // Flag for interpolation if required // if ( interp_lo[idim] == 1) { interpolate = 1; @@ -80,10 +74,10 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, } else { slice_ncells[idim] = ( slice_hi[idim] - slice_lo[idim] + 1 ) - / cr_ratio[idim]; + / slice_cr_ratio[idim]; - int refined_ncells = (slice_hi[idim] - slice_lo[idim]) + 1 ; - if ( cr_ratio[idim] > 1) { + int refined_ncells = slice_hi[idim] - slice_lo[idim] + 1 ; + if ( slice_cr_ratio[idim] > 1) { coarsen = true; // modify slice_grid_size if >= refines_cells // @@ -114,13 +108,8 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, Vector<DistributionMapping> sdmap(1); sdmap[0] = DistributionMapping{sba[0]}; - if ( slicetypeToBeConverted==1 ) { - smf.reset(new MultiFab(amrex::convert(sba[0],SliceType), - sdmap[0],ncomp,nghost)); - } - else { - smf.reset(new MultiFab(sba[0],sdmap[0],ncomp,nghost)); - } + smf.reset(new MultiFab(amrex::convert(sba[0],SliceType), sdmap[0], + ncomp, nghost)); // Copy data from domain to slice that has same cell size as that of // // the domain mf. src and dst have the same number of ghost cells // @@ -132,14 +121,14 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, ncomp, nghost, slice_lo, slice_hi, SliceType, real_box); } - if ( coarsen == false ) { -// amrex::Print() << " Cell sizes are equal. No averaging required for slice data " << "\n"; + + if (coarsen == false) { return smf; } else if ( coarsen == true ) { Vector<BoxArray> crse_ba(1); crse_ba[0] = sba[0]; - crse_ba[0].coarsen(cr_ratio); + crse_ba[0].coarsen(slice_cr_ratio); AMREX_ALWAYS_ASSERT(crse_ba[0].size() == sba[0].size()); @@ -165,39 +154,39 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, int scomp = 0; int dcomp = 0; - amrex::IntVect cctype(AMREX_D_DECL(0,0,0)); + IntVect cctype(AMREX_D_DECL(0,0,0)); if( SliceType==cctype ) { amrex::amrex_avgdown(Dst_bx, Dst_fabox, Src_fabox, dcomp, scomp, - ncomp, cr_ratio); + ncomp, slice_cr_ratio); } - amrex::IntVect ndtype(AMREX_D_DECL(1,1,1)); + IntVect ndtype(AMREX_D_DECL(1,1,1)); if( SliceType == ndtype ) { amrex::amrex_avgdown_nodes(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio); + scomp, ncomp, slice_cr_ratio); } if( SliceType == WarpX::Ex_nodal_flag ) { amrex::amrex_avgdown_edges(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio, 0); + scomp, ncomp, slice_cr_ratio, 0); } if( SliceType == WarpX::Ey_nodal_flag) { amrex::amrex_avgdown_edges(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio, 1); + scomp, ncomp, slice_cr_ratio, 1); } if( SliceType == WarpX::Ez_nodal_flag ) { amrex::amrex_avgdown_edges(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio, 2); + scomp, ncomp, slice_cr_ratio, 2); } if( SliceType == WarpX::Bx_nodal_flag) { amrex::amrex_avgdown_faces(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio, 0); + scomp, ncomp, slice_cr_ratio, 0); } if( SliceType == WarpX::By_nodal_flag ) { amrex::amrex_avgdown_faces(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio, 1); + scomp, ncomp, slice_cr_ratio, 1); } if( SliceType == WarpX::Bz_nodal_flag ) { amrex::amrex_avgdown_faces(Dst_bx, Dst_fabox, Src_fabox, dcomp, - scomp, ncomp, cr_ratio, 2); + scomp, ncomp, slice_cr_ratio, 2); } if ( mfi_dst.isValid() ) { @@ -219,7 +208,6 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, * for example, if the coarsening ratio is 3, 5 or 6, which is not an exponent of 2, * then the value of coarsening ratio is modified to the nearest exponent of 2. * The default value for coarsening ratio is 1. - * \param slice_realbox a Real box for defining the slice dimension. * slice_realbox.lo and slice_realbox.hi are set equal to the simulation domain lo and hi * if for the user-input for the slice lo and hi coordinates are outside the domain. * If the slice_realbox.lo and slice_realbox.hi coordinates do not align with the data @@ -229,28 +217,36 @@ std::unique_ptr<MultiFab> CreateSlice( const amrex::MultiFab& mf, * If slice_realbox.lo==slice_realbox.hi, then that dimension has only one cell and no * modifications are made to the value. If the lo and hi do not align with a data point, * then it is flagged for interpolation. + * \param real_box a Real box defined for the underlying domain. + * \param slice_realbox a Real box for defining the slice dimension. + * \param slice_cc_nd_box a Real box for defining the modified lo and hi of the slice + * such that the coordinates align with the underlying data points. + * If the dimension is reduced to have only one cell, the slice_realbox is not modified and * instead the values are interpolated to the coordinate from the nearest data points. + * \param slice_cr_ratio contains values of the coarsening ratio which may be modified + * if the input values do not satisfy coarsenability conditions. * \param slice_lo and slice_hi are the index values of the slice + * \param interp_lo are set to 0 or 1 if they are flagged for interpolation. * The slice shares the same index space as that of the simulation domain. */ -void CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, - RealBox &slice_realbox, IntVect &cr_ratio, IntVect slice_cr_ratio, - Vector<Geometry> dom_geom, IntVect const SliceType, IntVect &slice_lo, - IntVect &slice_hi, IntVect &interp_lo) +void +CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, + RealBox &slice_realbox, IntVect &slice_cr_ratio, + Vector<Geometry> dom_geom, IntVect const SliceType, + IntVect &slice_lo, IntVect &slice_hi, IntVect &interp_lo) { + IntVect slice_lo2(AMREX_D_DECL(0,0,0)); for ( int idim = 0; idim < AMREX_SPACEDIM; ++idim) - { - // Modify coarsening ratio if the input valie is not an exponnt of 2 for AMR // + { + // Modify coarsening ratio if the input value is not an exponent of 2 for AMR // if ( slice_cr_ratio[idim] > 0 ) { int log_cr_ratio = floor ( log2( slice_cr_ratio[idim])); slice_cr_ratio[idim] = exp2( log_cr_ratio ); } - - // Default coarsening ratio i 1 // - if ( slice_cr_ratio[idim] > 0 ) cr_ratio[idim] = slice_cr_ratio[idim]; - + + //// Default coarsening ratio is 1 // // Modify lo if input is out of bounds // if ( slice_realbox.lo(idim) < real_box.lo(idim) ) { slice_realbox.setLo( idim, real_box.lo(idim)); @@ -266,17 +262,16 @@ void CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, " Modified it in dimension " << idim << " to be aligned with the domain box\n"; } - + // Factor to ensure index values computation depending on index type // double fac = ( 1.0 - SliceType[idim] )*dom_geom[0].CellSize(idim)*0.5; - // if dimension is reduced to one cell length // if ( slice_realbox.hi(idim) - slice_realbox.lo(idim) <= 0) { slice_cc_nd_box.setLo( idim, slice_realbox.lo(idim) ); slice_cc_nd_box.setHi( idim, slice_realbox.hi(idim) ); - - if ( slice_cr_ratio[idim] > 1) cr_ratio[idim] = 1; + + if ( slice_cr_ratio[idim] > 1) slice_cr_ratio[idim] = 1; // check for interpolation -- compute index lo with floor and ceil if ( slice_cc_nd_box.lo(idim) - real_box.lo(idim) >= fac ) { @@ -295,7 +290,7 @@ void CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, - (real_box.lo(idim) ) ) / dom_geom[0].CellSize(idim) ); } - + // flag for interpolation -- if reduced dimension location // // does not align with data point // if ( slice_lo[idim] == slice_lo2[idim]) { @@ -306,60 +301,59 @@ void CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, else { interp_lo[idim] = 1; } - + // ncells = 1 if dimension is reduced // slice_hi[idim] = slice_lo[idim] + 1; } else { // moving realbox.lo and reabox.hi to nearest coarsenable grid point // - int index_lo = floor(((slice_realbox.lo(idim) - (real_box.lo(idim))) - / dom_geom[0].CellSize(idim) ) +fac * 1E-10); - int index_hi = ceil(((slice_realbox.hi(idim) - (real_box.lo(idim))) - / dom_geom[0].CellSize(idim) ) - fac * 1E-10); + int index_lo = floor(((slice_realbox.lo(idim) + 1E-10 + - (real_box.lo(idim))) / dom_geom[0].CellSize(idim))); + int index_hi = ceil(((slice_realbox.hi(idim) - 1E-10 + - (real_box.lo(idim))) / dom_geom[0].CellSize(idim))); bool modify_cr = true; - + while ( modify_cr == true) { int lo_new = index_lo; - int hi_new = index_hi; - int mod_lo = index_lo % cr_ratio[idim]; - int mod_hi = index_hi % cr_ratio[idim]; + int hi_new = index_hi; + int mod_lo = index_lo % slice_cr_ratio[idim]; + int mod_hi = index_hi % slice_cr_ratio[idim]; modify_cr = false; - + // To ensure that the index.lo is coarsenable // if ( mod_lo > 0) { lo_new = index_lo - mod_lo; } // To ensure that the index.hi is coarsenable // if ( mod_hi > 0) { - hi_new = index_hi + (cr_ratio[idim] - mod_hi); + hi_new = index_hi + (slice_cr_ratio[idim] - mod_hi); } - + // If modified index.hi is > baselinebox.hi, reduce cr ratio // // and provide more points that asked for // if ( (hi_new * dom_geom[0].CellSize(idim)) - > real_box.hi(idim) - real_box.lo(idim) ) + > real_box.hi(idim) - real_box.lo(idim) + dom_geom[0].CellSize(idim)*0.01 ) { - cr_ratio[idim] = cr_ratio[idim]/2; + slice_cr_ratio[idim] = slice_cr_ratio[idim]/2; modify_cr = true; } - + int ncells = (hi_new - lo_new); - // If refined cells is not an integer multiple of cr ratio // // then reduce coarsening ratio by factor of 2 // - - if ( ( ncells % cr_ratio[idim] ) != 0 ) { - cr_ratio[idim] = cr_ratio[idim]/2; + + if ( ( ncells % slice_cr_ratio[idim] ) != 0 ) { + slice_cr_ratio[idim] = slice_cr_ratio[idim]/2; modify_cr = true; } - + if ( modify_cr == false ) { index_lo = lo_new; index_hi = hi_new; } - + slice_lo[idim] = index_lo; slice_hi[idim] = index_hi - 1; // since default is cell-centered } @@ -371,7 +365,6 @@ void CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, slice_cc_nd_box.setHi( idim, slice_realbox.hi(idim) - fac ); } } - } @@ -380,10 +373,11 @@ void CheckSliceInput( const RealBox real_box, RealBox &slice_cc_nd_box, * \param interp_lo is an IntVect which is flagged as 1, if interpolation is required in the dimension. */ -void InterpolateSliceValues( MultiFab& smf, IntVect interp_lo, - RealBox slice_realbox, Vector<Geometry> geom, int ncomp, - int nghost, IntVect slice_lo, IntVect slice_hi, - IntVect SliceType, const RealBox real_box) +void +InterpolateSliceValues(MultiFab& smf, IntVect interp_lo, RealBox slice_realbox, + Vector<Geometry> geom, int ncomp, int nghost, + IntVect slice_lo, IntVect slice_hi, IntVect SliceType, + const RealBox real_box) { for (MFIter mfi(smf); mfi.isValid(); ++mfi) { @@ -403,9 +397,11 @@ void InterpolateSliceValues( MultiFab& smf, IntVect interp_lo, } -void InterpolateLo(const Box& bx, FArrayBox &fabox, IntVect slice_lo, - Vector<Geometry> geom, int idir, IntVect IndType, RealBox slice_realbox, - int srccomp, int ncomp, int nghost, const RealBox real_box ) +void +InterpolateLo(const Box& bx, FArrayBox &fabox, IntVect slice_lo, + Vector<Geometry> geom, int idir, IntVect IndType, + RealBox slice_realbox, int srccomp, int ncomp, + int nghost, const RealBox real_box ) { auto fabarr = fabox.array(); const auto lo = amrex::lbound(bx); diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index 18d89951d..8ae219ddf 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -65,6 +65,26 @@ WarpX::MoveWindow (bool move_j) RealBox new_box(new_lo, new_hi); Geometry::ProbDomain(new_box); + // Moving slice coordinates - lo and hi - with moving window // + // slice box is modified only if slice diagnostics is initialized in input // + if ( slice_plot_int > 0 ) + { + Real new_slice_lo[AMREX_SPACEDIM]; + Real new_slice_hi[AMREX_SPACEDIM]; + const Real* current_slice_lo = slice_realbox.lo(); + const Real* current_slice_hi = slice_realbox.hi(); + for ( int i = 0; i < AMREX_SPACEDIM; i++) { + new_slice_lo[i] = current_slice_lo[i]; + new_slice_hi[i] = current_slice_hi[i]; + } + int num_shift_base_slice = static_cast<int> ((moving_window_x - + current_slice_lo[dir]) / cdx[dir]); + new_slice_lo[dir] = current_slice_lo[dir] + num_shift_base_slice*cdx[dir]; + new_slice_hi[dir] = current_slice_hi[dir] + num_shift_base_slice*cdx[dir]; + slice_realbox.setLo(new_slice_lo); + slice_realbox.setHi(new_slice_hi); + } + int num_shift = num_shift_base; int num_shift_crse = num_shift; |