]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Also introduce MGDoFAccessor<0,1,spacedim>.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 18 Dec 2010 05:00:38 +0000 (05:00 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 18 Dec 2010 05:00:38 +0000 (05:00 +0000)
git-svn-id: https://svn.dealii.org/trunk@23000 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/dofs/dof_accessor.h
deal.II/include/deal.II/multigrid/mg_dof_accessor.h
deal.II/include/deal.II/multigrid/mg_dof_iterator_selector.h
deal.II/source/multigrid/mg_dof_accessor.cc
deal.II/source/multigrid/mg_dof_accessor.inst.in

index b74f47048ed55958e6fb72df9debe000fb34efbd..d33b447bacced43e9525ad78ae3e80562c6f7a7a 100644 (file)
@@ -859,7 +859,7 @@ class DoFAccessor<0,DH<1,spacedim> > : public TriaAccessor<0,1,spacedim>
                                      * do anything useful here and so
                                      * may not actually be called.
                                      */
-    DoFAccessor (const Triangulation<1,spacedim> * =  0,
+    DoFAccessor (const Triangulation<1,spacedim> *,
                 const int = 0,
                 const int = 0,
                 const DH<1,spacedim> * = 0);
index 42a14d8942b0b98f5ba713269afd5cc08ae80b6b..438add06a8aaa48df0d98f1cef7732438f758115 100644 (file)
@@ -338,6 +338,14 @@ class MGDoFAccessor : public internal::MGDoFAccessor::Inheritance<structdim,dim,
                                      * @}
                                      */
 
+                                    /**
+                                     * Return a handle on the
+                                     * DoFHandler object which we
+                                     * are using.
+                                     */
+    const MGDoFHandler<dim,spacedim> &
+    get_mg_dof_handler () const;
+
                                     /**
                                      * Return an iterator pointing to
                                      * the the parent.
@@ -413,6 +421,161 @@ class MGDoFAccessor : public internal::MGDoFAccessor::Inheritance<structdim,dim,
 };
 
 
+
+/**
+ * Specialization of the general MGDoFAccessor class template for the
+ * case of zero-dimensional objects (a vertex) that are the face of a
+ * one-dimensional cell in spacedim space dimensions. Since vertices
+ * function differently than general faces, this class does a few
+ * things differently than the general template, but the interface
+ * should look the same.
+ *
+ * @author Wolfgang Bangerth, 2010
+ */
+template <int spacedim>
+class MGDoFAccessor<0,1,spacedim> : public internal::MGDoFAccessor::Inheritance<0,1,spacedim>::BaseClass
+{
+  public:
+                                    /**
+                                     * Declare a typedef to the base
+                                     * class to make accessing some
+                                     * of the exception classes
+                                     * simpler.
+                                     */
+    typedef
+    typename internal::MGDoFAccessor::Inheritance<0,1,spacedim>::BaseClass
+    BaseClass;
+
+                                    /**
+                                     * A typedef necessary for the
+                                     * iterator classes.
+                                     */
+    typedef MGDoFHandler<1,spacedim> AccessorData;
+
+                                    /**
+                                     * The dof handler type used by
+                                     * this accessor.
+                                     */
+    typedef MGDoFHandler<1, spacedim> Container;
+
+                                    /**
+                                     * @name Constructors
+                                     */
+                                    /**
+                                     * @{
+                                     */
+
+                                    /**
+                                     * Default constructor. Provides
+                                     * an accessor that can't be
+                                     * used.
+                                     */
+    MGDoFAccessor ();
+
+                                    /**
+                                     * Constructor to be used if the
+                                     * object here refers to a vertex
+                                     * of a one-dimensional
+                                     * triangulation, i.e. a face of
+                                     * the triangulation.
+                                     *
+                                     * Since there is no mapping from
+                                     * vertices to cells, an accessor
+                                     * object for a point has no way
+                                     * to figure out whether it is at
+                                     * the boundary of the domain or
+                                     * not. Consequently, the second
+                                     * argument must be passed by the
+                                     * object that generates this
+                                     * accessor -- e.g. a 1d cell
+                                     * that can figure out whether
+                                     * its left or right vertex are
+                                     * at the boundary.
+                                     *
+                                     * The third argument is the
+                                     * global index of the vertex we
+                                     * point to.
+                                     *
+                                     * The fourth argument is a
+                                     * pointer to the MGDoFHandler
+                                     * object.
+                                     *
+                                     * This iterator can only be
+                                     * called for one-dimensional
+                                     * triangulations.
+                                     */
+    MGDoFAccessor (const Triangulation<1,spacedim> * tria,
+                  const typename TriaAccessor<0,1,spacedim>::VertexKind vertex_kind,
+                  const unsigned int    vertex_index,
+                  const MGDoFHandler<1,spacedim> * dof_handler);
+
+                                    /**
+                                     * Constructor. This constructor
+                                     * exists in order to maintain
+                                     * interface compatibility with
+                                     * the other accessor
+                                     * classes. However, it doesn't
+                                     * do anything useful here and so
+                                     * may not actually be called.
+                                     */
+    MGDoFAccessor (const Triangulation<1,spacedim> * =  0,
+                  const int = 0,
+                  const int = 0,
+                  const MGDoFHandler<1,spacedim> * = 0);
+
+                                    /**
+                                     * Conversion constructor. This
+                                     * constructor exists to make certain
+                                     * constructs simpler to write in
+                                     * dimension independent code. For
+                                     * example, it allows assigning a face
+                                     * iterator to a line iterator, an
+                                     * operation that is useful in 2d but
+                                     * doesn't make any sense in 3d. The
+                                     * constructor here exists for the
+                                     * purpose of making the code conform to
+                                     * C++ but it will unconditionally abort;
+                                     * in other words, assigning a face
+                                     * iterator to a line iterator is better
+                                     * put into an if-statement that checks
+                                     * that the dimension is two, and assign
+                                     * to a quad iterator in 3d (an operator
+                                     * that, without this constructor would
+                                     * be illegal if we happen to compile for
+                                     * 2d).
+                                     */
+    template <int structdim2, int dim2, int spacedim2>
+    MGDoFAccessor (const InvalidAccessor<structdim2,dim2,spacedim2> &);
+
+                                    /**
+                                     * Another conversion operator
+                                     * between objects that don't
+                                     * make sense, just like the
+                                     * previous one.
+                                     */
+    template <int structdim2, int dim2, int spacedim2>
+    MGDoFAccessor (const MGDoFAccessor<structdim2, dim2, spacedim2> &);
+
+
+                                    /**
+                                     * Another conversion operator
+                                     * between objects that don't
+                                     * make sense, just like the
+                                     * previous one.
+                                     */
+    template <int dim2, int spacedim2>
+    MGDoFAccessor (const MGDoFCellAccessor<dim2, spacedim2> &);
+
+  protected:
+                                    /**
+                                     * Store the address of the @p MGDoFHandler object
+                                     * to be accessed.
+                                     */
+    MGDoFHandler<1,spacedim> *mg_dof_handler;
+};
+
+
+
 /* -------------------------------------------------------------------------- */
 
 
index bcbdc8fbfc4f6b23803cb3bf88891e296232207e..716a864cbb49667376a15dd237c6475c86474f72 100644 (file)
@@ -45,8 +45,8 @@ namespace internal
     {
       public:
        typedef dealii::MGDoFCellAccessor<1,spacedim> CellAccessor;
-       typedef InvalidAccessor<0,1,spacedim> FaceAccessor;
-       
+       typedef dealii::MGDoFAccessor<0,1,spacedim> FaceAccessor;
+
        typedef TriaRawIterator   <CellAccessor> raw_line_iterator;
         typedef TriaIterator      <CellAccessor> line_iterator;
         typedef TriaActiveIterator<CellAccessor> active_line_iterator;
@@ -80,11 +80,11 @@ namespace internal
       public:
        typedef dealii::MGDoFCellAccessor<2,spacedim> CellAccessor;
        typedef dealii::MGDoFAccessor<1,2,spacedim> FaceAccessor;
-       
+
         typedef TriaRawIterator   <FaceAccessor> raw_line_iterator;
         typedef TriaIterator      <FaceAccessor> line_iterator;
         typedef TriaActiveIterator<FaceAccessor> active_line_iterator;
-    
+
         typedef TriaRawIterator   <CellAccessor> raw_quad_iterator;
         typedef TriaIterator      <CellAccessor> quad_iterator;
         typedef TriaActiveIterator<CellAccessor> active_quad_iterator;
@@ -99,7 +99,7 @@ namespace internal
 
         typedef raw_line_iterator    raw_face_iterator;
         typedef line_iterator        face_iterator;
-        typedef active_line_iterator active_face_iterator;    
+        typedef active_line_iterator active_face_iterator;
     };
 
 
@@ -133,7 +133,7 @@ namespace internal
 
         typedef raw_quad_iterator    raw_face_iterator;
         typedef quad_iterator        face_iterator;
-        typedef active_quad_iterator active_face_iterator;    
+        typedef active_quad_iterator active_face_iterator;
     };
   }
 }
