]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Introduce a way to pre-compute the number of lines, quads, etc instead of always...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 5 Aug 1999 14:00:40 +0000 (14:00 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 5 Aug 1999 14:00:40 +0000 (14:00 +0000)
git-svn-id: https://svn.dealii.org/trunk@1627 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/Todo
deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/source/grid/tria.all_dimensions.cc
deal.II/deal.II/source/grid/tria.cc

index 885ba452646a7477992dec3314291732ba8cbc2a..a502b6138a564bdf267f33ad703c7da307046905 100644 (file)
@@ -212,8 +212,6 @@ Remove the hack in Make_forward_declarations when using gcc2.95. See
   the cvs diffs of 1.6 to 1.7 for the reason.
 
 
-Introduce a number cache, where functions like n_active_quads get
-  their data from. This is worth the work, as can be seen from
-  profiling data, since it prevents us from cycling through all
-  elements each time and above that allows to write
-     for (; i<n_active_cells(); ++i)
\ No newline at end of file
+Move the TriaNumberCache into the Triangulation class once the
+  compiler supports this (gcc2.95 does not allow this at present)
+  
index 43767785a93bb45b3d0eb91d74c79e526b4d517f..5da8eb673aae1447eccef77aec1f9c0b4713d580 100644 (file)
@@ -290,6 +290,187 @@ class TriaDimensionInfo<3> {
 
 
 
+/*------------------------------------------------------------------------*/
+
+/**
+ * Cache class used to store the number of used and active elements
+ * (lines or quads etc) within the levels of a triangulation. This
+ * is only the declaration of the template, concrete instantiations
+ * are below.
+ *
+ * In the old days, whenever one wanted to access one of these
+ * numbers, one had to perform a loop over all lines, e.g., and count
+ * the elements until we hit the end iterator. This is time consuming
+ * and since access to the number of lines etc is a rather frequent
+ * operation, this was not an optimal solution.
+ *
+ * Note: these class should be made local to the triangulation class
+ * once the compiler supports that (gcc2.95 does not at present).
+ *
+ * @author Wolfgang Bangerth, 1999 */
+template <int dim>
+struct TriaNumberCache;
+
+/**
+ * Cache class used to store the number of used and active elements
+ * (lines or quads etc) within the levels of a triangulation. This
+ * specialization stores the numbers of lines.
+ *
+ * In the old days, whenever one wanted to access one of these
+ * numbers, one had to perform a loop over all lines, e.g., and count
+ * the elements until we hit the end iterator. This is time consuming
+ * and since access to the number of lines etc is a rather frequent
+ * operation, this was not an optimal solution.
+ *
+ * Note: these class should be made local to the triangulation class
+ * once the compiler supports that (gcc2.95 does not at present).
+ *
+ * @author Wolfgang Bangerth, 1999
+ */
+template <>
+struct TriaNumberCache<1> 
+{
+                                    /**
+                                     * Number of used lines in the whole
+                                     * triangulation.
+                                     */
+    unsigned int n_lines;
+
+                                    /**
+                                     * Array holding the number of used
+                                     * lines on each level.
+                                     */
+    vector<unsigned int> n_lines_level;
+    
+                                    /**
+                                     * Number of active lines in the
+                                     * whole triangulation.
+                                     */
+    unsigned int n_active_lines;
+
+                                    /**
+                                     * Array holding the number of active
+                                     * lines on each level.
+                                     */
+    vector<unsigned int> n_active_lines_level;
+
+                                    /**
+                                     * Constructor. Set values to zero
+                                     * by default.
+                                     */
+    TriaNumberCache ();
+};
+
+
+
+/**
+ * Cache class used to store the number of used and active elements
+ * (lines or quads etc) within the levels of a triangulation. This
+ * specialization stores the numbers of quads. Due to the inheritance
+ * from the base class #TriaNumberCache<1>#, the numbers of lines
+ * are also within this class.
+ *
+ * In the old days, whenever one wanted to access one of these
+ * numbers, one had to perform a loop over all lines, e.g., and count
+ * the elements until we hit the end iterator. This is time consuming
+ * and since access to the number of lines etc is a rather frequent
+ * operation, this was not an optimal solution.
+ *
+ * Note: these class should be made local to the triangulation class
+ * once the compiler supports that (gcc2.95 does not at present).
+ *
+ * @author Wolfgang Bangerth, 1999
+ */
+template <>
+struct TriaNumberCache<2> : public TriaNumberCache<1>
+{
+                                    /**
+                                     * Number of used quads in the whole
+                                     * triangulation.
+                                     */
+    unsigned int n_quads;
+
+                                    /**
+                                     * Array holding the number of used
+                                     * quads on each level.
+                                     */
+    vector<unsigned int> n_quads_level;
+    
+                                    /**
+                                     * Number of active quads in the
+                                     * whole triangulation.
+                                     */
+    unsigned int n_active_quads;
+
+                                    /**
+                                     * Array holding the number of active
+                                     * quads on each level.
+                                     */
+    vector<unsigned int> n_active_quads_level;
+
+                                    /**
+                                     * Constructor. Set values to zero
+                                     * by default.
+                                     */
+    TriaNumberCache ();
+};
+
+
+
+/**
+ * Cache class used to store the number of used and active elements
+ * (lines or quads etc) within the levels of a triangulation. This
+ * specialization stores the numbers of hexes. Due to the inheritance
+ * from the base class #TriaNumberCache<2>#, the numbers of lines
+ * and quads are also within this class.
+ *
+ * In the old days, whenever one wanted to access one of these
+ * numbers, one had to perform a loop over all lines, e.g., and count
+ * the elements until we hit the end iterator. This is time consuming
+ * and since access to the number of lines etc is a rather frequent
+ * operation, this was not an optimal solution.
+ *
+ * Note: these class should be made local to the triangulation class
+ * once the compiler supports that (gcc2.95 does not at present).
+ *
+ * @author Wolfgang Bangerth, 1999
+ */
+template <>
+struct TriaNumberCache<3> : public TriaNumberCache<2>
+{
+                                    /**
+                                     * Number of used hexes in the whole
+                                     * triangulation.
+                                     */
+    unsigned int n_hexes;
+
+                                    /**
+                                     * Array holding the number of used
+                                     * hexes on each level.
+                                     */
+    vector<unsigned int> n_hexes_level;
+    
+                                    /**
+                                     * Number of active hexes in the
+                                     * whole triangulation.
+                                     */
+    unsigned int n_active_hexes;
+
+                                    /**
+                                     * Array holding the number of active
+                                     * hexes on each level.
+                                     */
+    vector<unsigned int> n_active_hexes_level;
+
+                                    /**
+                                     * Constructor. Set values to zero
+                                     * by default.
+                                     */
+    TriaNumberCache ();
+};
+
+
+
 
 /*------------------------------------------------------------------------*/
 
@@ -2544,119 +2725,71 @@ class Triangulation
                                     /**
                                      *  Return total number of used lines,
                                      *  active or not.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_lines () const;
     
                                     /**
                                      *  Return total number of used lines,
                                      *  active or not on level #level#.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_lines (const unsigned int level) const;
     
                                     /**
                                      * Return total number of active lines.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_lines () const;
     
                                     /**
                                      *  Return total number of active lines,
                                      *  on level #level#.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_lines (const unsigned int level) const;
     
                                     /**
                                      *  Return total number of used quads,
                                      *  active or not.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_quads () const;
 
                                     /**
                                      *  Return total number of used quads,
                                      *  active or not on level #level#.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_quads (const unsigned int level) const;
     
                                     /**
                                      *  Return total number of active quads,
                                      *  active or not.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_quads () const;
     
                                     /**
                                      *  Return total number of active quads,
                                      *  active or not on level #level#.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_quads (const unsigned int level) const;
     
                                     /**
                                      *  Return total number of used hexahedra,
                                      *  active or not.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_hexs() const;
 
                                     /**
                                      *  Return total number of used hexahedra,
                                      *  active or not on level #level#.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_hexs(const unsigned int level) const;
     
                                     /**
                                      *  Return total number of active hexahedra,
                                      *  active or not.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_hexs() const;
     
                                     /**
                                      *  Return total number of active hexahedra,
                                      *  active or not on level #level#.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_hexs(const unsigned int level) const;
 
@@ -2665,18 +2798,6 @@ class Triangulation
                                      *  active or not.
                                      *  Maps to #n_lines()# in one space
                                      *  dimension and so on.
-                                     *
-                                     * Be a bit careful with this function,
-                                     * since it needs quite some computational
-                                     * ressources: the number of cells is
-                                     * computed by looping over all cells,
-                                     * rather than by returning a pre-stored
-                                     * value. It is therefore often better
-                                     * to rewrite lines like
-                                     * #for (i=0; i<tria.n_cells() const; ++i) ...#
-                                     * by
-                                     * #const unsigned int n=tria.n_cells();#
-                                     * #for (i=0; i<n; ++i) ...#.
                                      */
     unsigned int n_cells () const;
 
@@ -2685,10 +2806,6 @@ class Triangulation
                                      *  active or not, on level #level#.
                                      *  Maps to #n_lines(level)# in one space
                                      *  dimension and so on.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_cells (const unsigned int level) const;
 
@@ -2696,10 +2813,6 @@ class Triangulation
                                      *  Return total number of active cells.
                                      *  Maps to #n_active_lines()# in one space
                                      *  dimension and so on.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_cells () const;
 
@@ -2708,10 +2821,6 @@ class Triangulation
                                      * on level #level#.
                                      * Maps to #n_active_lines(level)# in one
                                      * space dimension and so on.
-                                     *
-                                     * Regarding the computational effort of
-                                     * this function, the same applies as
-                                     * for the #n_cells()# function.
                                      */
     unsigned int n_active_cells (const unsigned int level) const;
 
@@ -2965,6 +3074,35 @@ class Triangulation
                                      * coarsening.
                                      */
     void fix_coarsen_flags ();
+
+                                    /**
+                                     * Re-compute the number of lines, quads,
+                                     * etc. This function is called by
+                                     * #execute_{coarsening,refinement}# and
+                                     * by #create_triangulation# after the
+                                     * grid was changed.
+                                     *
+                                     * This function simply delegates to the
+                                     * functions below, which count
+                                     * only a certain class of objects.
+                                     */
+    void update_number_cache ();
+
+                                    /**
+                                     * Recompute the number of lines.
+                                     */
+    void update_number_cache_lines ();
+
+                                    /**
+                                     * Recompute the number of quads.
+                                     */
+    void update_number_cache_quads ();
+
+                                    /**
+                                     * Recompute the number of hexes.
+                                     */
+    void update_number_cache_hexes ();
+
     
                                     /**
                                      *  Array of pointers pointing to the
@@ -3001,6 +3139,24 @@ class Triangulation
                                      */
     MeshSmoothing                    smooth_grid;
 
+                                    /**
+                                     * Cache to hold the numbers of lines,
+                                     * quads, hexes, etc. These numbers
+                                     * are set at the end of the refinement
+                                     * and coarsening functions and enable
+                                     * faster access later on. In the old
+                                     * days, whenever one wanted to access
+                                     * one of these numbers, one had to
+                                     * perform a loop over all lines, e.g.,
+                                     * and count the elements until we hit
+                                     * the end iterator. This is time
+                                     * consuming and since access to the
+                                     * number of lines etc is a rather
+                                     * frequent operation, this was not
+                                     * an optimal solution.
+                                     */
+    TriaNumberCache<dim>             number_cache;
+    
                                     // Friendship includes local classes.
     friend class TriaAccessor<dim>;
     friend class TriaObjectAccessor<1, dim>;
index edf06a5945421eeada505b6f918934a1f5f1912a..5258724bd9e1a2f052530e1c8b4cc4be4220b3d1 100644 (file)
@@ -670,3 +670,32 @@ void TriangulationLevel<3>::block_write (ostream &out) const
 };
 
 
+
+
+
+TriaNumberCache<1>::TriaNumberCache () :
+               n_lines (0),
+               n_active_lines (0) 
+                                                // all other fields are
+                                                // default constructed
+{};
+
+
+
+
+TriaNumberCache<2>::TriaNumberCache () :
+               n_quads (0),
+               n_active_quads (0) 
+                                                // all other fields are
+                                                // default constructed
+{};
+
+
+
+
+TriaNumberCache<3>::TriaNumberCache () :
+               n_hexes (0),
+               n_active_hexes (0) 
+                                                // all other fields are
+                                                // default constructed
+{};
index 8fa9399927947acc6b5de67e6ac0ec43fe041022..98e3785f8aa1e79d64dac77ea6977343086114ce 100644 (file)
@@ -36,8 +36,6 @@ Triangulation<dim>::Triangulation (const MeshSmoothing smooth_grid) :
       boundary[i] = &straight_boundary;
       straight_boundary.subscribe();
     }
