diff options
Diffstat (limited to 'Carpet/CarpetRegrid/src/manualgridpoints.cc')
-rw-r--r-- | Carpet/CarpetRegrid/src/manualgridpoints.cc | 131 |
1 files changed, 131 insertions, 0 deletions
diff --git a/Carpet/CarpetRegrid/src/manualgridpoints.cc b/Carpet/CarpetRegrid/src/manualgridpoints.cc new file mode 100644 index 000000000..2a8bedb5e --- /dev/null +++ b/Carpet/CarpetRegrid/src/manualgridpoints.cc @@ -0,0 +1,131 @@ +#include <assert.h> + +#include <sstream> +#include <vector> + +#include "cctk.h" +#include "cctk_Parameters.h" + +#include "gh.hh" + +#include "carpet.hh" +#include "regrid.hh" + +extern "C" { + static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetRegrid/src/manualgridpoints.cc,v 1.5 2004/07/02 10:14:51 tradke Exp $"; + CCTK_FILEVERSION(Carpet_CarpetRegrid_manualgridpoints_cc); +} + + + +namespace CarpetRegrid { + + using namespace std; + using namespace Carpet; + + + + int ManualGridpoints (cGH const * const cctkGH, + gh<dim> const & hh, + gh<dim>::rexts & bbsss, + gh<dim>::rbnds & obss, + gh<dim>::rprocs & pss) + { + DECLARE_CCTK_PARAMETERS; + + if (refinement_levels > 4) { + CCTK_WARN (0, "Cannot currently specify manual refinement regions for more than 4 refinement levels"); + } + assert (refinement_levels >= 1 && refinement_levels <= 4); + + // do nothing if the levels already exist + if (reflevel == refinement_levels) return 0; + + assert (bbsss.size() >= 1); + + bbsss.resize (refinement_levels); + obss.resize (refinement_levels); + pss.resize (refinement_levels); + + vector<ivect> ilower(3), iupper(3); + ilower.at(0) = ivect (l1ixmin, l1iymin, l1izmin); + iupper.at(0) = ivect (l1ixmax, l1iymax, l1izmax); + ilower.at(1) = ivect (l2ixmin, l2iymin, l2izmin); + iupper.at(1) = ivect (l2ixmax, l2iymax, l2izmax); + ilower.at(2) = ivect (l3ixmin, l3iymin, l3izmin); + iupper.at(2) = ivect (l3ixmax, l3iymax, l3izmax); + + assert (! smart_outer_boundaries); + + for (size_t rl=1; rl<bbsss.size(); ++rl) { + + bbvect const ob (false); + + vector<ibbox> bbs; + gh<dim>::cbnds obs; + + ManualGridpoints_OneLevel + (cctkGH, hh, rl,refinement_levels, + ilower.at(rl-1), iupper.at(rl-1), ob, bbs, obs); + + // make multiprocessor aware + gh<dim>::cprocs ps; + SplitRegions (cctkGH, bbs, obs, ps); + + // make multigrid aware + vector<vector<ibbox> > bbss; + MakeMultigridBoxes (cctkGH, bbs, obs, bbss); + + bbsss.at(rl) = bbss; + obss.at(rl) = obs; + pss.at(rl) = ps; + + } // for rl + + return 1; + } + + + + void ManualGridpoints_OneLevel (const cGH * const cctkGH, + const gh<dim> & hh, + const int rl, + const int numrl, + const ivect ilower, + const ivect iupper, + const bbvect obound, + vector<ibbox> & bbs, + vector<bbvect> & obs) + { + const ivect rstr = hh.baseextent.stride(); + const ivect rlb = hh.baseextent.lower(); + const ivect rub = hh.baseextent.upper(); + + const int levfac = ipow(hh.reffact, rl); + assert (all (rstr % levfac == 0)); + const ivect str (rstr / levfac); + const ivect lb (ilower); + const ivect ub (iupper); + if (! all(lb>=rlb && ub<=rub)) { + ostringstream buf; + buf << "The refinement region boundaries for refinement level #" << rl << " are not within the main grid. Allowed are the grid point boundaries " << rlb << " - " << rub << "; specified were " << lb << " - " << ub << ends; + CCTK_WARN (0, buf.str().c_str()); + } + if (! all(lb<=ub)) { + ostringstream buf; + buf << "The refinement region boundaries for refinement level #" << rl << " have the upper boundary (" << ub << ") less than the lower boundary (" << lb << ")" << ends; + CCTK_WARN (0, buf.str().c_str()); + } + if (! all(lb%str==0 && ub%str==0)) { + CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, + "The refinement region boundaries for refinement level #%d are not a multiple of the stride for that level", rl); + } + assert (all(lb>=rlb && ub<=rub)); + assert (all(lb<=ub)); + assert (all(lb%str==0 && ub%str==0)); + + bbs.push_back (ibbox(lb, ub, str)); + obs.push_back (obound); + } + +} // namespace CarpetRegrid |