]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
remove specializations in TriaObjects
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 12 Sep 2012 15:26:51 +0000 (15:26 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 12 Sep 2012 15:26:51 +0000 (15:26 +0000)
git-svn-id: https://svn.dealii.org/trunk@26313 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/grid/tria_object.h
deal.II/include/deal.II/grid/tria_objects.h
deal.II/source/grid/tria.cc
deal.II/source/grid/tria_objects.cc
deal.II/source/grid/tria_objects.inst.in

index 8fa11c17aee5ffc1123e82b9a5092e6f6ba9308d..2cd849306e1b855c92d4cbbc0d4487a3e3ab9c37 100644 (file)
@@ -26,7 +26,7 @@ namespace internal
 
 /**
  * Class template for the <tt>structdim</tt>-dimensional cells constituting
- * a Triangulation of dimension <tt>structdim</tt> or lower dimensional
+ * a dealii::Triangulation of dimension <tt>structdim</tt> or lower dimensional
  * objects of higher dimensions.  They are characterized by the
  * (global) indices of their faces, which are cells of dimension
  * <tt>structdim-1</tt> or vertices if <tt>structdim=1</tt>.
index f64e7fd589b05764daeb4ae9b2f7dd1874f053cd..b9c60b5567d694033c82fd4dd02774f63eb4ea3d 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
-//TODO: See if we can unify the class hierarchy a bit using partial
-// specialization of classes here, e.g. declare a general class
-// TriaObjects<dim,G>.
-// To consider for this: in that case we would have to duplicate quite a few
-// things, e.g. TriaObjectsHex is derived from TriaObjects<Hexahedron> and
-// declares mainly additional data. This would have to be changed in case of a
-// partial specialization.
+//TODO: This should all be cleaned up. Currently, only a single
+//function in the library makes use of the odd specializations, and
+//this function is Triangulation::execute_refinement() in 3D. I
+//assume, that the other refinement functions would profit from using
+//next_free_single_object() and next_free_pair_object, but they seem
+//to get around it.
 
 //TODO: The TriaObjects class contains a std::vector<G>. This is only an
 //efficient storage scheme if G is relatively well packed, i.e. it's not a
@@ -35,7 +34,8 @@ DEAL_II_NAMESPACE_OPEN
 //actually the case.
 
 template <int dim, int spacedim> class Triangulation;
-
+template <class Accessor> class TriaRawIterator;
+template <int, int, int> class TriaAccessor;
 
 namespace internal
 {
@@ -51,7 +51,7 @@ namespace internal
  * Objects of these classes are included in the TriaLevel and TriaFaces
  * classes.
  *
- * @author Tobias Leicht, Guido Kanschat, 2006, 2007
+ * @author Tobias Leicht, Guido Kanschat, 2006, 2007, 2012
  */
 
     template <typename G>
@@ -124,7 +124,7 @@ namespace internal
                                           *  already been processed.
                                           *
                                           *  You can clear all used flags using
-                                          *  Triangulation::clear_user_flags().
+                                          *  dealii::Triangulation::clear_user_flags().
                                           */
         std::vector<bool> user_flags;
 
@@ -221,48 +221,52 @@ namespace internal
                                          /**
                                           * Return an iterator to the
                                           * next free slot for a
-                                          * single line. Only
-                                          * implemented for
-                                          * <code>G=TriaObject<1>
-                                          * </code>.
+                                          * single object. This
+                                          * function is only used by
+                                          * dealii::Triangulation::execute_refinement()
+                                          * in 3D.
+                                         *
+                                         * @warning Interestingly,
+                                         * this function is not used
+                                         * for 1D or 2D
+                                         * triangulations, where it
+                                         * seems the authors of the
+                                         * refinement function insist
+                                         * on reimplementing its
+                                         * contents.
+                                         *
+                                         * @todo This function is
+                                         * not instantiated for the
+                                         * codim-one case
                                           */
         template <int dim, int spacedim>
-        typename dealii::Triangulation<dim,spacedim>::raw_line_iterator
-        next_free_single_line (const dealii::Triangulation<dim,spacedim> &tria);
+       dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >
+        next_free_single_object (const dealii::Triangulation<dim,spacedim> &tria);
 
                                          /**
                                           * Return an iterator to the
                                           * next free slot for a pair
-                                          * of lines. Only implemented
-                                          * for <code>G=TriaObject<1>
-                                          * </code>.
-                                          */
-        template <int dim, int spacedim>
-        typename dealii::Triangulation<dim,spacedim>::raw_line_iterator
-        next_free_pair_line (const dealii::Triangulation<dim,spacedim> &tria);
-
-                                         /**
-                                          * Return an iterator to the
-                                          * next free slot for a
-                                          * single quad. Only
-                                          * implemented for
-                                          * <code>G=TriaObject@<2@>
-                                          * </code>.
+                                          * of objects. This
+                                          * function is only used by
+                                          * dealii::Triangulation::execute_refinement()
+                                          * in 3D.
+                                         *
+                                         * @warning Interestingly,
+                                         * this function is not used
+                                         * for 1D or 2D
+                                         * triangulations, where it
+                                         * seems the authors of the
+                                         * refinement function insist
+                                         * on reimplementing its
+                                         * contents.
+                                         *
+                                         * @todo This function is
+                                         * not instantiated for the
+                                         * codim-one case
                                           */
         template <int dim, int spacedim>
-        typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator
-        next_free_single_quad (const dealii::Triangulation<dim,spacedim> &tria);
-
-                                         /**
-                                          * Return an iterator to the
-                                          * next free slot for a pair
-                                          * of quads. Only implemented
-                                          * for <code>G=TriaObject@<2@>
-                                          * </code>.
-                                          */
-        template <int dim, int spacedim>
-        typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator
-        next_free_pair_quad (const dealii::Triangulation<dim,spacedim> &tria);
+       dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >
+        next_free_pair_object (const dealii::Triangulation<dim,spacedim> &tria);
 
                                          /**
                                           * Return an iterator to the
@@ -401,7 +405,7 @@ namespace internal
                         "but you can only ask for " << arg2 <<"_iterators.");
 
                                          /**
-                                          * Triangulation objects can
+                                          * dealii::Triangulation objects can
                                           * either access a user
                                           * pointer or a user
                                           * index. What you tried to
@@ -930,16 +934,6 @@ namespace internal
 
 // declaration of explicit specializations
 
-    template<>
-    void
-    TriaObjects<TriaObject<1> >::reserve_space (const unsigned int new_lines_in_pairs,
-                                      const unsigned int new_lines_single);
-
-    template<>
-    void
-    TriaObjects<TriaObject<2> >::reserve_space (const unsigned int new_quads_in_pairs,
-                                      const unsigned int new_quads_single);
-
     template<>
     void
     TriaObjects<TriaObject<2> >::monitor_memory (const unsigned int) const;
index 3743949c8d08a52a3d503f711e9c44aaef961876..0e0f18ad22a910a0ff9d1812170ec0e1ca3813a5 100644 (file)
@@ -5374,7 +5374,7 @@ namespace internal
                                                        // up the two child lines
                                                        // (++ takes care of the
                                                        // end of the vector)
-                      next_unused_line=triangulation.faces->lines.next_free_pair_line(triangulation);
+                      next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation);
                       Assert(next_unused_line.state() == IteratorState::valid,
                              ExcInternalError());
 
@@ -5506,7 +5506,7 @@ namespace internal
                                                          // the quad
                         typename Triangulation<dim,spacedim>::raw_line_iterator new_line;
 
-                        new_line=triangulation.faces->lines.next_free_single_line(triangulation);
+                        new_line=triangulation.faces->lines.next_free_single_object(triangulation);
                         Assert (new_line->used() == false,
                                 ExcCellShouldBeUnused());
 
@@ -5561,7 +5561,7 @@ namespace internal
                                                          // created quads.
                         typename Triangulation<dim,spacedim>::raw_quad_iterator new_quads[2];
 
-                        next_unused_quad=triangulation.faces->quads.next_free_pair_quad(triangulation);
+                        next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation);
                         new_quads[0] = next_unused_quad;
                         Assert (new_quads[0]->used() == false, ExcCellShouldBeUnused());
 
@@ -5691,7 +5691,7 @@ namespace internal
                                                                  // all info
                                 typename Triangulation<dim,spacedim>::raw_line_iterator new_child[2];
 
-                                new_child[0]=new_child[1]=triangulation.faces->lines.next_free_pair_line(triangulation);
+                                new_child[0]=new_child[1]=triangulation.faces->lines.next_free_pair_object(triangulation);
                                 ++new_child[1];
 
                                 new_child[0]->set_used_flag();
@@ -5960,7 +5960,7 @@ namespace internal
 
                                                                  // now search a slot for the two
                                                                  // child lines
-                                next_unused_line=triangulation.faces->lines.next_free_pair_line(triangulation);
+                                next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation);
 
                                                                  // set the child
                                                                  // pointer of the present
@@ -6096,7 +6096,7 @@ namespace internal
                                                                // anisotropically and the
                                                                // two lines end up as
                                                                // children of new line
-                              next_unused_line=triangulation.faces->lines.next_free_pair_line(triangulation);
+                              next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation);
 
                             new_lines[i] = next_unused_line;
                             ++next_unused_line;
@@ -6193,7 +6193,7 @@ namespace internal
                                                          // created quads.
                         typename Triangulation<dim,spacedim>::raw_quad_iterator new_quads[4];
 
-                        next_unused_quad=triangulation.faces->quads.next_free_pair_quad(triangulation);
+                        next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation);
 
                         new_quads[0] = next_unused_quad;
                         Assert (new_quads[0]->used() == false, ExcCellShouldBeUnused());
@@ -6202,7 +6202,7 @@ namespace internal
                         new_quads[1] = next_unused_quad;
                         Assert (new_quads[1]->used() == false, ExcCellShouldBeUnused());
 
-                        next_unused_quad=triangulation.faces->quads.next_free_pair_quad(triangulation);
+                        next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation);
                         new_quads[2] = next_unused_quad;
                         Assert (new_quads[2]->used() == false, ExcCellShouldBeUnused());
 
@@ -6366,7 +6366,7 @@ namespace internal
                         new_lines(n_new_lines);
                       for (unsigned int i=0; i<n_new_lines; ++i)
                         {
-                          new_lines[i] = triangulation.faces->lines.next_free_single_line(triangulation);
+                          new_lines[i] = triangulation.faces->lines.next_free_single_object(triangulation);
 
                           Assert (new_lines[i]->used() == false,
                                   ExcCellShouldBeUnused());
@@ -6385,7 +6385,7 @@ namespace internal
                         new_quads(n_new_quads);
                       for (unsigned int i=0; i<n_new_quads; ++i)
                         {
-                          new_quads[i] = triangulation.faces->quads.next_free_single_quad(triangulation);
+                          new_quads[i] = triangulation.faces->quads.next_free_single_object(triangulation);
 
                           Assert (new_quads[i]->used() == false,
                                   ExcCellShouldBeUnused());
@@ -6412,7 +6412,6 @@ namespace internal
                         {
                           if (i%2==0)
                             next_unused_hex=triangulation.levels[level+1]->cells.next_free_hex(triangulation,level+1);
-
                           else
                             ++next_unused_hex;
 
index 7ce3ccfd6909a448a432dac7b9d4e16ec95c4b08..ff30938e7eedcb11f0a810c703e5f17e39bda98e 100644 (file)
@@ -29,105 +29,10 @@ namespace internal
 {
   namespace Triangulation
   {
-    template<>
-    void
-    TriaObjects<TriaObject<1> >::reserve_space (const unsigned int new_lines_in_pairs,
-                                      const unsigned int new_lines_single)
-    {
-      Assert(new_lines_in_pairs%2==0, ExcInternalError());
-
-      next_free_single=0;
-      next_free_pair=0;
-      reverse_order_next_free_single=false;
-
-                                       // count the number of lines, of
-                                       // unused single lines and of
-                                       // unused pairs of lines
-      unsigned int n_lines=0;
-      unsigned int n_unused_pairs=0;
-      unsigned int n_unused_singles=0;
-      for (unsigned int i=0; i<used.size(); ++i)
-        {
-          if (used[i])
-            ++n_lines;
-          else if (i+1<used.size())
-            {
-              if (used[i+1])
-                {
-                  ++n_unused_singles;
-                  if (next_free_single==0)
-                    next_free_single=i;
-                }
-              else
-                {
-                  ++n_unused_pairs;
-                  if (next_free_pair==0)
-                    next_free_pair=i;
-                  ++i;
-                }
-            }
-          else
-            ++n_unused_singles;
-        }
-      Assert(n_lines+2*n_unused_pairs+n_unused_singles==used.size(),
-             ExcInternalError());
-
-                                       // how many single lines are needed in
-                                       // addition to n_unused_singles?
-      const int additional_single_lines=
-        new_lines_single-n_unused_singles;
-
-      unsigned int new_size=
-        used.size() + new_lines_in_pairs - 2*n_unused_pairs;
-      if (additional_single_lines>0)
-        new_size+=additional_single_lines;
-
-                                       // only allocate space if necessary
-      if (new_size>cells.size())
-        {
-          cells.reserve (new_size);
-          cells.insert (cells.end(),
-                        new_size-cells.size(),
-                        TriaObject<1> ());
-
-          used.reserve (new_size);
-          used.insert (used.end(),
-                       new_size-used.size(),
-                       false);
-
-          user_flags.reserve (new_size);
-          user_flags.insert (user_flags.end(),
-                             new_size-user_flags.size(),
-                             false);
-
-          children.reserve (new_size);
-          children.insert (children.end(),
-                           new_size-children.size(),
-                           -1);
-
-          boundary_or_material_id.reserve (new_size);
-          boundary_or_material_id.insert (boundary_or_material_id.end(),
-                              new_size-boundary_or_material_id.size(),
-                              BoundaryOrMaterialId());
-
-          user_data.reserve (new_size);
-          user_data.insert (user_data.end(),
-                            new_size-user_data.size(),
-                            UserData());
-        }
-
-      if (n_unused_singles==0)
-        {
-          next_free_single=new_size-1;
-          reverse_order_next_free_single=true;
-        }
-    }
-
-
-    template <>
+    template <class G>
     template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_line_iterator
-    TriaObjects<TriaObject<1> >::next_free_single_line (const dealii::Triangulation<dim,spacedim> &tria)
+    dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >
+    TriaObjects<G>::next_free_single_object (const dealii::Triangulation<dim,spacedim> &tria)
     {
                                        // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this.
 
@@ -166,19 +71,19 @@ namespace internal
           if (pos>0)
             next_free_single=pos-1;
           else
-                                             // no valid single line anymore
-            return tria.end_line();
+                                             // no valid single object anymore
+            return dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >(&tria, -1, -1);
         }
 
-      return typename dealii::Triangulation<dim,spacedim>::raw_line_iterator(&tria,0,pos);
+      return dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >(&tria, 0, pos);
     }
 
 
 
-    template <>
+    template <class G>
     template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_line_iterator
-    TriaObjects<TriaObject<1> >::next_free_pair_line (const dealii::Triangulation<dim,spacedim> &tria)
+    dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >
+    TriaObjects<G>::next_free_pair_object (const dealii::Triangulation<dim,spacedim> &tria)
     {
                                        // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this.
 
@@ -194,58 +99,35 @@ namespace internal
             }
       if (pos>=last)
                                          // no free slot
-        return tria.end_line();
+        return dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >(&tria, -1, -1);
       else
         next_free_pair=pos+2;
 
-      return typename dealii::Triangulation<dim,spacedim>::raw_line_iterator(&tria,0,pos);
-    }
-
-
-
-    template <>
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator
-    TriaObjects<TriaObject<1> >::next_free_single_quad (const dealii::Triangulation<dim,spacedim> &tria)
-    {
-      Assert(false, ExcWrongIterator("quad","line"));
-      return tria.end_quad();
-    }
-
-
-
-    template <>
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator
-    TriaObjects<TriaObject<1> >::next_free_pair_quad (const dealii::Triangulation<dim,spacedim> &tria)
-    {
-      Assert(false, ExcWrongIterator("quad","line"));
-      return tria.end_quad();
+      return dealii::TriaRawIterator<dealii::TriaAccessor<G::dimension,dim,spacedim> >(&tria, 0, pos);
     }
 
 
-
-    template<>
+    template<class G>
     void
-    TriaObjects<TriaObject<2> >::reserve_space (const unsigned int new_quads_in_pairs,
-                                      const unsigned int new_quads_single)
+    TriaObjects<G>::reserve_space (const unsigned int new_objects_in_pairs,
+                                  const unsigned int new_objects_single)
     {
-      Assert(new_quads_in_pairs%2==0, ExcInternalError());
+      Assert(new_objects_in_pairs%2==0, ExcInternalError());
 
       next_free_single=0;
       next_free_pair=0;
       reverse_order_next_free_single=false;
 
-                                       // count the number of lines, of
-                                       // unused single lines and of
-                                       // unused pairs of lines
-      unsigned int n_quads=0;
+                                       // count the number of objects, of
+                                       // unused single objects and of
+                                       // unused pairs of objects
+      unsigned int n_objects=0;
       unsigned int n_unused_pairs=0;
       unsigned int n_unused_singles=0;
       for (unsigned int i=0; i<used.size(); ++i)
         {
           if (used[i])
-            ++n_quads;
+            ++n_objects;
           else if (i+1<used.size())
             {
               if (used[i+1])
@@ -265,18 +147,18 @@ namespace internal
           else
             ++n_unused_singles;
         }
-      Assert(n_quads+2*n_unused_pairs+n_unused_singles==used.size(),
+      Assert(n_objects+2*n_unused_pairs+n_unused_singles==used.size(),
              ExcInternalError());
 
-                                       // how many single quads are needed in
-                                       // addition to n_unused_quads?
-      const int additional_single_quads=
-        new_quads_single-n_unused_singles;
+                                       // how many single objects are needed in
+                                       // addition to n_unused_objects?
+      const int additional_single_objects=
+        new_objects_single-n_unused_singles;
 
       unsigned int new_size=
-        used.size() + new_quads_in_pairs - 2*n_unused_pairs;
-      if (additional_single_quads>0)
-        new_size+=additional_single_quads;
+        used.size() + new_objects_in_pairs - 2*n_unused_pairs;
+      if (additional_single_objects>0)
+        new_size+=additional_single_objects;
 
                                        // only allocate space if necessary
       if (new_size>cells.size())
@@ -284,7 +166,7 @@ namespace internal
           cells.reserve (new_size);
           cells.insert (cells.end(),
                         new_size-cells.size(),
-                        TriaObject<2> ());
+                        G ());
 
           used.reserve (new_size);
           used.insert (used.end(),
@@ -296,16 +178,19 @@ namespace internal
                              new_size-user_flags.size(),
                              false);
 
-          children.reserve (2*new_size);
+         const unsigned int factor = GeometryInfo<G::dimension>::max_children_per_cell / 2;
+          children.reserve (factor*new_size);
           children.insert (children.end(),
-                           2*new_size-children.size(),
+                           factor*new_size-children.size(),
                            -1);
 
-          refinement_cases.reserve (new_size);
-          refinement_cases.insert (refinement_cases.end(),
-                               new_size - refinement_cases.size(),
-                               RefinementCase<2>::no_refinement);
-
+         if (G::dimension > 1)
+           {
+             refinement_cases.reserve (new_size);
+             refinement_cases.insert (refinement_cases.end(),
+                                      new_size - refinement_cases.size(),
+                                      RefinementCase<G::dimension>::no_refinement);
+           }
 
           boundary_or_material_id.reserve (new_size);
           boundary_or_material_id.insert (boundary_or_material_id.end(),
@@ -326,109 +211,6 @@ namespace internal
     }
 
 
-
-
-    template <>
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator
-    TriaObjects<TriaObject<2> >::next_free_single_quad (const dealii::Triangulation<dim,spacedim> &tria)
-    {
-                                       // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this.
-
-      int pos=next_free_single,
-         last=used.size()-1;
-      if (!reverse_order_next_free_single)
-        {
-                                           // first sweep forward, only use
-                                           // really single slots, do not use
-                                           // pair slots
-          for (; pos<last; ++pos)
-            if (!used[pos])
-              if (used[++pos])
-                {
-                                                   // this was a single slot
-                  pos-=1;
-                  break;
-                }
-          if (pos>=last)
-            {
-              reverse_order_next_free_single=true;
-              next_free_single=used.size()-1;
-              pos=used.size()-1;
-            }
-          else
-            next_free_single=pos+1;
-        }
-
-      if(reverse_order_next_free_single)
-        {
-                                           // second sweep, use all slots, even
-                                           // in pairs
-          for(;pos>=0;--pos)
-            if (!used[pos])
-              break;
-          if (pos>0)
-            next_free_single=pos-1;
-          else
-                                             // no valid single quad anymore
-            return tria.end_quad();
-        }
-
-      return typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator(&tria,0,pos);
-    }
-
-
-
-    template <>
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator
-    TriaObjects<TriaObject<2> >::next_free_pair_quad (const dealii::Triangulation<dim,spacedim> &tria)
-    {
-                                       // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this.
-
-      int pos=next_free_pair,
-         last=used.size()-1;
-      for (; pos<last; ++pos)
-        if (!used[pos])
-          if (!used[++pos])
-            {
-                                               // this was a pair slot
-              pos-=1;
-              break;
-            }
-      if (pos>=last)
-                                         // no free slot
-        return tria.end_quad();
-      else
-        next_free_pair=pos+2;
-
-      return typename dealii::Triangulation<dim,spacedim>::raw_quad_iterator(&tria,0,pos);
-    }
-
-
-    template <>
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_line_iterator
-    TriaObjects<TriaObject<2> >::next_free_single_line (const dealii::Triangulation<dim,spacedim> &tria)
-    {
-      Assert(false, ExcWrongIterator("line","quad"));
-      return tria.end_line();
-    }
-
-
-
-    template <>
-    template <int dim, int spacedim>
-    typename dealii::Triangulation<dim,spacedim>::raw_line_iterator
-    TriaObjects<TriaObject<2> >::next_free_pair_line (const dealii::Triangulation<dim,spacedim> &tria)
-    {
-      Assert(false, ExcWrongIterator("line","quad"));
-      return tria.end_line();
-    }
-
-
-
-
     template <>
     template <int dim, int spacedim>
     typename dealii::Triangulation<dim,spacedim>::raw_hex_iterator
@@ -536,9 +318,9 @@ namespace internal
       next_free_pair=0;
       reverse_order_next_free_single=false;
 
-                                       // count the number of lines, of unused
-                                       // single lines and of unused pairs of
-                                       // lines
+                                       // count the number of objects, of unused
+                                       // single objects and of unused pairs of
+                                       // objects
       unsigned int n_quads=0;
       unsigned int n_unused_pairs=0;
       unsigned int n_unused_singles=0;
index 76fc1259ad1191b908c3c2ae2035c7d84c010666..777e0ba6cb09e42e37fcd7058441690109b740b5 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2006, 2007, 2010 by the deal.II authors
+//    Copyright (C) 2006, 2007, 2010, 2012 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 for (deal_II_dimension : DIMENSIONS)
   {
 #if deal_II_dimension > 1
-    template dealii::Triangulation<deal_II_dimension>::raw_line_iterator
-    TriaObjects<TriaObject<1> >::next_free_single_line(const dealii::Triangulation<deal_II_dimension> &);
-    template dealii::Triangulation<deal_II_dimension>::raw_line_iterator
-    TriaObjects<TriaObject<1> >::next_free_pair_line(const dealii::Triangulation<deal_II_dimension> &);
-    template dealii::Triangulation<deal_II_dimension>::raw_quad_iterator
-    TriaObjects<TriaObject<1> >::next_free_single_quad(const dealii::Triangulation<deal_II_dimension> &);
-    template dealii::Triangulation<deal_II_dimension>::raw_quad_iterator
-    TriaObjects<TriaObject<1> >::next_free_pair_quad(const dealii::Triangulation<deal_II_dimension> &);
 
-    template dealii::Triangulation<deal_II_dimension>::raw_line_iterator
-    TriaObjects<TriaObject<2> >::next_free_single_line(const dealii::Triangulation<deal_II_dimension> &);
-    template dealii::Triangulation<deal_II_dimension>::raw_line_iterator
-    TriaObjects<TriaObject<2> >::next_free_pair_line(const dealii::Triangulation<deal_II_dimension> &);
-    template dealii::Triangulation<deal_II_dimension>::raw_quad_iterator
-    TriaObjects<TriaObject<2> >::next_free_single_quad(const dealii::Triangulation<deal_II_dimension> &);
-    template dealii::Triangulation<deal_II_dimension>::raw_quad_iterator
-    TriaObjects<TriaObject<2> >::next_free_pair_quad(const dealii::Triangulation<deal_II_dimension> &);
+    template dealii::TriaRawIterator<dealii::TriaAccessor<1,deal_II_dimension,deal_II_dimension> >
+    TriaObjects<TriaObject<1> >::next_free_single_object (const dealii::Triangulation<deal_II_dimension> &tria);
+    template dealii::TriaRawIterator<dealii::TriaAccessor<1,deal_II_dimension,deal_II_dimension> >
+    TriaObjects<TriaObject<1> >::next_free_pair_object (const dealii::Triangulation<deal_II_dimension> &tria);
+    template dealii::TriaRawIterator<dealii::TriaAccessor<2,deal_II_dimension,deal_II_dimension> >
+    TriaObjects<TriaObject<2> >::next_free_single_object (const dealii::Triangulation<deal_II_dimension> &tria);
+    template dealii::TriaRawIterator<dealii::TriaAccessor<2,deal_II_dimension,deal_II_dimension> >
+    TriaObjects<TriaObject<2> >::next_free_pair_object (const dealii::Triangulation<deal_II_dimension> &tria);
 #endif
 #if deal_II_dimension == 3
     template dealii::Triangulation<deal_II_dimension>::raw_hex_iterator

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.