-  
-  levels.push_back (new TriangulationLevel<dim>);
 };
 
 
@@ -161,8 +159,7 @@ template <int dim>
 void Triangulation<dim>::copy_triangulation (const Triangulation<dim> &old_tria)
 {
   Assert (vertices.size() == 0, ExcTriangulationNotEmpty());
-  Assert (n_cells () == 0, ExcTriangulationNotEmpty());
-  Assert (levels.size () == 1, ExcTriangulationNotEmpty());
+  Assert (levels.size () == 0, ExcTriangulationNotEmpty());
 
                                   // copy normal elements
   vertices      = old_tria.vertices;
@@ -175,14 +172,13 @@ void Triangulation<dim>::copy_triangulation (const Triangulation<dim> &old_tria)
       boundary[i]      = old_tria.boundary[i];
       boundary[i]->subscribe ();
     }
-                                  // delete only level previously existing,
-                                  // reserve space for new and copy
-  delete levels[0];
-  levels.clear ();
+
   levels.reserve (old_tria.levels.size());
   for (unsigned int level=0; level<old_tria.levels.size(); ++level)
     levels.push_back (new TriangulationLevel<dim>(*old_tria.levels[level]));
 
+  number_cache = old_tria.number_cache;
+  
                                   // note that we need not copy the
                                   // subscriptor!
 };
