]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Deprecate hp::DoFHandler.
authorMarc Fehling <mafehling.git@gmail.com>
Fri, 13 Nov 2020 02:43:09 +0000 (19:43 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Mon, 16 Nov 2020 03:20:54 +0000 (20:20 -0700)
Added missing deprecate flags in matrix_free.h.

Added DEAL_II_*_EXTRA_DIAGNOSTICS pragmas in matrix_free.h
to avoid deprecation warnings during CI.

include/deal.II/dofs/dof_handler.h
include/deal.II/hp/dof_handler.h
include/deal.II/matrix_free/matrix_free.h
source/hp/dof_handler.cc

index bd297f1b012a8e16d95a47f37647a311a64e53ae..91228ebebc1f47dfcabebc7ad11584775ee7dd60 100644 (file)
@@ -213,6 +213,101 @@ namespace parallel
  * using the same kind of finite element before re-loading data from the
  * serialization archive.
  *
+ *
+ * <h3>hp-adaptive finite element methods</h3>
+ *
+ * Instead of only using one particular FiniteElement on all cells, this class
+ * also allows for an enumeration of degrees of freedom on different finite
+ * elements on every cells. To this end, one assigns an
+ * <code>active_fe_index</code> to every cell that indicates which element
+ * within a collection of finite elements (represented by an object of type
+ * hp::FECollection) is the one that lives on this cell. The class then
+ * enumerates the degree of freedom associated with these finite elements on
+ * each cell of a triangulation and, if possible, identifies degrees of
+ * freedom at the interfaces of cells if they match. If neighboring cells
+ * have degrees of freedom along the common interface that do not immediate
+ * match (for example, if you have $Q_2$ and $Q_3$ elements meeting at a
+ * common face), then one needs to compute constraints to ensure that the
+ * resulting finite element space on the mesh remains conforming.
+ *
+ * The whole process of working with objects of this type is explained in
+ * step-27. Many of the algorithms this class implements are described in
+ * the
+ * @ref hp_paper "hp paper".
+ *
+ *
+ * <h3>Active FE indices and their behavior under mesh refinement</h3>
+ *
+ * The typical workflow for using this class is to create a mesh, assign an
+ * active FE index to every active cell, call DoFHandler::distribute_dofs(),
+ * and then assemble a linear system and solve a problem on this finite element
+ * space.
+ *
+ * Active FE indices will be automatically transferred during mesh adaptation
+ * from the old to the new mesh. Future FE indices are meant to determine the
+ * active FE index after mesh adaptation, and are used to prepare data on the
+ * old mesh for the new one. If no future FE index is specified, the finite
+ * element prevails.
+ *
+ * In particular, the following rules apply during adaptation:
+ * - Upon mesh refinement, child cells inherit the future FE index of
+ *   the parent.
+ * - When coarsening cells, the (now active) parent cell will be assigned
+ *   a future FE index that is determined from its (no longer active)
+ *   children, following the FiniteElementDomination logic: Out of the set of
+ *   elements previously assigned to the former children, we choose the one
+ *   dominated by all children for the parent cell. If none was found, we pick
+ *   the most dominant element in the whole collection that is dominated by
+ *   all former children. See hp::FECollection::find_dominated_fe_extended()
+ *   for further information on this topic.
+ *
+ * Strategies for automatic hp-adaptation which will set future FE indices based
+ * on criteria are available in the hp::Refinement namespace.
+ *
+ *
+ * <h3>Active FE indices and parallel meshes</h3>
+ *
+ * When this class is used with either a parallel::shared::Triangulation
+ * or a parallel::distributed::Triangulation, you can only set active
+ * FE indices on cells that are locally owned,
+ * using a call such as <code>cell-@>set_active_fe_index(...)</code>.
+ * On the other hand, setting the active FE index on ghost
+ * or artificial cells is not allowed.
+ *
+ * Ghost cells do acquire the information what element
+ * is active on them, however: whenever you call DoFHandler::distribute_dofs(),
+ * all processors that participate in the parallel mesh exchange information in
+ * such a way that the active FE index on ghost cells equals the active FE index
+ * that was set on that processor that owned that particular ghost cell.
+ * Consequently, one can <i>query</i> the @p active_fe_index on ghost
+ * cells, just not set it by hand.
+ *
+ * On artificial cells, no information is available about the
+ * @p active_fe_index used there. That's because we don't even know
+ * whether these cells exist at all, and even if they did, the
+ * current processor does not know anything specific about them.
+ * See
+ * @ref GlossArtificialCell "the glossary entry on artificial cells"
+ * for more information.
+ *
+ * During refinement and coarsening, information about the @p active_fe_index
+ * of each cell will be automatically transferred.
+ *
+ * However, using a parallel::distributed::Triangulation with a DoFHandler
+ * in hp-mode requires additional attention during serialization, since no
+ * information on active FE indices will be automatically transferred. This
+ * has to be done manually using the
+ * prepare_for_serialization_of_active_fe_indices() and
+ * deserialize_active_fe_indices() functions. The former has to be called
+ * before parallel::distributed::Triangulation::save() is invoked, and the
+ * latter needs to be run after parallel::distributed::Triangulation::load().
+ * If further data will be attached to the triangulation via the
+ * parallel::distributed::CellDataTransfer,
+ * parallel::distributed::SolutionTransfer, or Particles::ParticleHandler
+ * classes, all corresponding preparation and deserialization function calls
+ * need to happen in the same order. Consult the documentation of
+ * parallel::distributed::SolutionTransfer for more information.
+ *
  * @ingroup dofs
  */
 template <int dim, int spacedim = dim>
index cf383640b98e48d96d41085a9b7b4af31ec8472c..62de381d4359d60bb8b74139c19acd38a4631140 100644 (file)
@@ -119,45 +119,13 @@ namespace hp
    * @ingroup dofs
    * @ingroup hp
    *
-   * @note Task is delegated to the base class dealii::DoFHandler.
+   * @deprecated The basic dealii::DoFHandler is capable of hp-adaptation now.
    */
   template <int dim, int spacedim = dim>
-  class DoFHandler : public dealii::DoFHandler<dim, spacedim>
+  class DEAL_II_DEPRECATED DoFHandler : public dealii::DoFHandler<dim, spacedim>
   {
-  public:
-    /**
-     * Default Constructor.
-     *
-     * @deprecated Use the overload taking a `bool` instead.
-     */
-    DEAL_II_DEPRECATED
-    DoFHandler();
-
-    /**
-     * Constructor. Take @p tria as the triangulation to work on.
-     *
-     * @deprecated Use the overload taking a
-     * `const Triangulation<dim, spacedim>&` and a `bool` instead.
-     */
-    DEAL_II_DEPRECATED
-    DoFHandler(const Triangulation<dim, spacedim> &tria);
-
-    /**
-     * Copy constructor. DoFHandler objects are large and expensive.
-     * They should not be copied, in particular not by accident, but
-     * rather deliberately constructed. As a consequence, this constructor
-     * is explicitly removed from the interface of this class.
-     */
-    DoFHandler(const DoFHandler &) = delete;
-
-    /**
-     * Copy operator. DoFHandler objects are large and expensive.
-     * They should not be copied, in particular not by accident, but
-     * rather deliberately constructed. As a consequence, this operator
-     * is explicitly removed from the interface of this class.
-     */
-    DoFHandler &
-    operator=(const DoFHandler &) = delete;
+    using dealii::DoFHandler<dim, spacedim>::DoFHandler;
+    using dealii::DoFHandler<dim, spacedim>::operator=;
   };
 
 } // namespace hp
