]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Write necessary documentation. Avoid one #include file.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 30 Aug 2013 03:00:19 +0000 (03:00 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 30 Aug 2013 03:00:19 +0000 (03:00 +0000)
git-svn-id: https://svn.dealii.org/trunk@30532 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/hp/dof_levels.h

index 7bd350cdf6c49f69d662ad69798e710ecf08d488..154ebe9c9d0de3417cdcf80653f231b27e12de1b 100644 (file)
@@ -19,7 +19,7 @@
 
 
 #include <deal.II/base/config.h>
-#include <deal.II/hp/dof_objects.h>
+#include <deal.II/base/exceptions.h>
 
 #include <vector>
 
@@ -30,76 +30,13 @@ namespace internal
   namespace hp
   {
     /**
-     * Store the indices of the degrees of freedom that are located on
-     * objects of dimension @p structdim.
-     *
-     * The things we store here are very similar to what is stored in the
-     * internal::DoFHandler::DoFLevel class hierarchy (see there for more
-     * information, in particular on the layout of the class hierarchy,
-     * and the use of file names). There are two main
-     * differences, discussed in the following subsections. In addition to
-     * the data already stored by the internal::DoFHandler::DoFLevel
-     * classes, we also have to store which finite element each cell
-     * uses.
-     *
-     *
-     * <h4>Offset computations</h4>
-     *
-     * For hp methods, not all cells may use the same finite element, and
-     * it is consequently more complicated to determine where the DoF
-     * indices for a given line, quad, or hex are stored. As described in
-     * the documentation of the internal::DoFHandler::DoFLevel class, we
-     * can compute the location of the first line DoF, for example, by
-     * calculating the offset as <code>line_index *
-     * dof_handler.get_fe().dofs_per_line</code>. This of course doesn't
-     * work any more if different lines may have different numbers of
-     * degrees of freedom associated with them. Consequently, rather than
-     * using this simple multiplication, each of the lines.dofs, quads.dofs,
-     * and hexes.dofs arrays has an associated array lines.dof_offsets,
-     * quads.dof_offsets, and hexes.dof_offsets. The data corresponding to a
-     * line then starts at index <code>line_dof_offsets[line_index]</code>
-     * within the <code>line_dofs</code> array.
-     *
-     *
-     * <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
-     * 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.
-     *
-     * Consequently, for objects that have a lower dimensionality than
-     * cells, we have to store a map from the finite element index to the
-     * set of DoF indices associated. Since real sets are typically very
-     * inefficient to store, and since most of the time we expect the
-     * number of individual keys to be small (frequently, adjacent cells
-     * will have the same finite element, and only a single entry will
-     * exist in the map), what we do is instead to store a linked list. In
-     * this format, the first entry starting at position
-     * <code>lines.dofs[lines.dof_offsets[line_index]]</code> will denote
-     * the finite element index of the set of DoF indices following; after
-     * this set, we will store the finite element index of the second set
-     * followed by the corresponding DoF indices; and so on. Finally, when
-     * all finite element indices adjacent to this object have been
-     * covered, we write a -1 to indicate the end of the list.
-     *
-     * Access to this kind of data, as well as the distinction between
-     * cells and objects of lower dimensionality are encoded in the
-     * accessor functions, DoFLevel::set_dof_index() and
-     * 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.
-     *
-     *
-     * @ingroup hp
-     * @author Wolfgang Bangerth, 1998, 2006, Oliver Kayser-Herold 2003.
+     * This is the class that stores the degrees of freedom on cells in a hp
+     * hierarchy. Compared to faces and edges, the task here is simple since
+     * each cell can only have a single active finite element index. Consequently,
+     * all we need is one long array with DoF indices and one array of offsets
+     * where each cell's indices start within the array of indices. This is in
+     * contrast to the DoFObjects class where each face or edge may have more than
+     * one associated finite element with corresponding degrees of freedom.
      */
     template <int dim>
     class DoFLevel
@@ -312,7 +249,8 @@ namespace internal
     fe_index_is_active (const unsigned int                obj_index,
                         const unsigned int                fe_index) const
     {
-      Assert ((fe_index != dealii::hp::DoFHandler<1,1>::default_fe_index),
+      const unsigned int invalid_fe_index = numbers::invalid_unsigned_int;
+      Assert ((fe_index != invalid_fe_index),
               ExcMessage ("You need to specify a FE index when working "
                           "with hp DoFHandlers"));
 

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.