@@ -225,7 +221,7 @@ void Triangulation<1>::create_triangulation (const vector<Point<1> >    &v,
   const unsigned int dim=1;
   
   Assert (vertices.size() == 0, ExcTriangulationNotEmpty());
-  Assert (n_lines() == 0, ExcTriangulationNotEmpty());
+  Assert (levels.size() == 0, ExcTriangulationNotEmpty());
                                   // check that no forbidden arrays are used
   Assert (subcelldata.check_consistency(dim), ExcInternalError());
 
@@ -238,6 +234,7 @@ void Triangulation<1>::create_triangulation (const vector<Point<1> >    &v,
   vector<vector<int> > lines_at_vertex (v.size());
 
                                   // reserve enough space
+  levels.push_back (new TriangulationLevel<dim>);
   levels[0]->TriangulationLevel<0>::reserve_space (cells.size(), dim);
   levels[0]->TriangulationLevel<1>::reserve_space (cells.size());
   
@@ -310,6 +307,10 @@ void Triangulation<1>::create_triangulation (const vector<Point<1> >    &v,
                                        lines_at_vertex[line->vertex_index(vertex)][0]);
          line->set_neighbor (vertex, neighbor);
        };
+
+
+                                  // re-compute numbers of lines, etc
+  update_number_cache ();  
 };
 
 #endif
