From ef8e396d1487a452033add20ad33ef8c502c4272 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Thu, 12 Nov 2020 19:43:09 -0700 Subject: [PATCH] Deprecate hp::DoFHandler. 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 | 95 +++++++++++++++++++++++ include/deal.II/hp/dof_handler.h | 40 +--------- include/deal.II/matrix_free/matrix_free.h | 64 +++++++++------ source/hp/dof_handler.cc | 16 ---- 4 files changed, 138 insertions(+), 77 deletions(-) diff --git a/include/deal.II/dofs/dof_handler.h b/include/deal.II/dofs/dof_handler.h index bd297f1b01..91228ebebc 100644 --- a/include/deal.II/dofs/dof_handler.h +++ b/include/deal.II/dofs/dof_handler.h @@ -213,6 +213,101 @@ namespace parallel * using the same kind of finite element before re-loading data from the * serialization archive. * + * + *

hp-adaptive finite element methods

+ * + * 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 + * active_fe_index 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". + * + * + *

Active FE indices and their behavior under mesh refinement

+ * + * 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. + * + * + *

Active FE indices and parallel meshes

+ * + * 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 cell-@>set_active_fe_index(...). + * 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 query 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 diff --git a/include/deal.II/hp/dof_handler.h b/include/deal.II/hp/dof_handler.h index cf383640b9..62de381d43 100644 --- a/include/deal.II/hp/dof_handler.h +++ b/include/deal.II/hp/dof_handler.h @@ -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 - class DoFHandler : public dealii::DoFHandler + class DEAL_II_DEPRECATED DoFHandler : public dealii::DoFHandler { - 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&` and a `bool` instead. - */ - DEAL_II_DEPRECATED - DoFHandler(const Triangulation &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::DoFHandler; + using dealii::DoFHandler::operator=; }; } // namespace hp diff --git a/include/deal.II/matrix_free/matrix_free.h b/include/deal.II/matrix_free/matrix_free.h index 82bc2adab1..165e5a23ee 100644 --- a/include/deal.II/matrix_free/matrix_free.h +++ b/include/deal.II/matrix_free/matrix_free.h @@ -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 + template DEAL_II_DEPRECATED void reinit(const Mapping & mapping, - const std::vector *> & dof_handler_hp, + const std::vector & dof_handler, const std::vector *> &constraint, const std::vector & 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 + template DEAL_II_DEPRECATED void - reinit(const std::vector *> & dof_handler_hp, + reinit(const std::vector & dof_handler, const std::vector *> &constraint, const std::vector & 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 + template DEAL_II_DEPRECATED void reinit(const Mapping & mapping, - const std::vector *> & dof_handler_hp, + const std::vector & dof_handler, const std::vector *> &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 + template DEAL_II_DEPRECATED void - reinit(const std::vector *> & dof_handler_hp, + reinit(const std::vector & dof_handler, const std::vector *> &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::active_cell_iterator @@ -2922,19 +2924,22 @@ MatrixFree::reinit( template -template +template void MatrixFree::reinit( const Mapping & mapping, - const std::vector *> & dof_handler_hp, + const std::vector & dof_handler, const std::vector *> &constraint, const std::vector & quad, const AdditionalData & additional_data) { + static_assert(dim == DoFHandlerType::dimension, + "Dimension dim not equal to DoFHandlerType::dimension."); + std::vector *> 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::reinit( template -template +template void MatrixFree::reinit( - const std::vector *> & dof_handler_hp, + const std::vector & dof_handler, const std::vector *> &constraint, const std::vector & quad, const AdditionalData & additional_data) { + static_assert(dim == DoFHandlerType::dimension, + "Dimension dim not equal to DoFHandlerType::dimension."); + std::vector *> 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::reinit( template -template +template void MatrixFree::reinit( const Mapping & mapping, - const std::vector *> & dof_handler_hp, + const std::vector & dof_handler, const std::vector *> &constraint, const QuadratureType & quad, const AdditionalData & additional_data) { + static_assert(dim == DoFHandlerType::dimension, + "Dimension dim not equal to DoFHandlerType::dimension."); + std::vector *> 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::reinit( template -template +template void MatrixFree::reinit( - const std::vector *> & dof_handler_hp, + const std::vector & dof_handler, const std::vector *> &constraint, const QuadratureType & quad, const AdditionalData & additional_data) { + static_assert(dim == DoFHandlerType::dimension, + "Dimension dim not equal to DoFHandlerType::dimension."); + std::vector *> 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); diff --git a/source/hp/dof_handler.cc b/source/hp/dof_handler.cc index 43ee32b005..598e6476c7 100644 --- a/source/hp/dof_handler.cc +++ b/source/hp/dof_handler.cc @@ -17,22 +17,6 @@ DEAL_II_NAMESPACE_OPEN -namespace hp -{ - template - DoFHandler::DoFHandler() - : dealii::DoFHandler() - {} - - - - template - DoFHandler::DoFHandler( - const Triangulation &tria) - : dealii::DoFHandler(tria) - {} - -} // namespace hp /*-------------- Explicit Instantiations -------------------------------*/ #include "dof_handler.inst" -- 2.39.5