+/*------------------------------------------------------------------------*/
+
+/**
+ * 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 ();
+};
+
+
+
/*------------------------------------------------------------------------*/
/**
* 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;
* 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;
* 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;
* 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;
* 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;
* 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
*/
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>;
boundary[i] = &straight_boundary;
straight_boundary.subscribe();
}
-
- levels.push_back (new TriangulationLevel<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;
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!
};
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());
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());
lines_at_vertex[line->vertex_index(vertex)][0]);
line->set_neighbor (vertex, neighbor);
};
+
+
+ // re-compute numbers of lines, etc
+ update_number_cache ();
};
#endif
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());
};
// 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());
else
cell->set_neighbor (side,
adjacent_cells[cell->line(side)->index()][0]);
+
+
+ // re-compute numbers of lines, etc
+ update_number_cache ();
};
#endif
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());
// 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());
else
cell->set_neighbor (face,
adjacent_cells[cell->quad(face)->index()][0]);
+
+
+ // re-compute numbers of lines, etc
+ update_number_cache ();
};
#endif
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());
};
+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());
+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());
+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());
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];
};
};
};
+
+ // re-compute number of lines
+ update_number_cache ();
+
#ifdef DEBUG
for (unsigned int level=0; level<levels.size(); ++level)
levels[level]->monitor_memory (1);
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);
///////////////////////////////////////
// now refine marked quads
+ ///////////////////////////////////////
for (unsigned int level=0; level!=levels.size()-1; ++level)
{
// only active objects can be refined
};
};
-
///////////////////////////////////
// Now, finally, set up the new
// cells
+ ///////////////////////////////////
for (unsigned int level=0; level!=levels.size()-1; ++level)
{
// only active objects can be refined
};
};
+ // re-compute number of lines and quads
+ update_number_cache ();
+
#ifdef DEBUG
for (unsigned int level=0; level<levels.size(); ++level)
// 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
+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());
+ };
+};