]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove DoFHandlerType from hp::FEValues templates. 10559/head
authorDavid Wells <drwells@email.unc.edu>
Fri, 19 Jun 2020 18:36:36 +0000 (14:36 -0400)
committerDavid Wells <drwells@email.unc.edu>
Fri, 19 Jun 2020 18:36:36 +0000 (14:36 -0400)
include/deal.II/hp/fe_values.h
source/hp/fe_values.cc
source/hp/fe_values.inst.in

index c599df056e0b8039a7c75ce24b82249773360415..44fa1b968e98392c725d42dcac241856f5a3a18d 100644 (file)
@@ -18,6 +18,8 @@
 
 #include <deal.II/base/config.h>
 
+#include <deal.II/dofs/dof_handler.h>
+
 #include <deal.II/fe/fe.h>
 #include <deal.II/fe/fe_values.h>
 
@@ -318,7 +320,7 @@ namespace hp
      * passed to the constructor of this class with index given by
      * <code>cell-@>active_fe_index()</code>. Consequently, the
      * hp::FECollection argument given to this object should really be the
-     * same as that used in the construction of the hp::DoFHandler associated
+     * same as that used in the construction of the DoFHandler associated
      * with the present cell. On the other hand, if a value is given for this
      * argument, it overrides the choice of
      * <code>cell-@>active_fe_index()</code>.
@@ -347,19 +349,18 @@ namespace hp
      * all finite elements in an hp discretization), then this single mapping
      * is used unless a different value for this argument is specified.
      */
-    template <typename DoFHandlerType, bool lda>
+    template <bool lda>
     void
-    reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType, lda>> cell,
+    reinit(const TriaIterator<
+             DoFCellAccessor<dealii::DoFHandler<dim, spacedim>, lda>> cell,
            const unsigned int q_index       = numbers::invalid_unsigned_int,
            const unsigned int mapping_index = numbers::invalid_unsigned_int,
            const unsigned int fe_index      = numbers::invalid_unsigned_int);
 
     /**
-     * Like the previous function, but for non-hp iterators. The reason this
-     * (and the other non-hp iterator) function exists is so that one can use
-     * hp::FEValues not only for hp::DoFhandler objects, but for all sorts of
-     * DoFHandler objects, and triangulations not associated with DoFHandlers
-     * in general.
+     * Like the previous function, but for non-DoFHandler iterators. The reason
+     * this function exists is so that one can use hp::FEValues for
+     * Triangulation objects too.
      *
      * Since <code>cell-@>active_fe_index()</code> doesn't make sense for
      * triangulation iterators, this function chooses the zero-th finite
@@ -439,7 +440,7 @@ namespace hp
      * passed to the constructor of this class with index given by
      * <code>cell-@>active_fe_index()</code>. Consequently, the
      * hp::FECollection argument given to this object should really be the
-     * same as that used in the construction of the hp::DoFHandler associated
+     * same as that used in the construction of the DoFHandler associated
      * with the present cell. On the other hand, if a value is given for this
      * argument, it overrides the choice of
      * <code>cell-@>active_fe_index()</code>.
@@ -468,20 +469,19 @@ namespace hp
      * all finite elements in an hp discretization), then this single mapping
      * is used unless a different value for this argument is specified.
      */
-    template <typename DoFHandlerType, bool lda>
+    template <bool lda>
     void
-    reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType, lda>> cell,
-           const unsigned int                                       face_no,
+    reinit(const TriaIterator<
+             DoFCellAccessor<dealii::DoFHandler<dim, spacedim>, lda>> cell,
+           const unsigned int                                         face_no,
            const unsigned int q_index       = numbers::invalid_unsigned_int,
            const unsigned int mapping_index = numbers::invalid_unsigned_int,
            const unsigned int fe_index      = numbers::invalid_unsigned_int);
 
     /**
-     * Like the previous function, but for non-hp iterators. The reason this
-     * (and the other non-hp iterator) function exists is so that one can use
-     * hp::FEValues not only for hp::DoFhandler objects, but for all sorts of
-     * DoFHandler objects, and triangulations not associated with DoFHandlers
-     * in general.
+     * Like the previous function, but for non-DoFHandler iterators. The reason
+     * this function exists is so that one can use this class for
+     * Triangulation objects too.
      *
      * Since <code>cell-@>active_fe_index()</code> doesn't make sense for
      * triangulation iterators, this function chooses the zero-th finite
@@ -566,24 +566,23 @@ namespace hp
      * all finite elements in an hp discretization), then this single mapping
      * is used unless a different value for this argument is specified.
      */
-    template <typename DoFHandlerType, bool lda>
+    template <bool lda>
     void
