#include <base/exceptions.h>
#include <base/subscriptor.h>
#include <basic/forward-declarations.h>
+#include <lac/forward-declarations.h>
#include <vector>
/**
* Specialization of #TimeStepBase# which addresses some aspects of grid handling.
- * In particular, this class is though to make handling of grids available that
+ * In particular, this class is thought to make handling of grids available that
* are adaptively refined on each time step separately or with a loose coupling
* between time steps. It also takes care of deleting and rebuilding grids when
* memory resources are a point, through the #sleep# and #wake_up# functions
* declared in the base class.
*
- * @author Wolfgang Bangerth, 1999
+ * In ddition to that, it offers a function which do some rather hairy refinement
+ * rules for time dependent problems. trying to avoid to much change in the grids
+ * between subsequent time levels, while also trying to retain the freedom of
+ * refining each grid separately. There are lots of flags and numbers controlling
+ * this function, which might drastically change the behaviour of the function -- see
+ * the description of the flags for further information.
+ *
+ * @author Wolfgang Bangerth, 1999; large parts taken from the wave program, by Wolfgang Bangerth 1998
*/
template <int dim>
class TimeStepBase_Tria : public TimeStepBase
public:
// forward declaration
struct Flags;
-
+ struct RefinementFlags;
+ struct RefinementData;
+
/**
* Default constructor. Does nothing but
* throws an exception. We need to have
* the lifetime of the coarse grid
* is longer than it is needed by
* this object.
+ *
+ * You need to give the general flags
+ * structure to this function since it
+ * is needed anyway; the refinement
+ * flags can be omitted if you do
+ * not intend to call the refinement
+ * function of this class.
*/
TimeStepBase_Tria (const double time,
const Triangulation<dim> &coarse_grid,
- const Flags &flags);
+ const Flags &flags,
+ const RefinementFlags &refinement_flags = RefinementFlags());
/**
* Destructor. At present, this does
*/
virtual void sleep (const unsigned);
+ /**
+ * Do the refinement according to the
+ * flags passed to the constructor of this
+ * object and the data passed to this
+ * function. For a description of the
+ * working of this function refer to the
+ * general documentation of this class.
+ *
+ * In fact, this function does not
+ * actually refine or coarsen the
+ * triangulation, but only sets the
+ * respective flags. This is done because
+ * usually you will not need the grid
+ * immediately afterwards but only
+ * in the next sweep, so it suffices
+ * to store the flags and rebuild it
+ * the next time we need it. Also, it
+ * may be that the next time step
+ * would like to add or delete some
+ * flags, so we have to wait anyway
+ * with the use of this grid.
+ *
+ * The function is made virtual so
+ * that you can overload it if needed.
+ */
+ virtual void refine_grid (const RefinementData &data);
+
protected:
/**
*/
const Flags flags;
+ /**
+ * Flags controlling the refinement
+ * process; see the documentation of the
+ * respective structure for more
+ * information.
+ */
+ const RefinementFlags refinement_flags;
+
private:
/**
* Vectors holding the refinement and
* rebuilt using the saved flags.
*/
void restore_grid ();
+ void mirror_refinement_flags (typename Triangulation<dim>::cell_iterator &new_cell,
+ typename Triangulation<dim>::cell_iterator &old_cell);
+
+ void adapt_grids (Triangulation<dim> &tria1, Triangulation<dim> &tria2);
};
+
template <int dim>
struct TimeStepBase_Tria<dim>::Flags
{
* fields for a description of the
* meaning of the parameters.
*/
- Flags (const unsigned int max_refinement_level,
- const bool delete_and_rebuild_tria,
+ Flags (const bool delete_and_rebuild_tria,
const unsigned int wakeup_level_to_build_grid,
const unsigned int sleep_level_to_delete_grid);
- /**
- * Maximum level of a cell in the
- * triangulation of a time level. If it
- * is set to zero, then no limit is imposed
- * on the number of refinements a coarse
- * grid cell may undergo. Usually, this
- * field is used, if for some reason you
- * want to limit refinement in an
- * adaptive process, for example to avoid
- * overly large numbers of cells or to
- * compare with grids which have a certain
- * number of refinements.
- */
- const unsigned int max_refinement_level;
-
/**
* This flag determines whether the
* #sleep# and #wake_up# functions shall
+/**
+ * Terminology:
+ * Correction: change the number of cells on this grid according to a criterion
+ * that the number of cells may be only a certain fraction more or less
+ * then the number of cells on the previous grid.
+ * Adaption: flag some cells such that there are no to grave differences.
+ */
+template <int dim>
+struct TimeStepBase_Tria<dim>::RefinementFlags
+{
+ /**
+ * Constructor. The default values are
+ * chosen such that almost no restriction
+ * on the mesh refinement is imposed.
+ */
+ RefinementFlags (const unsigned int max_refinement_level = 0,
+ const double cell_number_corridor_top = (1<<dim),
+ const double cell_number_corridor_bottom = 1,
+ const unsigned int cell_number_correction_steps = 0);
+
+ /**
+ * Maximum level of a cell in the
+ * triangulation of a time level. If it
+ * is set to zero, then no limit is imposed
+ * on the number of refinements a coarse
+ * grid cell may undergo. Usually, this
+ * field is used, if for some reason you
+ * want to limit refinement in an
+ * adaptive process, for example to avoid
+ * overly large numbers of cells or to
+ * compare with grids which have a certain
+ * number of refinements.
+ */
+ const unsigned int max_refinement_level;
+
+ /**
+ * First sweep to perform cell number
+ * correction steps on; for sweeps
+ * before, cells are only flagged and
+ * no number-correction to previous grids
+ * is performed.
+ */
+ const unsigned int first_sweep_with_correction;
+
+
+ /**
+ * Apply cell number correction with the
+ * previous time level only if there are
+ * more than this number of cells.
+ */
+ const unsigned int min_cells_for_correction;
+
+ /**
+ * Fraction by which the number of cells
+ * on a time level may differ from the
+ * number on the previous time level
+ * (first: top deviation, second: bottom
+ * deviation).
+ */
+ const double cell_number_corridor_top, cell_number_corridor_bottom;
+
+ const vector<vector<pair<unsigned int,double> > > correction_relaxations;
+
+ /**
+ * Number of iterations to be performed
+ * to adjust the number of cells on a
+ * time level to those on the previous
+ * one.
+ */
+ const unsigned int cell_number_correction_steps;
+
+ /**
+ * Flag all cells which are flagged on this
+ * timestep for refinement on the previous
+ * one also. This is useful in case the
+ * error indicator was computed by
+ * integration over time-space cells, but
+ * are now associated to a grid on a
+ * discrete time level. Since the error
+ * contribution comes from both grids,
+ * however, it is appropriate to refine
+ * both grids.
+ *
+ * Since the previous grid does not mirror
+ * the flags to the one before it, this
+ * does not lead to an almost infinite
+ * growth of cell numbers. You should
+ * use this flag with cell number
+ * correction switched on only, however.
+ *
+ * Mirroring is done after cell number
+ * correction is done, but before grid
+ * adaption, so the cell number on
+ * this grid is not noticably influenced
+ * by the cells flagged additionally on
+ * the previous grid.
+ */
+ const bool mirror_flags_to_previous_grid;
+
+ const bool adapt_grids;
+
+ /**
+ * Exception
+ */
+ DeclException1 (ExcInvalidValue,
+ int,
+ << "The following value does not fulfil the requirements: " << arg1);
+};
+
+
+
+template <int dim>
+struct TimeStepBase_Tria<dim>::RefinementData
+{
+ /**
+ * Constructor
+ */
+ RefinementData (const Vector<float> &criteria,
+ const double refinement_threshold,
+ const double coarsening_threshold=0);
+
+ /**
+ * Vector of values indicating the
+ * error per cell or some other
+ * quantity used to determine which cells
+ * are to be refined and which are to be
+ * coarsened.
+ *
+ * The vector contains exactly one value
+ * per active cell, in the usual order
+ * in which active cells are sorted in
+ * loops over cells. It is assumed that
+ * all values in this vector have a
+ * nonnegative value.
+ */
+ const Vector<float> &criteria;
+
+ /**
+ * Threshold for refinement: cells having
+ * a larger value will be refined (at least
+ * in the first round; subsequent steps
+ * of the refinement process may flag
+ * other cells as well or remove the
+ * flag from cells with a criterion higher
+ * than this threshold).
+ */
+ const double refinement_threshold;
+
+ /**
+ * Same threshold for coarsening: cells
+ * with a smaller threshold will be
+ * coarsened if possible.
+ */
+ const double coarsening_threshold;
+
+ /**
+ * Exception
+ */
+ DeclException1 (ExcInvalidValue,
+ int,
+ << "The following value does not fulfil the requirements: " << arg1);
+};
+
+
/*---------------------------- time-dependent.h ---------------------------*/
/* end of #ifndef __time_dependent_H */
#include <grid/tria.h>
#include <grid/tria_accessor.h>
#include <grid/tria_iterator.h>
+#include <lac/vector.h>
-#include <function.h>
-
+#include <functional>
+#include <algorithm>
+#include <numeric>
TimeDependent::TimeSteppingData::TimeSteppingData (const unsigned int look_ahead,
/* ------------------------------------------------------------------------- */
-template <int dim>
-TimeStepBase_Tria<dim>::Flags::Flags () {
- Assert (false, ExcInternalError());
-};
-
-
-
-template <int dim>
-TimeStepBase_Tria<dim>::Flags::Flags (const unsigned int max_refinement_level,
- const bool delete_and_rebuild_tria,
- const unsigned int wakeup_level_to_build_grid,
- const unsigned int sleep_level_to_delete_grid):
- max_refinement_level (max_refinement_level),
- delete_and_rebuild_tria (delete_and_rebuild_tria),
- wakeup_level_to_build_grid (wakeup_level_to_build_grid),
- sleep_level_to_delete_grid (sleep_level_to_delete_grid)
-{
-// Assert (!delete_and_rebuild_tria || (wakeup_level_to_build_grid>=1),
-// ExcInvalidParameter(wakeup_level_to_build_grid));
-// Assert (!delete_and_rebuild_tria || (sleep_level_to_delete_grid>=1),
-// ExcInvalidParameter(sleep_level_to_delete_grid));
-};
-
template <int dim>
TimeStepBase_Tria<dim>::TimeStepBase_Tria (const double time,
const Triangulation<dim> &coarse_grid,
- const Flags &flags) :
+ const Flags &flags,
+ const RefinementFlags &refinement_flags) :
TimeStepBase (time),
tria(0),
coarse_grid (coarse_grid),
- flags (flags)
+ flags (flags),
+ refinement_flags (refinement_flags)
{
coarse_grid.subscribe();
};
// limit refinement depth if the user
// desired so
- if (flags.max_refinement_level != 0)
- {
- Triangulation<dim>::active_cell_iterator cell, endc;
- for (cell = tria->begin_active(),
- endc = tria->end();
- cell!=endc; ++cell)
- if (static_cast<unsigned int>(cell->level()) >=
- flags.max_refinement_level)
- cell->clear_refine_flag();
- };
+// if (flags.max_refinement_level != 0)
+// {
+// Triangulation<dim>::active_cell_iterator cell, endc;
+// for (cell = tria->begin_active(),
+// endc = tria->end();
+// cell!=endc; ++cell)
+// if (static_cast<unsigned int>(cell->level()) >=
+// flags.max_refinement_level)
+// cell->clear_refine_flag();
+// };
tria->execute_coarsening_and_refinement ();
};
};
+
+
+template <int dim>
+void TimeStepBase_Tria<dim>::refine_grid (const RefinementData &refinement_data)
+{
+ // copy the following two values since
+ // we may need modified values in the
+ // process of this function
+ double refinement_threshold = refinement_data.refinement_threshold,
+ coarsening_threshold = refinement_data.coarsening_threshold;
+
+ // prepare an array where the partial
+ // sums of the criteria are stored
+ // we need this if cell number correction
+ // is switched on
+ Vector<float> partial_sums;
+ // two pointers into this array denoting
+ // the position where the two thresholds
+ // are assumed (sorry for the stupid
+ // naming of these variables, they were
+ // chosen when the whole function had about
+ // 30 lines of code and I am now too lazy
+ // to change it...)
+ Vector<float>::const_iterator p=0, q=0;
+
+ if ((timestep_no != 0) &&
+ (sweep_no>=refinement_flags.first_sweep_with_correction) &&
+ (refinement_flags.cell_number_correction_steps > 0))
+ {
+ partial_sums = refinement_data.criteria;
+ sort (partial_sums.begin(), partial_sums.end(),
+ greater<float>());
+ partial_sum (partial_sums.begin(), partial_sums.end(),
+ partial_sums.begin());
+ p = upper_bound (partial_sums.begin(), partial_sums.end(),
+ coarsening_threshold);
+ q = lower_bound (partial_sums.begin(), partial_sums.end(),
+ refinement_threshold);
+ };
+
+
+
+ // actually flag cells the first time
+ tria->refine (refinement_data.criteria, refinement_threshold);
+ tria->coarsen (refinement_data.criteria, coarsening_threshold);
+
+ // store this number for the following
+ // since its computation is rather
+ // expensive and since it doesn't change
+ const unsigned int n_active_cells = tria->n_active_cells ();
+
+ // if not on first time level: try to
+ // adjust the number of resulting
+ // cells to those on the previous
+ // time level. Only do the cell number
+ // correction for higher sweeps and if
+ // there are sufficiently many cells
+ // already to avoid "grid stall" i.e.
+ // that the grid's evolution is hindered
+ // by the correction (this usually
+ // happens if there are very few cells,
+ // since then the number of cells touched
+ // by the correction step may exceed the
+ // number of cells which are flagged for
+ // refinement; in this case it often
+ // happens that the number of cells
+ // does not grow between sweeps, which
+ // clearly is not the wanted behaviour)
+ //
+ // however, if we do not do anything, we
+ // can get into trouble somewhen later.
+ // therefore, we also use the correction
+ // step for the first sweep or if the
+ // number of cells is between 100 and 300
+ // (unlike in the first version of the
+ // algorithm), but relax the conditions
+ // for the correction to allow deviations
+ // which are three times as high than
+ // allowed (sweep==1 || cell number<200)
+ // or twice as high (sweep==2 ||
+ // cell number<300). Also, since
+ // refinement never does any harm other
+ // than increased work, we allow for
+ // arbitrary growth of cell number if
+ // the estimated cell number is below
+ // 200.
+ //
+ // repeat this loop several times since
+ // the first estimate may not be totally
+ // correct
+ if ((timestep_no != 0) && (sweep_no>=refinement_flags.first_sweep_with_correction))
+ for (unsigned int loop=0;
+ loop<refinement_flags.cell_number_correction_steps; ++loop)
+ {
+ Triangulation<dim> *previous_tria
+ = dynamic_cast<const TimeStepBase_Tria<dim>*>(previous_timestep)->tria;
+
+ // do one adaption step if desired
+ // (there are more coming below then
+ // also)
+ if (refinement_flags.adapt_grids)
+ adapt_grids (*previous_tria, *tria);
+
+ // perform flagging of cells
+ // needed to regularize the
+ // triangulation
+ tria->prepare_coarsening ();
+ tria->prepare_refinement ();
+ previous_tria->prepare_coarsening ();
+ previous_tria->prepare_refinement ();
+
+
+
+ // now count the number of elements
+ // which will result on the previous
+ // grid after it will be refined. The
+ // number which will really result
+ // should be approximately that that we
+ // compute here, since we already
+ // performed most of the prepare*
+ // steps for the previous grid
+ //
+ // use a double value since for each
+ // four cells (in 2D) that we flagged
+ // for coarsening we result in one
+ // new. but since we loop over flagged
+ // cells, we have to subtract 3/4 of
+ // a cell for each flagged cell
+ double previous_cells = previous_tria->n_active_cells();
+ Triangulation<dim>::active_cell_iterator cell, endc;
+ cell = previous_tria->begin_active();
+ endc = previous_tria->end();
+ for (; cell!=endc; ++cell)
+ if (cell->refine_flag_set())
+ previous_cells += (GeometryInfo<dim>::children_per_cell-1);
+ else
+ if (cell->coarsen_flag_set())
+ previous_cells -= (GeometryInfo<dim>::children_per_cell-1) /
+ GeometryInfo<dim>::children_per_cell;
+
+ // #previous_cells# now gives the
+ // number of cells which would result
+ // from the flags on the previous grid
+ // if we refined it now. However, some
+ // more flags will be set when we adapt
+ // the previous grid with this one
+ // after the flags have been set for
+ // this time level; on the other hand,
+ // we don't account for this, since the
+ // number of cells on this time level
+ // will be changed afterwards by the
+ // same way, when it is adapted to the
+ // next time level
+
+ // now estimate the number of cells which
+ // will result on this level
+ double estimated_cells = n_active_cells;
+ cell = tria->begin_active();
+ endc = tria->end();
+ for (; cell!=endc; ++cell)
+ if (cell->refine_flag_set())
+ estimated_cells += (GeometryInfo<dim>::children_per_cell-1);
+ else
+ if (cell->coarsen_flag_set())
+ estimated_cells -= (GeometryInfo<dim>::children_per_cell-1) /
+ GeometryInfo<dim>::children_per_cell;
+
+ // calculate the allowed delta in
+ // cell numbers; be more lenient
+ // if there are few cells
+ double delta_up = refinement_flags.cell_number_corridor_top,
+ delta_down = refinement_flags.cell_number_corridor_bottom;
+
+ const vector<pair<unsigned int,double> > &relaxations
+ = (sweep_no > refinement_flags.correction_relaxations.size() ?
+ refinement_flags.correction_relaxations.back() :
+ refinement_flags.correction_relaxations[sweep_no]);
+ for (unsigned int r=0; r!=relaxations.size(); ++r)
+ if (n_active_cells < relaxations[r].first)
+ {
+ delta_up *= relaxations[r].second;
+ delta_down *= relaxations[r].second;
+ break;
+ };
+
+ // now, if the number of estimated
+ // cells exceeds the number of cells
+ // on the old time level by more than
+ // delta: cut the top threshold
+ //
+ // note that for each cell that
+ // we unflag we have to diminish the
+ // estimated number of cells by
+ // #children_per_cell#.
+ if (estimated_cells > previous_cells*(1.+delta_up))
+ {
+ // only limit the cell number
+ // if there will not be less
+ // than some number of cells
+ //
+ // also note that when using the
+ // dual estimator, the initial
+ // time level is not refined
+ // on its own, so we may not
+ // limit the number of the second
+ // time level on the basis of
+ // the initial one; since for
+ // the dual estimator, we
+ // mirror the refinement
+ // flags, the initial level
+ // will be passively refined
+ // later on.
+ if (estimated_cells>refinement_flags.min_cells_for_correction)
+ {
+ // number of cells by which the
+ // new grid is to be diminished
+ double delta_cells = estimated_cells -
+ previous_cells*(1.+delta_up);
+
+ for (unsigned int i=0; i<delta_cells;
+ i += (GeometryInfo<dim>::children_per_cell-1))
+ if (q!=partial_sums.begin())
+ --q;
+ else
+ break;
+ }
+ else
+ // too many cells, but we
+ // won't do anything about
+ // that
+ break;
+ }
+ else
+ // likewise: if the estimated number
+ // of cells is less than 90 per cent
+ // of those at the previous time level:
+ // raise threshold by refining
+ // additional cells. if we start to
+ // run into the area of cells
+ // which are to be coarsened, we
+ // raise the limit for these too
+ if (estimated_cells < previous_cells*(1.-delta_down))
+ {
+ // number of cells by which the
+ // new grid is to be enlarged
+ double delta_cells = previous_cells*(1.-delta_down)-estimated_cells;
+ // heuristics: usually, if we
+ // add #delta_cells# to the
+ // present state, we end up
+ // with much more than only
+ // (1-delta_down)*prev_cells
+ // because of the effect of
+ // regularization and because
+ // of adaption to the
+ // following grid. Therefore,
+ // if we are not in the last
+ // correction loop, we try not
+ // to add as many cells as seem
+ // necessary at first and hope
+ // to get closer to the limit
+ // this way. Only in the last
+ // loop do we have to take the
+ // full number to guarantee the
+ // wanted result.
+ //
+ // The value 0.9 is taken from
+ // practice, as the additional
+ // number of cells introduced
+ // by regularization is
+ // approximately 10 per cent
+ // of the flagged cells.
+ if (loop != refinement_flags.cell_number_correction_steps-1)
+ delta_cells *= 0.9;
+
+
+ for (unsigned int i=0; i<delta_cells;
+ i += (GeometryInfo<dim>::children_per_cell-1))
+ if (q!=p)
+ ++q;
+ else
+ if (p!=partial_sums.end())
+ ++p, ++q;
+ else
+ break;
+ }
+ else
+ // estimated cell number is ok,
+ // stop correction steps
+ break;
+
+
+ if (p==partial_sums.end())
+ coarsening_threshold = 0;
+ else
+ coarsening_threshold = *p - *(p-1);
+
+ if (q==partial_sums.begin())
+ refinement_threshold = *q;
+ else
+ refinement_threshold = *q - *(q-1);
+ if (coarsening_threshold>=refinement_threshold)
+ coarsening_threshold = 0.999*refinement_threshold;
+
+ // now that we have re-adjusted
+ // thresholds: clear all refine and
+ // coarsening flags and do it all
+ // over again
+ cell = tria->begin_active();
+ endc = tria->end();
+ for (; cell!=endc; ++cell)
+ {
+ cell->clear_refine_flag ();
+ cell->clear_coarsen_flag ();
+ };
+
+
+ // flag cells finally
+ tria->refine (refinement_data.criteria, refinement_threshold);
+ tria->coarsen (refinement_data.criteria, coarsening_threshold);
+ };
+
+
+
+ // if step number is greater than
+ // one: adapt this and the previous
+ // grid to each other. Don't do so
+ // for the initial grid because
+ // it is always taken to be the first
+ // grid and needs therefore no
+ // treatment of its own.
+ if ((timestep_no >= 1) && (refinement_flags.adapt_grids))
+ {
+ Triangulation<dim> *previous_tria
+ = dynamic_cast<const TimeStepBase_Tria<dim>*>(previous_timestep)->tria;
+
+
+ // if we used the dual estimator, we
+ // computed the error information on
+ // a time slab, rather than on a level
+ // of its own. we then mirror the
+ // refinement flags we determined for
+ // the present level to the previous
+ // one
+ //
+ // do this mirroring only, if cell number
+ // adjustment is on, since otherwise
+ // strange things may happen
+ if (refinement_flags.mirror_flags_to_previous_grid)
+ {
+ adapt_grids (*previous_tria, *tria);
+
+ Triangulation<dim>::cell_iterator old_cell, new_cell, endc;
+ old_cell = previous_tria->begin(0);
+ new_cell = tria->begin(0);
+ endc = tria->end(0);
+ for (; new_cell!=endc; ++new_cell, ++old_cell)
+ mirror_refinement_flags (new_cell, old_cell);
+ };
+
+ tria->prepare_coarsening ();
+ tria->prepare_refinement ();
+ previous_tria->prepare_coarsening ();
+ previous_tria->prepare_refinement ();
+
+ // adapt present and previous grids
+ // to each other: flag additional
+ // cells to avoid the previous grid
+ // to have cells refined twice more
+ // than the present one and vica versa.
+ adapt_grids (*tria, *tria);
+ tria->prepare_coarsening ();
+ tria->prepare_refinement ();
+ previous_tria->prepare_coarsening ();
+ previous_tria->prepare_refinement ();
+ };
+};
+
+
+
+
+
+
+template <int dim>
+TimeStepBase_Tria<dim>::Flags::Flags () {
+ Assert (false, ExcInternalError());
+};
+
+
+
+template <int dim>
+TimeStepBase_Tria<dim>::Flags::Flags (const bool delete_and_rebuild_tria,
+ const unsigned int wakeup_level_to_build_grid,
+ const unsigned int sleep_level_to_delete_grid):
+ delete_and_rebuild_tria (delete_and_rebuild_tria),
+ wakeup_level_to_build_grid (wakeup_level_to_build_grid),
+ sleep_level_to_delete_grid (sleep_level_to_delete_grid)
+{
+// Assert (!delete_and_rebuild_tria || (wakeup_level_to_build_grid>=1),
+// ExcInvalidParameter(wakeup_level_to_build_grid));
+// Assert (!delete_and_rebuild_tria || (sleep_level_to_delete_grid>=1),
+// ExcInvalidParameter(sleep_level_to_delete_grid));
+};
+
+
+
+template <int dim>
+TimeStepBase_Tria<dim>::RefinementFlags::
+RefinementFlags (const unsigned int max_refinement_level,
+ const double cell_number_corridor_top,
+ const double cell_number_corridor_bottom,
+ const unsigned int cell_number_correction_steps) :
+ max_refinement_level(max_refinement_level),
+ cell_number_corridor_top(cell_number_corridor_top),
+ cell_number_corridor_bottom(cell_number_corridor_bottom),
+ cell_number_correction_steps(cell_number_correction_steps)
+{
+ Assert (cell_number_corridor_top>=0, ExcInvalidValue (cell_number_corridor_top));
+ Assert (cell_number_corridor_bottom>=0, ExcInvalidValue (cell_number_corridor_bottom));
+ Assert (cell_number_corridor_bottom<=1, ExcInvalidValue (cell_number_corridor_bottom));
+};
+
+
+
+template <int dim>
+TimeStepBase_Tria<dim>::RefinementData::
+RefinementData (const Vector<float> &criteria,
+ const double _refinement_threshold,
+ const double _coarsening_threshold) :
+ criteria (criteria),
+ refinement_threshold(_refinement_threshold),
+ // in some rare cases it may happen that
+ // both thresholds are the same (e.g. if
+ // there are many cells with the same
+ // error indicator). That would mean that
+ // all cells will be flagged for
+ // refinement or coarsening, but some will
+ // be flagged for both, namely those for
+ // which the indicator equals the
+ // thresholds. This is forbidden, however.
+ //
+ // In some rare cases with very few cells
+ // we also could get integer round off
+ // errors and get problems with
+ // the top and bottom fractions.
+ //
+ // In these case we arbitrarily reduce the
+ // bottom threshold by one permille below
+ // the top threshold
+ coarsening_threshold((_coarsening_threshold == _refinement_threshold ?
+ _coarsening_threshold :
+ 0.999*_coarsening_threshold))
+{
+ Assert (*min_element(criteria.begin(), criteria.end()) >= 0,
+ ExcInvalidValue(*min_element(criteria.begin(), criteria.end())));
+
+ Assert (refinement_threshold >= 0, ExcInvalidValue(refinement_threshold));
+ Assert (coarsening_threshold >= 0, ExcInvalidValue(coarsening_threshold));
+ Assert (coarsening_threshold < refinement_threshold,
+ ExcInvalidValue (coarsening_threshold));
+};
+
+
+
+
+
// explicit instantiations
template class TimeStepBase_Tria<deal_II_dimension>;