index b5efa2acfc82dce646e06bf28d5363e6f10197ca..61007de2d5e28d25e820ff1bcb49ad4aee2aebc9 100644 (file)
@@ -65,6 +65,15 @@ MGDoFAccessor<structdim, dim, spacedim>::set_mg_dof_handler (MGDoFHandler<dim,sp
 
 
 
+template <int structdim, int dim, int spacedim>
+const MGDoFHandler<dim,spacedim> &
+MGDoFAccessor<structdim, dim, spacedim>::get_mg_dof_handler () const
+{
+  return *mg_dof_handler;
+}
+
+
+
 template <int structdim,int dim, int spacedim>
 unsigned int
 MGDoFAccessor<structdim, dim, spacedim>::mg_vertex_dof_index (const int level,
@@ -454,9 +463,78 @@ get_mg_dof_values (const int level,
 
 
 
+/*------------------------- Functions: MGDoFAccessor<0,1,spacedim> -----------------------*/
+
+template <int spacedim>
+inline
+MGDoFAccessor<0,1,spacedim>::MGDoFAccessor ()
+               :
+               BaseClass ()
+{
+  Assert (false, typename BaseClass::ExcInvalidObject());
+}
+
+
+
+template <int spacedim>
+inline
+MGDoFAccessor<0,1,spacedim>::
+MGDoFAccessor (const Triangulation<1,spacedim> *tria,
+              const typename TriaAccessor<0,1,spacedim>::VertexKind vertex_kind,
+              const unsigned int    vertex_index,
+              const MGDoFHandler<1,spacedim> * dof_handler)
+               :
+               BaseClass (tria,
+                          vertex_kind,
+                          vertex_index),
+                mg_dof_handler(const_cast<MGDoFHandler<1,spacedim>*>(dof_handler))
+{}
+
+
+
+template <int spacedim>
+inline
+MGDoFAccessor<0,1,spacedim>::
+MGDoFAccessor (const Triangulation<1,spacedim> *,
+              const int                 ,
+              const int                 ,
+              const MGDoFHandler<1,spacedim>     *)
+               :
+                mg_dof_handler(0)
+{
+  Assert (false,
+         ExcMessage ("This constructor can not be called for face iterators in 1d."));
+}
+
+
+
+template <int spacedim>
+template <int structdim2, int dim2, int spacedim2>
+MGDoFAccessor<0,1,spacedim>::MGDoFAccessor (const InvalidAccessor<structdim2,dim2,spacedim2> &)
+{
+  Assert (false, typename BaseClass::ExcInvalidObject());
+}
 
 
 
+template <int spacedim>
+template <int structdim2, int dim2, int spacedim2>
+inline
+MGDoFAccessor<0,1,spacedim>::MGDoFAccessor (const MGDoFAccessor<structdim2, dim2, spacedim2> &)
+{
+  Assert (false, typename BaseClass::ExcInvalidObject());
+}
+
+
+
+template <int spacedim>
+template <int dim2, int spacedim2>
+inline
+MGDoFAccessor<0,1,spacedim>::MGDoFAccessor (const MGDoFCellAccessor<dim2, spacedim2> &)
+{
+  Assert (false, typename BaseClass::ExcInvalidObject());
+}
+
 
 
 /*------------------------- Functions: MGDoFCellAccessor -----------------------*/
@@ -538,25 +616,65 @@ MGDoFCellAccessor<dim, spacedim>::parent () const
 }
 
 
+namespace internal
+{
+  namespace MGDoFCellAccessor
+  {
+    template <int spacedim>
+    inline
+    typename internal::MGDoFHandler::Iterators<1,spacedim>::face_iterator
+    get_face (const dealii::MGDoFCellAccessor<1,spacedim> &cell,
+             const unsigned int i,
+             const internal::int2type<1>)
+    {
+      dealii::MGDoFAccessor<0,1,spacedim>
+       a (&cell.get_triangulation(),
+          ((i == 0) && cell.at_boundary(0)
+           ?
+           dealii::TriaAccessor<0, 1, spacedim>::left_vertex
+           :
+           ((i == 1) && cell.at_boundary(1)
+            ?
+            dealii::TriaAccessor<0, 1, spacedim>::right_vertex
+            :
+            dealii::TriaAccessor<0, 1, spacedim>::interior_vertex)),
+          cell.vertex_index(i),
+          &cell.get_mg_dof_handler());
+      return typename internal::MGDoFHandler::Iterators<1,spacedim>::face_iterator(a);
+    }
+
+
+    template <int spacedim>
+    inline
+    typename internal::MGDoFHandler::Iterators<2,spacedim>::face_iterator
+    get_face (const dealii::MGDoFCellAccessor<2,spacedim> &cell,
+             const unsigned int i,
+             const internal::int2type<2>)
+    {
+      return cell.line(i);
+    }
+
+
+    template <int spacedim>
+    inline
+    typename internal::MGDoFHandler::Iterators<3,spacedim>::face_iterator
+    get_face (const dealii::MGDoFCellAccessor<3,spacedim> &cell,
+             const unsigned int i,
+             const internal::int2type<3>)
+    {
+      return cell.quad(i);
+    }
+  }
+}
+
+
 template <int dim, int spacedim>
 typename internal::MGDoFHandler::Iterators<dim,spacedim>::face_iterator
 MGDoFCellAccessor<dim,spacedim>::face (const unsigned int i) const
 {
-  switch (dim)
-    {
-      case 1:
-           Assert (false, ExcImpossibleInDim(1));
-           return
-             typename internal::MGDoFHandler::Iterators<dim,spacedim>::face_iterator();
-      case 2:
-           return this->line(i);
-      case 3:
-           return this->quad(i);
-      default:
-           Assert (false, ExcNotImplemented());
-           return
-             typename internal::MGDoFHandler::Iterators<dim,spacedim>::face_iterator();
-    }
+  Assert (i<GeometryInfo<dim>::faces_per_cell, ExcIndexRange (i, 0, GeometryInfo<dim>::faces_per_cell));
+
+  return internal::MGDoFCellAccessor::get_face (*this, i, internal::int2type<dim>());
 }
 
 
index 4888f59abbf1faac24219cfa12722c2bfc5ce82d..4fdcec74de5d83acfa926e7b568f9da73bc5d06b 100644 (file)
@@ -82,6 +82,9 @@ for (deal_II_dimension : DIMENSIONS)
 
 
 #if deal_II_dimension == 1
+  template class MGDoFAccessor<0, 1, 1>;
+  template class MGDoFAccessor<0, 1, 2>;
+
   template class MGDoFAccessor<1, 1, 1>;
   template class MGDoFAccessor<1, 1, 2>;
 

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.