diff --git a/Docs/sphinx_documentation/source/RuntimeParameters.rst b/Docs/sphinx_documentation/source/RuntimeParameters.rst index 4e9f419680..5840445624 100644 --- a/Docs/sphinx_documentation/source/RuntimeParameters.rst +++ b/Docs/sphinx_documentation/source/RuntimeParameters.rst @@ -787,6 +787,15 @@ Embedded Boundary .. tip:: Because AMReX currently does not support multi-cut cells, it would be a runtime error if multi-cut cells are left unfixed. +.. py:data:: eb2.plt_multiple_cuts + :type: bool + :value: false + + If this parameter is set to true and multi-cut cells are present, + we produce plot files to see the location of the multi-cut issue. + In the plot file, a cell value of 2 indicates regular, 0 means covered, + and values greater than 2 indicate multi-cut cells. + .. py:data:: eb2.maxiter :type: int :value: 32 diff --git a/Src/EB/AMReX_EB2_2D_C.cpp b/Src/EB/AMReX_EB2_2D_C.cpp index b99b5559c7..e420265583 100644 --- a/Src/EB/AMReX_EB2_2D_C.cpp +++ b/Src/EB/AMReX_EB2_2D_C.cpp @@ -203,7 +203,8 @@ int build_faces (Box const& bx, Array4 const& cell, Array4 const& fcx, Array4 const& fcy, GpuArray const& dx, GpuArray const& problo, - bool cover_multiple_cuts, int& nsmallfaces) noexcept + bool cover_multiple_cuts, int& nsmallfaces, + bool plt_multiple_cuts, Array4 const& mcx) noexcept { #ifdef AMREX_USE_FLOAT constexpr Real small = 1.e-5_rt; @@ -311,6 +312,10 @@ int build_faces (Box const& bx, Array4 const& cell, if (ncuts > 2) { Gpu::Atomic::Add(dp,1); } + if (plt_multiple_cuts) + { + mcx(i,j,k) = ncuts; + } } } }); @@ -342,7 +347,10 @@ int build_faces (Box const& bx, Array4 const& cell, nsmallfaces += *(hp+1); if (*hp > 0 && !cover_multiple_cuts) { - amrex::Abort("amrex::EB2::build_faces: more than 2 cuts not supported"); + if (!plt_multiple_cuts) + { + amrex::Abort("amrex::EB2::build_faces: more than 2 cuts not supported"); + } } return *hp; diff --git a/Src/EB/AMReX_EB2_3D_C.cpp b/Src/EB/AMReX_EB2_3D_C.cpp index 2d02e53bdc..ce6dd8ec5a 100644 --- a/Src/EB/AMReX_EB2_3D_C.cpp +++ b/Src/EB/AMReX_EB2_3D_C.cpp @@ -370,7 +370,8 @@ int build_faces (Box const& bx, Array4 const& cell, Array4 const& m2z, GpuArray const& dx, GpuArray const& problo, - bool cover_multiple_cuts) noexcept + bool cover_multiple_cuts, bool plt_multiple_cuts, + Array4 const& mcx, Array4 const& mcy, Array4 const& mcz) noexcept { Gpu::Buffer nmulticuts = {0}; int* hp = nmulticuts.hostData(); @@ -491,6 +492,11 @@ int build_faces (Box const& bx, Array4 const& cell, } else if (apx(i,j,k) == 1.0_rt) { fx(i,j,k) = Type::regular; } + + if (plt_multiple_cuts){ + mcx(i,j,k) = ncuts; + } + } }); @@ -599,6 +605,10 @@ int build_faces (Box const& bx, Array4 const& cell, } else if (apy(i,j,k) == 1.0_rt) { fy(i,j,k) = Type::regular; } + + if (plt_multiple_cuts){ + mcy(i,j,k) = ncuts; + } } }); @@ -707,6 +717,10 @@ int build_faces (Box const& bx, Array4 const& cell, } else if (apz(i,j,k) == 1.0_rt) { fz(i,j,k) = Type::regular; } + + if (plt_multiple_cuts){ + mcz(i,j,k) = ncuts; + } } }); @@ -768,7 +782,10 @@ int build_faces (Box const& bx, Array4 const& cell, } }); } else { - amrex::Abort("amrex::EB2::build_faces: more than 2 cuts not supported"); + if (!plt_multiple_cuts) + { + amrex::Abort("amrex::EB2::build_faces: more than 2 cuts not supported"); + } } } @@ -787,7 +804,7 @@ void build_cells (Box const& bx, Array4 const& cell, Array4 const& barea, Array4 const& bcent, Array4 const& bnorm, Array4 const& ctmp, Array4 const& levset, Real small_volfrac, Geometry const& geom, - bool extend_domain_face, bool cover_multiple_cuts, + bool extend_domain_face, bool cover_multiple_cuts, bool plt_multiple_cuts, int& nsmallcells, int& nmulticuts) noexcept { Gpu::Buffer n_smallcell_multicuts = {0,0}; @@ -932,7 +949,10 @@ void build_cells (Box const& bx, Array4 const& cell, if (nsmallcells > 0 || nmulticuts > 0) { if (!cover_multiple_cuts && nmulticuts > 0) { - amrex::Abort("amrex::EB2::build_cells: multi-cuts not supported"); + if (!plt_multiple_cuts) + { + amrex::Abort("amrex::EB2::build_cells: multi-cuts not supported"); + } } return; } else { diff --git a/Src/EB/AMReX_EB2_C.H b/Src/EB/AMReX_EB2_C.H index 49db3a270b..aeaac8d231 100644 --- a/Src/EB/AMReX_EB2_C.H +++ b/Src/EB/AMReX_EB2_C.H @@ -25,7 +25,8 @@ int build_faces (Box const& bx, Array4 const& cell, Array4 const& fcx, Array4 const& fcy, GpuArray const& dx, GpuArray const& problo, - bool cover_multiple_cuts, int& nsmallfaces) noexcept; + bool cover_multiple_cuts, int& nsmallfaces, + bool plt_multiple_cuts, Array4 const& mcx) noexcept; void build_cells (Box const& bx, Array4 const& cell, Array4 const& fx, Array4 const& fy, @@ -55,7 +56,8 @@ int build_faces (Box const& bx, Array4 const& cell, Array4 const& m2z, GpuArray const& dx, GpuArray const& problo, - bool cover_multiple_cuts) noexcept; + bool cover_multiple_cuts, bool plt_multiple_cuts, + Array4 const& mcx, Array4 const& mcy, Array4 const& mcz) noexcept; void build_cells (Box const& bx, Array4 const& cell, Array4 const& fx, Array4 const& fy, @@ -69,7 +71,7 @@ void build_cells (Box const& bx, Array4 const& cell, Array4 const& barea, Array4 const& bcent, Array4 const& bnorm, Array4 const& ctmp, Array4 const& levset, Real small_volfrac, Geometry const& geom, - bool extend_domain_face, bool cover_multiple_cuts, + bool extend_domain_face, bool cover_multiple_cuts, bool plt_multiple_cuts, int& nsmallcells, int& nmulticuts) noexcept; void set_connection_flags(Box const& bx, Box const& bxg1, diff --git a/Src/EB/AMReX_EB2_Level.H b/Src/EB/AMReX_EB2_Level.H index 7e30e51a6b..2d9470bfa6 100644 --- a/Src/EB/AMReX_EB2_Level.H +++ b/Src/EB/AMReX_EB2_Level.H @@ -14,6 +14,7 @@ #include #include #include +#include #ifdef AMREX_USE_OMP #include @@ -92,6 +93,9 @@ protected: Array m_areafrac; Array m_facecent; Array m_edgecent; + MultiFab m_multicut_fcx; + MultiFab m_multicut_fcy; + MultiFab m_multicut_fcz; iMultiFab m_cutcellmask; bool m_allregular = false; bool m_ok = false; @@ -167,12 +171,14 @@ GShopLevel::define_fine (G const& gshop, const Geometry& geom, Real small_volfrac = 1.e-14; #endif bool cover_multiple_cuts = false; + bool plt_multiple_cuts = false; int maxiter = 32; { ParmParse pp("eb2"); pp.queryAdd("small_volfrac", small_volfrac); pp.queryAdd("cover_multiple_cuts", cover_multiple_cuts); pp.queryAdd("maxiter", maxiter); + pp.queryAdd("plt_multiple_cuts", plt_multiple_cuts); } maxiter = std::min(100000, maxiter); @@ -290,6 +296,9 @@ GShopLevel::define_fine (G const& gshop, const Geometry& geom, m_bndryarea.define(m_grids, m_dmap, 1, ng, mf_info); m_bndrycent.define(m_grids, m_dmap, AMREX_SPACEDIM, ng, mf_info); m_bndrynorm.define(m_grids, m_dmap, AMREX_SPACEDIM, ng, mf_info); + m_multicut_fcx.define(m_grids, m_dmap, 1, ng, mf_info); + m_multicut_fcy.define(m_grids, m_dmap, 1, ng, mf_info); + m_multicut_fcz.define(m_grids, m_dmap, 1, ng, mf_info); for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) { m_areafrac[idim].define(amrex::convert(m_grids, IntVect::TheDimensionVector(idim)), m_dmap, 1, ng, mf_info); @@ -413,9 +422,14 @@ GShopLevel::define_fine (G const& gshop, const Geometry& geom, Array4 const& ym2 = M2[1].array(); Array4 const& zm2 = M2[2].array(); + Array4 const& mcx = m_multicut_fcx.array(mfi); + Array4 const& mcy = m_multicut_fcy.array(mfi); + Array4 const& mcz = m_multicut_fcz.array(mfi); + nmc = build_faces(vbx, cfg, ftx, fty, ftz, xdg, ydg, zdg, lst, xip, yip, zip, apx, apy, apz, fcx, fcy, fcz, - xm2, ym2, zm2, dx, problo, cover_multiple_cuts); + xm2, ym2, zm2, dx, problo, cover_multiple_cuts, + plt_multiple_cuts, mcx, mcy, mcz); cellflagtmp.resize(m_cellflag[mfi].box()); Elixir cellflagtmp_eli = cellflagtmp.elixir(); @@ -424,7 +438,8 @@ GShopLevel::define_fine (G const& gshop, const Geometry& geom, build_cells(vbx, cfg, ftx, fty, ftz, apx, apy, apz, fcx, fcy, fcz, xm2, ym2, zm2, dx, vfr, ctr, bar, bct, bnm, cfgtmp, lst, - small_volfrac, geom, extend_domain_face, cover_multiple_cuts, + small_volfrac, geom, extend_domain_face, + cover_multiple_cuts, plt_multiple_cuts, nsm, nmc); // Because it is used in a synchronous reduction kernel in @@ -466,8 +481,11 @@ GShopLevel::define_fine (G const& gshop, const Geometry& geom, clst, geom); } + Array4 const& mcx = m_multicut_fcx.array(mfi); + nmc = build_faces(vbx, cfg, ftx, fty, lst, xip, yip, apx, apy, fcx, fcy, - dx, problo, cover_multiple_cuts, nsm); + dx, problo, cover_multiple_cuts, nsm, + plt_multiple_cuts, mcx); build_cells(vbx, cfg, ftx, fty, apx, apy, dx, vfr, ctr, bar, bct, bnm, lst, small_volfrac, geom, extend_domain_face, @@ -479,6 +497,34 @@ GShopLevel::define_fine (G const& gshop, const Geometry& geom, } ParallelAllReduce::Sum({nsmallcells,nmulticuts}, ParallelContext::CommunicatorSub()); + + if (plt_multiple_cuts && nmulticuts > 0) + { + amrex::Print() << "Total number of multicuts = " << nmulticuts << "\n"; + amrex::Print() << "Plotting multicut locations..." << "\n"; + + #if (AMREX_SPACEDIM == 3) + if (m_multicut_fcx.max(0) > 2) + { + WriteSingleLevelPlotfile("plt.multicut.x", m_multicut_fcx, {"multicut_fcx"}, geom, 0.0, 0); + } + if (m_multicut_fcy.max(0) > 2) + { + WriteSingleLevelPlotfile("plt.multicut.y", m_multicut_fcy, {"multicut_fcy"}, geom, 0.0, 0); + } + if (m_multicut_fcz.max(0) > 2) + { + WriteSingleLevelPlotfile("plt.multicut.z", m_multicut_fcz, {"multicut_fcz"}, geom, 0.0, 0); + } + #elif (AMREX_SPACEDIM == 2) + if (m_multicut_fcx.max(0) > 2) + { + WriteSingleLevelPlotfile("plt.multicut", m_multicut_fcx, {"multicut_fcx"}, geom, 0.0, 0); + } + #endif + amrex::Abort("amrex::EB2:: more than 2 cuts not supported"); + } + if (nsmallcells == 0 && nmulticuts == 0) { break; } else {