]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update the documentation. Rename template arguments.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 30 Aug 2013 12:57:12 +0000 (12:57 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 30 Aug 2013 12:57:12 +0000 (12:57 +0000)
git-svn-id: https://svn.dealii.org/trunk@30537 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/hp/dof_objects.h

index 1a035b80bfeaa43ea73f2832c1f93019f28a452e..59ec451d53289626ffb5e694c2780ce4d987f0b0 100644 (file)
@@ -32,7 +32,9 @@ namespace internal
 
     /**
      * Store the indices of the degrees of freedom which are located on
-     * objects of dimension @p dim.
+     * objects of dimension @p structdim < dim, i.e., for faces or edges
+     * of cells. This is opposed to the internal::hp::DoFLevels class
+     * that stores the DoF indices on cells.
      *
      * The things we store here is very similar to what is stored in the
      * internal::DoFHandler::DoFObjects classes (see there for more
@@ -58,16 +60,12 @@ namespace internal
      *
      * <h4>Multiple data sets per object</h4>
      *
-     * If an object corresponds to a cell, the global dof indices of this
-     * cell are stored at the location indicated above in sequential
-     * order.
-     *
-     * However, if two adjacent cells use different finite elements, then
+     * If two adjacent cells use different finite elements, then
      * the face that they share needs to store DoF indices for both
      * involved finite elements. While faces therefore have to have at
-     * most two sets of DoF indices, it is easy to see that vertices for
-     * example can have as many sets of DoF indices associated with them
-     * as there are adjacent cells, and the same holds for lines in 3d.
+     * most two sets of DoF indices, it is easy to see that edges and
+     * vertices can have as many sets of DoF indices associated with them
+     * as there are adjacent cells.
      *
      * Consequently, for objects that have a lower dimensionality than
      * cells, we have to store a map from the finite element index to the
@@ -87,7 +85,7 @@ namespace internal
      * Access to this kind of data, as well as the distinction between
      * cells and objects of lower dimensionality are encoded in the
      * accessor functions, DoFObjects::set_dof_index() and
-     * DoFLevel::get_dof_index() They are able to traverse this
+     * DoFLevel::get_dof_index(). They are able to traverse this
      * list and pick out or set a DoF index given the finite element index
      * and its location within the set of DoFs corresponding to this
      * finite element.
@@ -96,7 +94,7 @@ namespace internal
      * @ingroup hp
      * @author Tobias Leicht, 2006
      */
-    template <int dim>
+    template <int structdim>
     class DoFObjects
     {
     public:
@@ -142,9 +140,9 @@ namespace internal
        * class template for more
        * information.
        */
-      template <int dimm, int spacedim>
+      template <int dim, int spacedim>
       void
-      set_dof_index (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+      set_dof_index (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                      const unsigned int               obj_index,
                      const unsigned int               fe_index,
                      const unsigned int               local_index,
@@ -174,9 +172,9 @@ namespace internal
        * class template for more
        * information.
        */
-      template <int dimm, int spacedim>
+      template <int dim, int spacedim>
       types::global_dof_index
-      get_dof_index (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+      get_dof_index (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                      const unsigned int               obj_index,
                      const unsigned int               fe_index,
                      const unsigned int               local_index,
@@ -205,9 +203,9 @@ namespace internal
        * been distributed and zero
        * is returned.
        */
-      template <int dimm, int spacedim>
+      template <int dim, int spacedim>
       unsigned int
-      n_active_fe_indices (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+      n_active_fe_indices (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                            const unsigned int               obj_index) const;
 
       /**
@@ -215,9 +213,9 @@ namespace internal
        * n-th active finite element
        * on this object.
        */
-      template <int dimm, int spacedim>
+      template <int dim, int spacedim>
       types::global_dof_index
-      nth_active_fe_index (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+      nth_active_fe_index (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                            const unsigned int               obj_level,
                            const unsigned int               obj_index,
                            const unsigned int               n) const;
@@ -228,9 +226,9 @@ namespace internal
        * used on the present
        * object or not.
        */
-      template <int dimm, int spacedim>
+      template <int dim, int spacedim>
       bool
-      fe_index_is_active (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+      fe_index_is_active (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                           const unsigned int               obj_index,
                           const unsigned int               fe_index,
                           const unsigned int               obj_level) const;
@@ -246,18 +244,18 @@ namespace internal
 
     // --------------------- inline and template functions ------------------
 
-    template <int dim>
-    template <int dimm, int spacedim>
+    template <int structdim>
+    template <int dim, int spacedim>
     inline
     types::global_dof_index
-    DoFObjects<dim>::
-    get_dof_index (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+    DoFObjects<structdim>::
+    get_dof_index (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                    const unsigned int                obj_index,
                    const unsigned int                fe_index,
                    const unsigned int                local_index,
                    const unsigned int                obj_level) const
     {
-      Assert ((fe_index != dealii::hp::DoFHandler<dimm,spacedim>::default_fe_index),
+      Assert ((fe_index != dealii::hp::DoFHandler<dim,spacedim>::default_fe_index),
               ExcMessage ("You need to specify a FE index when working "
                           "with hp DoFHandlers"));
       Assert (&dof_handler != 0,
@@ -268,10 +266,10 @@ namespace internal
       Assert (fe_index < dof_handler.get_fe().size(),
               ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
       Assert (local_index <
-              dof_handler.get_fe()[fe_index].template n_dofs_per_object<dim>(),
+              dof_handler.get_fe()[fe_index].template n_dofs_per_object<structdim>(),
               ExcIndexRange(local_index, 0,
                             dof_handler.get_fe()[fe_index]
-                            .template n_dofs_per_object<dim>()));
+                            .template n_dofs_per_object<structdim>()));
       Assert (obj_index < dof_offsets.size(),
               ExcIndexRange (obj_index, 0, dof_offsets.size()));
 
@@ -283,7 +281,7 @@ namespace internal
                           "information for an object on which no such "
                           "information is available"));
 
-      Assert (dim<dimm, ExcMessage ("This object can not be used for cells."));
+      Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
 
       // there may be multiple finite elements associated with
       // this object. hop along the list of index sets until we
@@ -301,25 +299,25 @@ namespace internal
          else
            pointer += static_cast<types::global_dof_index>(
              dof_handler.get_fe()[*pointer]
-             .template n_dofs_per_object<dim>() + 1);
+             .template n_dofs_per_object<structdim>() + 1);
        }
     }
 
 
 
-    template <int dim>
-    template <int dimm, int spacedim>
+    template <int structdim>
+    template <int dim, int spacedim>
     inline
     void
-    DoFObjects<dim>::
-    set_dof_index (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+    DoFObjects<structdim>::
+    set_dof_index (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                    const unsigned int                obj_index,
                    const unsigned int                fe_index,
                    const unsigned int                local_index,
                    const types::global_dof_index     global_index,
                    const unsigned int                obj_level)
     {
-      Assert ((fe_index != dealii::hp::DoFHandler<dimm,spacedim>::default_fe_index),
+      Assert ((fe_index != dealii::hp::DoFHandler<dim,spacedim>::default_fe_index),
               ExcMessage ("You need to specify a FE index when working "
                           "with hp DoFHandlers"));
       Assert (&dof_handler != 0,
@@ -330,10 +328,10 @@ namespace internal
       Assert (fe_index < dof_handler.get_fe().size(),
               ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
       Assert (local_index <
-              dof_handler.get_fe()[fe_index].template n_dofs_per_object<dim>(),
+              dof_handler.get_fe()[fe_index].template n_dofs_per_object<structdim>(),
               ExcIndexRange(local_index, 0,
                             dof_handler.get_fe()[fe_index]
-                            .template n_dofs_per_object<dim>()));
+                            .template n_dofs_per_object<structdim>()));
       Assert (obj_index < dof_offsets.size(),
               ExcIndexRange (obj_index, 0, dof_offsets.size()));
 
@@ -345,7 +343,7 @@ namespace internal
                           "information for an object on which no such "
                           "information is available"));
 
-      Assert (dim<dimm, ExcMessage ("This object can not be used for cells."));
+      Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
 
       // there may be multiple finite elements associated with
       // this object.  hop along the list of index sets until we
@@ -365,21 +363,20 @@ namespace internal
            }
          else
            pointer += dof_handler.get_fe()[*pointer]
-                      .template n_dofs_per_object<dim>() + 1;
+                      .template n_dofs_per_object<structdim>() + 1;
        }
     }
 
 
 
-    template <int dim>
-    template <int dimm, int spacedim>
+    template <int structdim>
+    template <int dim, int spacedim>
     inline
     unsigned int
-    DoFObjects<dim>::
-    n_active_fe_indices (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+    DoFObjects<structdim>::
+    n_active_fe_indices (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                          const unsigned int                obj_index) const
     {
-      Assert (dim <= dimm, ExcInternalError());
       Assert (&dof_handler != 0,
               ExcMessage ("No DoFHandler is specified for this iterator"));
       Assert (&dof_handler.get_fe() != 0,
@@ -394,7 +391,7 @@ namespace internal
       if (dof_offsets[obj_index] == numbers::invalid_dof_index)
         return 0;
 
-      Assert (dim<dimm, ExcMessage ("This object can not be used for cells."));
+      Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
 
       // there may be multiple finite elements associated with this
       // object. hop along the list of index sets until we find the
@@ -413,24 +410,23 @@ namespace internal
            {
              ++counter;
              pointer += dof_handler.get_fe()[*pointer]
-                        .template n_dofs_per_object<dim>() + 1;
+                        .template n_dofs_per_object<structdim>() + 1;
            }
        }
     }
 
 
 
-    template <int dim>
-    template <int dimm, int spacedim>
+    template <int structdim>
+    template <int dim, int spacedim>
     inline
     types::global_dof_index
-    DoFObjects<dim>::
-    nth_active_fe_index (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+    DoFObjects<structdim>::
+    nth_active_fe_index (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                          const unsigned int                obj_level,
                          const unsigned int                obj_index,
                          const unsigned int                n) const
     {
-      Assert (dim <= dimm, ExcInternalError());
       Assert (&dof_handler != 0,
               ExcMessage ("No DoFHandler is specified for this iterator"));
       Assert (&dof_handler.get_fe() != 0,
@@ -447,7 +443,7 @@ namespace internal
                           "information for an object on which no such "
                           "information is available"));
 
-      Assert (dim<dimm, ExcMessage ("This object can not be used for cells."));
+      Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
 
       Assert (n < n_active_fe_indices(dof_handler, obj_index),
              ExcIndexRange (n, 0,
@@ -476,18 +472,18 @@ namespace internal
 
          ++counter;
          pointer += dof_handler.get_fe()[fe_index]
-                    .template n_dofs_per_object<dim>() + 1;
+                    .template n_dofs_per_object<structdim>() + 1;
        }
     }
 
 
 
-    template <int dim>
-    template <int dimm, int spacedim>
+    template <int structdim>
+    template <int dim, int spacedim>
     inline
     bool
-    DoFObjects<dim>::
-    fe_index_is_active (const dealii::hp::DoFHandler<dimm,spacedim> &dof_handler,
+    DoFObjects<structdim>::
+    fe_index_is_active (const dealii::hp::DoFHandler<dim,spacedim> &dof_handler,
                         const unsigned int                obj_index,
                         const unsigned int                fe_index,
                         const unsigned int                obj_level) const
@@ -499,7 +495,7 @@ namespace internal
                           "this DoFHandler"));
       Assert (obj_index < dof_offsets.size(),
               ExcIndexRange (obj_index, 0, static_cast<unsigned int>(dof_offsets.size())));
-      Assert ((fe_index != dealii::hp::DoFHandler<dimm,spacedim>::default_fe_index),
+      Assert ((fe_index != dealii::hp::DoFHandler<dim,spacedim>::default_fe_index),
               ExcMessage ("You need to specify a FE index when working "
                           "with hp DoFHandlers"));
       Assert (fe_index < dof_handler.get_fe().size(),
@@ -513,7 +509,7 @@ namespace internal
                           "information for an object on which no such "
                           "information is available"));
 
-      Assert (dim<dimm, ExcMessage ("This object can not be used for cells."));
+      Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
 
       // there may be multiple finite elements associated with
       // this object. hop along the list of index sets until we
@@ -532,7 +528,7 @@ namespace internal
          else
            pointer += static_cast<types::global_dof_index>(
              dof_handler.get_fe()[*pointer]
-             .template n_dofs_per_object<dim>()+1);
+             .template n_dofs_per_object<structdim>()+1);
        }
     }
 

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.