index 82bc2adab177a100e3f986db16c83406b436ed14..165e5a23eecd7ff47d4c347259fcc61f5b4aa420 100644 (file)
@@ -631,14 +631,14 @@ public:
          const AdditionalData &additional_data = AdditionalData());
 
   /**
-   * Initializes the data structures. Same as above, but  using hp::DoFHandlers.
+   * Initializes the data structures. Same as above, but  using DoFHandlerType.
    *
    * @deprecated Use the overload taking a DoFHandler object instead.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename DoFHandlerType>
   DEAL_II_DEPRECATED void
   reinit(const Mapping<dim> &                                   mapping,
-         const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+         const std::vector<const DoFHandlerType *> &            dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const std::vector<QuadratureType> &                    quad,
          const AdditionalData &additional_data = AdditionalData());
@@ -657,13 +657,13 @@ public:
          const AdditionalData &additional_data = AdditionalData());
 
   /**
-   * Initializes the data structures. Same as above, but  using hp::DoFHandlers.
+   * Initializes the data structures. Same as above, but  using DoFHandlerType.
    *
    * @deprecated Use the overload taking a DoFHandler object instead.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename DoFHandlerType>
   DEAL_II_DEPRECATED void
-  reinit(const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+  reinit(const std::vector<const DoFHandlerType *> &            dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const std::vector<QuadratureType> &                    quad,
          const AdditionalData &additional_data = AdditionalData());
@@ -684,14 +684,14 @@ public:
          const AdditionalData &additional_data = AdditionalData());
 
   /**
-   * Initializes the data structures. Same as above, but  using hp::DoFHandlers.
+   * Initializes the data structures. Same as above, but  using DoFHandlerType.
    *
    * @deprecated Use the overload taking a DoFHandler object instead.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename DoFHandlerType>
   DEAL_II_DEPRECATED void
   reinit(const Mapping<dim> &                                   mapping,
-         const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+         const std::vector<const DoFHandlerType *> &            dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const QuadratureType &                                 quad,
          const AdditionalData &additional_data = AdditionalData());
@@ -710,13 +710,13 @@ public:
          const AdditionalData &additional_data = AdditionalData());
 
   /**
-   * Initializes the data structures. Same as above, but  using hp::DoFHandlers.
+   * Initializes the data structures. Same as above, but  using DoFHandlerType.
    *
    * @deprecated Use the overload taking a DoFHandler object instead.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename DoFHandlerType>
   DEAL_II_DEPRECATED void
-  reinit(const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+  reinit(const std::vector<const DoFHandlerType *> &            dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const QuadratureType &                                 quad,
          const AdditionalData &additional_data = AdditionalData());
@@ -1749,6 +1749,8 @@ public:
                     const unsigned int fe_component = 0) const;
 
   /**
+   * @copydoc MatrixFree::get_cell_iterator()
+   *
    * @deprecated Use get_cell_iterator() instead.
    */
   DEAL_II_DEPRECATED typename DoFHandler<dim>::active_cell_iterator
