]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
replace classes for Line, Quad and Hexahedron by a single template
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 27 May 2007 19:57:42 +0000 (19:57 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 27 May 2007 19:57:42 +0000 (19:57 +0000)
git-svn-id: https://svn.dealii.org/trunk@14716 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria_accessor.h
deal.II/deal.II/include/grid/tria_accessor.templates.h
deal.II/deal.II/include/grid/tria_faces.h
deal.II/deal.II/include/grid/tria_hex.h [deleted file]
deal.II/deal.II/include/grid/tria_levels.h
deal.II/deal.II/include/grid/tria_line.h [deleted file]
deal.II/deal.II/include/grid/tria_object.h [new file with mode: 0644]
deal.II/deal.II/include/grid/tria_objects.h
deal.II/deal.II/include/grid/tria_quad.h [deleted file]
deal.II/deal.II/source/grid/tria.cc
deal.II/deal.II/source/grid/tria_accessor.cc

index 26bfa2b015f251c350c60689d91d76279b272c15..325f7693828a9dd5cb88624281c5eb043b3aaba6 100644 (file)
@@ -39,9 +39,7 @@ namespace internal
 {
   namespace Triangulation
   {
-    class Line;
-    class Quad;
-    class Hexahedron;
+    template <int dim> class TriaObject;
     template <typename G> class TriaObjects;
   }
 }
@@ -515,21 +513,21 @@ class TriaObjectAccessor :  public TriaAccessor<celldim,dim>
                                      * <tt>celldim==1</tt>.
                                      */
 
-    void set (const internal::Triangulation::Line&) const;
+    void set (const internal::Triangulation::TriaObject<1>&) const;
     
                                     /**
                                      * Copy the data of the given
                                      * quad. Only implemented for
                                      * <tt>celldim==2</tt>.
                                      */
-    void set (const internal::Triangulation::Quad&) const;
+    void set (const internal::Triangulation::TriaObject<2>&) const;
     
                                     /**
                                      * Copy the data of the given
                                      * hex. Only implemented for
                                      * <tt>celldim==3</tt>.
                                      */
-    void set (const internal::Triangulation::Hexahedron&) const;
+    void set (const internal::Triangulation::TriaObject<3>&) const;
 
                                     /**
                                      *  Index of vertex. The convention
@@ -1239,7 +1237,7 @@ class TriaObjectAccessor<1, dim> :  public TriaAccessor<1,dim>
                                      *  Copy the data of the given
                                      *  line.
                                      */
-    void set (const internal::Triangulation::Line &l) const;
+    void set (const internal::Triangulation::TriaObject<1> &l) const;
 
                                     /**
                                      *  Return the index of vertex
@@ -1792,7 +1790,7 @@ class TriaObjectAccessor<1, dim> :  public TriaAccessor<1,dim>
                                      * for all dim without
                                      * specialization.
                                      */
-    internal::Triangulation::TriaObjects<internal::Triangulation::Line> & lines() const;
+    internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<1> > & lines() const;
 
   private:
     
@@ -1843,7 +1841,7 @@ class TriaObjectAccessor<2, dim> :  public TriaAccessor<2,dim>
                                     /**
                                      *  Copy the data of the given quad.
                                      */
-    void set (const internal::Triangulation::Quad &q) const;
+    void set (const internal::Triangulation::TriaObject<2> &q) const;
 
                                     /**
                                      * Return index of a vertex of a
@@ -2446,7 +2444,7 @@ class TriaObjectAccessor<2, dim> :  public TriaAccessor<2,dim>
                                      * implemented for all dimensions
                                      * without specialization.
                                      */
-    internal::Triangulation::TriaObjects<internal::Triangulation::Quad> & quads() const;
+    internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<2> > & quads() const;
 
   private:
     
@@ -2497,7 +2495,7 @@ class TriaObjectAccessor<3, dim> :  public TriaAccessor<3,dim>
                                      *  Copy the data of the given
                                      *  hex.
                                      */
-    void set (const internal::Triangulation::Hexahedron &h) const;
+    void set (const internal::Triangulation::TriaObject<3> &h) const;
     
                                     /**
                                      *  Return index of a vertex of a hex in the internal structures of Triangulation.
@@ -3549,8 +3547,8 @@ class CellAccessor :  public TriaObjectAccessor<dim,dim>
 template <> Point<2> TriaObjectAccessor<2, 2>::barycenter () const;
 template <> Point<3> TriaObjectAccessor<2, 3>::barycenter () const;
 template <> Point<3> TriaObjectAccessor<3, 3>::barycenter () const;
-template <> internal::Triangulation::TriaObjects<internal::Triangulation::Line> &TriaObjectAccessor<1, 1>::lines() const;
-template <> internal::Triangulation::TriaObjects<internal::Triangulation::Quad> &TriaObjectAccessor<2, 2>::quads() const;
+template <> internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<1> > &TriaObjectAccessor<1, 1>::lines() const;
+template <> internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<2> > &TriaObjectAccessor<2, 2>::quads() const;
 template <> bool CellAccessor<1>::at_boundary () const;
 template <> unsigned char CellAccessor<1>::material_id () const;
 template <> void CellAccessor<1>::set_material_id (const unsigned char mat_id) const;
index 0e03942d1da0417f029b7c65d8837fcdab9444bb..12e354bfa368725911e8829c8c85a8086d150c58 100644 (file)
@@ -616,7 +616,7 @@ unsigned int
 TriaObjectAccessor<2,dim>::line_index (const unsigned int i) const
 {
   Assert (i<4, ExcIndexRange(i,0,4));
-  return quads().cells[this->present_index].line(i);
+  return quads().cells[this->present_index].face(i);
 }
 
 
@@ -989,7 +989,7 @@ TriaObjectAccessor<3,3>::quad_index (const unsigned int i) const
 {
   Assert (i<6, ExcIndexRange(i,0,6));
 
-  return this->tria->levels[this->present_level]->hexes.cells[this->present_index].quad(i);
+  return this->tria->levels[this->present_level]->hexes.cells[this->present_index].face(i);
 }
 
 
index 8981bdbb098301bf9fbce3e4c156cc0a8d8158bd..897aeb7d7c4dff5edcd3287a42eb1919e30489f9 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2006 by the deal.II authors
+//    Copyright (C) 2006, 2007 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -14,9 +14,7 @@
 #define __deal2__tria_faces_h
 
 #include <base/config.h>
-#include <grid/tria_line.h>
-#include <grid/tria_quad.h>
-#include <grid/tria_hex.h>
+#include <grid/tria_object.h>
 #include <grid/tria_objects.h>
 
 DEAL_II_NAMESPACE_OPEN
diff --git a/deal.II/deal.II/include/grid/tria_hex.h b/deal.II/deal.II/include/grid/tria_hex.h
deleted file mode 100644 (file)
index d26e690..0000000
+++ /dev/null
@@ -1,157 +0,0 @@
-//---------------------------------------------------------------------------
-//    $Id$
-//    Version: $Name$
-//
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
-//
-//    This file is subject to QPL and may not be  distributed
-//    without copyright and license information. Please refer
-//    to the file deal.II/doc/license.html for the  text  and
-//    further information on this license.
-//
-//---------------------------------------------------------------------------
-#ifndef __deal2__tria_hex_h
-#define __deal2__tria_hex_h
-
-
-#include <base/config.h>
-#include <base/exceptions.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-namespace internal
-{
-  namespace Triangulation
-  {
-    
-/**
- *   @p Hexahedrons denote the fundamental entities of triangulations
- *   in three dimensions. They are characterized by the (global)
- *   indices of the six bounding quadrilaterals.
- *
- *   A heaxhedron itself has one index, as far as the topological part
- *   handled in the triangulation is concerned: the index in the level
- *   it belongs to. The level index is implicitly given by the
- *   position in the <tt>hexes.hexes</tt> list attached to the information
- *   of each level of the triangulation.
- *
- *   For conventions on the numbering of faces, lines and vertices
- *   within a hexahedron, refer to the documentation of the
- *   Triangulation class.
- *
- * @ingroup grid
- *   @author Wolfgang Bangerth, 1998
- */
-    class Hexahedron
-    {
-      public:
-
-                                         /**
-                                          *  Construct a Hex with quad
-                                          *  indices @p i0 through
-                                          *  @p i5. By default, indices
-                                          *  are set to -1, i.e. an
-                                          *  invalid value.
-                                          *
-                                          *  No convention is set as of
-                                          *  now on the order of quads
-                                          */
-        Hexahedron (const int i0 = -1,
-                    const int i1 = -1,
-                    const int i2 = -1,
-                    const int i3 = -1,
-                    const int i4 = -1,
-                    const int i5 = -1);
-    
-                                         /**
-                                          *  Return the index of quad
-                                          *  @p i=0 through 5.
-                                          */
-        int quad (const int i) const;
-    
-                                         /**
-                                          *  Set the index of the @p ith
-                                          *  quad to @p index. @p i=0
-                                          *  through 5.
-                                          */
-        void set_quad (const int i, const int index);
-    
-                                         /**
-                                          * Determine an estimate for the
-                                          * memory consumption (in bytes)
-                                          * of this object.
-                                          */
-        static unsigned int memory_consumption ();
-
-                                         /**
-                                          *  Exception
-                                          */ 
-        DeclException1 (ExcRange,
-                        int,
-                        << "Indices for the quad number must be 0 through 5, "
-                        << "but you gave " << arg1); 
-      protected:
-                                         /**
-                                          * Indices of the six quads that
-                                          * bound this hexahedron. The
-                                          * level index is implicitly
-                                          * given by the level of this
-                                          * hexahedron, so this is only
-                                          * the index within this level.
-                                          */
-        int quads[6];
-    };
-
-
-/*----------------------------- Inline Function: Hexahedron ------------------------*/
-
-
-    inline
-    Hexahedron::Hexahedron (const int i0, const int i1,
-                            const int i2, const int i3,
-                            const int i4, const int i5)
-    {
-      quads[0] = i0;
-      quads[1] = i1;
-      quads[2] = i2;
-      quads[3] = i3;
-      quads[4] = i4;
-      quads[5] = i5;
-    }
-
-
-
-    inline
-    int Hexahedron::quad (const int i) const
-    {
-      Assert ((i>=0) && (i<6),
-              ExcRange(i));
-      return quads[i];
-    }
-
-
-
-    inline
-    void Hexahedron::set_quad (const int i, const int index)
-    {
-      Assert ((i>=0) && (i<6),
-              ExcRange(i));
-      quads[i] = index;
-    }
-
-
-
-    inline
-    unsigned int
-    Hexahedron::memory_consumption ()
-    {
-      return sizeof(Hexahedron);
-    }
-    
-  }
-  
-}
-
-DEAL_II_NAMESPACE_CLOSE
-
-#endif
index df0154c0cc8b044bdfe3b07f89b989aae065f44e..8d6405561d71fbb193376e4d3f97eba8214627dd 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -16,9 +16,7 @@
 
 #include <base/config.h>
 #include <vector>
