From 4a06c16aa75a23c86623d003b05d0fc5f84cbe98 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Thu, 2 Jun 2022 14:14:18 -0600 Subject: [PATCH] Move invalid_fe_index to numbers. --- include/deal.II/base/types.h | 5 ++ include/deal.II/dofs/dof_accessor.h | 60 +++++++++---------- include/deal.II/dofs/dof_accessor.templates.h | 22 +++---- include/deal.II/dofs/dof_handler.h | 18 ++++-- source/dofs/dof_accessor_get.cc | 6 +- source/dofs/dof_accessor_set.cc | 6 +- source/dofs/dof_handler.cc | 15 ++--- source/hp/fe_collection.cc | 5 ++ 8 files changed, 74 insertions(+), 63 deletions(-) diff --git a/include/deal.II/base/types.h b/include/deal.II/base/types.h index 0b8efbc481..5b7c96710d 100644 --- a/include/deal.II/base/types.h +++ b/include/deal.II/base/types.h @@ -215,6 +215,11 @@ namespace numbers const types::global_dof_index invalid_size_type = static_cast(-1); + /** + * An invalid value for active and future fe indices. + */ + const types::fe_index invalid_fe_index = static_cast(-1); + /** * An invalid value for indices of degrees of freedom. */ diff --git a/include/deal.II/dofs/dof_accessor.h b/include/deal.II/dofs/dof_accessor.h index 9de534f05f..19eae0834a 100644 --- a/include/deal.II/dofs/dof_accessor.h +++ b/include/deal.II/dofs/dof_accessor.h @@ -444,9 +444,9 @@ public: * cell-@>active_fe_index as last argument. */ void - get_dof_indices(std::vector &dof_indices, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + get_dof_indices( + std::vector &dof_indices, + const unsigned int fe_index = numbers::invalid_fe_index) const; /** * Return the global multilevel indices of the degrees of freedom that live @@ -455,19 +455,18 @@ public: * level this line lives on. */ void - get_mg_dof_indices(const int level, - std::vector &dof_indices, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + get_mg_dof_indices( + const int level, + std::vector &dof_indices, + const unsigned int fe_index = numbers::invalid_fe_index) const; /** * Set the level DoF indices that are returned by get_mg_dof_indices. */ void - set_mg_dof_indices( - const int level, - const std::vector &dof_indices, - const unsigned int fe_index = DoFHandler::invalid_fe_index); + set_mg_dof_indices(const int level, + const std::vector &dof_indices, + const unsigned int fe_index = numbers::invalid_fe_index); /** * Global DoF index of the i degree associated with the @p vertexth @@ -491,10 +490,10 @@ public: * this is interpreted as equal to `cell->active_fe_index()`. */ types::global_dof_index - vertex_dof_index(const unsigned int vertex, - const unsigned int i, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + vertex_dof_index( + const unsigned int vertex, + const unsigned int i, + const unsigned int fe_index = numbers::invalid_fe_index) const; /** * Return the global DoF index of the ith degree of freedom @@ -502,11 +501,11 @@ public: * see vertex_dof_index(). */ types::global_dof_index - mg_vertex_dof_index(const int level, - const unsigned int vertex, - const unsigned int i, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + mg_vertex_dof_index( + const int level, + const unsigned int vertex, + const unsigned int i, + const unsigned int fe_index = numbers::invalid_fe_index) const; /** * Index of the ith degree of freedom of this object. @@ -537,8 +536,7 @@ public: */ types::global_dof_index dof_index(const unsigned int i, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + const unsigned int fe_index = numbers::invalid_fe_index) const; /** * Return the dof_index on the given level. Also see dof_index. @@ -715,8 +713,7 @@ protected: void set_dof_index(const unsigned int i, const types::global_dof_index index, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + const unsigned int fe_index = numbers::invalid_fe_index) const; void set_mg_dof_index(const int level, @@ -724,12 +721,12 @@ protected: const types::global_dof_index index) const; void - set_mg_vertex_dof_index(const int level, - const unsigned int vertex, - const unsigned int i, - const types::global_dof_index index, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + set_mg_vertex_dof_index( + const int level, + const unsigned int vertex, + const unsigned int i, + const types::global_dof_index index, + const unsigned int fe_index = numbers::invalid_fe_index) const; // Iterator classes need to be friends because they need to access // operator== and operator!=. @@ -1302,8 +1299,7 @@ public: void set_dof_index(const unsigned int i, const types::global_dof_index index, - const unsigned int fe_index = - DoFHandler::invalid_fe_index) const; + const unsigned int fe_index = numbers::invalid_fe_index) const; }; diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index cbad821178..3515d1ca24 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -199,7 +199,7 @@ namespace internal { /** * Convert an FE index that might contain the right value but also - * invalid_fe_index to a right value if needed/possible. + * numbers::invalid_fe_index to a right value if needed/possible. */ template unsigned int @@ -212,7 +212,7 @@ namespace internal // 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 == numbers::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!")); @@ -226,15 +226,14 @@ namespace internal // 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)) + if ((dim == structdim) && (fe_index == numbers::invalid_fe_index)) { AssertDimension(cell.n_active_fe_indices(), 1); return cell.nth_active_fe_index(0); } - Assert((fe_index != DoFHandler::invalid_fe_index), + Assert((fe_index != numbers::invalid_fe_index), ExcMessage( "You need to specify a FE index if hp support is used!")); @@ -1610,7 +1609,7 @@ DoFAccessor::vertex_dof_index( { const unsigned int fe_index = (((this->dof_handler->hp_capability_enabled == false) && - (fe_index_ == DoFHandler::invalid_fe_index)) ? + (fe_index_ == numbers::invalid_fe_index)) ? // No hp enabled, and the argument is at its default value -> we // can translate to the default active fe index DoFHandler::default_fe_index : @@ -1619,8 +1618,7 @@ DoFAccessor::vertex_dof_index( // 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: - ((dim == structdim) && - (fe_index_ == DoFHandler::invalid_fe_index) ? + ((dim == structdim) && (fe_index_ == numbers::invalid_fe_index) ? this->nth_active_fe_index(0) : fe_index_)); @@ -2213,6 +2211,8 @@ namespace internal Assert(static_cast(accessor.level()) < accessor.dof_handler->hp_cell_future_fe_indices.size(), ExcMessage("DoFHandler not initialized")); + Assert(i != numbers::invalid_fe_index, + ExcMessage("Invalid finite element index.")); accessor.dof_handler ->hp_cell_active_fe_indices[accessor.level()] @@ -2277,6 +2277,8 @@ namespace internal Assert(static_cast(accessor.level()) < accessor.dof_handler->hp_cell_future_fe_indices.size(), ExcMessage("DoFHandler not initialized")); + Assert(i != numbers::invalid_fe_index, + ExcMessage("Invalid finite element index.")); accessor.dof_handler ->hp_cell_future_fe_indices[accessor.level()] @@ -2308,7 +2310,7 @@ namespace internal return accessor.dof_handler ->hp_cell_future_fe_indices[accessor.level()] [accessor.present_index] != - DoFHandler::invalid_active_fe_index; + numbers::invalid_fe_index; } @@ -2336,7 +2338,7 @@ namespace internal accessor.dof_handler ->hp_cell_future_fe_indices[accessor.level()] [accessor.present_index] = - DoFHandler::invalid_active_fe_index; + numbers::invalid_fe_index; } }; } // namespace DoFCellAccessorImplementation diff --git a/include/deal.II/dofs/dof_handler.h b/include/deal.II/dofs/dof_handler.h index 11af8d9d7a..786a46129c 100644 --- a/include/deal.II/dofs/dof_handler.h +++ b/include/deal.II/dofs/dof_handler.h @@ -519,8 +519,11 @@ public: /** * Invalid index of the finite element to be used on a given cell. + * + * @deprecated Use numbers::invalid_fe_index instead. */ - static const unsigned int invalid_fe_index = numbers::invalid_unsigned_int; + static const unsigned int invalid_fe_index DEAL_II_DEPRECATED = + numbers::invalid_fe_index; /** * The type in which we store the active FE index. @@ -537,9 +540,11 @@ public: /** * Invalid active FE index which will be used as a default value to determine * whether a future FE index has been set or not. + * + * @deprecated Use numbers::invalid_fe_index instead. */ - static const types::fe_index invalid_active_fe_index = - static_cast(-1); + static const types::fe_index invalid_active_fe_index DEAL_II_DEPRECATED = + numbers::invalid_fe_index; /** * Standard constructor, not initializing any data. After constructing an @@ -583,7 +588,8 @@ public: /** * Go through the triangulation and return a vector of active FE indices of - * all locally relevant cells. + * all locally relevant cells. Artificial cells will have the value + * numbers::invalid_fe_index assigned. */ std::vector get_active_fe_indices() const; @@ -602,7 +608,7 @@ public: /** * Go through the triangulation and set the future FE indices of all * locally owned cells to the values given in @p future_fe_indices. - * Cells corresponding to invalid_active_fe_index will be skipped. + * Cells corresponding to numbers::invalid_fe_index will be skipped. */ void set_future_fe_indices(const std::vector &future_fe_indices); @@ -610,7 +616,7 @@ public: /** * Go through the triangulation and return a vector of future FE indices of * all locally owned cells. If no future FE index has been set on a cell, - * its value will be invalid_active_fe_index. + * its value will be numbers::invalid_fe_index. */ std::vector get_future_fe_indices() const; diff --git a/source/dofs/dof_accessor_get.cc b/source/dofs/dof_accessor_get.cc index 0909f6805e..cd7861068f 100644 --- a/source/dofs/dof_accessor_get.cc +++ b/source/dofs/dof_accessor_get.cc @@ -50,7 +50,7 @@ DoFCellAccessor::get_interpolated_dof_values( { const unsigned int fe_index = (this->dof_handler->hp_capability_enabled == false && - fe_index_ == DoFHandler::invalid_fe_index) ? + fe_index_ == numbers::invalid_fe_index) ? DoFHandler::default_fe_index : fe_index_; @@ -64,7 +64,7 @@ DoFCellAccessor::get_interpolated_dof_values( // active cells, you either don't specify an fe_index, // or that you specify the correct one (fe_index == this->active_fe_index()) || - (fe_index == DoFHandler::invalid_fe_index)) + (fe_index == numbers::invalid_fe_index)) this->get_dof_values(values, interpolated_values); else { @@ -102,7 +102,7 @@ DoFCellAccessor::get_interpolated_dof_values( // space to this cell's (unknown) FE space unless an explicit // fe_index is given Assert((this->dof_handler->hp_capability_enabled == false) || - (fe_index != DoFHandler::invalid_fe_index), + (fe_index != numbers::invalid_fe_index), ExcMessage( "You cannot call this function on non-active cells " "of DoFHandler objects unless you provide an explicit " diff --git a/source/dofs/dof_accessor_set.cc b/source/dofs/dof_accessor_set.cc index 6ff52a32cd..feea52deff 100644 --- a/source/dofs/dof_accessor_set.cc +++ b/source/dofs/dof_accessor_set.cc @@ -186,7 +186,7 @@ namespace internal { const unsigned int fe_index = (cell.get_dof_handler().has_hp_capabilities() == false && - fe_index_ == DoFHandler::invalid_fe_index) ? + fe_index_ == numbers::invalid_fe_index) ? DoFHandler::default_fe_index : fe_index_; @@ -197,7 +197,7 @@ namespace internal // active cells, you either don't specify an fe_index, // or that you specify the correct one (fe_index == cell.active_fe_index()) || - (fe_index == DoFHandler::invalid_fe_index)) + (fe_index == numbers::invalid_fe_index)) // simply set the values on this cell processor(cell, local_values, values); else @@ -228,7 +228,7 @@ namespace internal // otherwise distribute them to the children { Assert((cell.get_dof_handler().has_hp_capabilities() == false) || - (fe_index != DoFHandler::invalid_fe_index), + (fe_index != numbers::invalid_fe_index), ExcMessage( "You cannot call this function on non-active cells " "of DoFHandler objects unless you provide an explicit " diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index 12dca976a7..846a0f60c6 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -46,10 +46,6 @@ template const unsigned int DoFHandler::default_fe_index; -template -const types::fe_index DoFHandler::invalid_active_fe_index; - - namespace internal { template @@ -2638,7 +2634,7 @@ std::vector DoFHandler::get_active_fe_indices() const { std::vector active_fe_indices( - this->get_triangulation().n_active_cells(), invalid_active_fe_index); + this->get_triangulation().n_active_cells(), numbers::invalid_fe_index); // we could try to extract the values directly, since they are // stored as protected data of this object, but for simplicity we @@ -2678,7 +2674,8 @@ DoFHandler::set_future_fe_indices( // tests which we would have to duplicate ourselves otherwise for (const auto &cell : this->active_cell_iterators()) if (cell->is_locally_owned() && - future_fe_indices[cell->active_cell_index()] != invalid_active_fe_index) + future_fe_indices[cell->active_cell_index()] != + numbers::invalid_fe_index) cell->set_future_fe_index(future_fe_indices[cell->active_cell_index()]); } @@ -2689,7 +2686,7 @@ std::vector DoFHandler::get_future_fe_indices() const { std::vector future_fe_indices( - this->get_triangulation().n_active_cells(), invalid_active_fe_index); + this->get_triangulation().n_active_cells(), numbers::invalid_fe_index); // we could try to extract the values directly, since they are // stored as protected data of this object, but for simplicity we @@ -2820,7 +2817,7 @@ DoFHandler::create_active_fe_table() this->hp_cell_active_fe_indices[level].resize( this->tria->n_raw_cells(level), 0); this->hp_cell_future_fe_indices[level].resize( - this->tria->n_raw_cells(level), invalid_active_fe_index); + this->tria->n_raw_cells(level), numbers::invalid_fe_index); } else { @@ -2879,7 +2876,7 @@ DoFHandler::update_active_fe_table() // We have used future FE indices to update all active FE indices // before refinement happened, thus we are safe to clear them now. this->hp_cell_future_fe_indices[i].assign(this->tria->n_raw_cells(i), - invalid_active_fe_index); + numbers::invalid_fe_index); } } diff --git a/source/hp/fe_collection.cc b/source/hp/fe_collection.cc index c840c63da2..dfee75622a 100644 --- a/source/hp/fe_collection.cc +++ b/source/hp/fe_collection.cc @@ -72,6 +72,11 @@ namespace hp ExcMessage("All elements inside a collection need to have the " "same number of vector components!")); + Assert(this->size() <= std::numeric_limits::max() && + this->size() != numbers::invalid_fe_index, + ExcMessage( + "You reached the maximum possible number of finite elements.")); + Collection>::push_back(new_fe.clone()); } -- 2.39.5