@@ -2922,19 +2924,22 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename DoFHandlerType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
   const Mapping<dim> &                                   mapping,
-  const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+  const std::vector<const DoFHandlerType *> &            dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const std::vector<QuadratureType> &                    quad,
   const AdditionalData &                                 additional_data)
 {
+  static_assert(dim == DoFHandlerType::dimension,
+                "Dimension dim not equal to DoFHandlerType::dimension.");
+
   std::vector<const DoFHandler<dim> *> dof_handlers;
 
-  for (const auto dof_handler : dof_handler_hp)
-    dof_handlers.push_back(dof_handler);
+  for (const auto dh : dof_handler)
+    dof_handlers.push_back(dh);
 
   this->reinit(mapping, dof_handlers, constraint, quad, additional_data);
 }
@@ -2967,17 +2972,20 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename DoFHandlerType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+  const std::vector<const DoFHandlerType *> &            dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const std::vector<QuadratureType> &                    quad,
   const AdditionalData &                                 additional_data)
 {
+  static_assert(dim == DoFHandlerType::dimension,
+                "Dimension dim not equal to DoFHandlerType::dimension.");
+
   std::vector<const DoFHandler<dim> *> dof_handlers;
 
-  for (const auto dof_handler : dof_handler_hp)
+  for (const auto dh : dof_handler)
     dof_handlers.push_back(dof_handler);
 
   this->reinit(dof_handlers, constraint, quad, additional_data);
@@ -3039,18 +3047,21 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename DoFHandlerType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
   const Mapping<dim> &                                   mapping,
-  const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+  const std::vector<const DoFHandlerType *> &            dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const QuadratureType &                                 quad,
   const AdditionalData &                                 additional_data)
 {
+  static_assert(dim == DoFHandlerType::dimension,
+                "Dimension dim not equal to DoFHandlerType::dimension.");
+
   std::vector<const DoFHandler<dim> *> dof_handlers;
 
-  for (const auto dof_handler : dof_handler_hp)
+  for (const auto dh : dof_handler)
     dof_handlers.push_back(dof_handler);
 
   this->reinit(mapping, dof_handlers, constraint, quad, additional_data);
@@ -3059,17 +3070,20 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename DoFHandlerType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const std::vector<const hp::DoFHandler<dim> *> &       dof_handler_hp,
+  const std::vector<const DoFHandlerType *> &            dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const QuadratureType &                                 quad,
   const AdditionalData &                                 additional_data)
 {
+  static_assert(dim == DoFHandlerType::dimension,
+                "Dimension dim not equal to DoFHandlerType::dimension.");
+
   std::vector<const DoFHandler<dim> *> dof_handlers;
 
-  for (const auto dof_handler : dof_handler_hp)
+  for (const auto dh : dof_handler)
     dof_handlers.push_back(dof_handler);
 
   this->reinit(dof_handlers, constraint, quad, additional_data);
index 43ee32b005adde7c6cde0fa5dbb8c568c3f0ef43..598e6476c70bfe7acd1edd346ed0181b4e760998 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
-namespace hp
-{
-  template <int dim, int spacedim>
-  DoFHandler<dim, spacedim>::DoFHandler()
-    : dealii::DoFHandler<dim, spacedim>()
-  {}
-
-
-
-  template <int dim, int spacedim>
-  DoFHandler<dim, spacedim>::DoFHandler(
-    const Triangulation<dim, spacedim> &tria)
-    : dealii::DoFHandler<dim, spacedim>(tria)
-  {}
-
-} // namespace hp
 
 /*-------------- Explicit Instantiations -------------------------------*/
 #include "dof_handler.inst"

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.