-    reinit(const TriaIterator<DoFCellAccessor<DoFHandlerType, lda>> cell,
-           const unsigned int                                       face_no,
-           const unsigned int                                       subface_no,
+    reinit(const TriaIterator<
+             DoFCellAccessor<dealii::DoFHandler<dim, spacedim>, lda>> cell,
+           const unsigned int                                         face_no,
+           const unsigned int subface_no,
            const unsigned int q_index       = numbers::invalid_unsigned_int,
            const unsigned int mapping_index = numbers::invalid_unsigned_int,
            const unsigned int fe_index      = numbers::invalid_unsigned_int);
 
     /**
-     * Like the previous function, but for non-hp iterators. The reason this
-     * (and the other non-hp iterator) function exists is so that one can use
-     * hp::FEValues not only for hp::DoFhandler objects, but for all sorts of
-     * DoFHandler objects, and triangulations not associated with DoFHandlers
-     * in general.
+     * Like the previous function, but for non-DoFHandler iterators. The reason
+     * this function exists is so that one can use this class for
+     * Triangulation objects too.
      *
      * Since <code>cell-@>active_fe_index()</code> doesn't make sense for
-     * triangulation iterators, this function chooses the zero-th finite
+     * Triangulation iterators, this function chooses the zero-th finite
      * element, mapping, and quadrature object from the relevant constructions
      * passed to the constructor of this object. The only exception is if you
      * specify a value different from the default value for any of these last
index 8bc477b42e1d2c1604842453c3b04b6ccc02cf3e..233b474ab5e0d3538dcc6925e12999a231b8a2ce 100644 (file)
@@ -223,13 +223,14 @@ namespace hp
 
 
   template <int dim, int spacedim>
-  template <typename DoFHandlerType, bool lda>
+  template <bool lda>
   void
   FEValues<dim, spacedim>::reinit(
-    const TriaIterator<DoFCellAccessor<DoFHandlerType, lda>> cell,
-    const unsigned int                                       q_index,
-    const unsigned int                                       mapping_index,
-    const unsigned int                                       fe_index)
+    const TriaIterator<DoFCellAccessor<dealii::DoFHandler<dim, spacedim>, lda>>
+                       cell,
+    const unsigned int q_index,
+    const unsigned int mapping_index,
+    const unsigned int fe_index)
   {
     // determine which indices we
     // should actually use
@@ -334,14 +335,15 @@ namespace hp
 
 
   template <int dim, int spacedim>
-  template <typename DoFHandlerType, bool lda>
+  template <bool lda>
   void
   FEFaceValues<dim, spacedim>::reinit(
-    const TriaIterator<DoFCellAccessor<DoFHandlerType, lda>> cell,
-    const unsigned int                                       face_no,
-    const unsigned int                                       q_index,
-    const unsigned int                                       mapping_index,
-    const unsigned int                                       fe_index)
+    const TriaIterator<DoFCellAccessor<dealii::DoFHandler<dim, spacedim>, lda>>
+                       cell,
+    const unsigned int face_no,
+    const unsigned int q_index,
+    const unsigned int mapping_index,
+    const unsigned int fe_index)
   {
     // determine which indices we
     // should actually use
@@ -447,15 +449,16 @@ namespace hp
 
 
   template <int dim, int spacedim>
-  template <typename DoFHandlerType, bool lda>
+  template <bool lda>
   void
   FESubfaceValues<dim, spacedim>::reinit(
-    const TriaIterator<DoFCellAccessor<DoFHandlerType, lda>> cell,
-    const unsigned int                                       face_no,
-    const unsigned int                                       subface_no,
-    const unsigned int                                       q_index,
-    const unsigned int                                       mapping_index,
-    const unsigned int                                       fe_index)
+    const TriaIterator<DoFCellAccessor<dealii::DoFHandler<dim, spacedim>, lda>>
+                       cell,
+    const unsigned int face_no,
+    const unsigned int subface_no,
+    const unsigned int q_index,
+    const unsigned int mapping_index,
+    const unsigned int fe_index)
   {
     // determine which indices we
     // should actually use
index 5b38d728245d1b88caff04bde2171b4e874e6549..dac83a8e8014c92bae276de7c006e7ee87e3b1a7 100644 (file)
@@ -82,8 +82,7 @@ for (deal_II_dimension : DIMENSIONS)
 #endif
   }
 
-for (dof_handler : DOFHANDLER_TEMPLATE; deal_II_dimension : DIMENSIONS;
-     deal_II_space_dimension : SPACE_DIMENSIONS;
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS;
      lda : BOOL)
   {
     namespace hp
@@ -93,7 +92,7 @@ for (dof_handler : DOFHANDLER_TEMPLATE; deal_II_dimension : DIMENSIONS;
       template void
       FEValues<deal_II_dimension, deal_II_space_dimension>::reinit(
         TriaIterator<DoFCellAccessor<
-          dealii::dof_handler<deal_II_dimension, deal_II_space_dimension>,
+          dealii::DoFHandler<deal_II_dimension, deal_II_space_dimension>,
           lda>>,
         unsigned int,
         unsigned int,
@@ -101,7 +100,7 @@ for (dof_handler : DOFHANDLER_TEMPLATE; deal_II_dimension : DIMENSIONS;
       template void
       FEFaceValues<deal_II_dimension, deal_II_space_dimension>::reinit(
         TriaIterator<DoFCellAccessor<
-          dealii::dof_handler<deal_II_dimension, deal_II_space_dimension>,
+          dealii::DoFHandler<deal_II_dimension, deal_II_space_dimension>,
           lda>>,
         unsigned int,
         unsigned int,
@@ -110,7 +109,7 @@ for (dof_handler : DOFHANDLER_TEMPLATE; deal_II_dimension : DIMENSIONS;
       template void
       FESubfaceValues<deal_II_dimension, deal_II_space_dimension>::reinit(
         TriaIterator<DoFCellAccessor<
-          dealii::dof_handler<deal_II_dimension, deal_II_space_dimension>,
+          dealii::DoFHandler<deal_II_dimension, deal_II_space_dimension>,
           lda>>,
         unsigned int,
         unsigned int,

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.