diff options
author | schnetter <> | 2004-04-18 11:29:00 +0000 |
---|---|---|
committer | schnetter <> | 2004-04-18 11:29:00 +0000 |
commit | 4513dabacb47c14dac8eb02274a3de502c537572 (patch) | |
tree | 601f8a3509c6af7fb2720bfc4a95f87726dd69b7 /Carpet/CarpetLib | |
parent | 013e3244bab2b0ec447d0526ffdafbda5553e951 (diff) |
Remove the parameters Carpet::prolongate_initial_data; this is now
Remove the parameters Carpet::prolongate_initial_data; this is now
always done.
Remove arguments initialise_from and do_prolongate from Regrid().
Regridding is now done in level mode instead of meta mode.
Furthermore, CarpetRegrid is called in singlemape mode.
darcs-hash:20040418112943-07bb3-2e392df1737ab75f3f0d553bb53bde2ed41f8773.gz
Diffstat (limited to 'Carpet/CarpetLib')
-rw-r--r-- | Carpet/CarpetLib/src/dh.cc | 8 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/dh.hh | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/gf.cc | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/ggf.cc | 122 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/ggf.hh | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/gh.cc | 12 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/gh.hh | 6 |
7 files changed, 89 insertions, 71 deletions
diff --git a/Carpet/CarpetLib/src/dh.cc b/Carpet/CarpetLib/src/dh.cc index 1cd9772e7..f7b7b88e5 100644 --- a/Carpet/CarpetLib/src/dh.cc +++ b/Carpet/CarpetLib/src/dh.cc @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/dh.cc,v 1.50 2004/04/07 16:56:42 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/dh.cc,v 1.51 2004/04/18 13:29:43 schnetter Exp $ #include <assert.h> @@ -31,7 +31,7 @@ dh<D>::dh (gh<D>& h, assert (buffer_width>=0); h.add(this); CHECKPOINT; - recompose (0, true); + recompose (); } // Destructors @@ -51,7 +51,7 @@ int dh<D>::prolongation_stencil_size () const { // Modifiers template<int D> -void dh<D>::recompose (const int initialise_from, const bool do_prolongate) { +void dh<D>::recompose () { DECLARE_CCTK_PARAMETERS; CHECKPOINT; @@ -574,7 +574,7 @@ void dh<D>::recompose (const int initialise_from, const bool do_prolongate) { for (typename list<ggf<D>*>::iterator f=gfs.begin(); f!=gfs.end(); ++f) { - (*f)->recompose (initialise_from, do_prolongate); + (*f)->recompose (); } } diff --git a/Carpet/CarpetLib/src/dh.hh b/Carpet/CarpetLib/src/dh.hh index d525a0be4..c78d31291 100644 --- a/Carpet/CarpetLib/src/dh.hh +++ b/Carpet/CarpetLib/src/dh.hh @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/dh.hh,v 1.17 2004/03/03 16:20:19 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/dh.hh,v 1.18 2004/04/18 13:29:43 schnetter Exp $ #ifndef DH_HH #define DH_HH @@ -115,7 +115,7 @@ public: int prolongation_stencil_size () const; // Modifiers - void recompose (const int initialise_from, const bool do_prolongate); + void recompose (); // Grid function management void add (ggf<D>* f); diff --git a/Carpet/CarpetLib/src/gf.cc b/Carpet/CarpetLib/src/gf.cc index 462c179b8..36b601850 100644 --- a/Carpet/CarpetLib/src/gf.cc +++ b/Carpet/CarpetLib/src/gf.cc @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/gf.cc,v 1.17 2004/03/23 19:30:14 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/gf.cc,v 1.18 2004/04/18 13:29:43 schnetter Exp $ #include <assert.h> @@ -23,7 +23,7 @@ gf<T,D>::gf (const int varindex, const operator_type transport_operator, t, d, tmin, tmax, prolongation_order_time, vectorlength, vectorindex, vectorleader) { - this->recompose (0, true); + this->recompose (); } // Destructors diff --git a/Carpet/CarpetLib/src/ggf.cc b/Carpet/CarpetLib/src/ggf.cc index 38657d367..6b78c43dd 100644 --- a/Carpet/CarpetLib/src/ggf.cc +++ b/Carpet/CarpetLib/src/ggf.cc @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/ggf.cc,v 1.36 2004/03/23 19:30:14 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/ggf.cc,v 1.37 2004/04/18 13:28:25 schnetter Exp $ #include <assert.h> #include <stdlib.h> @@ -55,7 +55,7 @@ bool ggf<D>::operator== (const ggf<D>& f) const { // Modifiers template<int D> -void ggf<D>::recompose (const int initialise_from, const bool do_prolongate) { +void ggf<D>::recompose () { // TODO: restructure storage only when needed @@ -66,7 +66,7 @@ void ggf<D>::recompose (const int initialise_from, const bool do_prolongate) { storage.resize(tmax-tmin+1); for (int tl=tmin; tl<=tmax; ++tl) { storage.at(tl-tmin).resize(h.reflevels()); - for (int rl=initialise_from; rl<h.reflevels(); ++rl) { + for (int rl=0; rl<h.reflevels(); ++rl) { storage.at(tl-tmin).at(rl).resize(h.components(rl)); for (int c=0; c<h.components(rl); ++c) { storage.at(tl-tmin).at(rl).at(c).resize(h.mglevels(rl,c)); @@ -78,53 +78,73 @@ void ggf<D>::recompose (const int initialise_from, const bool do_prolongate) { } // for tl // Initialise the new storage - for (int rl=initialise_from; rl<h.reflevels(); ++rl) { - for (int tl=tmin; tl<=tmax; ++tl) { - for (int c=0; c<h.components(rl); ++c) { - for (int ml=0; ml<h.mglevels(rl,c); ++ml) { - - storage.at(tl-tmin).at(rl).at(c).at(ml) = typed_data(tl,rl,c,ml); - - // Allocate storage - storage.at(tl-tmin).at(rl).at(c).at(ml)->allocate - (d.boxes.at(rl).at(c).at(ml).exterior, h.proc(rl,c)); - - if (do_prolongate) { - // Initialise from coarser level, if possible - // TODO: init only un-copied regions - if (rl>0) { - for (comm_state<D> state; !state.done(); state.step()) { - const CCTK_REAL time = t.time(tl,rl,ml); - ref_prolongate (state,tl,rl,c,ml,time); - } - } // if rl - } // if do_prolongate - - // Copy from old storage, if possible - // todo: copy only from interior regions? - if (rl<(int)oldstorage.at(tl-tmin).size()) { - for (comm_state<D> state; !state.done(); state.step()) { + // TODO: overlap this initialisation for all variables + for (int rl=0; rl<h.reflevels(); ++rl) { + bool firsttime=true; + for (comm_state<D> state; !state.done(); state.step(), firsttime=false) { + for (int tl=tmin; tl<=tmax; ++tl) { + for (int c=0; c<h.components(rl); ++c) { + for (int ml=0; ml<h.mglevels(rl,c); ++ml) { + + if (firsttime) { + storage.at(tl-tmin).at(rl).at(c).at(ml) = typed_data(tl,rl,c,ml); + + // Allocate storage + storage.at(tl-tmin).at(rl).at(c).at(ml)->allocate + (d.boxes.at(rl).at(c).at(ml).exterior, h.proc(rl,c)); + } + + // Find out which regions need to be prolongated + // TODO: do this once in the dh instead of for each variable here + ibset work (d.boxes.at(rl).at(c).at(ml).exterior); + + // Copy from old storage, if possible + // TODO: copy only from interior regions? + if (rl<(int)oldstorage.at(tl-tmin).size()) { for (int cc=0; cc<(int)oldstorage.at(tl-tmin).at(rl).size(); ++cc) { if (ml<(int)oldstorage.at(tl-tmin).at(rl).at(cc).size()) { - const ibbox ovlp = - (d.boxes.at(rl).at(c).at(ml).exterior - & oldstorage.at(tl-tmin).at(rl).at(cc).at(ml)->extent()); - storage.at(tl-tmin).at(rl).at(c).at(ml)->copy_from - (state, oldstorage.at(tl-tmin).at(rl).at(cc).at(ml), ovlp); + // TODO: prefer same processor, etc., see dh.cc + ibset ovlp = work & oldstorage.at(tl-tmin).at(rl).at(cc).at(ml)->extent(); + ovlp.normalize(); + work -= ovlp; + for (typename ibset::const_iterator r=ovlp.begin(); r!=ovlp.end(); ++r) { + storage.at(tl-tmin).at(rl).at(c).at(ml)->copy_from (state, oldstorage.at(tl-tmin).at(rl).at(cc).at(ml), *r); + } } // if ml } // for cc - } // for step - } // if rl - - } // for ml - } // for c - - } // for tl + } // if rl + work.normalize(); + + // Initialise from coarser level, if possible + if (rl>0) { + if (transport_operator != op_none) { + assert (tmax-tmin+1 >= prolongation_order_time+1); + vector<int> tls(prolongation_order_time+1); + for (int i=0; i<=prolongation_order_time; ++i) { + tls.at(i) = tmax - i; + } + vector<const gdata<D>*> gsrcs(tls.size()); + vector<CCTK_REAL> times(tls.size()); + for (int i=0; i<(int)gsrcs.size(); ++i) { + gsrcs.at(i) = storage.at(tls.at(i)-tmin).at(rl-1).at(c).at(ml); + times.at(i) = t.time(tls.at(i),rl-1,ml); + } + const CCTK_REAL time = t.time(tl,rl,ml); + for (typename ibset::const_iterator r=work.begin(); r!=work.end(); ++r) { + storage.at(tl-tmin).at(rl).at(c).at(ml)->interpolate_from (state, gsrcs, times, *r, time, d.prolongation_order_space, prolongation_order_time); + } + } // if transport_operator + } // if rl + + } // for ml + } // for c + } // for tl + } // for step } // for rl // Delete old storage for (int tl=tmin; tl<=tmax; ++tl) { - for (int rl=initialise_from; rl<(int)oldstorage.at(tl-tmin).size(); ++rl) { + for (int rl=0; rl<(int)oldstorage.at(tl-tmin).size(); ++rl) { for (int c=0; c<(int)oldstorage.at(tl-tmin).at(rl).size(); ++c) { for (int ml=0; ml<(int)oldstorage.at(tl-tmin).at(rl).at(c).size(); ++ml) { delete oldstorage.at(tl-tmin).at(rl).at(c).at(ml); @@ -134,21 +154,19 @@ void ggf<D>::recompose (const int initialise_from, const bool do_prolongate) { } // for tl for (int tl=tmin; tl<=tmax; ++tl) { - for (int rl=initialise_from; rl<h.reflevels(); ++rl) { + for (int rl=0; rl<h.reflevels(); ++rl) { // Set boundaries for (int c=0; c<h.components(rl); ++c) { for (int ml=0; ml<h.mglevels(rl,c); ++ml) { - if (do_prolongate) { - // TODO: assert that reflevel 0 boundaries are copied - if (rl>0) { - for (comm_state<D> state; !state.done(); state.step()) { - const CCTK_REAL time = t.time(tl,rl,ml); - ref_bnd_prolongate (state,tl,rl,c,ml,time); - } - } // if rl - } + // TODO: assert that reflevel 0 boundaries are copied + if (rl>0) { + for (comm_state<D> state; !state.done(); state.step()) { + const CCTK_REAL time = t.time(tl,rl,ml); + ref_bnd_prolongate (state,tl,rl,c,ml,time); + } + } // if rl for (comm_state<D> state; !state.done(); state.step()) { sync (state,tl,rl,c,ml); diff --git a/Carpet/CarpetLib/src/ggf.hh b/Carpet/CarpetLib/src/ggf.hh index 64998c8a2..d228a863e 100644 --- a/Carpet/CarpetLib/src/ggf.hh +++ b/Carpet/CarpetLib/src/ggf.hh @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/ggf.hh,v 1.20 2004/03/23 12:40:27 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/ggf.hh,v 1.21 2004/04/18 13:29:43 schnetter Exp $ #ifndef GGF_HH #define GGF_HH @@ -87,7 +87,7 @@ public: // Modifiers - void recompose (const int initialise_from, const bool do_prolongate); + void recompose (); // Cycle the time levels by rotating the data sets void cycle (int rl, int c, int ml); diff --git a/Carpet/CarpetLib/src/gh.cc b/Carpet/CarpetLib/src/gh.cc index 172da013d..b68be34d8 100644 --- a/Carpet/CarpetLib/src/gh.cc +++ b/Carpet/CarpetLib/src/gh.cc @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/gh.cc,v 1.26 2004/04/07 16:59:47 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/gh.cc,v 1.28 2004/04/18 13:29:43 schnetter Exp $ #include <assert.h> #include <stdlib.h> @@ -36,9 +36,7 @@ gh<D>::~gh () { } template<int D> void gh<D>::recompose (const rexts& exts, const rbnds& outer_bounds, - const rprocs& procs, - const int initialise_from, - const bool do_prolongate) + const rprocs& procs) { DECLARE_CCTK_PARAMETERS; @@ -96,7 +94,11 @@ void gh<D>::recompose (const rexts& exts, // Check base grid extent if (reflevels()>0) { for (int c=0; c<components(0); ++c) { + // TODO: put the check back in, taking outer boundaries into + // account +#if 0 assert (extents.at(0).at(c).at(0).is_contained_in(baseextent)); +#endif } } @@ -169,7 +171,7 @@ void gh<D>::recompose (const rexts& exts, } for (typename list<dh<D>*>::iterator d=dhs.begin(); d!=dhs.end(); ++d) { - (*d)->recompose (initialise_from, do_prolongate); + (*d)->recompose (); } } diff --git a/Carpet/CarpetLib/src/gh.hh b/Carpet/CarpetLib/src/gh.hh index babc85cfb..69f9f31c0 100644 --- a/Carpet/CarpetLib/src/gh.hh +++ b/Carpet/CarpetLib/src/gh.hh @@ -1,4 +1,4 @@ -// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/gh.hh,v 1.16 2004/03/23 19:30:14 schnetter Exp $ +// $Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetLib/src/gh.hh,v 1.17 2004/04/18 13:29:43 schnetter Exp $ #ifndef GH_HH #define GH_HH @@ -85,9 +85,7 @@ public: // Modifiers void recompose (const rexts& exts, const rbnds& outer_bounds, - const rprocs& procs, - const int initialise_from, - const bool do_prolongate); + const rprocs& procs); // Accessors int reflevels () const { |