@@ -325,8 +326,7 @@ void Triangulation<2>::create_triangulation (const vector<Point<2> >    &v,
   const unsigned int dim=2;
 
   Assert (vertices.size() == 0, ExcTriangulationNotEmpty());
-  Assert (n_lines() == 0, ExcTriangulationNotEmpty());
-  Assert (n_quads() == 0, ExcTriangulationNotEmpty());
+  Assert (levels.size() == 0, ExcTriangulationNotEmpty());
                                   // check that no forbidden arrays are used
   Assert (subcelldata.check_consistency(dim), ExcInternalError());
 
@@ -460,6 +460,7 @@ void Triangulation<2>::create_triangulation (const vector<Point<2> >    &v,
     };
        
                                   // reserve enough space
+  levels.push_back (new TriangulationLevel<dim>);
   levels[0]->TriangulationLevel<0>::reserve_space (cells.size(), dim);
   levels[0]->TriangulationLevel<1>::reserve_space (needed_lines.size());
   levels[0]->TriangulationLevel<2>::reserve_space (cells.size());
@@ -601,6 +602,10 @@ void Triangulation<2>::create_triangulation (const vector<Point<2> >    &v,
       else
        cell->set_neighbor (side,
                            adjacent_cells[cell->line(side)->index()][0]);
+
+
+                                  // re-compute numbers of lines, etc
+  update_number_cache ();  
 };
 
 #endif
@@ -647,9 +652,7 @@ void Triangulation<3>::create_triangulation (const vector<Point<3> >    &v,
   const unsigned int dim=3;
 
   Assert (vertices.size() == 0, ExcTriangulationNotEmpty());
-  Assert (n_lines() == 0, ExcTriangulationNotEmpty());
-  Assert (n_quads() == 0, ExcTriangulationNotEmpty());
-  Assert (n_hexs() == 0, ExcTriangulationNotEmpty());
+  Assert (levels.size() == 0, ExcTriangulationNotEmpty());
                                   // check that no forbidden arrays are used
   Assert (subcelldata.check_consistency(dim), ExcInternalError());
 
@@ -761,6 +764,7 @@ void Triangulation<3>::create_triangulation (const vector<Point<3> >    &v,
                                   // actually set up data structures
                                   // for the lines
                                   // reserve enough space
+  levels.push_back (new TriangulationLevel<dim>);
   levels[0]->TriangulationLevel<0>::reserve_space (cells.size(), dim);
   levels[0]->TriangulationLevel<1>::reserve_space (needed_lines.size());
 
@@ -1104,6 +1108,10 @@ void Triangulation<3>::create_triangulation (const vector<Point<3> >    &v,
       else
        cell->set_neighbor (face,
                            adjacent_cells[cell->quad(face)->index()][0]);
+
+
+                                  // re-compute numbers of lines, etc
+  update_number_cache ();  
 };
 
 #endif
@@ -3120,64 +3128,48 @@ unsigned int Triangulation<3>::n_active_cells (const unsigned int level) const {
 
 template <int dim>
 unsigned int Triangulation<dim>::n_lines () const {
-  unsigned int n=0;
-  for (unsigned int l=0; l<levels.size(); ++l)
-    n += n_lines (l);
-  return n;
+  return number_cache.n_lines;
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_lines (const unsigned int level) const {
-  if (levels[level]->lines.lines.size() == 0)
-    return 0;
-
-                                  // only evaluate begin_/end_line if there
-                                  // are lines.
-  line_iterator line = begin_line (level),
-               endc = (level == levels.size()-1 ?
-                       line_iterator(end_line()) :
-                       begin_line (level+1));
-  unsigned int n=0;
-  for (; line!=endc; ++line)
-    ++n;
-  return n;
+  Assert (level < number_cache.n_lines_level.size(),
+         ExcIndexRange (level, 0, number_cache.n_lines_level.size()));
+  
+  return number_cache.n_lines_level[level];
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_active_lines () const {
-  unsigned int n=0;
-  for (unsigned int l=0; l<levels.size(); ++l)
-    n += n_active_lines (l);
-  return n;
+  return number_cache.n_active_lines;
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_active_lines (const unsigned int level) const {
-  if (levels[level]->lines.lines.size() == 0)
-    return 0;
-
-                                  // only evaluate begin_/end_line if there
-                                  // are lines.
-  active_line_iterator line = begin_active_line (level),
-                      endc = (level == levels.size()-1 ?
-                              active_line_iterator(end_line()) :
-                              begin_active_line (level+1));
-  unsigned int n=0;
-  for (; line!=endc; ++line)
-    ++n;
-  return n;
+  Assert (level < number_cache.n_lines_level.size(),
+         ExcIndexRange (level, 0, number_cache.n_lines_level.size()));
+  
+  return number_cache.n_active_lines_level[level];
 };
 
 
 
 #if deal_II_dimension == 1
 
+template <>
+unsigned int Triangulation<1>::n_quads () const {
+  Assert (false, ExcFunctionNotUseful());
+  return 0;
+};
+
+
+
 template <>
 unsigned int Triangulation<1>::n_quads (const unsigned int) const {
   Assert (false, ExcFunctionNotUseful());
@@ -3193,69 +3185,60 @@ unsigned int Triangulation<1>::n_active_quads (const unsigned int) const {
 };
 
 
+template <>
+unsigned int Triangulation<1>::n_active_quads () const {
+  Assert (false, ExcFunctionNotUseful());
+  return 0;
+};
+
+
 #endif
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_quads () const {
-  unsigned int n=0;
-  for (unsigned int l=0; l<levels.size(); ++l)
-    n += n_quads (l);
-  return n;
+  return number_cache.n_quads;
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_quads (const unsigned int level) const {
-  if (levels[level]->quads.quads.size() == 0)
-    return 0;
-
-                                  // only evaluate begin_/end_quad if there
-                                  // are quads.
-  quad_iterator quad = begin_quad (level),
-               endc = (level == levels.size()-1 ?
-                       quad_iterator(end_quad()) :
-                       begin_quad (level+1));
-  unsigned int n=0;
-  for (; quad!=endc; ++quad)
-    ++n;
-  return n;
+  Assert (level < number_cache.n_quads_level.size(),
+         ExcIndexRange (level, 0, number_cache.n_quads_level.size()));
+  
+  return number_cache.n_quads_level[level];
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_active_quads () const {
-  unsigned int n=0;
-  for (unsigned int l=0; l<levels.size(); ++l)
-    n += n_active_quads (l);
-  return n;
+  return number_cache.n_active_quads;
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_active_quads (const unsigned int level) const {
-  if (levels[level]->quads.quads.size() == 0)
-    return 0;
-
-                                  // only evaluate begin_/end_quad if there
-                                  // are quads.
-  active_quad_iterator quad = begin_active_quad (level),
-                      endc = (level == levels.size()-1 ?
-                              active_quad_iterator(end_quad()) :
-                              begin_active_quad (level+1));
-  unsigned int n=0;
-  for (; quad!=endc; ++quad)
-    ++n;
-  return n;
+  Assert (level < number_cache.n_quads_level.size(),
+         ExcIndexRange (level, 0, number_cache.n_quads_level.size()));
+  
+  return number_cache.n_active_quads_level[level];
 };
 
 
 
 #if deal_II_dimension < 3
 
+template <>
+unsigned int Triangulation<1>::n_hexs () const {
+  Assert (false, ExcFunctionNotUseful());
+  return 0;
+};
+
+
+
 template <>
 unsigned int Triangulation<1>::n_hexs (const unsigned int) const {
   Assert (false, ExcFunctionNotUseful());
@@ -3264,6 +3247,14 @@ unsigned int Triangulation<1>::n_hexs (const unsigned int) const {
 
 
 
+template <>
+unsigned int Triangulation<1>::n_active_hexs () const {
+  Assert (false, ExcFunctionNotUseful());
+  return 0;
+};
+
+
+
 template <>
 unsigned int Triangulation<1>::n_active_hexs (const unsigned int) const {
   Assert (false, ExcFunctionNotUseful());
@@ -3272,6 +3263,14 @@ unsigned int Triangulation<1>::n_active_hexs (const unsigned int) const {
 
 
 
+template <>
+unsigned int Triangulation<2>::n_hexs () const {
+  Assert (false, ExcFunctionNotUseful());
+  return 0;
+};
+
+
+
 template <>
 unsigned int Triangulation<2>::n_hexs (const unsigned int) const {
   Assert (false, ExcFunctionNotUseful());
@@ -3281,70 +3280,52 @@ unsigned int Triangulation<2>::n_hexs (const unsigned int) const {
 
 
 template <>
-unsigned int Triangulation<2>::n_active_hexs (const unsigned int) const {
+unsigned int Triangulation<2>::n_active_hexs () const {
   Assert (false, ExcFunctionNotUseful());
   return 0;
 };
 
-#endif
 
 
+template <>
+unsigned int Triangulation<2>::n_active_hexs (const unsigned int) const {
+  Assert (false, ExcFunctionNotUseful());
+  return 0;
+};
+
+#endif
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_hexs () const {
-  unsigned int n=0;
-  for (unsigned int l=0; l<levels.size(); ++l)
-    n += n_hexs (l);
-  return n;
+  return number_cache.n_hexes;
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_hexs (const unsigned int level) const {
-  if (levels[level]->hexes.hexes.size() == 0)
-    return 0;
-
-                                  // only evaluate begin_/end_hex if there
-                                  // are hexs.
-  hex_iterator hex = begin_hex (level),
-             endc = (level == levels.size()-1 ?
-                     hex_iterator(end_hex()) :
-                     begin_hex (level+1));
-  unsigned int n=0;
-  for (; hex!=endc; ++hex)
-    ++n;
-  return n;
+  Assert (level < number_cache.n_hexes_level.size(),
+         ExcIndexRange (level, 0, number_cache.n_hexes_level.size()));
+  
+  return number_cache.n_hexes_level[level];
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_active_hexs () const {
-  unsigned int n=0;
-  for (unsigned int l=0; l<levels.size(); ++l)
-    n += n_active_hexs (l);
-  return n;
+  return number_cache.n_active_hexes;
 };
 
 
 
 template <int dim>
 unsigned int Triangulation<dim>::n_active_hexs (const unsigned int level) const {
-  if (levels[level]->hexes.hexes.size() == 0)
-    return 0;
-
-                                  // only evaluate begin_/end_hex if there
-                                  // are hexs.
-  active_hex_iterator hex = begin_active_hex (level),
-                    endc = (level == levels.size()-1 ?
-                            active_hex_iterator(end_hex()) :
-                            begin_active_hex (level+1));
-  unsigned int n=0;
-  for (; hex!=endc; ++hex)
-    ++n;
-  return n;
+  Assert (level < number_cache.n_hexes_level.size(),
+         ExcIndexRange (level, 0, number_cache.n_hexes_level.size()));
+  
+  return number_cache.n_active_hexes_level[level];
 };
 
 
@@ -3823,6 +3804,10 @@ void Triangulation<1>::execute_refinement () {
          };      
     };
 
+
+                                  // re-compute number of lines
+  update_number_cache ();
+  
 #ifdef DEBUG
   for (unsigned int level=0; level<levels.size(); ++level) 
     levels[level]->monitor_memory (1);
@@ -4445,6 +4430,11 @@ void Triangulation<2>::execute_refinement () {
            subcells[3]->set_material_id (cell->material_id());
          };
     };
+
+                                  // re-compute number of lines and quads
+  update_number_cache ();
+
+
 #ifdef DEBUG
   for (unsigned int level=0; level<levels.size(); ++level) 
     levels[level]->monitor_memory (2);
@@ -4744,6 +4734,7 @@ void Triangulation<3>::execute_refinement () {
 
                                   ///////////////////////////////////////
                                   // now refine marked quads
+                                  ///////////////////////////////////////
   for (unsigned int level=0; level!=levels.size()-1; ++level)
     {
                                       // only active objects can be refined
@@ -4915,10 +4906,10 @@ void Triangulation<3>::execute_refinement () {
          };
     };
 
-
                                   ///////////////////////////////////
                                   // Now, finally, set up the new
                                   // cells
+                                  ///////////////////////////////////
   for (unsigned int level=0; level!=levels.size()-1; ++level)
     {
                                       // only active objects can be refined
@@ -5665,6 +5656,9 @@ void Triangulation<3>::execute_refinement () {
          };
     };
 
+                                  // re-compute number of lines and quads
+  update_number_cache ();
+
 
 #ifdef DEBUG
   for (unsigned int level=0; level<levels.size(); ++level) 
@@ -5754,6 +5748,10 @@ void Triangulation<dim>::execute_coarsening () {
                                         // use a separate function, since this
                                         // is dimension specific
        delete_children (cell);
+
+                                  // re-compute number of lines and quads
+  update_number_cache ();
+
                                   // in principle no user flags should be
                                   // set any more at this point
 #if DEBUG
@@ -7063,6 +7061,212 @@ void Triangulation<dim>::read_bool_vector (const unsigned int  magic_number1,
 
 
 
+template <int dim>
+void Triangulation<dim>::update_number_cache_lines () 
+{
+                                  ///////////////////////////////////
+                                  // update the number of lines
+                                  // on the different levels in
+                                  // the cache
+  number_cache.n_lines_level.resize (levels.size());
+  number_cache.n_lines = 0;
+  for (unsigned int level=0; level<levels.size(); ++level)
+    {
+                                      // count lines on this level
+      number_cache.n_lines_level[level] = 0;
+      if (levels[level]->lines.lines.size() != 0) 
+       {
+         line_iterator line = begin_line (level),
+                       endc = (level == levels.size()-1 ?
+                               line_iterator(end_line()) :
+                               begin_line (level+1));
+         for (; line!=endc; ++line)
+           ++number_cache.n_lines_level[level];
+       };
+      
+                                      // update total number of lines
+      number_cache.n_lines += number_cache.n_lines_level[level];
+    };
+
+                                  // do the update for the number
+                                  // of active lines as well
+  number_cache.n_active_lines_level.resize (levels.size());
+  number_cache.n_active_lines = 0;
+  for (unsigned int level=0; level<levels.size(); ++level)
+    {
+                                      // count lines on this level
+      number_cache.n_active_lines_level[level] = 0;
+      if (levels[level]->lines.lines.size() != 0) 
+       {
+         active_line_iterator line = begin_active_line (level),
+                              endc = end_active_line (level);
+         for (; line!=endc; ++line)
+           ++number_cache.n_active_lines_level[level];
+       };
+      
+                                      // update total number of lines
+      number_cache.n_active_lines += number_cache.n_active_lines_level[level];
+    };
+};
+
+
+
+#if deal_II_dimension == 1
+
+template <>
+void Triangulation<1>::update_number_cache_quads () 
+{
+  Assert (false, ExcInternalError());
+};
+
+#endif
+
+#if deal_II_dimension >= 2
+
+template <int dim>
+void Triangulation<dim>::update_number_cache_quads () 
+{
+                                  ///////////////////////////////////
+                                  // update the number of quads
+                                  // on the different levels in
+                                  // the cache
+  number_cache.n_quads_level.resize (levels.size());
+  number_cache.n_quads = 0;
+  for (unsigned int level=0; level<levels.size(); ++level)
+    {
+                                      // count quads on this level
+      number_cache.n_quads_level[level] = 0;
+      if (levels[level]->quads.quads.size() != 0) 
+       {
+         quad_iterator quad = begin_quad (level),
+                       endc = (level == levels.size()-1 ?
+                               quad_iterator(end_quad()) :
+                               begin_quad (level+1));
+         for (; quad!=endc; ++quad)
+           ++number_cache.n_quads_level[level];
+       };
+      
+                                      // update total number of quads
+      number_cache.n_quads += number_cache.n_quads_level[level];
+    };
+
+                                  // do the update for the number
+                                  // of active quads as well
+  number_cache.n_active_quads_level.resize (levels.size());
+  number_cache.n_active_quads = 0;
+  for (unsigned int level=0; level<levels.size(); ++level)
+    {
+                                      // count quads on this level
+      number_cache.n_active_quads_level[level] = 0;
+      if (levels[level]->quads.quads.size() != 0) 
+       {
+         active_quad_iterator quad = begin_active_quad (level),
+                              endc = end_active_quad (level);
+         for (; quad!=endc; ++quad)
+           ++number_cache.n_active_quads_level[level];
+       };
+      
+                                      // update total number of quads
+      number_cache.n_active_quads += number_cache.n_active_quads_level[level];
+    };
+};
+
+#endif
+
+
+#if deal_II_dimension == 1
+
+template <>
+void Triangulation<1>::update_number_cache_hexes () 
+{
+  Assert (false, ExcInternalError());
+};
+
+#endif
+
+
+#if deal_II_dimension == 2
+
+template <>
+void Triangulation<2>::update_number_cache_hexes () 
+{
+  Assert (false, ExcInternalError());
+};
+
+#endif
+
+#if deal_II_dimension >= 3
+
+template <int dim>
+void Triangulation<dim>::update_number_cache_hexes () 
+{
+                                  ///////////////////////////////////
+                                  // update the number of hexes
+                                  // on the different levels in
+                                  // the cache
+  number_cache.n_hexes_level.resize (levels.size());
+  number_cache.n_hexes = 0;
+  for (unsigned int level=0; level<levels.size(); ++level)
+    {
+                                      // count hexes on this level
+      number_cache.n_hexes_level[level] = 0;
+      if (levels[level]->hexes.hexes.size() != 0) 
+       {
+         hex_iterator hex = begin_hex (level),
+                     endc = (level == levels.size()-1 ?
+                             hex_iterator(end_hex()) :
+                             begin_hex (level+1));
+         for (; hex!=endc; ++hex)
+           ++number_cache.n_hexes_level[level];
+       };
+      
+                                      // update total number of hexes
+      number_cache.n_hexes += number_cache.n_hexes_level[level];
+    };
+
+                                  // do the update for the number
+                                  // of active hexes as well
+  number_cache.n_active_hexes_level.resize (levels.size());
+  number_cache.n_active_hexes = 0;
+  for (unsigned int level=0; level<levels.size(); ++level)
+    {
+                                      // count hexes on this level
+      number_cache.n_active_hexes_level[level] = 0;
+      if (levels[level]->hexes.hexes.size() != 0) 
+       {
+         active_hex_iterator hex = begin_active_hex (level),
+                            endc = end_active_hex (level);
+         for (; hex!=endc; ++hex)
+           ++number_cache.n_active_hexes_level[level];
+       };
+      
+                                      // update total number of hexes
+      number_cache.n_active_hexes += number_cache.n_active_hexes_level[level];
+    };
+};
+
+#endif
+
+
+template <int dim>
+void Triangulation<dim>::update_number_cache () 
+{
+                                  // play a little bit with switch
+                                  // statements without 'break's
+  switch (dim)
+    {
+      case 3:
+           update_number_cache_hexes();
+      case 2:
+           update_number_cache_quads();
+      case 1:
+           update_number_cache_lines();
+           break;
+           
+      default:
+           Assert (false, ExcInternalError());
+    };
+};
 
 
 

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.