virtual void execute_coarsening_and_refinement ();
/**
- * Restore the grid according to the saved
- * data. For this, the coarse grid is
- * copied and the grid is stepwise
- * rebuilt using the saved flags.
+ * Restore the grid according to
+ * the saved data. For this, the
+ * coarse grid is copied and the
+ * grid is stepwise rebuilt using
+ * the saved flags.
*
- * Note that this function will result in
- * an error if the underlying triangulation
- * is not empty, i.e. it will only succeed
- * if this object is newly created or
- * @p{clear()} was called on it before.
+ * Note that this function will
+ * result in an error if the
+ * underlying triangulation is
+ * not empty, i.e. it will only
+ * succeed if this object is
+ * newly created or @p{clear()}
+ * was called on it before.
*/
void restore ();
+ /**
+ * Differential restore. Performs
+ * the @p{step}th local
+ * refinement and coarsening step.
+ * Step 0 stands for the copying
+ * of the coarse grid.
+ */
+ void restore (const unsigned int step);
+
+ /**
+ * Returns the number of
+ * refinement and coarsening
+ * steps. This is given by the
+ * size of the @p{refine_flags}
+ * vector.
+ */
+ unsigned int n_refinement_steps () const;
+
/**
* Overload this function to use
* @p{tria} as a new coarse grid. The
*/
virtual void read_flags(std::istream &in);
+ /**
+ * Clears all flags. Retains the
+ * same coarse grid.
+ */
+ virtual void clear_flags();
+
/**
* Determine an estimate for the
* memory consumption (in bytes)
* Exception.
*/
DeclException0 (ExcTriaNotEmpty);
+ /**
+ * Exception.
+ */
+ DeclException0 (ExcFlagsNotCleared);
private:
/**
template <int dim>
void
-PersistentTriangulation<dim>::restore () {
- // copy the old triangulation.
- // this will yield an error if the
- // underlying triangulation was not
- // empty
- Triangulation<dim>::copy_triangulation (*coarse_grid);
-
- // for each of the previous refinement
- // sweeps
- for (unsigned int i=0; i<refine_flags.size(); ++i)
+PersistentTriangulation<dim>::restore ()
+{
+ // for each of the previous
+ // refinement sweeps
+ for (unsigned int i=0; i<refine_flags.size()+1; ++i)
+ restore(i);
+};
+
+
+template <int dim>
+void
+PersistentTriangulation<dim>::restore (const unsigned int step) {
+
+ if (step==0)
+ // copy the old triangulation.
+ // this will yield an error if
+ // the underlying triangulation
+ // was not empty
+ Triangulation<dim>::copy_triangulation (*coarse_grid);
+ else
+ // for each of the previous
+ // refinement sweeps
{
- // get flags
- load_refine_flags (refine_flags[i]);
- load_coarsen_flags (coarsen_flags[i]);
+ Assert(step<refine_flags.size()+1,
+ ExcDimensionMismatch(step, refine_flags.size()+1));
+
+ load_refine_flags (refine_flags[step-1]);
+ load_coarsen_flags (coarsen_flags[step-1]);
Triangulation<dim>::execute_coarsening_and_refinement ();
};
};
+
+template <int dim>
+unsigned int
+PersistentTriangulation<dim>::n_refinement_steps() const
+{
+ return refine_flags.size();
+}
+
+
template <int dim>
void
PersistentTriangulation<dim>::copy_triangulation (const Triangulation<dim> &old_grid)
PersistentTriangulation<dim>::read_flags(std::istream &in)
{
Assert(refine_flags.size()==0 && coarsen_flags.size()==0,
- ExcTriaNotEmpty());
-
+ ExcFlagsNotCleared());
AssertThrow (in, ExcIO());
-
+
unsigned int magic_number;
in >> magic_number;
AssertThrow(magic_number==mn_persistent_tria_flags_begin,
}
+template <int dim>
+void
+PersistentTriangulation<dim>::clear_flags()
+{
+ refine_flags.clear();
+ coarsen_flags.clear();
+}
+
template <int dim>
unsigned int