]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New functions restore(uint step), clear_flags, n_refinement_steps. This allows to...
authorRalf Hartmann <Ralf.Hartmann@dlr.de>
Tue, 20 Mar 2001 10:21:02 +0000 (10:21 +0000)
committerRalf Hartmann <Ralf.Hartmann@dlr.de>
Tue, 20 Mar 2001 10:21:02 +0000 (10:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@4247 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/persistent_tria.h
deal.II/deal.II/source/grid/persistent_tria.cc

index c0cd5d3f33c7233a92f2940f4dd2b6190e82c887..ba4fc4f2f29615f30b0c168b4d7b0b0f8515f219 100644 (file)
@@ -140,19 +140,40 @@ class PersistentTriangulation : public Triangulation<dim>
     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
@@ -196,6 +217,12 @@ class PersistentTriangulation : public Triangulation<dim>
                                      */
     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)
@@ -211,6 +238,10 @@ class PersistentTriangulation : public Triangulation<dim>
                                      * Exception.
                                      */
     DeclException0 (ExcTriaNotEmpty);
+                                    /**
+                                     * Exception.
+                                     */
+    DeclException0 (ExcFlagsNotCleared);
     
   private:
                                     /**
index 409f79def2de7a4bcc52de9f5a089eacdf1c7012..7e3dcc2ddcd1137843ce0d122ce98ee7b8a5be8c 100644 (file)
@@ -69,26 +69,49 @@ PersistentTriangulation<dim>::execute_coarsening_and_refinement ()
 
 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) 
@@ -139,10 +162,9 @@ void
 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,
@@ -168,6 +190,14 @@ PersistentTriangulation<dim>::read_flags(std::istream &in)
 }
 
 
+template <int dim>
+void
+PersistentTriangulation<dim>::clear_flags()
+{
+  refine_flags.clear();
+  coarsen_flags.clear();
+}
+
 
 template <int dim>
 unsigned int

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.