-#include <grid/tria_line.h>
-#include <grid/tria_quad.h>
-#include <grid/tria_hex.h>
+#include <grid/tria_object.h>
 #include <base/point.h>
 #include <grid/tria_objects.h>
 
diff --git a/deal.II/deal.II/include/grid/tria_line.h b/deal.II/deal.II/include/grid/tria_line.h
deleted file mode 100644 (file)
index b528ed2..0000000
+++ /dev/null
@@ -1,131 +0,0 @@
-//---------------------------------------------------------------------------
-//    $Id$
-//    Version: $Name$
-//
-//    Copyright (C) 1998, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
-//
-//    This file is subject to QPL and may not be  distributed
-//    without copyright and license information. Please refer
-//    to the file deal.II/doc/license.html for the  text  and
-//    further information on this license.
-//
-//---------------------------------------------------------------------------
-#ifndef __deal2__tria_line_h
-#define __deal2__tria_line_h
-
-
-#include <base/config.h>
-#include <base/exceptions.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-namespace internal
-{
-  namespace Triangulation
-  {
-    
-/**
- *   Lines denote the boundaries of quads and the edges of hexaeders. They are
- *   characterized by the (global) indices of the endpoints.
- *
- *   A line itself has one index, as far as the topological part handled in
- *   the triangulation is concerned: the index in the level
- *   it belongs to. The level index is implicitly given by the position
- *   in the <tt>lines.lines</tt> list attached to the information of each level.
- *
- * @ingroup grid
- *   @author Wolfgang Bangerth, 1998
- */
-    class Line
-    {
-      public:
-                                         /**
-                                          *  Construct a line with end point
-                                          *  indices @p i0 and @p i1. By default,
-                                          *  indices are set to -1, i.e. an
-                                          *  invalid value.
-                                          */
-        Line (const int i0 = -1,
-              const int i1 = -1);
-    
-                                         /**
-                                          *  Return the index of end point @p i=0
-                                          *  or 1.
-                                          */
-        int vertex (const int i) const;
-    
-                                         /**
-                                          *  Set the index of the @p ith vertex to
-                                          *  @p index. @p i=0 or 1.
-                                          */
-        void set_vertex (const int i, const int index);
-    
-                                         /**
-                                          * Determine an estimate for the
-                                          * memory consumption (in bytes)
-                                          * of this object.
-                                          */
-        static unsigned int memory_consumption ();
-
-                                         /**
-                                          *  Exception
-                                          */
-        DeclException1 (ExcRange,
-                        int,
-                        << "Indices for the point number must be 0 or 1, "
-                        << "but you gave " << arg1);
-      
-      protected:
-                                         /**
-                                          *  Global indices of the two end points.
-                                          */
-        int end_points[2];
-    };
-
-
-/*----------------------------- Inline Function: Line ------------------------*/
-
-
-    inline                               // truly in-line here!
-    Line::Line (const int i0, const int i1)
-    {
-      end_points[0] = i0;
-      end_points[1] = i1;
-    }
-
-
-
-    inline
-    int Line::vertex (const int i) const
-    {
-      Assert ((i==0) || (i==1),
-              ExcRange(i));
-      return end_points[i];
-    }
-
-
-
-    inline
-    void Line::set_vertex (const int i, const int index)
-    {
-      Assert ((i==0) || (i==1),
-              ExcRange(i));
-      end_points[i] = index;
-    }
-
-
-
-    inline
-    unsigned int
-    Line::memory_consumption ()
-    {
-      return sizeof(Line);
-    }
-    
-  }
-  
-}
-
-DEAL_II_NAMESPACE_CLOSE
-
-#endif
diff --git a/deal.II/deal.II/include/grid/tria_object.h b/deal.II/deal.II/include/grid/tria_object.h
new file mode 100644 (file)
index 0000000..5e40775
--- /dev/null
@@ -0,0 +1,202 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__tria_object_h
+#define __deal2__tria_object_h
+
+
+#include <base/config.h>
+#include <base/exceptions.h>
+#include <base/geometry_info.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace internal
+{
+  namespace Triangulation
+  {
+    
+/**
+ * Class template for the <tt>dim</tt>-dimensional cells constituting
+ * a Triangulation of dimension <tt>dim</tt> or lower dimensional
+ * objects of higher dimensions.  They are characterized by the
+ * (global) indices of their faces, which are cells of dimension
+ * <tt>dim-1</tt> or vertices if <tt>dim=1</tt>.
+ *
+ * @author Guido Kanschat, 2007
+ */
+    template <int dim>
+    class TriaObject
+    {
+      public:
+                                        /**
+                                         * Default constructor,
+                                         * setting all face indices
+                                         * to invalid values.
+                                         */                              
+       TriaObject ();
+       
+                                        /**
+                                         * Constructor for a line
+                                         * object with the numbers of
+                                         * its two end points.
+                                         *
+                                         * Throws an exception if
+                                         * dimension is not one.
+                                         */
+       TriaObject (const int i0, const int i1);
+
+                                        /**
+                                         * Constructor for a quadrilateral
+                                         * object with the numbers of
+                                         * its four lines.
+                                         *
+                                         * Throws an exception if
+                                         * dimension is not two.
+                                         */
+       TriaObject (const int i0, const int i1,
+                   const int i2, const int i3);
+
+                                        /**
+                                         * Constructor for a hexahedron
+                                         * object with the numbers of
+                                         * its six quadrilaterals.
+                                         *
+                                         * Throws an exception if
+                                         * dimension is not two.
+                                         */
+       TriaObject (const int i0, const int i1,
+                    const int i2, const int i3,
+                    const int i4, const int i5);
+    
+
+                                        /**
+                                         * Return the index of the
+                                         * ith face object.
+                                         */
+       int face (const unsigned int i) const;
+
+                                        /**
+                                         * Set the index of the ith
+                                         * face object.
+                                         */
+       void set_face (const unsigned int i, const int index);
+
+                                        /**
+                                          * Determine an estimate for the
+                                          * memory consumption (in bytes)
+                                          * of this object.
+                                          */
+        static unsigned int memory_consumption ();
+       
+      protected:
+                                         /**
+                                          *  Global indices of the two end points.
+                                          */
+        int faces[GeometryInfo<dim>::faces_per_cell];
+    };
+
+                                    /// Legacy typedef
+    typedef TriaObject<1> Line;
+
+                                    /// Legacy typedef
+    typedef TriaObject<2> Quad;
+
+                                    /// Legacy typedef
+    typedef TriaObject<3> Hexahedron;
+
+//----------------------------------------------------------------------//
+
+    template <int dim>
+    inline
+    TriaObject<dim>::TriaObject ()
+    {
+      for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+       faces[i] = -1;
+    }
+
+
+    template <int dim>
+    inline
+    TriaObject<dim>::TriaObject (const int i0, const int i1)
+    {
+      Assert (dim==1, ExcImpossibleInDim(dim));
+      faces[0] = i0;
+      faces[1] = i1;
+    }
+
+
+    template <int dim>
+    inline
+    TriaObject<dim>::TriaObject (const int i0, const int i1, const int i2, const int i3)
+    {
+      Assert (dim==2, ExcImpossibleInDim(dim));
+      faces[0] = i0;
+      faces[1] = i1;
+      faces[2] = i2;
+      faces[3] = i3;
+    }
+
+
+    template <int dim>
+    inline
+    TriaObject<dim>::TriaObject (const int i0, const int i1, const int i2,
+                                const int i3, const int i4, const int i5)
+    {
+      Assert (dim==3, ExcImpossibleInDim(dim));
+      faces[0] = i0;
+      faces[1] = i1;
+      faces[2] = i2;
+      faces[3] = i3;
+      faces[4] = i4;
+      faces[5] = i5;
+    }
+
+
+    template <int dim>
+    inline
+    int TriaObject<dim>::face (const unsigned int i) const
+    {
+      Assert (i<GeometryInfo<dim>::faces_per_cell,
+              ExcIndexRange(i,0,GeometryInfo<dim>::faces_per_cell));
+      return faces[i];
+    }
+    
+    
+    
+    template <int dim>
+    inline
+    void TriaObject<dim>::set_face (const unsigned int i, const int index)
+    {
+      Assert (i<GeometryInfo<dim>::faces_per_cell,
+              ExcIndexRange(i,0,GeometryInfo<dim>::faces_per_cell));
+      faces[i] = index;
+    }
+    
+    
+    
+    template <int dim>
+    inline
+    unsigned int
+    TriaObject<dim>::memory_consumption ()
+    {
+      return sizeof(TriaObject<dim>);
+    }
+    
+    
+  }
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index c8b001f9ecf6de7361267216288e8a08eb1db9c6..defdce2a0c3470ae65911745fbb5d0e9f17564dd 100644 (file)
@@ -16,9 +16,7 @@
 #include <base/config.h>
 #include <base/exceptions.h>
 #include <base/geometry_info.h>
-#include <grid/tria_line.h>
-#include <grid/tria_quad.h>
-#include <grid/tria_hex.h>
+#include <grid/tria_object.h>
 #include <vector>
 
 DEAL_II_NAMESPACE_OPEN
diff --git a/deal.II/deal.II/include/grid/tria_quad.h b/deal.II/deal.II/include/grid/tria_quad.h
deleted file mode 100644 (file)
index 9dbf31b..0000000
+++ /dev/null
@@ -1,142 +0,0 @@
-//---------------------------------------------------------------------------
-//    $Id$
-//    Version: $Name$
-//
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
-//
-//    This file is subject to QPL and may not be  distributed
-//    without copyright and license information. Please refer
-//    to the file deal.II/doc/license.html for the  text  and
-//    further information on this license.
-//
-//---------------------------------------------------------------------------
-#ifndef __deal2__tria_quad_h
-#define __deal2__tria_quad_h
-
-
-#include <base/config.h>
-#include <base/exceptions.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-namespace internal
-{
-  namespace Triangulation
-  {
-    
-/**
- *   @p Quads denote the fundamental entities of triangulations in two dimensions
- *   and the boundaries of hexaeders in three dimensions. They are
- *   characterized by the (global) indices of the four bounding lines.
- *
- *   A quad itself has one index, as far as the topological part handled in
- *   the triangulation is concerned: the index in the level
- *   it belongs to. The level index is implicitly given by the position
- *   in the <tt>quads.quads</tt> list attached to the information of each level
- *   of the triangulation.
- *
- * @ingroup grid
- *   @author Wolfgang Bangerth, 1998
- */
-    class Quad
-    {
-      public:
-
-                                         /**
-                                          *  Construct a Quad with line
-                                          *  indices @p i0 through @p i3. By default,
-                                          *  indices are set to -1, i.e. an
-                                          *  invalid value.
-                                          *
-                                          * By convention, the four lines
-                                          * must be numbered as follows
-                                          *  .--3--.
-                                          *  |     |
-                                          *  0     1
-                                          *  |     |
-                                          *  .--2--.
-                                          */
-        Quad (const int i0 = -1,
-              const int i1 = -1,
-              const int i2 = -1,
-              const int i3 = -1);
-    
-                                         /**
-                                          *  Return the index of line @p i=0
-                                          *  through 3.
-                                          */
-        int line (const int i) const;
-    
-                                         /**
-                                          *  Set the index of the @p ith line to
-                                          *  @p index. @p i=0 through 3.
-                                          */
-        void set_line (const int i, const int index);
-    
-                                         /**
-                                          * Determine an estimate for the
-                                          * memory consumption (in bytes)
-                                          * of this object.
-                                          */
-        static unsigned int memory_consumption ();
-
-                                         /**
-                                          *  Exception
-                                          */ 
-        DeclException1 (ExcRange,
-                        int,
-                        << "Indices for the line number must be 0, 1, 2 or 3, "
-                        << "but you gave " << arg1); 
-      protected:
-        int lines[4];
-    };
-
-
-/*----------------------------- Inline Function: Quad ------------------------*/
-
-
-    inline
-    Quad::Quad (const int i0, const int i1, const int i2, const int i3)
-    {
-      lines[0] = i0;
-      lines[1] = i1;
-      lines[2] = i2;
-      lines[3] = i3;
-    }
-
-
-
-    inline
-    int Quad::line (const int i) const
-    {
-      Assert ((i>=0) && (i<4),
-              ExcRange(i));
-      return lines[i];
-    }
-
-
-
-    inline
-    void Quad::set_line (const int i, const int index)
-    {
-      Assert ((i>=0) && (i<4),
-              ExcRange(i));
-      lines[i] = index;
-    }
-
-
-
-    inline
-    unsigned int
-    Quad::memory_consumption ()
-    {
-      return sizeof(Quad);
-    }
-    
-  }
-  
-}
-
-DEAL_II_NAMESPACE_CLOSE
-
-#endif
index 8ccd9087f63397e8c22321d7507b682f5af00c44..f8bf1d9ce71c1b5f9762bca9ed253d9775cccb67 100644 (file)
@@ -778,16 +778,16 @@ struct QuadComparator
                                         // the repeated equality test
                                         // of the previous lines, but
                                         // I don't care at present
-       if ((q1.line(0) < q2.line(0))          ||
-           ((q1.line(0) == q2.line(0)) &&
-            (q1.line(1) <  q2.line(1)))       ||
-           ((q1.line(0) == q2.line(0)) &&
-            (q1.line(1) == q2.line(1)) &&
-            (q1.line(2) <  q2.line(2)))       ||
-           ((q1.line(0) == q2.line(0)) &&
-            (q1.line(1) == q2.line(1)) &&
-            (q1.line(2) == q2.line(2)) &&
-            (q1.line(3) <  q2.line(3))))
+       if ((q1.face(0) < q2.face(0))          ||
+           ((q1.face(0) == q2.face(0)) &&
+            (q1.face(1) <  q2.face(1)))       ||
+           ((q1.face(0) == q2.face(0)) &&
+            (q1.face(1) == q2.face(1)) &&
+            (q1.face(2) <  q2.face(2)))       ||
+           ((q1.face(0) == q2.face(0)) &&
+            (q1.face(1) == q2.face(1)) &&
+            (q1.face(2) == q2.face(2)) &&
+            (q1.face(3) <  q2.face(3))))
          return true;
        else
          return false;
@@ -1168,21 +1168,21 @@ Triangulation<3>::create_triangulation (const std::vector<Point<3> >    &v,
                                            // new one and instead
                                            // later set the
                                            // face_orientation flag
-          const internal::Triangulation::Quad
-           test_quad_1(quad.line(2), quad.line(3),
-                       quad.line(0), quad.line(1)),//face_orientation=false, face_flip=false, face_rotation=false
-           test_quad_2(quad.line(0), quad.line(1),
-                       quad.line(3), quad.line(2)),//face_orientation=false, face_flip=false, face_rotation=true
-           test_quad_3(quad.line(3), quad.line(2),
-                       quad.line(1), quad.line(0)),//face_orientation=false, face_flip=true,  face_rotation=false
-           test_quad_4(quad.line(1), quad.line(0),
-                       quad.line(2), quad.line(3)),//face_orientation=false, face_flip=true,  face_rotation=true
-           test_quad_5(quad.line(2), quad.line(3),
-                       quad.line(1), quad.line(0)),//face_orientation=true,  face_flip=false, face_rotation=true
-           test_quad_6(quad.line(1), quad.line(0),
-                       quad.line(3), quad.line(2)),//face_orientation=true,  face_flip=true,  face_rotation=false
-           test_quad_7(quad.line(3), quad.line(2),
-                       quad.line(0), quad.line(1));//face_orientation=true,  face_flip=true,  face_rotation=true
+          const internal::Triangulation::TriaObject<2>
+           test_quad_1(quad.face(2), quad.face(3),
+                       quad.face(0), quad.face(1)),//face_orientation=false, face_flip=false, face_rotation=false
+           test_quad_2(quad.face(0), quad.face(1),
+                       quad.face(3), quad.face(2)),//face_orientation=false, face_flip=false, face_rotation=true
+           test_quad_3(quad.face(3), quad.face(2),
+                       quad.face(1), quad.face(0)),//face_orientation=false, face_flip=true,  face_rotation=false
+           test_quad_4(quad.face(1), quad.face(0),
+                       quad.face(2), quad.face(3)),//face_orientation=false, face_flip=true,  face_rotation=true
+           test_quad_5(quad.face(2), quad.face(3),
+                       quad.face(1), quad.face(0)),//face_orientation=true,  face_flip=false, face_rotation=true
+           test_quad_6(quad.face(1), quad.face(0),
+                       quad.face(3), quad.face(2)),//face_orientation=true,  face_flip=true,  face_rotation=false
+           test_quad_7(quad.face(3), quad.face(2),
+                       quad.face(0), quad.face(1));//face_orientation=true,  face_flip=true,  face_rotation=true
           if (needed_quads.find (test_quad_1) == needed_quads.end() &&
              needed_quads.find (test_quad_2) == needed_quads.end() &&
              needed_quads.find (test_quad_3) == needed_quads.end() &&
@@ -1313,21 +1313,21 @@ Triangulation<3>::create_triangulation (const std::vector<Point<3> >    &v,
                                                   // then. construct all
                                                   // possibilities and check
                                                   // them one after the other
-                 const internal::Triangulation::Quad
-                   test_quad_1(quad.line(2), quad.line(3),
-                               quad.line(0), quad.line(1)),//face_orientation=false, face_flip=false, face_rotation=false
-                   test_quad_2(quad.line(0), quad.line(1),
-                               quad.line(3), quad.line(2)),//face_orientation=false, face_flip=false, face_rotation=true
-                   test_quad_3(quad.line(3), quad.line(2),
-                               quad.line(1), quad.line(0)),//face_orientation=false, face_flip=true,  face_rotation=false
-                   test_quad_4(quad.line(1), quad.line(0),
-                               quad.line(2), quad.line(3)),//face_orientation=false, face_flip=true,  face_rotation=true
-                   test_quad_5(quad.line(2), quad.line(3),
-                               quad.line(1), quad.line(0)),//face_orientation=true,  face_flip=false, face_rotation=true
-                   test_quad_6(quad.line(1), quad.line(0),
-                               quad.line(3), quad.line(2)),//face_orientation=true,  face_flip=true,  face_rotation=false
-                   test_quad_7(quad.line(3), quad.line(2),
-                               quad.line(0), quad.line(1));//face_orientation=true,  face_flip=true,  face_rotation=true
+                 const internal::Triangulation::TriaObject<2>
+                   test_quad_1(quad.face(2), quad.face(3),
+                               quad.face(0), quad.face(1)),//face_orientation=false, face_flip=false, face_rotation=false
+                   test_quad_2(quad.face(0), quad.face(1),
+                               quad.face(3), quad.face(2)),//face_orientation=false, face_flip=false, face_rotation=true
+                   test_quad_3(quad.face(3), quad.face(2),
+                               quad.face(1), quad.face(0)),//face_orientation=false, face_flip=true,  face_rotation=false
+                   test_quad_4(quad.face(1), quad.face(0),
+                               quad.face(2), quad.face(3)),//face_orientation=false, face_flip=true,  face_rotation=true
+                   test_quad_5(quad.face(2), quad.face(3),
+                               quad.face(1), quad.face(0)),//face_orientation=true,  face_flip=false, face_rotation=true
+                   test_quad_6(quad.face(1), quad.face(0),
+                               quad.face(3), quad.face(2)),//face_orientation=true,  face_flip=true,  face_rotation=false
+                   test_quad_7(quad.face(3), quad.face(2),
+                               quad.face(0), quad.face(1));//face_orientation=true,  face_flip=true,  face_rotation=true
                  if (needed_quads.find (test_quad_1) != needed_quads.end())
                    {
                      face_iterator[face] = needed_quads[test_quad_1].first;
@@ -1711,8 +1711,8 @@ Triangulation<3>::create_triangulation (const std::vector<Point<3> >    &v,
                                           // lexicographic ordering)
          for (unsigned int i=0; i<4; ++i)
            {
-             quad_compare_1.set_line(i,       line_counterclock[lex2cclock[i]]->index());
-             quad_compare_2.set_line((i+2)%4, line_counterclock[lex2cclock[i]]->index());
+             quad_compare_1.set_face(i,       line_counterclock[lex2cclock[i]]->index());
+             quad_compare_2.set_face((i+2)%4, line_counterclock[lex2cclock[i]]->index());
            }
              
          ++n_rotations;
index ec1b369fb604d60bb9c1d7211cd327f2473e887e..95b071aa09c0337df3621f7adde6628bec5cd22b 100644 (file)
@@ -50,7 +50,7 @@ template <int dim>
 int TriaObjectAccessor<1, dim>::vertex_index (const unsigned int i) const
 {
   Assert (i<2, ExcIndexRange(i,0,2));
-  return lines().cells[this->present_index].vertex (i);
+  return lines().cells[this->present_index].face (i);
 }
 
 

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.