From: Peter Munch Date: Mon, 21 Feb 2022 10:20:36 +0000 (+0100) Subject: Improve assert message for invalid FE indices in DoFAccessor X-Git-Tag: v9.4.0-rc1~454^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F13430%2Fhead;p=dealii.git Improve assert message for invalid FE indices in DoFAccessor --- diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index 62c99a95d7..865d206a7f 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -197,6 +197,51 @@ namespace internal { namespace DoFAccessorImplementation { + /** + * Convert an FE index that might contain the right value but also + * invalid_fe_index to a right value if needed/possible. + */ + template + unsigned int + get_fe_index_or_default( + const DoFAccessor &cell, + const unsigned int fe_index) + { + if (cell.get_dof_handler().has_hp_capabilities() == false) + { + // No hp enabled, and the argument is at its default value -> we + // can translate to the default active fe index + Assert( + (fe_index == DoFHandler::invalid_fe_index) || + (fe_index == DoFHandler::default_fe_index), + ExcMessage( + "It is not possible to specify a FE index if no hp support is used!")); + + return DoFHandler::default_fe_index; + } + else + { + // Otherwise: If anything other than the default is provided by + // the caller, then we should take just that. As an exception, if + // we are on a cell (rather than a face/edge/vertex), then we know + // that there is only one active fe index on this cell and we can + // use that: + if ((dim == structdim) && + (fe_index == DoFHandler::invalid_fe_index)) + { + AssertDimension(cell.n_active_fe_indices(), 1); + + return cell.nth_active_fe_index(0); + } + + Assert((fe_index != DoFHandler::invalid_fe_index), + ExcMessage( + "You need to specify a FE index if hp support is used!")); + + return fe_index; + } + } + /** * A class like the one with same name in tria.cc. See there for more * information. @@ -884,11 +929,9 @@ namespace internal } else { - const unsigned int fe_index = - (accessor.get_dof_handler().hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default( + accessor, fe_index_); unsigned int index = 0; @@ -1054,11 +1097,9 @@ namespace internal const unsigned int fe_index_, const DoFProcessor &dof_processor) const { - const unsigned int fe_index = - (accessor.get_dof_handler().hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default( + accessor, fe_index_); process_object( accessor.get_dof_handler(), @@ -1092,11 +1133,9 @@ namespace internal const unsigned int fe_index_, const DoFProcessor &dof_processor) const { - const unsigned int fe_index = - (accessor.get_dof_handler().hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default( + accessor, fe_index_); process_object(accessor.get_dof_handler(), accessor.level(), @@ -1407,11 +1446,9 @@ DoFAccessor::dof_index( const unsigned int i, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); // access the respective DoF return dealii::internal::DoFAccessorImplementation::Implementation:: @@ -1442,11 +1479,9 @@ DoFAccessor::set_dof_index( const types::global_dof_index index, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); // access the respective DoF dealii::internal::DoFAccessorImplementation::Implementation::set_dof_index( @@ -1562,11 +1597,9 @@ DoFAccessor::mg_vertex_dof_index( const unsigned int i, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); (void)fe_index; Assert(this->dof_handler != nullptr, ExcInvalidObject()); Assert(this->dof_handler->mg_vertex_dofs.size() > 0, @@ -1594,11 +1627,9 @@ DoFAccessor:: const types::global_dof_index index, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); (void)fe_index; Assert(this->dof_handler != nullptr, ExcInvalidObject()); AssertIndexRange(vertex, this->n_vertices()); @@ -1646,11 +1677,9 @@ DoFAccessor::get_dof_indices( { Assert(this->dof_handler != nullptr, ExcInvalidObject()); - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); Assert(static_cast(this->level()) < this->dof_handler->object_dof_indices.size(), @@ -1689,11 +1718,9 @@ DoFAccessor::get_mg_dof_indices( ExcMessage("Multigrid DoF indices can only be accessed after " "DoFHandler::distribute_mg_dofs() has been called!")); - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); internal::DoFAccessorImplementation::Implementation::get_mg_dof_indices( *this, level, dof_indices, fe_index); @@ -1709,11 +1736,9 @@ DoFAccessor::set_mg_dof_indices( { Assert(this->dof_handler != nullptr, ExcInvalidObject()); - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? - DoFHandler::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); internal::DoFAccessorImplementation::Implementation::set_mg_dof_indices( *this, level, dof_indices, fe_index); @@ -1905,11 +1930,9 @@ DoFAccessor<0, 1, spacedim, level_dof_access>::get_dof_indices( std::vector &dof_indices, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler<1, spacedim>::invalid_fe_index) ? - DoFHandler<1, spacedim>::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); for (unsigned int i = 0; i < dof_indices.size(); ++i) dof_indices[i] = dealii::internal::DoFAccessorImplementation:: @@ -1930,11 +1953,9 @@ DoFAccessor<0, 1, spacedim, level_dof_access>::get_mg_dof_indices( std::vector &dof_indices, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler<1, spacedim>::invalid_fe_index) ? - DoFHandler<1, spacedim>::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); (void)fe_index; AssertDimension(fe_index, (DoFHandler<1, spacedim>::default_fe_index)); @@ -1953,11 +1974,9 @@ DoFAccessor<0, 1, spacedim, level_dof_access>::vertex_dof_index( const unsigned int i, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler<1, spacedim>::invalid_fe_index) ? - DoFHandler<1, spacedim>::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); (void)vertex; AssertIndexRange(vertex, 1); @@ -1978,11 +1997,9 @@ DoFAccessor<0, 1, spacedim, level_dof_access>::dof_index( const unsigned int i, const unsigned int fe_index_) const { - const unsigned int fe_index = - (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler<1, spacedim>::invalid_fe_index) ? - DoFHandler<1, spacedim>::default_fe_index : - fe_index_; + const auto fe_index = + internal::DoFAccessorImplementation::get_fe_index_or_default(*this, + fe_index_); return dealii::internal::DoFAccessorImplementation::Implementation:: get_dof_index(*this->dof_handler,