]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
implement user_index functions
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 23 May 2007 20:50:10 +0000 (20:50 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 23 May 2007 20:50:10 +0000 (20:50 +0000)
git-svn-id: https://svn.dealii.org/trunk@14687 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/include/grid/tria_accessor.h
deal.II/deal.II/source/grid/tria.cc
deal.II/deal.II/source/grid/tria_accessor.cc

index 8684358f454198dac5a7a7364c7838e432f5d230..0c92e0b52084b7e00a359da280846858a1f19226 100644 (file)
@@ -978,17 +978,17 @@ namespace internal
  *   execute_coarsening_and_refinement() function uses the face user
  *   flags.
  *
- *   There is another set of user data, namely a <tt>void *</tt>, for
- *   each line, quad, etc. You can access these user pointers through
- *   the functions <tt>user_pointer()</tt>,
- *   <tt>clear_user_pointer()</tt> and <tt>set_user_pointer(p)</tt> in
+ *   There is another set of user data, which can be either an
+ *   <tt>unsigned int</tt> or a <tt>void *</tt>, for
+ *   each line, quad, etc. You can access these through
+ *   the functions listed under <tt>User data</tt> in
  *   the accessor classes. These pointers are not used nor changed in
  *   many places of the library, and those classes and functions that
  *   do use them should document this clearly; the most prominent user
  *   of these pointers is the SolutionTransfer class which uses the
  *   cell->user_pointers.
  *
- *   The value of these user pointers is @p NULL by default. Note that
+ *   The value of these user indices or pointers is @p NULL by default. Note that
  *   the pointers are not inherited to children upon
  *   refinement. Still, after a remeshing they are available on all
  *   cells, where they were set on the previous mesh (unless, of
@@ -1745,7 +1745,7 @@ class Triangulation : public Subscriptor
 
 
                                     /**
-                                     *  @name User flag handling
+                                     *  @name User data
                                      */
                                     /*@{*/
                                     /**
@@ -1866,12 +1866,34 @@ class Triangulation : public Subscriptor
                                      * Load the user flags located on hexs.
                                      */
     void load_user_flags_hex (const std::vector<bool> &v);
-                                    /*@}*/
 
                                     /**
+                                     * Clear all user pointers and
+                                     * indices and allow the use of
+                                     * both for next access.
+                                     */
+    void clear_user_data ();
+    
+                                    /**
+                                     * @deprecated User
+                                     * clear_user_data() instead.
+                                     *
                                      *  Clear all user pointers.
                                      */
     void clear_user_pointers ();
+    
+                                    /**
+                                     * Save all user indices. The
+                                     * output vector is resized if
+                                     * necessary.
+                                     */
+    void save_user_indices (std::vector<unsigned int> &v) const;
+
+                                    /**
+                                     * Read the information stored by
+                                     * save_user_indices().
+                                     */
+    void load_user_indices (const std::vector<unsigned int> &v);
 
                                     /**
                                      * Save all user pointers. The
@@ -1882,12 +1904,50 @@ class Triangulation : public Subscriptor
 
                                     /**
                                      * Read the information stored by
-                                     * @p save_user_pointers.
+                                     * save_user_pointers().
                                      */
     void load_user_pointers (const std::vector<void *> &v);
 
                                     /**
-                                     * Save the user pointers on
+                                     * Save the user indices on
+                                     * lines. The output vector is
+                                     * resized if necessary.
+                                     */
+    void save_user_indices_line (std::vector<unsigned int> &v) const;
+
+                                    /**
+                                     * Load the user indices located
+                                     * on lines.
+                                     */
+    void load_user_indices_line (const std::vector<unsigned int> &v);
+
+                                    /**
+                                     * Save the user indices on
+                                     * quads. The output vector is
+                                     * resized if necessary.
+                                     */
+    void save_user_indices_quad (std::vector<unsigned int> &v) const;
+
+                                    /**
+                                     * Load the user indices located
+                                     * on quads.
+                                     */
+    void load_user_indices_quad (const std::vector<unsigned int> &v);
+
+                                    /**
+                                     * Save the user indices on
+                                     * hexes. The output vector is
+                                     * resized if necessary.
+                                     */
+    void save_user_indices_hex (std::vector<unsigned int> &v) const;
+
+                                    /**
+                                     * Load the user indices located
+                                     * on hexs.
+                                     */
+    void load_user_indices_hex (const std::vector<unsigned int> &v);
+                                    /**
+                                     * Save the user indices on
                                      * lines. The output vector is
                                      * resized if necessary.
                                      */
@@ -1924,7 +1984,8 @@ class Triangulation : public Subscriptor
                                      * on hexs.
                                      */
     void load_user_pointers_hex (const std::vector<void *> &v);
-                                    /* ------------------------------------ */
+                                    /*@}*/
+    
     
                                     /**
                                      *  @name Cell iterator functions
@@ -1936,8 +1997,8 @@ class Triangulation : public Subscriptor
                                      *  has no cells, a past-the-end iterator
                                      *  is returned.
                                      *
-                                     *  This function calls @p begin_raw_line
-                                     *  in 1D and @p begin_raw_quad in 2D.
+                                     *  This function calls begin_raw_line()
+                                     *  in 1D and begin_raw_quad() in 2D.
                                      */
     raw_cell_iterator    begin_raw   (const unsigned int level = 0) const;
 
index 1832356a21106c7aeaf7a339cfe917243d2fa716..afe9f67b19f959488ee4d7a7b2d8ee21ea0f7a70 100644 (file)
@@ -605,6 +605,10 @@ class TriaObjectAccessor :  public TriaAccessor<celldim,dim>
                                      */
     bool used () const;
 
+                                    /**
+                                     * @name User data
+                                     */
+                                    /*@{*/
                                     /**
                                      *  Set the @p used flag. You
                                      *  should know quite exactly
@@ -706,7 +710,7 @@ class TriaObjectAccessor :  public TriaAccessor<celldim,dim>
                                      * the triangulation.
                                      */
     void recursively_set_user_pointer (void *p) const;
-
+    
                                     /**
                                      * Clear the user pointer of this
                                      * object and all of its
@@ -716,7 +720,55 @@ class TriaObjectAccessor :  public TriaAccessor<celldim,dim>
                                      * function.
                                      */
     void recursively_clear_user_pointer () const;
+    
+                                    /**
+                                     * Set the user index
+                                     * to @p p.
+                                     */
+    void set_user_index (unsigned int p) const;
+    
+                                    /**
+                                     * Reset the user index to 0.
+                                     */
+    void clear_user_index () const;
+    
+                                    /**
+                                     * Access the value of the user
+                                     * index.
+                                     */
+    unsigned int user_index () const;
+
+                                    /**
+                                     * Set the user index of this
+                                     * object and all its children.
+                                     *
+                                     * Note that the user index is
+                                     * not inherited under mesh
+                                     * refinement, so after mesh
+                                     * refinement there might be
+                                     * cells or faces that don't have
+                                     * the expected user indices. In
+                                     * this case, simply loop over
+                                     * all the elements of the
+                                     * coarsest level that has this
+                                     * information, and use this
+                                     * function to recursively set
+                                     * the user index of all finer
+                                     * levels of the triangulation.
+                                     */
+    void recursively_set_user_index (unsigned int p) const;
 
+                                    /**
+                                     * Clear the user index of this
+                                     * object and all of its
+                                     * descendants. The same holds as
+                                     * said for the
+                                     * recursively_set_user_index()
+                                     * function.
+                                     */
+    void recursively_clear_user_index () const;
+                                    /*@}*/
+    
                                     /**
                                      *  Pointer to the @p ith
                                      *  child.
@@ -1335,6 +1387,53 @@ class TriaObjectAccessor<1, dim> :  public TriaAccessor<1,dim>
                                      */
     void * user_pointer () const;
     
+                                    /**
+                                     * Set the user index
+                                     * to @p p.
+                                     */
+    void set_user_index (unsigned int p) const;
+    
+                                    /**
+                                     * Reset the user index to 0.
+                                     */
+    void clear_user_index () const;
+    
+                                    /**
+                                     * Access the value of the user
+                                     * index.
+                                     */
+    unsigned int user_index () const;
+
+                                    /**
+                                     * Set the user index of this
+                                     * object and all its children.
+                                     *
+                                     * Note that the user index is
+                                     * not inherited under mesh
+                                     * refinement, so after mesh
+                                     * refinement there might be
+                                     * cells or faces that don't have
+                                     * the expected user indices. In
+                                     * this case, simply loop over
+                                     * all the elements of the
+                                     * coarsest level that has this
+                                     * information, and use this
+                                     * function to recursively set
+                                     * the user index of all finer
+                                     * levels of the triangulation.
+                                     */
+    void recursively_set_user_index (unsigned int p) const;
+
+                                    /**
+                                     * Clear the user index of this
+                                     * object and all of its
+                                     * descendants. The same holds as
+                                     * said for the
+                                     * recursively_set_user_index()
+                                     * function.
+                                     */
+    void recursively_clear_user_index () const;
+
                                     /**
                                      *  Return a pointer to the
                                      *  @p ith child.
@@ -1920,6 +2019,52 @@ class TriaObjectAccessor<2, dim> :  public TriaAccessor<2,dim>
                                      */
     void recursively_clear_user_pointer () const;
 
+                                    /**
+                                     * Set the user index
+                                     * to @p p.
+                                     */
+    void set_user_index (unsigned int p) const;
+    
+                                    /**
+                                     * Reset the user index to 0.
+                                     */
+    void clear_user_index () const;
+    
+                                    /**
+                                     * Access the value of the user
+                                     * index.
+                                     */
+    unsigned int user_index () const;
+
+                                    /**
+                                     * Set the user index of this
+                                     * object and all its children.
+                                     *
+                                     * Note that the user index is
+                                     * not inherited under mesh
+                                     * refinement, so after mesh
+                                     * refinement there might be
+                                     * cells or faces that don't have
+                                     * the expected user indices. In
+                                     * this case, simply loop over
+                                     * all the elements of the
+                                     * coarsest level that has this
+                                     * information, and use this
+                                     * function to recursively set
+                                     * the user index of all finer
+                                     * levels of the triangulation.
+                                     */
+    void recursively_set_user_index (unsigned int p) const;
+
+                                    /**
+                                     * Clear the user index of this
+                                     * object and all of its
+                                     * descendants. The same holds as
+                                     * said for the
+                                     * recursively_set_user_index()
+                                     * function.
+                                     */
+    void recursively_clear_user_index () const;
                                     /**
                                      *  Return a pointer to the @p ith
                                      *  child.
@@ -2541,6 +2686,53 @@ class TriaObjectAccessor<3, dim> :  public TriaAccessor<3,dim>
                                      */
     void recursively_clear_user_pointer () const;
 
+                                    /**
+                                     * Set the user index
+                                     * to @p p.
+                                     */
+    void set_user_index (unsigned int p) const;
+    
+                                    /**
+                                     * Reset the user index to 0.
+                                     */
+    void clear_user_index () const;
+    
+                                    /**
+                                     * Access the value of the user
+                                     * index.
+                                     */
+    unsigned int user_index () const;
+
+                                    /**
+                                     * Set the user index of this
+                                     * object and all its children.
+                                     *
+                                     * Note that the user index is
+                                     * not inherited under mesh
+                                     * refinement, so after mesh
+                                     * refinement there might be
+                                     * cells or faces that don't have
+                                     * the expected user indices. In
+                                     * this case, simply loop over
+                                     * all the elements of the
+                                     * coarsest level that has this
+                                     * information, and use this
+                                     * function to recursively set
+                                     * the user index of all finer
+                                     * levels of the triangulation.
+                                     */
+    void recursively_set_user_index (unsigned int p) const;
+
+                                    /**
+                                     * Clear the user index of this
+                                     * object and all of its
+                                     * descendants. The same holds as
+                                     * said for the
+                                     * recursively_set_user_index()
+                                     * function.
+                                     */
+    void recursively_clear_user_index () const;
+    
                                     /**
                                      *  Return a pointer to the
                                      *  @p ith child.
index 776a9515023320c78434409cd4680ea89a511c37..15c156928ae8f39f71c44fb025fa4098169c7d99 100644 (file)
@@ -2124,7 +2124,7 @@ void Triangulation<dim>::load_coarsen_flags (const std::vector<bool> &v)
 #if deal_II_dimension == 1
 
 template <>
-void Triangulation<1>::clear_user_pointers ()
+void Triangulation<1>::clear_user_data ()
 {
   for (unsigned int level=0;level<levels.size();++level)
     levels[level]->lines.clear_user_data();
@@ -2132,6 +2132,14 @@ void Triangulation<1>::clear_user_pointers ()
 
 
 
+template <>
+void Triangulation<1>::clear_user_pointers ()
+{
+  clear_user_data
+}
+
+
+
 template <>
 void Triangulation<1>::clear_user_flags ()
 {
@@ -2156,7 +2164,7 @@ void Triangulation<1>::clear_user_flags_hex ()
 #if deal_II_dimension == 2
 
 template <>
-void Triangulation<2>::clear_user_pointers ()
+void Triangulation<2>::clear_user_data ()
 {
   for (unsigned int level=0;level<levels.size();++level)
     levels[level]->quads.clear_user_data();
@@ -2165,6 +2173,14 @@ void Triangulation<2>::clear_user_pointers ()
 
 
 
+template <>
+void Triangulation<2>::clear_user_pointers ()
+{
+  clear_user_data();
+}
+
+
+
 template <>
 void Triangulation<2>::clear_user_flags ()
 {
@@ -2185,7 +2201,7 @@ void Triangulation<2>::clear_user_flags_hex ()
 #if deal_II_dimension == 3
 
 template <>
-void Triangulation<3>::clear_user_pointers ()
+void Triangulation<3>::clear_user_data ()
 {
   for (unsigned int level=0;level<levels.size();++level)
     levels[level]->hexes.clear_user_data();
@@ -2195,6 +2211,14 @@ void Triangulation<3>::clear_user_pointers ()
 
 
 
+template <>
+void Triangulation<3>::clear_user_pointers ()
+{
+  clear_user_data();
+}
+
+
+
 template <>
 void Triangulation<3>::clear_user_flags ()
 {
@@ -2597,6 +2621,209 @@ void Triangulation<dim>::load_user_flags_hex (const std::vector<bool> &v)
 
 
 
+template <int dim>
+void Triangulation<dim>::save_user_indices (std::vector<unsigned int> &v) const
+{
+                                  // clear vector and append all the
+                                  // stuff later on
+  v.clear ();
+
+  std::vector<unsigned int> tmp;
+
+  save_user_indices_line (tmp);
+  v.insert (v.end(), tmp.begin(), tmp.end());
+
+  if (dim >= 2)
+    {
+      save_user_indices_quad (tmp);
+      v.insert (v.end(), tmp.begin(), tmp.end());
+    }
+  
+  if (dim >= 3)
+    {
+      save_user_indices_hex (tmp);
+      v.insert (v.end(), tmp.begin(), tmp.end());
+    }      
+
+  if (dim >= 4)
+    Assert (false, ExcNotImplemented());
+}
+
+
+
+template <int dim>
+void Triangulation<dim>::load_user_indices (const std::vector<unsigned int> &v)
+{
+  Assert (v.size() == n_lines()+n_quads()+n_hexs(), ExcInternalError());
+  std::vector<unsigned int> tmp;
+
+                                  // first extract the indices
+                                  // belonging to lines
+  tmp.insert (tmp.end(),
+             v.begin(), v.begin()+n_lines());
+                                  // and set the lines
+  load_user_indices_line (tmp);
+
+  if (dim >= 2)
+    {
+      tmp.clear ();
+      tmp.insert (tmp.end(),
+                 v.begin()+n_lines(), v.begin()+n_lines()+n_quads());
+      load_user_indices_quad (tmp);
+    }
+  
+  if (dim >= 3)
+    {
+      tmp.clear ();
+      tmp.insert (tmp.end(),
+                 v.begin()+n_lines()+n_quads(), v.begin()+n_lines()+n_quads()+n_hexs());
+      load_user_indices_hex (tmp);
+    }      
+
+  if (dim >= 4)
+    Assert (false, ExcNotImplemented());
+}
+
+
+
+template <int dim>
+void Triangulation<dim>::save_user_indices_line (std::vector<unsigned int> &v) const
+{
+  v.resize (n_lines(), 0);
+  std::vector<unsigned int>::iterator  i = v.begin();
+  line_iterator line = begin_line(),
+               endl = end_line();
+  for (; line!=endl; ++line, ++i)
+    *i = line->user_index();
+}
+
+
+
+template <int dim>
+void Triangulation<dim>::load_user_indices_line (const std::vector<unsigned int> &v)
+{
+  Assert (v.size() == n_lines(), ExcGridReadError());
+  
+  line_iterator line = begin_line(),
+               endl = end_line();
+  std::vector<unsigned int>::const_iterator i = v.begin();
+  for (; line!=endl; ++line, ++i)
+    line->set_user_index(*i);
+}
+
+#if deal_II_dimension == 1
+
+
+template <>
+void Triangulation<1>::save_user_indices_quad (std::vector<unsigned int> &) const
+{
+  Assert (false, ExcImpossibleInDim(1));
+}
+
+
+
+template <>
+void Triangulation<1>::load_user_indices_quad (const std::vector<unsigned int> &)
+{
+  Assert (false, ExcImpossibleInDim(1));
+}
+
+
+
+template <>
+void Triangulation<1>::save_user_indices_hex (std::vector<unsigned int> &) const
+{
+  Assert (false, ExcImpossibleInDim(1));
+}
+
+
+
+template <>
+void Triangulation<1>::load_user_indices_hex (const std::vector<unsigned int> &)
+{
+  Assert (false, ExcImpossibleInDim(1));
+}
+
+#endif
+
+
+template <int dim>
+void Triangulation<dim>::save_user_indices_quad (std::vector<unsigned int> &v) const
+{
+  v.resize (n_quads(), 0);
+  std::vector<unsigned int>::iterator  i = v.begin();
+  quad_iterator quad = begin_quad(),
+               endq = end_quad();
+  for (; quad!=endq; ++quad, ++i)
+    *i = quad->user_index();
+}
+
+
+
+template <int dim>
+void Triangulation<dim>::load_user_indices_quad (const std::vector<unsigned int> &v)
+{
+  Assert (v.size() == n_quads(), ExcGridReadError());
+  
+  quad_iterator quad = begin_quad(),
+               endq = end_quad();
+  std::vector<unsigned int>::const_iterator i = v.begin();
+  for (; quad!=endq; ++quad, ++i)
+    quad->set_user_index(*i);
+}
+
+
+#if deal_II_dimension == 2
+
+
+
+template <>
+void Triangulation<2>::save_user_indices_hex (std::vector<unsigned int> &) const
+{
+  Assert (false, ExcImpossibleInDim(2));
+}
+
+
+
+template <>
+void Triangulation<2>::load_user_indices_hex (const std::vector<unsigned int> &)
+{
+  Assert (false, ExcImpossibleInDim(2));
+}
+
+
+#endif
+
+
+template <int dim>
+void Triangulation<dim>::save_user_indices_hex (std::vector<unsigned int> &v) const
+{
+  v.resize (n_hexs(), 0);
+  std::vector<unsigned int>::iterator  i = v.begin();
+  hex_iterator hex = begin_hex(),
+             endh = end_hex();
+  for (; hex!=endh; ++hex, ++i)
+    *i = hex->user_index();
+}
+
+
+
+template <int dim>
+void Triangulation<dim>::load_user_indices_hex (const std::vector<unsigned int> &v)
+{
+  Assert (v.size() == n_hexs(), ExcGridReadError());
+  
+  hex_iterator hex = begin_hex(),
+             endh = end_hex();
+  std::vector<unsigned int>::const_iterator i = v.begin();
+  for (; hex!=endh; ++hex, ++i)
+    hex->set_user_index(*i);
+}
+
+
+
+//----------------------------------------------------------------------//
+
 template <int dim>
 void Triangulation<dim>::save_user_pointers (std::vector<void *> &v) const
 {
index ac51722d2be4ec716673c7b11f7ba29a29eb3282..ec1b369fb604d60bb9c1d7211cd327f2473e887e 100644 (file)
@@ -161,6 +161,59 @@ TriaObjectAccessor<1, dim>::recursively_clear_user_pointer () const
 
 
 
+template <int dim>
+void TriaObjectAccessor<1, dim>::set_user_index (unsigned int p) const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  lines().user_index(this->present_index) = p;
+}
+
+
+
+template <int dim>
+void TriaObjectAccessor<1, dim>::clear_user_index () const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  lines().user_index(this->present_index) = 0;
+}
+
+
+
+template <int dim>
+unsigned int TriaObjectAccessor<1, dim>::user_index () const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  return lines().user_index(this->present_index);
+}
+
+
+
+template <int dim>
+void
+TriaObjectAccessor<1, dim>::recursively_set_user_index (unsigned int p) const
+{
+  set_user_index (p);
+
+  if (has_children())
+    for (unsigned int c=0; c<2; ++c)
+      child(c)->recursively_set_user_index (p);
+}
+
+
+
+template <int dim>
+void
+TriaObjectAccessor<1, dim>::recursively_clear_user_index () const
+{
+  clear_user_index ();
+
+  if (has_children())
+    for (unsigned int c=0; c<2; ++c)
+      child(c)->recursively_clear_user_index ();
+}
+
+
+
 template <int dim>
 void TriaObjectAccessor<1, dim>::set_children (const int index) const
 {
@@ -405,6 +458,59 @@ TriaObjectAccessor<2, dim>::recursively_clear_user_pointer () const
 
 
 
+template <int dim>
+void TriaObjectAccessor<2, dim>::set_user_index (unsigned int p) const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  quads().user_index(this->present_index) = p;
+}
+
+
+
+template <int dim>
+void TriaObjectAccessor<2, dim>::clear_user_index () const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  quads().user_index(this->present_index) = 0;
+}
+
+
+
+template <int dim>
+unsigned int  TriaObjectAccessor<2, dim>::user_index () const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  return quads().user_index(this->present_index);
+}
+
+
+
+template <int dim>
+void
+TriaObjectAccessor<2, dim>::recursively_set_user_index (unsigned int p) const
+{
+  set_user_index (p);
+
+  if (has_children())
+    for (unsigned int c=0; c<n_children(); ++c)
+      child(c)->recursively_set_user_index (p);
+}
+
+
+
+template <int dim>
+void
+TriaObjectAccessor<2, dim>::recursively_clear_user_index () const
+{
+  clear_user_index ();
+
+  if (has_children())
+    for (unsigned int c=0; c<n_children(); ++c)
+      child(c)->recursively_clear_user_index ();
+}
+
+
+
 template <int dim>
 void TriaObjectAccessor<2, dim>::set_children (const int index) const
 {
@@ -847,6 +953,30 @@ void * TriaObjectAccessor<3, 3>::user_pointer () const
   return this->tria->levels[this->present_level]->hexes.user_pointer(this->present_index);
 }
 
+
+template <>
+void TriaObjectAccessor<3, 3>::set_user_index (unsigned int p) const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  this->tria->levels[this->present_level]->hexes.user_index(this->present_index) = p;
+}
+
+
+template <>
+void TriaObjectAccessor<3, 3>::clear_user_index () const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  this->tria->levels[this->present_level]->hexes.user_index(this->present_index) = 0;
+}
+
+
+template <>
+unsigned int  TriaObjectAccessor<3, 3>::user_index () const
+{
+  Assert (used(), TriaAccessorExceptions::ExcCellNotUsed());
+  return this->tria->levels[this->present_level]->hexes.user_index(this->present_index);
+}
+
 #endif
 
 
@@ -875,6 +1005,31 @@ TriaObjectAccessor<3, dim>::recursively_clear_user_pointer () const
 }
 
 
+template <int dim>
+void
+TriaObjectAccessor<3, dim>::recursively_set_user_index (unsigned int p) const
+{
+  set_user_index (p);
+
+  if (has_children())
+    for (unsigned int c=0; c<n_children(); ++c)
+      child(c)->recursively_set_user_index (p);
+}
+
+
+
+template <int dim>
+void
+TriaObjectAccessor<3, dim>::recursively_clear_user_index () const
+{
+  clear_user_index ();
+
+  if (has_children())
+    for (unsigned int c=0; c<n_children(); ++c)
+      child(c)->recursively_clear_user_index ();
+}
+
+
 #if deal_II_dimension == 3
 
 template <>

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.