From 86f674a35266e2db2087cb962551ce7750856b22 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Tue, 18 Jan 2022 07:59:35 +0100 Subject: [PATCH] FEEval: do not use cell as indication if reinit() was called --- include/deal.II/matrix_free/fe_evaluation.h | 31 +++++++++++++------ .../deal.II/matrix_free/fe_evaluation_data.h | 17 +++++++--- 2 files changed, 34 insertions(+), 14 deletions(-) diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index 86e54949e2..13aa4761c9 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -5008,7 +5008,7 @@ FEEvaluationBase:: submit_value(const Tensor<1, n_components_, VectorizedArrayType> val_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -5046,7 +5046,7 @@ FEEvaluationBase:: const Tensor<1, n_components_, Tensor<1, dim, VectorizedArrayType>> grad_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -5158,7 +5158,7 @@ FEEvaluationBase:: hessian_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -5309,7 +5309,7 @@ inline Tensor<1, n_components_, VectorizedArrayType> FEEvaluationBase:: integrate_value() const { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); # ifdef DEBUG Assert(this->values_quad_submitted == true, internal::ExcAccessToUninitializedField()); @@ -5646,7 +5646,7 @@ FEEvaluationAccess::submit_value( const VectorizedArrayType val_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -5699,7 +5699,7 @@ FEEvaluationAccess:: submit_gradient(const Tensor<1, dim, VectorizedArrayType> grad_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -6016,7 +6016,7 @@ FEEvaluationAccess:: submit_divergence(const VectorizedArrayType div_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -6078,7 +6078,7 @@ FEEvaluationAccess:: // could have used base class operator, but that involves some overhead // which is inefficient. it is nice to have the symmetric tensor because // that saves some operations - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); Assert(this->J_value != nullptr, internal::ExcMatrixFreeAccessToUninitializedMappingField( @@ -6372,7 +6372,7 @@ FEEvaluationAccess<1, 1, Number, is_face, VectorizedArrayType>::submit_value( const VectorizedArrayType val_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); # ifdef DEBUG this->values_quad_submitted = true; @@ -6421,7 +6421,7 @@ FEEvaluationAccess<1, 1, Number, is_face, VectorizedArrayType>::submit_gradient( const VectorizedArrayType grad_in, const unsigned int q_point) { - Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(this->is_reinitialized, ExcNotInitialized()); AssertIndexRange(q_point, this->n_quadrature_points); # ifdef DEBUG this->gradients_quad_submitted = true; @@ -6907,6 +6907,7 @@ FEEvaluationcell_ids[i] = cell_index * VectorizedArrayType::size() + i; # ifdef DEBUG + this->is_reinitialized = true; this->dof_values_initialized = false; this->values_quad_initialized = false; this->gradients_quad_initialized = false; @@ -6945,6 +6946,10 @@ FEEvaluationget_mg_dof_indices(this->local_dof_indices); else cell->get_dof_indices(this->local_dof_indices); + +# ifdef DEBUG + this->is_reinitialized = true; +# endif } @@ -6971,6 +6976,10 @@ FEEvaluationmapped_geometry.get() != 0, ExcNotInitialized()); this->mapped_geometry->reinit(cell); + +# ifdef DEBUG + this->is_reinitialized = true; +# endif } @@ -7640,6 +7649,7 @@ FEFaceEvaluationnormals_times_jacobians[!this->is_interior_face][offsets]; # ifdef DEBUG + this->is_reinitialized = true; this->dof_values_initialized = false; this->values_quad_initialized = false; this->gradients_quad_initialized = false; @@ -7776,6 +7786,7 @@ FEFaceEvaluationis_interior_face][offsets]; # ifdef DEBUG + this->is_reinitialized = true; this->dof_values_initialized = false; this->values_quad_initialized = false; this->gradients_quad_initialized = false; diff --git a/include/deal.II/matrix_free/fe_evaluation_data.h b/include/deal.II/matrix_free/fe_evaluation_data.h index de3cb37535..ee815a8eba 100644 --- a/include/deal.II/matrix_free/fe_evaluation_data.h +++ b/include/deal.II/matrix_free/fe_evaluation_data.h @@ -472,7 +472,7 @@ public: get_cell_ids() const { // implemented inline to avoid compilation problems on Windows - Assert(cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(is_reinitialized, ExcNotInitialized()); return cell_ids; } @@ -484,7 +484,7 @@ public: get_cell_or_face_batch_id() const { // implemented inline to avoid compilation problems on Windows - Assert(cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(is_reinitialized, ExcNotInitialized()); return cell; } @@ -496,7 +496,7 @@ public: get_cell_or_face_ids() const { // implemented inline to avoid compilation problems on Windows - Assert(cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(is_reinitialized, ExcNotInitialized()); if (!is_face || dof_access_index == internal::MatrixFreeFunctions::DoFInfo::dof_access_cell) return cell_ids; @@ -705,6 +705,12 @@ protected: */ Number *hessians_quad; + /** + * Debug information to track whether the reinit() function of derived classes + * was called. + */ + bool is_reinitialized; + /** * Debug information to track whether dof values have been initialized * before accessed. Used to control exceptions when uninitialized data is @@ -895,6 +901,7 @@ inline FEEvaluationData::FEEvaluationData( , quadrature_weights( descriptor != nullptr ? descriptor->quadrature_weights.begin() : nullptr) # ifdef DEBUG + , is_reinitialized(false) , dof_values_initialized(false) , values_quad_initialized(false) , gradients_quad_initialized(false) @@ -944,6 +951,7 @@ inline FEEvaluationData::FEEvaluationData( , is_interior_face(true) , dof_access_index(internal::MatrixFreeFunctions::DoFInfo::dof_access_cell) , mapped_geometry(mapped_geometry) + , is_reinitialized(false) { mapping_data = &mapped_geometry->get_data_storage(); jacobian = mapped_geometry->get_data_storage().jacobians[0].begin(); @@ -974,6 +982,7 @@ FEEvaluationData::operator=(const FEEvaluationData &other) quadrature_weights = other.quadrature_weights; # ifdef DEBUG + is_reinitialized = false; dof_values_initialized = false; values_quad_initialized = false; gradients_quad_initialized = false; @@ -1242,7 +1251,7 @@ template inline internal::MatrixFreeFunctions::GeometryType FEEvaluationData::get_cell_type() const { - Assert(cell != numbers::invalid_unsigned_int, ExcNotInitialized()); + Assert(is_reinitialized, ExcNotInitialized()); return cell_type; } -- 2.39.5