From a2c6077906dcf45133dbe29d1664663934783bc4 Mon Sep 17 00:00:00 2001 From: Ian Hawke Date: Fri, 11 Feb 2005 13:57:00 +0000 Subject: AMR boundaries Add support for using boundaries with AMR. Seems to work with symmetry boundaries but gets a bit confused if you try refining outer boundaries. Requires that you setup the domain using CoordBase. darcs-hash:20050211135704-58c7f-6ed6ce1449de17599dd8b23abc909f5ba313d048.gz --- CarpetDev/CarpetAdaptiveRegrid/par/AMR1.par | 22 ++- CarpetDev/CarpetAdaptiveRegrid/par/AMR2.par | 76 ++++++++++ CarpetDev/CarpetAdaptiveRegrid/src/CAR.cc | 165 ++++++++++++++++++++- .../CarpetAdaptiveRegrid/src/CAR_Paramcheck.cc | 6 +- 4 files changed, 257 insertions(+), 12 deletions(-) create mode 100644 CarpetDev/CarpetAdaptiveRegrid/par/AMR2.par (limited to 'CarpetDev/CarpetAdaptiveRegrid') diff --git a/CarpetDev/CarpetAdaptiveRegrid/par/AMR1.par b/CarpetDev/CarpetAdaptiveRegrid/par/AMR1.par index 359c35d88..588d5de82 100644 --- a/CarpetDev/CarpetAdaptiveRegrid/par/AMR1.par +++ b/CarpetDev/CarpetAdaptiveRegrid/par/AMR1.par @@ -13,15 +13,25 @@ mol::initial_data_is_crap = "yes" ##carpet::veryverbose = "yes" grid::domain = "full" -grid::type = "byspacing" +grid::type = "coordbase" grid::avoid_origin = "no" -driver::global_nx = 51 -driver::global_ny = 51 -driver::global_nz = 51 -grid::dxyz = 0.02 + +coordbase::xmin = -0.5 +coordbase::ymin = -0.5 +coordbase::zmin = -0.5 +coordbase::xmax = 0.5 +coordbase::ymax = 0.5 +coordbase::zmax = 0.5 + +coordbase::dx = 0.02 +coordbase::dy = 0.02 +coordbase::dz = 0.02 + driver::ghost_size = 1 time::dtfac = 0.5 +carpet::domain_from_coordbase = yes + carpet::max_refinement_levels = 2 carpet::init_3_timelevels = yes carpetadaptiveregrid::max_error = 1.e-6 @@ -33,6 +43,8 @@ carpetadaptiveregrid::regrid_every = 1 #carpetadaptiveregrid::verbose = yes #carpetadaptiveregrid::veryverbose = yes +#carpet::verbose = yes +#carpet::veryverbose = yes cactus::terminate = "time" cactus::cctk_final_time = 0.5 diff --git a/CarpetDev/CarpetAdaptiveRegrid/par/AMR2.par b/CarpetDev/CarpetAdaptiveRegrid/par/AMR2.par new file mode 100644 index 000000000..e53527333 --- /dev/null +++ b/CarpetDev/CarpetAdaptiveRegrid/par/AMR2.par @@ -0,0 +1,76 @@ + +ActiveThorns = "coordbase SymBase NaNChecker carpetReduce CartGrid3D carpet carpetlib carpetadaptiveregrid Boundary IOBasic IOUtil carpetIOASCII IDWaveMoL carpetSlab WaveMoL Time MoL reflectionsymmetry" + + +IDWaveMoL::initial_data = "gaussian" + +wavemol::bound = "radiation" + +mol::initial_data_is_crap = "yes" + +##carpet::adaptive_stepsize = "yes" + +##carpet::veryverbose = "yes" + +grid::domain = "full" +grid::type = "coordbase" +grid::avoid_origin = "no" + +coordbase::xmin = -0.5 +coordbase::ymin = -0.5 +coordbase::zmin = 0.0 +coordbase::xmax = 0.5 +coordbase::ymax = 0.5 +coordbase::zmax = 0.5 + +coordbase::dx = 0.02 +coordbase::dy = 0.02 +coordbase::dz = 0.02 + +reflectionsymmetry::reflection_z = yes +reflectionsymmetry::avoid_origin_z = no +coordbase::boundary_size_z_lower = 2 +coordbase::boundary_shiftout_z_lower = 1 + +carpet::prolongation_order_space = 3 +carpet::prolongation_order_time = 2 +driver::ghost_size = 2 +time::dtfac = 0.5 + +carpet::domain_from_coordbase = yes + +carpet::max_refinement_levels = 2 +carpet::init_3_timelevels = yes +carpetadaptiveregrid::max_error = 1.e-6 +carpetadaptiveregrid::pad = 2 +carpetadaptiveregrid::min_width = 6 +carpetadaptiveregrid::error_var = "mol::errorestimate[0]" +#carpetadaptiveregrid::error_var = "wavemol::phi" +carpetadaptiveregrid::regrid_every = 1 + +#carpetadaptiveregrid::verbose = yes +#carpetadaptiveregrid::veryverbose = yes +#carpet::verbose = yes +#carpet::veryverbose = yes + +cactus::terminate = "time" +cactus::cctk_final_time = 0.5 + +iobasic::outScalar_every = 1 +iobasic::outScalar_vars = "wavemol::phi mol::errorestimate" + +iobasic::outInfo_every = 1 +iobasic::outInfo_vars = "wavemol::phi mol::errorestimate[0]" + +ioascii::out1D_every = 1 +ioascii::out1D_vars = "wavemol::scalarevolvemol_scalar mol::errorestimate" + +IO::out_dir = $. + +#mol::ode_method = "Generic" +#mol::mol_intermediate_steps = 4 +#mol::mol_num_scratch_levels = 3 +mol::ode_method = "RK45" +mol::mol_intermediate_steps = 6 +mol::mol_num_scratch_levels = 6 +mol::adaptive_stepsize = "no" diff --git a/CarpetDev/CarpetAdaptiveRegrid/src/CAR.cc b/CarpetDev/CarpetAdaptiveRegrid/src/CAR.cc index 5f198ba26..ce59cb533 100644 --- a/CarpetDev/CarpetAdaptiveRegrid/src/CAR.cc +++ b/CarpetDev/CarpetAdaptiveRegrid/src/CAR.cc @@ -44,6 +44,11 @@ namespace CarpetAdaptiveRegrid { const CCTK_INT& min_width, const CCTK_REAL& min_fraction, CCTK_INT& didit); } + + ivect pos2int (const cGH* const cctkGH, const gh& hh, + const rvect & rpos, const int rl); + rvect int2pos (const cGH* const cctkGH, const gh& hh, + const ivect & ipos, const int rl); CCTK_INT CarpetAdaptiveRegrid_Regrid (CCTK_POINTER_TO_CONST const cctkGH_, CCTK_POINTER const bbsss_, @@ -163,6 +168,16 @@ namespace CarpetAdaptiveRegrid { ivect imin = (bb.lower() - baseext.lower())/bb.stride(), imax = (bb.upper() - baseext.lower())/bb.stride(); + rvect physical_min, physical_max; + rvect interior_min, interior_max; + rvect exterior_min, exterior_max; + rvect base_spacing; + int ierr = GetDomainSpecification + (dim, &physical_min[0], &physical_max[0], + &interior_min[0], &interior_max[0], + &exterior_min[0], &exterior_max[0], &base_spacing[0]); + assert (!ierr); + BEGIN_LOCAL_COMPONENT_LOOP(cctkGH, CCTK_GF) { const CCTK_REAL *error_var_ptr = @@ -430,14 +445,104 @@ namespace CarpetAdaptiveRegrid { vector obs; while (! final.empty()) { ibbox bb = final.top(); final.pop(); - ivect lo = bb.lower(); - ivect hi = bb.upper(); - ivect str = bb.stride() / reffact; - newbbs.push_back (ibbox(lo,hi,str)); + ivect ilo = bb.lower(); + ivect ihi = bb.upper(); + rvect lo = int2pos(cctkGH, hh, ilo, reflevel); + rvect hi = int2pos(cctkGH, hh, ihi, reflevel); + rvect str = base_spacing * + ipow((CCTK_REAL)mgfact, basemglevel) / ipow(reffact, reflevel); + rbbox newbbcoord(lo, hi, str); + if (veryverbose) { + ostringstream buf; + buf << "Dealing with boundaries. Coord box is:" + << endl << newbbcoord; + CCTK_INFO(buf.str().c_str()); + } + // FIXME: Set the correct ob here. bbvect ob(false); + for (int d=0; d physical_max[d]) { + up[d] = exterior_max[d]; + } + rvect str = newbbcoord.stride(); + + // Set the ob if outside the physical domain + + ob[d][0] = + abs(lo[d] - exterior_min[d]) < 1.0e-6 * spacing[d]; + ob[d][1] = + abs(up[d] - exterior_max[d]) < 1.0e-6 * spacing[d]; + + newbbcoord = rbbox(lo, up, str); + } + if (veryverbose) { + ostringstream buf; + buf << "Done dealing with boundaries. Coord box is:" + << endl << newbbcoord << endl + << "obox is:" << endl << ob; + CCTK_INFO(buf.str().c_str()); + } + + // Convert back to integer coordinates + + ilo = pos2int(cctkGH, hh, newbbcoord.lower(), reflevel); + ihi = pos2int(cctkGH, hh, newbbcoord.upper(), reflevel); + ivect istr = bb.stride() / reffact; + + // Check that the width is sufficient + // This can only be too small if the domain was clipped + for (int d=0; d < dim; ++d) { + if (ihi[d] - ilo[d] < min_width * istr[d]) { + if (ob[d][0]) { + if (ob[d][1]) { + CCTK_WARN(0, "The domain is too small?!"); + } + ihi[d] = ilo[d] + min_width * istr[d]; + } + else if (ob[d][1]) { + if (ob[d][0]) { + CCTK_WARN(0, "The domain is too small?!"); + } + ilo[d] = ihi[d] - min_width * istr[d]; + } + else { + CCTK_WARN(0, "The grid is unclipped and too small?!"); + } + } + } + + ibbox newbb(ilo, ihi, istr); + + if (veryverbose) { + ostringstream buf; + buf << "After dealing with boundaries. Final box is:" + << endl << newbb; + CCTK_INFO(buf.str().c_str()); + } + + newbbs.push_back (newbb); obs.push_back(ob); } // FIXME: check if the newbbs is really different from the current bbs @@ -487,4 +592,56 @@ namespace CarpetAdaptiveRegrid { return do_recompose; } + ivect pos2int (const cGH* const cctkGH, const gh& hh, + const rvect & rpos, const int rl) + { + rvect global_lower, global_upper; + for (int d=0; d