From: Peter Munch <peterrmuench@gmail.com> Date: Sat, 7 May 2022 21:25:13 +0000 (+0200) Subject: Improve asserts in MF/FEEval regarding FENothing X-Git-Tag: v9.4.0-rc1~261^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F13690%2Fhead;p=dealii.git Improve asserts in MF/FEEval regarding FENothing --- diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index 5c80c2f994..120e9dbb31 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -6599,6 +6599,19 @@ FEEvaluation<dim, (void)dof_no; (void)first_selected_component; + Assert( + this->data->dofs_per_component_on_cell > 0, + ExcMessage( + "There is nothing useful you can do with an FEEvaluation object with " + "FE_Nothing, i.e., without DoFs! If you have passed to " + "MatrixFree::reinit() a collection of finite elements also containing " + "FE_Nothing, please check - before creating FEEvaluation - the category " + "of the current range by calling either " + "MatrixFree::get_cell_range_category(range) or " + "MatrixFree::get_face_range_category(range). The returned category " + "is the index of the active FE, which you can use to exclude " + "FE_Nothing.")); + # ifdef DEBUG // print error message when the dimensions do not match. Propose a possible // fix diff --git a/include/deal.II/matrix_free/matrix_free.h b/include/deal.II/matrix_free/matrix_free.h index f977f5d904..da75ab3aaa 100644 --- a/include/deal.II/matrix_free/matrix_free.h +++ b/include/deal.II/matrix_free/matrix_free.h @@ -1884,6 +1884,24 @@ public: get_face_quadrature(const unsigned int quad_index = 0, const unsigned int hp_active_fe_index = 0) const; + /** + * Return the category the current batch range of cells was assigned to. + * Categories run between the given values in the field + * AdditionalData::cell_vectorization_category for non-hp-DoFHandler types + * and return the active FE index in the hp-adaptive case. + */ + unsigned int + get_cell_range_category( + const std::pair<unsigned int, unsigned int> cell_batch_range) const; + + /** + * Return the category of the cells on the two sides of the current batch + * range of faces. + */ + std::pair<unsigned int, unsigned int> + get_face_range_category( + const std::pair<unsigned int, unsigned int> face_batch_range) const; + /** * Return the category the current batch of cells was assigned to. Categories * run between the given values in the field @@ -1894,7 +1912,7 @@ public: get_cell_category(const unsigned int cell_batch_index) const; /** - * Return the category on the cells on the two sides of the current batch of + * Return the category of the cells on the two sides of the current batch of * faces. */ std::pair<unsigned int, unsigned int> @@ -2780,6 +2798,44 @@ MatrixFree<dim, Number, VectorizedArrayType>::get_face_quadrature( +template <int dim, typename Number, typename VectorizedArrayType> +inline unsigned int +MatrixFree<dim, Number, VectorizedArrayType>::get_cell_range_category( + const std::pair<unsigned int, unsigned int> range) const +{ + const auto result = get_cell_category(range.first); + +# ifdef DEBUG + for (unsigned int i = range.first; i < range.second; ++i) + Assert(result == get_cell_category(i), + ExcMessage( + "The cell batches of the range have different categories!")); +# endif + + return result; +} + + + +template <int dim, typename Number, typename VectorizedArrayType> +inline std::pair<unsigned int, unsigned int> +MatrixFree<dim, Number, VectorizedArrayType>::get_face_range_category( + const std::pair<unsigned int, unsigned int> range) const +{ + const auto result = get_face_category(range.first); + +# ifdef DEBUG + for (unsigned int i = range.first; i < range.second; ++i) + Assert(result == get_face_category(i), + ExcMessage( + "The face batches of the range have different categories!")); +# endif + + return result; +} + + + template <int dim, typename Number, typename VectorizedArrayType> inline unsigned int MatrixFree<dim, Number, VectorizedArrayType>::get_cell_category( diff --git a/include/deal.II/matrix_free/util.h b/include/deal.II/matrix_free/util.h index 323f2093e0..b0b01f6fe1 100644 --- a/include/deal.II/matrix_free/util.h +++ b/include/deal.II/matrix_free/util.h @@ -110,6 +110,13 @@ namespace internal inline std::pair<Quadrature<dim - 1>, Quadrature<dim - 1>> get_unique_face_quadratures(const Quadrature<dim> &quad) { + AssertThrow( + quad.size() > 0, + ExcMessage( + "There is nothing useful you can do with a MatrixFree/FEEvaluation " + "object when using a quadrature formula with zero " + "quadrature points!")); + if (dim == 2 || dim == 3) { for (unsigned int i = 1; i <= 4; ++i)