From 04dff66c4ada038cbb2a6205a4d4b12caa14d309 Mon Sep 17 00:00:00 2001 From: schnetter <> Date: Mon, 17 Dec 2001 12:34:00 +0000 Subject: Fixed test cases for axial symmetry (changed parameter files.) Fixed test cases for axial symmetry (changed parameter files.) Extented Carpet to handle multiple grid components per processor. Currently works on one processor only. darcs-hash:20011217123400-07bb3-47796190d04061ea23798717f6714978bf7e1f55.gz --- Carpet/Carpet/options/carpet-lilypond-gcc3 | 5 +- Carpet/Carpet/src/Recompose.cc | 106 +++++++++++++++-------------- Carpet/Carpet/src/SetupGH.cc | 12 ++-- Carpet/Carpet/src/carpet_public.hh | 7 +- Carpet/Carpet/src/helpers.cc | 68 ++++++++---------- Carpet/Carpet/src/variables.cc | 5 +- Carpet/CarpetIOASCII/param.ccl | 6 +- Carpet/CarpetIOASCII/src/ioascii.cc | 11 ++- 8 files changed, 117 insertions(+), 103 deletions(-) (limited to 'Carpet') diff --git a/Carpet/Carpet/options/carpet-lilypond-gcc3 b/Carpet/Carpet/options/carpet-lilypond-gcc3 index 5900a6881..6bc774255 100644 --- a/Carpet/Carpet/options/carpet-lilypond-gcc3 +++ b/Carpet/Carpet/options/carpet-lilypond-gcc3 @@ -1,4 +1,4 @@ -# $Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/options/Attic/carpet-lilypond-gcc3,v 1.1 2001/12/07 18:24:15 schnetter Exp $ +# $Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/options/Attic/carpet-lilypond-gcc3,v 1.2 2001/12/17 13:34:00 schnetter Exp $ #CPP = /home/eschnett/gcc/bin/cpp -traditional -DTMPL_EXPLICIT -DCARPET_REAL CPP = cpp -traditional -DTMPL_EXPLICIT -DCARPET_REAL @@ -7,7 +7,8 @@ CXX = /home/eschnett/gcc/bin/g++ F77 = /home/eschnett/gcc/bin/g77 F90 = /home/eschnett/gcc/bin/g77 -CXXFLAGS = -DTMPL_EXPLICIT -DCARPET_REAL -ftemplate-depth-30 -Wl,-rpath,/home/eschnett/gcc/lib +CXXFLAGS = -DTMPL_EXPLICIT -DCARPET_REAL -ftemplate-depth-30 +LDFLAGS = -Wl,-rpath,/home/eschnett/gcc/lib C_OPTIMISE_FLAGS = -O2 -funroll-loops -march=pentiumpro CXX_OPTIMISE_FLAGS = -O2 -funroll-loops -march=pentiumpro diff --git a/Carpet/Carpet/src/Recompose.cc b/Carpet/Carpet/src/Recompose.cc index 1fded5cb8..2122d179f 100644 --- a/Carpet/Carpet/src/Recompose.cc +++ b/Carpet/Carpet/src/Recompose.cc @@ -15,7 +15,7 @@ #include "carpet.hh" -static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/Recompose.cc,v 1.12 2001/12/14 16:39:07 schnetter Exp $"; +static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/Recompose.cc,v 1.13 2001/12/17 13:34:01 schnetter Exp $"; @@ -31,8 +31,10 @@ namespace Carpet { - static void SplitRegions_AlongZ (const cGH* cgh, gh::cexts& bbss); - static void SplitRegions_AsSpecified (const cGH* cgh, gh::cexts& bbss); + static void SplitRegions_AlongZ (const cGH* cgh, + vector >& bbs); + static void SplitRegions_AsSpecified (const cGH* cgh, + vector >& bbs); static void MakeProcessors_RoundRobin (const cGH* cgh, const gh::rexts& bbss, @@ -40,32 +42,34 @@ namespace Carpet { - static bool RegionsOK (const gh::rexts& bbsss, - const gh::rprocs& pss); + static void CheckRegions (const gh::rexts& bbsss, + const gh::rprocs& pss); static void Adapt (const cGH* cgh, int reflevels, gh* hh); static void Output (const cGH* cgh, const gh* hh, const char* descr); - bool RegionsOK (const gh::rexts& bbsss, const gh::rprocs& pss) + void CheckRegions (const gh::rexts& bbsss, const gh::rprocs& pss) { // At least one level - if (bbsss.size() == 0) return false; + assert (bbsss.size() > 0); for (int rl=0; rl<(int)bbsss.size(); ++rl) { // No empty levels - if (bbsss[rl].size() == 0) return false; + assert (bbsss[rl].size() > 0); for (int c=0; c<(int)bbsss[rl].size(); ++c) { - // At least one multigrid level - if (bbsss[rl][c].size() == 0) return false; + // Just one multigrid level + assert (bbsss[rl][c].size() == 1); + for (int ml=0; ml<(int)bbsss[rl][c].size(); ++ml) { + assert (all(bbsss[rl][c][ml].stride() + == floor(pow((double)reffact,maxreflevels-rl-1)+0.5))); + } } } - if (pss.size() != bbsss.size()) return false; + assert (pss.size() == bbsss.size()); for (int rl=0; rl<(int)bbsss.size(); ++rl) { - if (pss[rl].size() != bbsss[rl].size()) return false; + assert (pss[rl].size() == bbsss[rl].size()); } - - return true; } @@ -74,7 +78,7 @@ namespace Carpet { const gh::rprocs& pss) { // Check the regions - assert (RegionsOK(bbsss,pss)); + CheckRegions (bbsss, pss); // Save the region information for the next regridding next_bbsss = bbsss; next_pss = pss; @@ -85,8 +89,6 @@ namespace Carpet { void Recompose (const cGH* cgh) { - DECLARE_CCTK_PARAMETERS; - assert (component == -1); Waypoint ("%*sRecompose", 2*reflevel, ""); @@ -94,7 +96,7 @@ namespace Carpet { if (!do_recompose) return; // Check the regions - assert (RegionsOK(next_bbsss, next_pss)); + CheckRegions (next_bbsss, next_pss); // Recompose hh->recompose (next_bbsss, next_pss); @@ -125,15 +127,13 @@ namespace Carpet { - // This routine is a leftover. It determines "automatically" how - // scalars and arrays should be refined. The user really should have - // a possibility to define how arrays are to be refined. + // This routine is a leftover. It determines "automatically" how + // scalars and arrays should be refined. The user really should + // have a possibility to define how arrays are to be refined. static void Adapt (const cGH* cgh, const int reflevels, gh* hh) { - DECLARE_CCTK_PARAMETERS; - - const int nprocs = CCTK_nProcs(cgh); - const int mglevels = 1; // for now + const int nprocs = CCTK_nProcs(cgh); + const int mglevels = 1; // for now vector > > bbss(reflevels); // note: what this routine calls "ub" is "ub+str" elsewhere vect rstr = hh->baseextent.stride(); @@ -174,8 +174,14 @@ namespace Carpet { assert (cub[dim-1] <= rub[dim-1]); bbs[c] = bbox(clb, cub-cstr, cstr); } - bbss[rl] = bbs; +// bbss[rl] = bbs; + bbss[rl].clear(); + assert (Carpet::hh->components(rl) % nprocs == 0); + for (int cc=0; cc<(int)Carpet::hh->components(rl)/nprocs; ++cc) { + bbss[rl].insert(bbss[rl].end(), bbs.begin(), bbs.end()); + } } + vector > > > bbsss = hh->make_multigrid_boxes(bbss, mglevels); @@ -239,14 +245,14 @@ namespace Carpet { - void SplitRegions (const cGH* cgh, gh::cexts& bbss) + void SplitRegions (const cGH* cgh, vector >& bbs) { DECLARE_CCTK_PARAMETERS; if (CCTK_EQUALS (processor_topology, "automatic")) { - SplitRegions_AlongZ (cgh, bbss); + SplitRegions_AlongZ (cgh, bbs); } else if (CCTK_EQUALS (processor_topology, "manual")) { - SplitRegions_AsSpecified (cgh, bbss); + SplitRegions_AsSpecified (cgh, bbs); } else { abort(); } @@ -254,22 +260,22 @@ namespace Carpet { - void SplitRegions_AlongZ (const cGH* cgh, gh::cexts& bbss) + void SplitRegions_AlongZ (const cGH* cgh, vector >& bbs) { // Something to do? - if (bbss.size() == 0) return; + if (bbs.size() == 0) return; + + const int nprocs = CCTK_nProcs(cgh); - const int nprocs = CCTK_nProcs(cgh); - const int mglevels = 1; // arbitrary value + if (nprocs==1) return; - assert (bbss.size() == 1); - assert (bbss[0].size() == 1); + assert (bbs.size() == 1); - const vect rstr = bbss[0][0].stride(); - const vect rlb = bbss[0][0].lower(); - const vect rub = bbss[0][0].upper() + rstr; + const vect rstr = bbs[0].stride(); + const vect rlb = bbs[0].lower(); + const vect rub = bbs[0].upper() + rstr; - vector > bbs(nprocs); + bbs.resize(nprocs); for (int c=0; c cstr = rstr; vect clb = rlb; @@ -285,27 +291,26 @@ namespace Carpet { assert (cub[dim-1] <= rub[dim-1]); bbs[c] = bbox(clb, cub-cstr, cstr); } - bbss = hh->make_reflevel_multigrid_boxes(bbs, mglevels); } - void SplitRegions_AsSpecified (const cGH* cgh, gh::cexts& bbss) + void SplitRegions_AsSpecified (const cGH* cgh, vector >& bbs) { DECLARE_CCTK_PARAMETERS; // Something to do? - if (bbss.size() == 0) return; + if (bbs.size() == 0) return; + + const int nprocs = CCTK_nProcs(cgh); - const int nprocs = CCTK_nProcs(cgh); - const int mglevels = 1; // arbitrary value + if (nprocs==1) return; - assert (bbss.size() == 1); - assert (bbss[0].size() == 1); + assert (bbs.size() == 1); - const vect rstr = bbss[0][0].stride(); - const vect rlb = bbss[0][0].lower(); - const vect rub = bbss[0][0].upper() + rstr; + const vect rstr = bbs[0].stride(); + const vect rlb = bbs[0].lower(); + const vect rub = bbs[0].upper() + rstr; const vect nprocs_dir (processor_topology_3d_x, processor_topology_3d_y, @@ -317,7 +322,7 @@ namespace Carpet { } assert (prod(nprocs_dir) == nprocs); - vector > bbs(nprocs); + bbs.resize(nprocs); assert (dim==3); for (int k=0; kmake_reflevel_multigrid_boxes(bbs, mglevels); } diff --git a/Carpet/Carpet/src/SetupGH.cc b/Carpet/Carpet/src/SetupGH.cc index 1efb5e3b3..431514da7 100644 --- a/Carpet/Carpet/src/SetupGH.cc +++ b/Carpet/Carpet/src/SetupGH.cc @@ -11,7 +11,7 @@ #include "carpet.hh" -static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/SetupGH.cc,v 1.16 2001/12/14 16:39:08 schnetter Exp $"; +static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/SetupGH.cc,v 1.17 2001/12/17 13:34:01 schnetter Exp $"; @@ -39,7 +39,8 @@ namespace Carpet { // Refinement information maxreflevels = max_refinement_levels; - maxreflevelfact = floor(pow((double)refinement_factor, maxreflevels-1) + 0.5); + reffact = refinement_factor; + maxreflevelfact = floor(pow((double)reffact, maxreflevels-1) + 0.5); // Ghost zones vect lghosts, ughosts; @@ -243,13 +244,14 @@ namespace Carpet { vector > bbs(1); bbs[0] = hh->baseextent; + SplitRegions (cgh, bbs); + vector > > bbss(1); bbss[0] = bbs; - SplitRegions (cgh, bbss); - gh::rexts bbsss; - bbsss = hh->make_multigrid_boxes(bbss, 1); + const int mglevels = 1; + bbsss = hh->make_multigrid_boxes(bbss, mglevels); gh::rprocs pss; MakeProcessors (cgh, bbsss, pss); diff --git a/Carpet/Carpet/src/carpet_public.hh b/Carpet/Carpet/src/carpet_public.hh index bcb483926..d5deb2b17 100644 --- a/Carpet/Carpet/src/carpet_public.hh +++ b/Carpet/Carpet/src/carpet_public.hh @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/carpet_public.hh,v 1.9 2001/12/14 16:39:09 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/carpet_public.hh,v 1.10 2001/12/17 13:34:01 schnetter Exp $ // It is assumed that the number of components of all arrays is equal // to the number of components of the grid functions, and that their @@ -35,6 +35,9 @@ namespace Carpet { // Maximum number of refinement levels extern int maxreflevels; + // Refinement factor + extern int reffact; + // Refinement factor on finest grid extern int maxreflevelfact; @@ -125,7 +128,7 @@ namespace Carpet { void RegisterRecomposeRegions (const gh::rexts& bbsss, const gh::rprocs& pss); - void SplitRegions (const cGH* cgh, gh::cexts& bbss); + void SplitRegions (const cGH* cgh, vector >& bbs); void MakeProcessors (const cGH* cgh, const gh::rexts& bbsss, gh::rprocs& pss); diff --git a/Carpet/Carpet/src/helpers.cc b/Carpet/Carpet/src/helpers.cc index fb36ba70c..fcc294cb6 100644 --- a/Carpet/Carpet/src/helpers.cc +++ b/Carpet/Carpet/src/helpers.cc @@ -11,7 +11,7 @@ #include "carpet.hh" -static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/helpers.cc,v 1.13 2001/12/14 16:39:10 schnetter Exp $"; +static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/helpers.cc,v 1.14 2001/12/17 13:34:01 schnetter Exp $"; @@ -294,37 +294,35 @@ namespace Carpet { // Local mode -- a component is active // Set Cactus parameters - for (int d=0; dboxes.size()); + assert (component < (int)dd->boxes[reflevel].size()); + assert (mglevel < (int)dd->boxes[reflevel][component].size()); const bbox& ext = dd->boxes[reflevel][component][mglevel].exterior; - cgh->cctk_lsh[d] = (ext.shape() / ext.stride())[d]; - cgh->cctk_lbnd[d] = (ext.lower() / ext.stride())[d]; - cgh->cctk_ubnd[d] = (ext.upper() / ext.stride())[d]; - assert (cgh->cctk_lsh[d]>=0 && cgh->cctk_lsh[d]<=cgh->cctk_gsh[d]); - assert (cgh->cctk_lbnd[d]>=0 && cgh->cctk_ubnd[d]cctk_gsh[d]); - assert (cgh->cctk_lbnd[d]<=cgh->cctk_ubnd[d]+1); -#if 0 - // No outer boundaries on the finer grids - cgh->cctk_bbox[2*d ] - = reflevel==0 && cgh->cctk_lbnd[d] == 0; - cgh->cctk_bbox[2*d+1] - = reflevel==0 && cgh->cctk_ubnd[d] == cgh->cctk_gsh[d]-1; -#else - // Do allow outer boundaries on the finer grids (but this is - // generally inconsistent -- c. f. periodicity) - const bbox& base = hh->baseextent; - cgh->cctk_bbox[2*d ] = (ext.lower() < base.lower())[d]; - cgh->cctk_bbox[2*d+1] = (ext.upper() > base.upper())[d]; -#endif - for (int stg=0; stgcctk_lssh[CCTK_LSSH_IDX(stg,d)] = cgh->cctk_lsh[d]; + for (int d=0; dcctk_lsh[d] = (ext.shape() / ext.stride())[d]; + cgh->cctk_lbnd[d] = (ext.lower() / ext.stride())[d]; + cgh->cctk_ubnd[d] = (ext.upper() / ext.stride())[d]; + assert (cgh->cctk_lsh[d]>=0 && cgh->cctk_lsh[d]<=cgh->cctk_gsh[d]); + assert (cgh->cctk_lbnd[d]>=0 && cgh->cctk_ubnd[d]cctk_gsh[d]); + assert (cgh->cctk_lbnd[d]<=cgh->cctk_ubnd[d]+1); + // Do allow outer boundaries on the finer grids + cgh->cctk_bbox[2*d ] = cgh->cctk_lbnd[d] == 0; + cgh->cctk_bbox[2*d+1] = cgh->cctk_ubnd[d] == cgh->cctk_gsh[d]-1; + for (int stg=0; stgcctk_lssh[CCTK_LSSH_IDX(stg,d)] = cgh->cctk_lsh[d]; + } } } for (int group=0; groupboxes.size()); + assert (component < (int)arrdata[group].dd->boxes[reflevel].size()); + assert (mglevel < (int)arrdata[group].dd->boxes[reflevel][component].size()); + const bbox& ext + = arrdata[group].dd->boxes[reflevel][component][mglevel].exterior; for (int d=0; d& ext - = arrdata[group].dd->boxes[reflevel][component][mglevel].exterior; ((int*)arrdata[group].info.lsh)[d] = (ext.shape() / ext.stride())[d]; ((int*)arrdata[group].info.lbnd)[d] @@ -336,21 +334,11 @@ namespace Carpet { assert (arrdata[group].info.lbnd[d]>=0 && arrdata[group].info.ubnd[d]& base = arrdata[group].hh->baseextent; + // Do allow outer boundaries on the finer grids ((int*)arrdata[group].info.bbox)[2*d ] - = (ext.lower() < base.lower())[d]; + = arrdata[group].info.lbnd[d] == 0; ((int*)arrdata[group].info.bbox)[2*d+1] - = (ext.upper() > base.upper())[d]; -#endif + = arrdata[group].info.ubnd[d] == arrdata[group].info.gsh[d]-1; } } @@ -359,7 +347,7 @@ namespace Carpet { const int group = CCTK_GroupIndexFromVarI(n); assert (group>=0); - const int var = n - CCTK_FirstVarIndexI(group); + const int var = n - CCTK_FirstVarIndexI(group); assert (var>=0); const int num_tl = CCTK_NumTimeLevelsFromVarI(n); assert (num_tl>0); diff --git a/Carpet/Carpet/src/variables.cc b/Carpet/Carpet/src/variables.cc index 91ee10bb6..b1d4e5a0d 100644 --- a/Carpet/Carpet/src/variables.cc +++ b/Carpet/Carpet/src/variables.cc @@ -5,7 +5,7 @@ #include "carpet.hh" -static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/variables.cc,v 1.2 2001/07/09 09:00:15 schnetter Exp $"; +static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/Carpet/src/variables.cc,v 1.3 2001/12/17 13:34:02 schnetter Exp $"; @@ -21,6 +21,9 @@ namespace Carpet { // Maximum number of refinement levels int maxreflevels; + // Refinement factor + int reffact; + // Refinement factor on finest grid int maxreflevelfact; diff --git a/Carpet/CarpetIOASCII/param.ccl b/Carpet/CarpetIOASCII/param.ccl index 2bf869427..f30df8d85 100644 --- a/Carpet/CarpetIOASCII/param.ccl +++ b/Carpet/CarpetIOASCII/param.ccl @@ -1,5 +1,5 @@ # Parameter definitions for thorn CarpetIOASCII -# $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetIOASCII/param.ccl,v 1.3 2001/12/14 17:59:54 schnetter Exp $ +# $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetIOASCII/param.ccl,v 1.4 2001/12/17 13:34:02 schnetter Exp $ @@ -45,6 +45,10 @@ BOOLEAN separate_grids "Separate grid levels in the output file by additional em { } "yes" +BOOLEAN separate_components "Separate grid components in the output file by additional empty lines" +{ +} "no" + CCTK_STRING outdir1D "Name of 1D ASCII output directory, overrides outdir" STEERABLE = ALWAYS diff --git a/Carpet/CarpetIOASCII/src/ioascii.cc b/Carpet/CarpetIOASCII/src/ioascii.cc index ce063b4f2..6db6e7f05 100644 --- a/Carpet/CarpetIOASCII/src/ioascii.cc +++ b/Carpet/CarpetIOASCII/src/ioascii.cc @@ -25,7 +25,7 @@ #include "ioascii.hh" -static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetIOASCII/src/ioascii.cc,v 1.24 2001/12/14 17:59:54 schnetter Exp $"; +static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetIOASCII/src/ioascii.cc,v 1.25 2001/12/17 13:34:03 schnetter Exp $"; @@ -337,6 +337,15 @@ int CarpetIOASCII data->write_ascii (file, cgh->cctk_iteration, offset1, dirs, tl, reflevel, component, mglevel); + // Append EOL after every component + if (CCTK_MyProc(cgh)==0) { + if (separate_components) { + assert (file.good()); + file << endl; + } + } + assert (file.good()); + } END_COMPONENT_LOOP(cgh); // Append EOL after every complete set of components -- cgit v1.2.3