From 41fa88691d0200428ddf65472acb27b3955d9bd5 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Mon, 2 Nov 2020 18:56:29 -0700 Subject: [PATCH] Rework DoFHandler connections to underlying Triangulation. --- include/deal.II/dofs/dof_accessor.templates.h | 2 +- include/deal.II/dofs/dof_handler.h | 71 ++- source/dofs/dof_handler.cc | 416 ++++++++---------- source/hp/fe_collection.cc | 48 +- 4 files changed, 230 insertions(+), 307 deletions(-) diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index 40f1038ad0..54a5023b9f 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -2943,7 +2943,7 @@ DoFCellAccessor:: Assert(!future_fe_indices_children.empty(), ExcInternalError()); const unsigned int future_fe_index = - this->dof_handler->get_fe_collection().find_dominated_fe_extended( + this->dof_handler->fe_collection.find_dominated_fe_extended( future_fe_indices_children, /*codim=*/0); Assert(future_fe_index != numbers::invalid_unsigned_int, diff --git a/include/deal.II/dofs/dof_handler.h b/include/deal.II/dofs/dof_handler.h index 6cb329aa62..bd297f1b01 100644 --- a/include/deal.II/dofs/dof_handler.h +++ b/include/deal.II/dofs/dof_handler.h @@ -1587,35 +1587,33 @@ private: connect_to_triangulation_signals(); /** - * Create default tables for the active_fe_indices. - * They are initialized with a zero - * indicator, meaning that fe[0] is going to be used by default. This - * method is called before refinement and while setting the finite elements - * via set_fe(). It ensures each cell has a valid active_fe_index. + * Create default tables for the active and future fe_indices. + * + * Active indices are initialized with a zero indicator, meaning that fe[0] is + * going to be used by default. Future indices are initialized with an invalid + * indicator, meaning that no p-adaptation is scheduled by default. + * + * This method is called upon construction and whenever the underlying + * triangulation gets created. This ensures that each cell has a valid active + * and future fe_index. */ void create_active_fe_table(); /** - * A function that will be triggered through a triangulation - * signal just before the triangulation is modified. + * Update tables for active and future fe_indices. * - * The function that stores the active_fe_flags of all cells that will - * be refined or coarsened before the refinement happens, so that - * they can be set again after refinement. - */ - void - pre_refinement_action(); - - /** - * A function that will be triggered through a triangulation - * signal just after the triangulation is modified. + * Whenever the underlying triangulation changes (either by adaptation or + * deserialization), active and future fe index tables will be adjusted to the + * current structure of the triangulation. Missing values of active and future + * indices will be initialized with their defaults (see + * create_active_fe_table()). * - * The function that restores the active_fe_flags of all cells that - * were refined. + * This method is called post refinement and post deserialization. This + * ensures that each cell has a valid active and future fe_index. */ void - post_refinement_action(); + update_active_fe_table(); /** * A function that will be triggered through a triangulation @@ -1627,18 +1625,7 @@ private: * they can be set again after refinement. */ void - pre_active_fe_index_transfer(); - - /** - * A function that will be triggered through a triangulation - * signal just before the associated parallel::distributed::Triangulation is - * modified. - * - * The function that stores all active_fe_indices on locally owned cells for - * distribution over all participating processors. - */ - void - pre_distributed_active_fe_index_transfer(); + pre_transfer_action(); /** * A function that will be triggered through a triangulation @@ -1649,29 +1636,29 @@ private: * were refined or coarsened. */ void - post_active_fe_index_transfer(); + post_transfer_action(); /** * A function that will be triggered through a triangulation - * signal just after the associated parallel::distributed::Triangulation is + * signal just before the associated parallel::distributed::Triangulation is * modified. * - * The function that restores all active_fe_indices on locally owned cells - * that have been communicated. + * The function that stores all active_fe_indices on locally owned cells for + * distribution over all participating processors. */ void - post_distributed_active_fe_index_transfer(); + pre_distributed_transfer_action(); /** * A function that will be triggered through a triangulation - * signal just after the associated parallel::distributed::Triangulation has - * been saved. + * signal just after the associated parallel::distributed::Triangulation is + * modified. * - * The function frees all memory related to the transfer of - * active_fe_indices. + * The function that restores all active_fe_indices on locally owned cells + * that have been communicated. */ void - post_distributed_serialization_of_active_fe_indices(); + post_distributed_transfer_action(); // Make accessor objects friends. diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index bf02ac40c6..d2e755e4e1 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -3041,16 +3041,11 @@ void DoFHandler::connect_to_triangulation_signals() { // connect functions to signals of the underlying triangulation - this->tria_listeners.push_back(this->tria->signals.pre_refinement.connect( - [this]() { this->pre_refinement_action(); })); - this->tria_listeners.push_back(this->tria->signals.post_refinement.connect( - [this]() { this->post_refinement_action(); })); this->tria_listeners.push_back(this->tria->signals.create.connect( - [this]() { this->post_refinement_action(); })); + [this]() { this->create_active_fe_table(); })); - // decide whether we need a sequential or a parallel shared/distributed - // policy and attach corresponding callback functions dealing with the - // transfer of active_fe_indices + // attach corresponding callback functions dealing with the transfer of + // active_fe_indices depending on the type of triangulation if (dynamic_cast< const dealii::parallel::DistributedTriangulationBase *>( &this->get_triangulation())) @@ -3063,24 +3058,26 @@ DoFHandler::connect_to_triangulation_signals() })); this->tria_listeners.push_back( this->tria->signals.pre_distributed_repartition.connect( - [this]() { this->pre_distributed_active_fe_index_transfer(); })); + [this]() { this->pre_distributed_transfer_action(); })); this->tria_listeners.push_back( this->tria->signals.post_distributed_repartition.connect( - [this] { this->post_distributed_active_fe_index_transfer(); })); + [this] { this->post_distributed_transfer_action(); })); // refinement signals this->tria_listeners.push_back( this->tria->signals.pre_distributed_refinement.connect( - [this]() { this->pre_distributed_active_fe_index_transfer(); })); + [this]() { this->pre_distributed_transfer_action(); })); this->tria_listeners.push_back( this->tria->signals.post_distributed_refinement.connect( - [this]() { this->post_distributed_active_fe_index_transfer(); })); + [this]() { this->post_distributed_transfer_action(); })); // serialization signals this->tria_listeners.push_back( - this->tria->signals.post_distributed_save.connect([this]() { - this->post_distributed_serialization_of_active_fe_indices(); - })); + this->tria->signals.post_distributed_save.connect( + [this]() { this->active_fe_index_transfer.reset(); })); + this->tria_listeners.push_back( + this->tria->signals.post_distributed_load.connect( + [this]() { this->update_active_fe_table(); })); } else if (dynamic_cast< const dealii::parallel::shared::Triangulation *>( @@ -3095,19 +3092,19 @@ DoFHandler::connect_to_triangulation_signals() // refinement signals this->tria_listeners.push_back(this->tria->signals.pre_refinement.connect( - [this] { this->pre_active_fe_index_transfer(); })); + [this] { this->pre_transfer_action(); })); this->tria_listeners.push_back( this->tria->signals.post_refinement.connect( - [this] { this->post_active_fe_index_transfer(); })); + [this] { this->post_transfer_action(); })); } else { // refinement signals this->tria_listeners.push_back(this->tria->signals.pre_refinement.connect( - [this] { this->pre_active_fe_index_transfer(); })); + [this] { this->pre_transfer_action(); })); this->tria_listeners.push_back( this->tria->signals.post_refinement.connect( - [this] { this->post_active_fe_index_transfer(); })); + [this] { this->post_transfer_action(); })); } } @@ -3169,16 +3166,7 @@ DoFHandler::create_active_fe_table() template void -DoFHandler::pre_refinement_action() -{ - create_active_fe_table(); -} - - - -template -void -DoFHandler::post_refinement_action() +DoFHandler::update_active_fe_table() { // // Normally only one level is added, but if this Triangulation // // is created by copy_triangulation, it can be more than one level. @@ -3226,27 +3214,21 @@ DoFHandler::post_refinement_action() template void -DoFHandler::pre_active_fe_index_transfer() +DoFHandler::pre_transfer_action() { - // Finite elements need to be assigned to each cell by calling - // distribute_dofs() first to make this functionality available. - if (this->fe_collection.size() > 0) - { - Assert(this->active_fe_index_transfer == nullptr, ExcInternalError()); + Assert(this->active_fe_index_transfer == nullptr, ExcInternalError()); - this->active_fe_index_transfer = - std::make_unique(); + this->active_fe_index_transfer = std::make_unique(); - dealii::internal::hp::DoFHandlerImplementation::Implementation:: - collect_fe_indices_on_cells_to_be_refined(*this); - } + dealii::internal::hp::DoFHandlerImplementation::Implementation:: + collect_fe_indices_on_cells_to_be_refined(*this); } template void -DoFHandler::pre_distributed_active_fe_index_transfer() +DoFHandler::pre_distributed_transfer_action() { #ifndef DEAL_II_WITH_P4EST Assert(false, @@ -3261,57 +3243,51 @@ DoFHandler::pre_distributed_active_fe_index_transfer() &this->get_triangulation()) != nullptr), ExcNotImplemented()); - // Finite elements need to be assigned to each cell by calling - // distribute_dofs() first to make this functionality available. - if (fe_collection.size() > 0) - { - Assert(active_fe_index_transfer == nullptr, ExcInternalError()); - - active_fe_index_transfer = std::make_unique(); - - // If we work on a p::d::Triangulation, we have to transfer all - // active_fe_indices since ownership of cells may change. We will - // use our p::d::CellDataTransfer member to achieve this. Further, - // we prepare the values in such a way that they will correspond to - // the active_fe_indices on the new mesh. - - // Gather all current future_fe_indices. - active_fe_index_transfer->active_fe_indices.resize( - get_triangulation().n_active_cells(), numbers::invalid_unsigned_int); - - for (const auto &cell : active_cell_iterators()) - if (cell->is_locally_owned()) - active_fe_index_transfer - ->active_fe_indices[cell->active_cell_index()] = - cell->future_fe_index(); - - // Create transfer object and attach to it. - const auto *distributed_tria = dynamic_cast< - const parallel::distributed::Triangulation *>( - &this->get_triangulation()); - - active_fe_index_transfer->cell_data_transfer = std::make_unique< - parallel::distributed:: - CellDataTransfer>>( - *distributed_tria, - /*transfer_variable_size_data=*/false, - /*refinement_strategy=*/ - &dealii::AdaptationStrategies::Refinement:: - preserve, - /*coarsening_strategy=*/ - [this]( - const typename Triangulation::cell_iterator &parent, - const std::vector &children_fe_indices) - -> unsigned int { - return dealii::internal::hp::DoFHandlerImplementation:: - Implementation::determine_fe_from_children( - parent, children_fe_indices, fe_collection); - }); - - active_fe_index_transfer->cell_data_transfer - ->prepare_for_coarsening_and_refinement( - active_fe_index_transfer->active_fe_indices); - } + Assert(active_fe_index_transfer == nullptr, ExcInternalError()); + + active_fe_index_transfer = std::make_unique(); + + // If we work on a p::d::Triangulation, we have to transfer all + // active_fe_indices since ownership of cells may change. We will + // use our p::d::CellDataTransfer member to achieve this. Further, + // we prepare the values in such a way that they will correspond to + // the active_fe_indices on the new mesh. + + // Gather all current future_fe_indices. + active_fe_index_transfer->active_fe_indices.resize( + get_triangulation().n_active_cells(), numbers::invalid_unsigned_int); + + for (const auto &cell : active_cell_iterators()) + if (cell->is_locally_owned()) + active_fe_index_transfer->active_fe_indices[cell->active_cell_index()] = + cell->future_fe_index(); + + // Create transfer object and attach to it. + const auto *distributed_tria = + dynamic_cast *>( + &this->get_triangulation()); + + active_fe_index_transfer->cell_data_transfer = std::make_unique< + parallel::distributed:: + CellDataTransfer>>( + *distributed_tria, + /*transfer_variable_size_data=*/false, + /*refinement_strategy=*/ + &dealii::AdaptationStrategies::Refinement:: + preserve, + /*coarsening_strategy=*/ + [this](const typename Triangulation::cell_iterator &parent, + const std::vector &children_fe_indices) + -> unsigned int { + return dealii::internal::hp::DoFHandlerImplementation::Implementation:: + determine_fe_from_children(parent, + children_fe_indices, + fe_collection); + }); + + active_fe_index_transfer->cell_data_transfer + ->prepare_for_coarsening_and_refinement( + active_fe_index_transfer->active_fe_indices); #endif } @@ -3319,61 +3295,54 @@ DoFHandler::pre_distributed_active_fe_index_transfer() template void -DoFHandler::post_active_fe_index_transfer() +DoFHandler::post_transfer_action() { - // Finite elements need to be assigned to each cell by calling - // distribute_dofs() first to make this functionality available. - if (this->fe_collection.size() > 0) - { - Assert(this->active_fe_index_transfer != nullptr, ExcInternalError()); + update_active_fe_table(); - dealii::internal::hp::DoFHandlerImplementation::Implementation:: - distribute_fe_indices_on_refined_cells(*this); + Assert(this->active_fe_index_transfer != nullptr, ExcInternalError()); - // We have to distribute the information about active_fe_indices - // of all cells (including the artificial ones) on all processors, - // if a parallel::shared::Triangulation has been used. - dealii::internal::hp::DoFHandlerImplementation::Implementation:: - communicate_active_fe_indices(*this); + dealii::internal::hp::DoFHandlerImplementation::Implementation:: + distribute_fe_indices_on_refined_cells(*this); - // Free memory. - this->active_fe_index_transfer.reset(); - } + // We have to distribute the information about active_fe_indices + // of all cells (including the artificial ones) on all processors, + // if a parallel::shared::Triangulation has been used. + dealii::internal::hp::DoFHandlerImplementation::Implementation:: + communicate_active_fe_indices(*this); + + // Free memory. + this->active_fe_index_transfer.reset(); } template void -DoFHandler::post_distributed_active_fe_index_transfer() +DoFHandler::post_distributed_transfer_action() { #ifndef DEAL_II_WITH_P4EST Assert(false, ExcInternalError()); #else - // Finite elements need to be assigned to each cell by calling - // distribute_dofs() first to make this functionality available. - if (this->fe_collection.size() > 0) - { - Assert(this->active_fe_index_transfer != nullptr, ExcInternalError()); + update_active_fe_table(); - // Unpack active_fe_indices. - this->active_fe_index_transfer->active_fe_indices.resize( - this->get_triangulation().n_active_cells(), - numbers::invalid_unsigned_int); - this->active_fe_index_transfer->cell_data_transfer->unpack( - this->active_fe_index_transfer->active_fe_indices); + Assert(this->active_fe_index_transfer != nullptr, ExcInternalError()); - // Update all locally owned active_fe_indices. - this->set_active_fe_indices( - this->active_fe_index_transfer->active_fe_indices); + // Unpack active_fe_indices. + this->active_fe_index_transfer->active_fe_indices.resize( + this->get_triangulation().n_active_cells(), numbers::invalid_unsigned_int); + this->active_fe_index_transfer->cell_data_transfer->unpack( + this->active_fe_index_transfer->active_fe_indices); - // Update active_fe_indices on ghost cells. - dealii::internal::hp::DoFHandlerImplementation::Implementation:: - communicate_active_fe_indices(*this); + // Update all locally owned active_fe_indices. + this->set_active_fe_indices( + this->active_fe_index_transfer->active_fe_indices); - // Free memory. - this->active_fe_index_transfer.reset(); - } + // Update active_fe_indices on ghost cells. + dealii::internal::hp::DoFHandlerImplementation::Implementation:: + communicate_active_fe_indices(*this); + + // Free memory. + this->active_fe_index_transfer.reset(); #endif } @@ -3396,70 +3365,42 @@ DoFHandler::prepare_for_serialization_of_active_fe_indices() &this->get_triangulation()) != nullptr), ExcNotImplemented()); - // Finite elements need to be assigned to each cell by calling - // distribute_dofs() first to make this functionality available. - if (fe_collection.size() > 0) - { - Assert(active_fe_index_transfer == nullptr, ExcInternalError()); - - active_fe_index_transfer = std::make_unique(); - - // Create transfer object and attach to it. - const auto *distributed_tria = dynamic_cast< - const parallel::distributed::Triangulation *>( - &this->get_triangulation()); - - active_fe_index_transfer->cell_data_transfer = std::make_unique< - parallel::distributed:: - CellDataTransfer>>( - *distributed_tria, - /*transfer_variable_size_data=*/false, - /*refinement_strategy=*/ - &dealii::AdaptationStrategies::Refinement:: - preserve, - /*coarsening_strategy=*/ - [this]( - const typename Triangulation::cell_iterator &parent, - const std::vector &children_fe_indices) - -> unsigned int { - return dealii::internal::hp::DoFHandlerImplementation:: - Implementation::determine_fe_from_children( - parent, children_fe_indices, fe_collection); - }); - - // If we work on a p::d::Triangulation, we have to transfer all - // active fe indices since ownership of cells may change. - - // Gather all current active_fe_indices - get_active_fe_indices(active_fe_index_transfer->active_fe_indices); - - // Attach to transfer object - active_fe_index_transfer->cell_data_transfer->prepare_for_serialization( - active_fe_index_transfer->active_fe_indices); - } -#endif -} - - - -template -void -DoFHandler::post_distributed_serialization_of_active_fe_indices() -{ -#ifndef DEAL_II_WITH_P4EST - Assert(false, - ExcMessage( - "You are attempting to use a functionality that is only available " - "if deal.II was configured to use p4est, but cmake did not find a " - "valid p4est library.")); -#else - if (this->fe_collection.size() > 0) - { - Assert(this->active_fe_index_transfer != nullptr, ExcInternalError()); - - // Free memory. - this->active_fe_index_transfer.reset(); - } + Assert(active_fe_index_transfer == nullptr, ExcInternalError()); + + active_fe_index_transfer = std::make_unique(); + + // Create transfer object and attach to it. + const auto *distributed_tria = + dynamic_cast *>( + &this->get_triangulation()); + + active_fe_index_transfer->cell_data_transfer = std::make_unique< + parallel::distributed:: + CellDataTransfer>>( + *distributed_tria, + /*transfer_variable_size_data=*/false, + /*refinement_strategy=*/ + &dealii::AdaptationStrategies::Refinement:: + preserve, + /*coarsening_strategy=*/ + [this](const typename Triangulation::cell_iterator &parent, + const std::vector &children_fe_indices) + -> unsigned int { + return dealii::internal::hp::DoFHandlerImplementation::Implementation:: + determine_fe_from_children(parent, + children_fe_indices, + fe_collection); + }); + + // If we work on a p::d::Triangulation, we have to transfer all + // active fe indices since ownership of cells may change. + + // Gather all current active_fe_indices + get_active_fe_indices(active_fe_index_transfer->active_fe_indices); + + // Attach to transfer object + active_fe_index_transfer->cell_data_transfer->prepare_for_serialization( + active_fe_index_transfer->active_fe_indices); #endif } @@ -3482,53 +3423,48 @@ DoFHandler::deserialize_active_fe_indices() &this->get_triangulation()) != nullptr), ExcNotImplemented()); - // Finite elements need to be assigned to each cell by calling - // distribute_dofs() first to make this functionality available. - if (fe_collection.size() > 0) - { - Assert(active_fe_index_transfer == nullptr, ExcInternalError()); - - active_fe_index_transfer = std::make_unique(); - - // Create transfer object and attach to it. - const auto *distributed_tria = dynamic_cast< - const parallel::distributed::Triangulation *>( - &this->get_triangulation()); - - active_fe_index_transfer->cell_data_transfer = std::make_unique< - parallel::distributed:: - CellDataTransfer>>( - *distributed_tria, - /*transfer_variable_size_data=*/false, - /*refinement_strategy=*/ - &dealii::AdaptationStrategies::Refinement:: - preserve, - /*coarsening_strategy=*/ - [this]( - const typename Triangulation::cell_iterator &parent, - const std::vector &children_fe_indices) - -> unsigned int { - return dealii::internal::hp::DoFHandlerImplementation:: - Implementation::determine_fe_from_children( - parent, children_fe_indices, fe_collection); - }); - - // Unpack active_fe_indices. - active_fe_index_transfer->active_fe_indices.resize( - get_triangulation().n_active_cells(), numbers::invalid_unsigned_int); - active_fe_index_transfer->cell_data_transfer->deserialize( - active_fe_index_transfer->active_fe_indices); - - // Update all locally owned active_fe_indices. - set_active_fe_indices(active_fe_index_transfer->active_fe_indices); - - // Update active_fe_indices on ghost cells. - dealii::internal::hp::DoFHandlerImplementation::Implementation:: - communicate_active_fe_indices(*this); - - // Free memory. - active_fe_index_transfer.reset(); - } + Assert(active_fe_index_transfer == nullptr, ExcInternalError()); + + active_fe_index_transfer = std::make_unique(); + + // Create transfer object and attach to it. + const auto *distributed_tria = + dynamic_cast *>( + &this->get_triangulation()); + + active_fe_index_transfer->cell_data_transfer = std::make_unique< + parallel::distributed:: + CellDataTransfer>>( + *distributed_tria, + /*transfer_variable_size_data=*/false, + /*refinement_strategy=*/ + &dealii::AdaptationStrategies::Refinement:: + preserve, + /*coarsening_strategy=*/ + [this](const typename Triangulation::cell_iterator &parent, + const std::vector &children_fe_indices) + -> unsigned int { + return dealii::internal::hp::DoFHandlerImplementation::Implementation:: + determine_fe_from_children(parent, + children_fe_indices, + fe_collection); + }); + + // Unpack active_fe_indices. + active_fe_index_transfer->active_fe_indices.resize( + get_triangulation().n_active_cells(), numbers::invalid_unsigned_int); + active_fe_index_transfer->cell_data_transfer->deserialize( + active_fe_index_transfer->active_fe_indices); + + // Update all locally owned active_fe_indices. + set_active_fe_indices(active_fe_index_transfer->active_fe_indices); + + // Update active_fe_indices on ghost cells. + dealii::internal::hp::DoFHandlerImplementation::Implementation:: + communicate_active_fe_indices(*this); + + // Free memory. + active_fe_index_transfer.reset(); #endif } diff --git a/source/hp/fe_collection.cc b/source/hp/fe_collection.cc index 9217b9337f..c3f84a30fb 100644 --- a/source/hp/fe_collection.cc +++ b/source/hp/fe_collection.cc @@ -28,13 +28,13 @@ namespace hp const std::set &fes, const unsigned int codim) const { +#ifdef DEBUG // Validate user inputs. Assert(codim <= dim, ExcImpossibleInDim(dim)); + Assert(size() > 0, ExcEmptyObject()); for (const auto &fe : fes) - { - (void)fe; - AssertIndexRange(fe, finite_elements.size()); - } + AssertIndexRange(fe, finite_elements.size()); +#endif // Check if any element of this FECollection is able to dominate all // elements of @p fes. If one was found, we add it to the set of @@ -68,13 +68,13 @@ namespace hp const std::set &fes, const unsigned int codim) const { +#ifdef DEBUG // Validate user inputs. Assert(codim <= dim, ExcImpossibleInDim(dim)); + Assert(size() > 0, ExcEmptyObject()); for (const auto &fe : fes) - { - (void)fe; - AssertIndexRange(fe, finite_elements.size()); - } + AssertIndexRange(fe, finite_elements.size()); +#endif // Check if any element of this FECollection is dominated by all // elements of @p fes. If one was found, we add it to the set of @@ -108,19 +108,19 @@ namespace hp const std::set &fes, const unsigned int codim) const { - // Validate user inputs. - Assert(codim <= dim, ExcImpossibleInDim(dim)); - for (const auto &fe : fes) - { - (void)fe; - AssertIndexRange(fe, finite_elements.size()); - } - // If the set of elements contains only a single element, // then this very element is considered to be the dominating one. if (fes.size() == 1) return *fes.begin(); +#ifdef DEBUG + // Validate user inputs. + Assert(codim <= dim, ExcImpossibleInDim(dim)); + Assert(size() > 0, ExcEmptyObject()); + for (const auto &fe : fes) + AssertIndexRange(fe, finite_elements.size()); +#endif + // There may also be others, in which case we'll check if any of these // elements is able to dominate all others. If one was found, we stop // looking further and return the dominating element. @@ -154,19 +154,19 @@ namespace hp const std::set &fes, const unsigned int codim) const { - // Validate user inputs. - Assert(codim <= dim, ExcImpossibleInDim(dim)); - for (const auto &fe : fes) - { - (void)fe; - AssertIndexRange(fe, finite_elements.size()); - } - // If the set of elements contains only a single element, // then this very element is considered to be the dominated one. if (fes.size() == 1) return *fes.begin(); +#ifdef DEBUG + // Validate user inputs. + Assert(codim <= dim, ExcImpossibleInDim(dim)); + Assert(size() > 0, ExcEmptyObject()); + for (const auto &fe : fes) + AssertIndexRange(fe, finite_elements.size()); +#endif + // There may also be others, in which case we'll check if any of these // elements is dominated by all others. If one was found, we stop // looking further and return the dominated element. -- 2.39.5