From d4edcbd243266d1d6bf29f5af037d4b973af15ea Mon Sep 17 00:00:00 2001 From: "Aaron M. Lattanzi" <103702284+AMLattanzi@users.noreply.github.com> Date: Thu, 20 Jul 2023 10:49:22 -0700 Subject: [PATCH] Towards Dynamic Mesh Refinement (#1172) * Allow only filling of external ghost cells and enable dynamic mesh refinement with relaxation/set regions. * Delete backtrace file. * Clean up. * Relocate crse register after it has been fillpatched. --------- Co-authored-by: Aaron Lattanzi --- Docs/sphinx_doc/MeshRefinement.rst | 11 +- Exec/RegTests/DynamicRefinement/inputs | 9 +- Source/BoundaryConditions/ERF_FillPatch.cpp | 4 +- Source/BoundaryConditions/ERF_FillPatcher.H | 40 +++--- Source/BoundaryConditions/ERF_FillPatcher.cpp | 79 +++++++++--- Source/ERF.H | 8 +- Source/ERF.cpp | 119 +++++++++++------- Source/TimeIntegration/ERF_TimeStep.cpp | 7 ++ Source/TimeIntegration/ERF_advance_dycore.cpp | 2 +- Source/TimeIntegration/TI_slow_rhs_fun.H | 12 +- 10 files changed, 197 insertions(+), 94 deletions(-) diff --git a/Docs/sphinx_doc/MeshRefinement.rst b/Docs/sphinx_doc/MeshRefinement.rst index ec9097fa4..cd0e36f72 100644 --- a/Docs/sphinx_doc/MeshRefinement.rst +++ b/Docs/sphinx_doc/MeshRefinement.rst @@ -135,9 +135,9 @@ and face-baced normal momenta on the coarse-fine interface, the coarse data is c interpolated to the fine mesh. The interpolated data is utilized to specify ghost cell data (outside of the valid fine region) as well as specify and relax data inside the lateral boundaries of the fine region. More specifically, a user may specify the total width of the interior -Dirichlet and relaxation region with ``erf.wrfbdy_width = `` (yellow + blue) +Dirichlet and relaxation region with ``erf.cf_width = `` (yellow + blue) and analogously the width of the interior Dirichlet region may be specified with -``erf.wrfbdy_set_width = `` (yellow). +``erf.cf_set_width = `` (yellow). .. |wrfbdy| image:: figures/wrfbdy_BCs.png :width: 600 @@ -169,7 +169,12 @@ the RHS (:math:`F`) is given by the following: where :math:`G` is the RHS of the NS equations, :math:`\psi^{\prime}` is the predicted update without relaxation, :math:`\psi^{FP}` is the fine patch data obtained from space-time interpolation of the -coarse mesh, and :math:`n` is the minimum number of grid point from a lateral boundary. +coarse mesh, and :math:`n` is the minimum number of grid point from a lateral boundary. Finally, we +note that time dependent Dirichlet data, provided via an external file, may be enforced on the +lateral boundary conditions of the domain (coarsest mesh). For such cases, the relaxation region width +at the domain edges may be specified with ``erf.wrfbdy_width = `` (yellow + blue) while the +interior Dirichlet region may be specified with ``erf.wrfbdy_set_width = `` (yellow). + By two-way coupling, we mean that in additional to specifying ghost cell data (outside of the valid fine region), the fine mesh communicates data back to the coarse mesh in two ways: diff --git a/Exec/RegTests/DynamicRefinement/inputs b/Exec/RegTests/DynamicRefinement/inputs index 8f777969b..ead24d4e8 100644 --- a/Exec/RegTests/DynamicRefinement/inputs +++ b/Exec/RegTests/DynamicRefinement/inputs @@ -78,9 +78,8 @@ erf.hi_scal2.max_level = 2 erf.hi_scal2.field_name = scalar erf.hi_scal2.value_greater = 10. -amr.n_error_buf = 4 -erf.regrid_int = 2 -erf.wrfbdy_width = 4 -erf.wrfbdy_set_width = 1 +amr.n_error_buf = 4 +erf.regrid_int = 2 +erf.cf_width = 4 +erf.cf_set_width = 1 -erf.use_NumDiff = 0.9 diff --git a/Source/BoundaryConditions/ERF_FillPatch.cpp b/Source/BoundaryConditions/ERF_FillPatch.cpp index 44dade13e..d2619abed 100644 --- a/Source/BoundaryConditions/ERF_FillPatch.cpp +++ b/Source/BoundaryConditions/ERF_FillPatch.cpp @@ -73,7 +73,7 @@ ERF::FillPatch (int lev, Real time, const Vector& mfs) } // var_idx // Coarse-Fine set region - if (lev>0 && coupling_type=="OneWay") { + if (lev>0 && coupling_type=="OneWay" && cf_set_width>0) { FPr_c[lev-1].fill(*mfs[Vars::cons], time, null_bc, domain_bcs_type, true); FPr_u[lev-1].fill(*mfs[Vars::xvel], time, null_bc, domain_bcs_type, true); FPr_v[lev-1].fill(*mfs[Vars::yvel], time, null_bc, domain_bcs_type, true); @@ -227,7 +227,7 @@ ERF::FillIntermediatePatch (int lev, Real time, } // var_idx // Coarse-Fine set region - if (lev>0 && coupling_type=="OneWay") { + if (lev>0 && coupling_type=="OneWay" && cf_set_width>0) { FPr_c[lev-1].fill(*mfs[Vars::cons], time, null_bc, domain_bcs_type, true); FPr_u[lev-1].fill(*mfs[Vars::xvel], time, null_bc, domain_bcs_type, true); FPr_v[lev-1].fill(*mfs[Vars::yvel], time, null_bc, domain_bcs_type, true); diff --git a/Source/BoundaryConditions/ERF_FillPatcher.H b/Source/BoundaryConditions/ERF_FillPatcher.H index 5874bf349..23c78fe14 100644 --- a/Source/BoundaryConditions/ERF_FillPatcher.H +++ b/Source/BoundaryConditions/ERF_FillPatcher.H @@ -13,6 +13,20 @@ public: amrex::Geometry const& cgeom, int nghost, int nghost_subset, int ncomp, amrex::InterpBase* interp); + ~ERFFillPatcher ( ) + { + delete m_cf_fine_data; + delete m_cf_fine_subset_data; + delete m_cf_crse_data[0]; + delete m_cf_crse_data[1]; + } + + void Define (amrex::BoxArray const& fba, amrex::DistributionMapping fdm, + amrex::Geometry const& fgeom, + amrex::BoxArray cba, amrex::DistributionMapping cdm, + amrex::Geometry const& cgeom, + int nghost, int nghost_subset, int ncomp, amrex::InterpBase* interp); + void registerCoarseData (amrex::Vector const& crse_data, amrex::Vector const& crse_time); @@ -33,9 +47,9 @@ private: int m_ncomp; amrex::InterpBase* m_interp; amrex::IntVect m_ratio; - amrex::Vector m_cf_crse_data; - amrex::MultiFab m_cf_fine_data; - amrex::MultiFab m_cf_fine_subset_data; + amrex::Vector m_cf_crse_data; + amrex::MultiFab* m_cf_fine_data; + amrex::MultiFab* m_cf_fine_subset_data; amrex::Vector m_crse_times; amrex::Real m_dt_crse; }; @@ -61,34 +75,32 @@ ERFFillPatcher::fill (amrex::MultiFab& mf, amrex::Real time, amrex::Real fac_old = 1.0 - fac_new; // Boundary condition operator - cbc(m_cf_crse_data[0], 0, m_ncomp, amrex::IntVect(0), time, 0); + cbc(*(m_cf_crse_data[0]), 0, m_ncomp, amrex::IntVect(0), time, 0); // Coarse MF to hold time interpolated data - amrex::MultiFab crse_data_time_interp(m_cf_crse_data[0].boxArray(), m_cf_crse_data[0].DistributionMap(), + amrex::MultiFab crse_data_time_interp(m_cf_crse_data[0]->boxArray(), m_cf_crse_data[0]->DistributionMap(), m_ncomp, amrex::IntVect{0}); // Time interpolate the coarse data amrex::MultiFab::LinComb(crse_data_time_interp, - fac_old, m_cf_crse_data[0], 0, - fac_new, m_cf_crse_data[1], 0, + fac_old, *(m_cf_crse_data[0]), 0, + fac_new, *(m_cf_crse_data[1]), 0, 0, m_ncomp, 0); // Ensure fine domain box is correct index type - amrex::Box fdest_dom = amrex::convert(m_fgeom.Domain(),m_cf_fine_data.boxArray().ixType()); + amrex::Box fdest_dom = amrex::convert(m_fgeom.Domain(),m_cf_fine_data->boxArray().ixType()); // Spatially interpolate the time-interpolated coarse data - amrex::FillPatchInterp(m_cf_fine_data, 0, crse_data_time_interp, 0, m_ncomp, amrex::IntVect(0), + amrex::FillPatchInterp(*m_cf_fine_data, 0, crse_data_time_interp, 0, m_ncomp, amrex::IntVect(0), m_cgeom, m_fgeom, fdest_dom, m_ratio, m_interp, bcs, 0); // Fill whole region or subset? if (fill_subset) { //amrex::MultiFab::Copy(m_cf_fine_subset_data, m_cf_fine_data, 0, 0, m_ncomp, amrex::IntVect{0}); - m_cf_fine_subset_data.ParallelCopy(m_cf_fine_data, 0, 0, m_ncomp, amrex::IntVect{0}, amrex::IntVect{0}); - mf.ParallelCopy(m_cf_fine_subset_data, 0, 0, m_ncomp, amrex::IntVect{0}, mf.nGrowVect()); + m_cf_fine_subset_data->ParallelCopy(*m_cf_fine_data, 0, 0, m_ncomp, amrex::IntVect{0}, amrex::IntVect{0}); + mf.ParallelCopy(*m_cf_fine_subset_data, 0, 0, m_ncomp, amrex::IntVect{0}, mf.nGrowVect()); } else { - mf.ParallelCopy(m_cf_fine_data, 0, 0, m_ncomp, amrex::IntVect{0}, mf.nGrowVect()); + mf.ParallelCopy(*m_cf_fine_data, 0, 0, m_ncomp, amrex::IntVect{0}, mf.nGrowVect()); } - - } #endif diff --git a/Source/BoundaryConditions/ERF_FillPatcher.cpp b/Source/BoundaryConditions/ERF_FillPatcher.cpp index 38e011317..e028ffc32 100644 --- a/Source/BoundaryConditions/ERF_FillPatcher.cpp +++ b/Source/BoundaryConditions/ERF_FillPatcher.cpp @@ -23,23 +23,64 @@ ERFFillPatcher::ERFFillPatcher (BoxArray const& fba, DistributionMapping fdm, int nghost, int nghost_subset, int ncomp, InterpBase* interp) : m_fba(fba), - m_cba(std::move(cba)), - m_fdm(std::move(fdm)), - m_cdm(std::move(cdm)), + m_cba(cba), + m_fdm(fdm), + m_cdm(cdm), m_fgeom(fgeom), - m_cgeom(cgeom), - m_nghost(nghost), - m_nghost_subset(nghost_subset), - m_ncomp(ncomp), - m_interp(interp) + m_cgeom(cgeom) { - AMREX_ALWAYS_ASSERT(nghost < 0); - AMREX_ALWAYS_ASSERT(nghost_subset < 0); - AMREX_ALWAYS_ASSERT(nghost < nghost_subset); AMREX_ALWAYS_ASSERT(fba.ixType() == cba.ixType()); // Vector to hold times for coarse data m_crse_times.resize(2); + m_cf_crse_data.resize(2); + + // Init MF patches + m_cf_fine_data = nullptr; m_cf_fine_subset_data = nullptr; + m_cf_crse_data[0] = nullptr; m_cf_crse_data[1] = nullptr; + + // Define the coarse and fine MFs + Define(fba, fdm, fgeom, cba, cdm, cgeom, + nghost, nghost_subset, ncomp, interp); +} + + +/* + * Redefine the coarse and fine patch MultiFabs. + * + * @param[in] fba BoxArray of data to be filled at fine level + * @param[in] fdm DistributionMapping of data to be filled at fine level + * @param[in] fgeom container of geometry infomation at fine level + * @param[in] cba BoxArray of data to be filled at coarse level + * @param[in] cdm DistributionMapping of data to be filled at coarse level + * @param[in] cgeom container of geometry infomation at coarse level + * @param[in] nghost number of ghost cells to be filled + * @param[in] ncomp number of components to be filled + * @param[in] interp interpolation operator to be used + */ +void ERFFillPatcher::Define (BoxArray const& fba, DistributionMapping fdm, + Geometry const& fgeom, + BoxArray cba, DistributionMapping cdm, + Geometry const& cgeom, + int nghost, int nghost_subset, + int ncomp, InterpBase* interp) +{ + AMREX_ALWAYS_ASSERT(nghost < 0); + AMREX_ALWAYS_ASSERT(nghost_subset <= 0); + AMREX_ALWAYS_ASSERT(nghost <= nghost_subset); + + // Set data memebers + m_fba = fba; m_cba = cba; + m_fdm = fdm; m_cdm = cdm; + m_fgeom = fgeom; m_cgeom = cgeom; + m_nghost = nghost; m_nghost_subset = nghost_subset; + m_ncomp = ncomp; m_interp = interp; + + // Delete old MFs if they exist + if (m_cf_fine_data) { + delete m_cf_fine_data; delete m_cf_fine_subset_data; + delete m_cf_crse_data[0]; delete m_cf_crse_data[1]; + } // Index type for the BL/BA IndexType m_ixt = fba.ixType(); @@ -108,16 +149,16 @@ ERFFillPatcher::ERFFillPatcher (BoxArray const& fba, DistributionMapping fdm, DistributionMapping gcf_dm(gcf_fba); DistributionMapping cf_dm_s(cf_fba_s); + // Fine patch to hold the time-interpolated state - m_cf_fine_data.define(gcf_fba, gcf_dm, m_ncomp, 0); + m_cf_fine_data = new MultiFab (gcf_fba, gcf_dm, m_ncomp, 0); // Fine subset patch to hold the time-interpolated state - //m_cf_fine_subset_data.define(cf_fba_s, cf_dm_s, m_ncomp, 0); - m_cf_fine_subset_data.define(cf_fba_s, gcf_dm, m_ncomp, 0); + m_cf_fine_subset_data = new MultiFab (cf_fba_s, gcf_dm, m_ncomp, 0); // Two coarse patches to hold the data to be interpolated - m_cf_crse_data.emplace_back(gcf_cba, gcf_dm, m_ncomp, 0); - m_cf_crse_data.emplace_back(gcf_cba, gcf_dm, m_ncomp, 0); + m_cf_crse_data[0] = new MultiFab (gcf_cba, gcf_dm, m_ncomp, 0); + m_cf_crse_data[1] = new MultiFab (gcf_cba, gcf_dm, m_ncomp, 0); } @@ -139,10 +180,10 @@ void ERFFillPatcher::registerCoarseData (Vector const& crse_dat // m_cf_crse_data into ghost cells in the z-dir. So we need // to include ghost cells for crse_data when doing the copy IntVect src_ng = crse_data[0]->nGrowVect(); - IntVect dst_ng = m_cf_crse_data[0].nGrowVect(); - m_cf_crse_data[0].ParallelCopy(*crse_data[0], 0, 0, m_ncomp, + IntVect dst_ng = m_cf_crse_data[0]->nGrowVect(); + m_cf_crse_data[0]->ParallelCopy(*crse_data[0], 0, 0, m_ncomp, src_ng, dst_ng, m_cgeom.periodicity()); // old data - m_cf_crse_data[1].ParallelCopy(*crse_data[1], 0, 0, m_ncomp, + m_cf_crse_data[1]->ParallelCopy(*crse_data[1], 0, 0, m_ncomp, src_ng, dst_ng, m_cgeom.periodicity()); // new data m_crse_times[0] = crse_time[0]; // time of "old" coarse data diff --git a/Source/ERF.H b/Source/ERF.H index 18cfed79a..22b6d4140 100644 --- a/Source/ERF.H +++ b/Source/ERF.H @@ -306,7 +306,11 @@ private: void define_grids_to_evolve (int lev, const amrex::BoxArray& ba); // NOLINT - void Define_ERFFillPatchers (); + void Construct_ERFFillPatchers (int lev); + + void Define_ERFFillPatchers (int lev); + + void Register_ERFFillPatchers (int lev); void init1DArrays (); @@ -508,6 +512,8 @@ private: #endif // Fillpatcher classes for coarse-fine boundaries + int cf_width{0}; + int cf_set_width{0}; amrex::Vector FPr_c; amrex::Vector FPr_u; amrex::Vector FPr_v; diff --git a/Source/ERF.cpp b/Source/ERF.cpp index a90bcb379..47594ed24 100644 --- a/Source/ERF.cpp +++ b/Source/ERF.cpp @@ -500,12 +500,6 @@ ERF::InitData () m_most->update_fluxes(lev); } - // NOTE: we must set up the FillPatcher object before calling - // WritePlotFile because WritePlotFile calls FillPatch - if (coupling_type=="OneWay") { - Define_ERFFillPatchers(); - } - if (solverChoice.use_rayleigh_damping) { initRayleigh(); @@ -562,6 +556,13 @@ ERF::InitData () // Fill ghost cells/faces for (int lev = 0; lev <= finest_level; ++lev) { + // NOTE: we must set up the FillPatcher object before calling + // FillPatch at a fine level + if (coupling_type=="OneWay" && cf_width>0 && lev>0) { + Construct_ERFFillPatchers(lev); + Register_ERFFillPatchers(lev); + } + auto& lev_new = vars_new[lev]; FillPatch(lev, t_new[lev], @@ -912,12 +913,20 @@ ERF::ReadParameters () // Specify whether ingest boundary planes of data pp.query("input_bndry_planes", input_bndry_planes); - // Query the set and total widths for interior ghost cells + // Query the set and total widths for wrfbdy interior ghost cells pp.query("wrfbdy_width", wrfbdy_width); pp.query("wrfbdy_set_width", wrfbdy_set_width); AMREX_ALWAYS_ASSERT(wrfbdy_width >= 0); AMREX_ALWAYS_ASSERT(wrfbdy_set_width >= 0); AMREX_ALWAYS_ASSERT(wrfbdy_width >= wrfbdy_set_width); + + // Query the set and total widths for crse-fine interior ghost cells + pp.query("cf_width", cf_width); + pp.query("cf_set_width", cf_set_width); + AMREX_ALWAYS_ASSERT(cf_width >= 0); + AMREX_ALWAYS_ASSERT(cf_set_width >= 0); + AMREX_ALWAYS_ASSERT(cf_width >= cf_set_width); + amrex::Print() << "READ CF params: " << cf_width << ' ' << cf_set_width << "\n"; } #ifdef ERF_USE_MULTIBLOCK @@ -1093,10 +1102,11 @@ ERF::AverageDownTo (int crse_lev) // NOLINT void ERF::define_grids_to_evolve (int lev, const BoxArray& ba) // NOLINT { - int width = wrfbdy_set_width; + Box domain(geom[lev].Domain()); if (lev == 0 && ( init_type == "real" || init_type == "metgrid" ) ) { + int width = wrfbdy_set_width; Box shrunk_domain(domain); shrunk_domain.grow(0,-width); shrunk_domain.grow(1,-width); @@ -1105,6 +1115,7 @@ ERF::define_grids_to_evolve (int lev, const BoxArray& ba) // NOLINT for (int i = 0; i < N; ++i) bl.push_back(ba[i] & shrunk_domain); grids_to_evolve[lev].define(std::move(bl)); } else if (lev == 1 && regrid_int < 0) { + int width = cf_set_width; Box shrunk_domain(boxes_at_level[lev][0]); shrunk_domain.grow(0,-width); shrunk_domain.grow(1,-width); @@ -1115,6 +1126,7 @@ ERF::define_grids_to_evolve (int lev, const BoxArray& ba) // NOLINT #if 0 if (num_boxes_at_level[lev] > 1) { for (int i = 1; i < num_boxes_at_level[lev]; i++) { + int width = cf_set_width; Box shrunk_domain(boxes_at_level[lev][i]); shrunk_domain.grow(0,-width); shrunk_domain.grow(1,-width); @@ -1132,43 +1144,64 @@ ERF::define_grids_to_evolve (int lev, const BoxArray& ba) // NOLINT } void -ERF::Define_ERFFillPatchers () +ERF::Construct_ERFFillPatchers (int lev) { - if (finest_level==0) return; - - // Construct the FillPatcher Objects for level>0 - for (int lev = 1; lev <= finest_level; ++lev) { - auto& fine_new = vars_new[lev]; - auto& crse_new = vars_new[lev-1]; - auto& ba_fine = fine_new[Vars::cons].boxArray(); - auto& ba_crse = crse_new[Vars::cons].boxArray(); - auto& dm_fine = fine_new[Vars::cons].DistributionMap(); - auto& dm_crse = crse_new[Vars::cons].DistributionMap(); - - // NOTE: crse-fine set/relaxation only done on Rho/RhoTheta - FPr_c.emplace_back(ba_fine, dm_fine, geom[lev] , - ba_crse, dm_crse, geom[lev-1], - -wrfbdy_width, -wrfbdy_set_width, 2, &cell_cons_interp); - FPr_u.emplace_back(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] , - convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1], - -wrfbdy_width, -wrfbdy_set_width, 1, &face_linear_interp); - FPr_v.emplace_back(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] , - convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1], - -wrfbdy_width, -wrfbdy_set_width, 1, &face_linear_interp); - FPr_w.emplace_back(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] , - convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1], - -wrfbdy_width, -wrfbdy_set_width, 1, &face_linear_interp); - } + auto& fine_new = vars_new[lev]; + auto& crse_new = vars_new[lev-1]; + auto& ba_fine = fine_new[Vars::cons].boxArray(); + auto& ba_crse = crse_new[Vars::cons].boxArray(); + auto& dm_fine = fine_new[Vars::cons].DistributionMap(); + auto& dm_crse = crse_new[Vars::cons].DistributionMap(); + + // NOTE: crse-fine set/relaxation only done on Rho/RhoTheta + FPr_c.emplace_back(ba_fine, dm_fine, geom[lev] , + ba_crse, dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 2, &cell_cons_interp); + FPr_u.emplace_back(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] , + convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 1, &face_linear_interp); + FPr_v.emplace_back(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] , + convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 1, &face_linear_interp); + FPr_w.emplace_back(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] , + convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 1, &face_linear_interp); +} - // Register the coarse data for each FPr object - for (int lev = 0; lev < finest_level; ++lev) { - auto& lev_new = vars_new[lev]; - auto& lev_old = vars_old[lev]; - FPr_c[lev].registerCoarseData({&lev_old[Vars::cons], &lev_new[Vars::cons]}, {t_old[lev], t_new[lev]}); - FPr_u[lev].registerCoarseData({&lev_old[Vars::xvel], &lev_new[Vars::xvel]}, {t_old[lev], t_new[lev]}); - FPr_v[lev].registerCoarseData({&lev_old[Vars::yvel], &lev_new[Vars::yvel]}, {t_old[lev], t_new[lev]}); - FPr_w[lev].registerCoarseData({&lev_old[Vars::zvel], &lev_new[Vars::zvel]}, {t_old[lev], t_new[lev]}); - } +void +ERF::Define_ERFFillPatchers (int lev) +{ + auto& fine_new = vars_new[lev]; + auto& crse_new = vars_new[lev-1]; + auto& ba_fine = fine_new[Vars::cons].boxArray(); + auto& ba_crse = crse_new[Vars::cons].boxArray(); + auto& dm_fine = fine_new[Vars::cons].DistributionMap(); + auto& dm_crse = crse_new[Vars::cons].DistributionMap(); + + // NOTE: crse-fine set/relaxation only done on Rho/RhoTheta + FPr_c[lev-1].Define(ba_fine, dm_fine, geom[lev] , + ba_crse, dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 2, &cell_cons_interp); + FPr_u[lev-1].Define(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] , + convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 1, &face_linear_interp); + FPr_v[lev-1].Define(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] , + convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 1, &face_linear_interp); + FPr_w[lev-1].Define(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] , + convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1], + -cf_width, -cf_set_width, 1, &face_linear_interp); +} + +void +ERF::Register_ERFFillPatchers (int lev) +{ + auto& lev_new = vars_new[lev-1]; + auto& lev_old = vars_old[lev-1]; + FPr_c[lev-1].registerCoarseData({&lev_old[Vars::cons], &lev_new[Vars::cons]}, {t_old[lev-1], t_new[lev-1]}); + FPr_u[lev-1].registerCoarseData({&lev_old[Vars::xvel], &lev_new[Vars::xvel]}, {t_old[lev-1], t_new[lev-1]}); + FPr_v[lev-1].registerCoarseData({&lev_old[Vars::yvel], &lev_new[Vars::yvel]}, {t_old[lev-1], t_new[lev-1]}); + FPr_w[lev-1].registerCoarseData({&lev_old[Vars::zvel], &lev_new[Vars::zvel]}, {t_old[lev-1], t_new[lev-1]}); } #ifdef ERF_USE_MULTIBLOCK diff --git a/Source/TimeIntegration/ERF_TimeStep.cpp b/Source/TimeIntegration/ERF_TimeStep.cpp index 074e37b43..2c26ab0d2 100644 --- a/Source/TimeIntegration/ERF_TimeStep.cpp +++ b/Source/TimeIntegration/ERF_TimeStep.cpp @@ -34,6 +34,13 @@ ERF::timeStep (int lev, Real time, int iteration) int old_finest = finest_level; regrid(lev, time); + // NOTE: Def & Reg index lev backwards (so we add 1 here) + // Redefine & register the ERFFillpatcher objects + if (coupling_type=="OneWay" && cf_width>0) { + Define_ERFFillPatchers(lev+1); + Register_ERFFillPatchers(lev+1); + } + // mark that we have regridded this level already for (int k = lev; k <= finest_level; ++k) { last_regrid_step[k] = istep[k]; diff --git a/Source/TimeIntegration/ERF_advance_dycore.cpp b/Source/TimeIntegration/ERF_advance_dycore.cpp index c9d0ae490..858bb60f8 100644 --- a/Source/TimeIntegration/ERF_advance_dycore.cpp +++ b/Source/TimeIntegration/ERF_advance_dycore.cpp @@ -308,7 +308,7 @@ void ERF::advance_dycore(int level, mri_integrator.advance(state_old, state_new, old_time, dt_advance); // Register coarse data for coarse-fine fill - if (level0) { FPr_c[level].registerCoarseData({&cons_old, &cons_new}, {old_time, old_time + dt_advance}); FPr_u[level].registerCoarseData({&xvel_old, &xvel_new}, {old_time, old_time + dt_advance}); FPr_v[level].registerCoarseData({&yvel_old, &yvel_new}, {old_time, old_time + dt_advance}); diff --git a/Source/TimeIntegration/TI_slow_rhs_fun.H b/Source/TimeIntegration/TI_slow_rhs_fun.H index ecc69a078..c281c51b0 100644 --- a/Source/TimeIntegration/TI_slow_rhs_fun.H +++ b/Source/TimeIntegration/TI_slow_rhs_fun.H @@ -203,9 +203,9 @@ #endif // Compute RHS for fine interior ghost - if (level>0 && coupling_type=="OneWay" && nrk==0) { + if (level>0 && coupling_type=="OneWay" && cf_width>0 && nrk==0) { fine_compute_interior_ghost_RHS(new_stage_time, slow_dt, - wrfbdy_width-1, wrfbdy_set_width, + cf_width-1, cf_set_width, &FPr_c[level-1], &FPr_u[level-1], &FPr_v[level-1], &FPr_w[level-1], boxes_at_level[level], domain_bcs_type, S_rhs, S_data); @@ -276,7 +276,7 @@ z_phys_nd_src[level], detJ_cc[level], detJ_cc_new[level], mapfac_m[level], mapfac_u[level], mapfac_v[level] #if defined(ERF_USE_NETCDF) && (defined(ERF_USE_MOISTURE) || defined(ERF_USE_WARM_NO_PRECIP)) - ,moist_relax, moist_zero, bdy_time_interval, start_bdy_time, new_stage_time, + ,moist_relax, moist_zero, bdy_time_interval, start_bdy_time, new_stage_time, wrfbdy_width-1, wrfbdy_set_width, bdy_data_xlo, bdy_data_xhi, bdy_data_ylo, bdy_data_yhi #endif @@ -291,7 +291,7 @@ z_phys_nd[level], detJ_cc[level], detJ_cc[level], mapfac_m[level], mapfac_u[level], mapfac_v[level] #if defined(ERF_USE_NETCDF) && (defined(ERF_USE_MOISTURE) || defined(ERF_USE_WARM_NO_PRECIP)) - ,moist_relax, moist_zero, bdy_time_interval, start_bdy_time, new_stage_time, + ,moist_relax, moist_zero, bdy_time_interval, start_bdy_time, new_stage_time, wrfbdy_width-1, wrfbdy_set_width, bdy_data_xlo, bdy_data_xhi, bdy_data_ylo, bdy_data_yhi #endif @@ -349,9 +349,9 @@ #endif // Compute RHS for fine interior ghost - if (level>0 && coupling_type=="OneWay" && nrk==0) { + if (level>0 && coupling_type=="OneWay" && cf_width>0 && nrk==0) { fine_compute_interior_ghost_RHS(new_stage_time, slow_dt, - wrfbdy_width-1, wrfbdy_set_width, + cf_width-1, cf_set_width, &FPr_c[level-1], &FPr_u[level-1], &FPr_v[level-1], &FPr_w[level-1], boxes_at_level[level], domain_bcs_type, S_rhs, S_data);