From f44a7fc6fc731e308e723bb1ee648163552a68b3 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Tue, 14 Dec 2021 11:18:56 +0100 Subject: [PATCH] Simplify call to detection for gather_evaluate --- .../matrix_free/evaluation_template_factory.h | 6 +- .../evaluation_template_factory.templates.h | 12 -- include/deal.II/matrix_free/fe_evaluation.h | 182 ++++++------------ 3 files changed, 57 insertions(+), 143 deletions(-) diff --git a/include/deal.II/matrix_free/evaluation_template_factory.h b/include/deal.II/matrix_free/evaluation_template_factory.h index 22c3591a66..d372e330f9 100644 --- a/include/deal.II/matrix_free/evaluation_template_factory.h +++ b/include/deal.II/matrix_free/evaluation_template_factory.h @@ -72,12 +72,10 @@ namespace internal const EvaluationFlags::EvaluationFlags integration_flag, Number * values_dofs, FEEvaluationData & eval); - - static bool - fast_evaluation_supported(const unsigned int given_degree, - const unsigned int n_q_points_1d); }; + + template struct FEFaceEvaluationGatherFactory { diff --git a/include/deal.II/matrix_free/evaluation_template_factory.templates.h b/include/deal.II/matrix_free/evaluation_template_factory.templates.h index 5fd88c1002..a644970b0a 100644 --- a/include/deal.II/matrix_free/evaluation_template_factory.templates.h +++ b/include/deal.II/matrix_free/evaluation_template_factory.templates.h @@ -184,18 +184,6 @@ namespace internal - template - bool - FEFaceEvaluationFactory::fast_evaluation_supported( - const unsigned int given_degree, - const unsigned int n_q_points_1d) - { - return instantiation_helper_run<1, FastEvaluationSupported>(given_degree, - n_q_points_1d); - } - - - template void CellwiseInverseMassFactory::apply( diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index 2a860a4a8c..a21bf2b821 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -2802,13 +2802,16 @@ namespace internal namespace internal { + // Extract all internal data pointers and indices in a single function that + // get passed on to the constructor of FEEvaluationData, avoiding to look + // things up multiple times template inline typename FEEvaluationData:: InitializationData - get_shape_info_and_indices( + extract_initialization_data( const MatrixFree &data, const unsigned int dof_no, const unsigned int first_selected_component, @@ -2862,6 +2865,7 @@ namespace internal } // namespace internal + /*----------------------- FEEvaluationBase ----------------------------------*/ template ( - internal::get_shape_info_and_indices(data_in, - dof_no, - first_selected_component, - quad_no, - fe_degree, - n_q_points, - active_fe_index, - active_quad_index, - face_type), + internal::extract_initialization_data(data_in, + dof_no, + first_selected_component, + quad_no, + fe_degree, + n_q_points, + active_fe_index, + active_quad_index, + face_type), is_interior_face, quad_no, first_selected_component) @@ -7147,6 +7151,7 @@ namespace internal * Implementation for standard vectors (that have the begin() methods). */ template ().begin()), Number *>::value, VectorType>::type * = nullptr> - bool - try_gather_evaluate_inplace( - EvaluatorType & phi, - const VectorType & input_vector, - const EvaluationFlags::EvaluationFlags evaluation_flag) + VectorizedArrayType * + check_vector_access_inplace(const EvaluatorType &phi, VectorType &vector) { const unsigned int cell = phi.get_cell_or_face_batch_id(); const auto & dof_info = phi.get_dof_info(); - using VectorizedArrayType = - typename std::remove_reference::type; - // If the index storage is interleaved and contiguous and the vector storage - // has the correct alignment, we can directly pass the pointer into the - // vector to the evaluate() call, without reading the vector entries into a - // separate data field. This saves some operations. + // If the index storage is interleaved and contiguous and the vector + // storage has the correct alignment, we can directly pass the pointer + // into the vector to the evaluate() and integrate() calls, without + // reading the vector entries into a separate data field. This saves some + // operations. if (std::is_same::value && dof_info.index_storage_variants [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell][cell] == internal::MatrixFreeFunctions::DoFInfo::IndexStorageVariants:: interleaved_contiguous && reinterpret_cast( - input_vector.begin() + + vector.begin() + dof_info.dof_indices_contiguous [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell] [cell * VectorizedArrayType::size()]) % sizeof(VectorizedArrayType) == 0) { - const VectorizedArrayType *vec_values = - reinterpret_cast( - input_vector.begin() + - dof_info.dof_indices_contiguous - [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell] - [cell * VectorizedArrayType::size()] + - dof_info.component_dof_indices_offset - [phi.get_active_fe_index()] - [phi.get_first_selected_component()] * - VectorizedArrayType::size()); - - phi.evaluate(vec_values, evaluation_flag); - - return true; - } - - return false; - } - - /** - * Implementation for block vectors. - */ - template ::value || - !std::is_same().begin()), - Number *>::value, - VectorType>::type * = nullptr> - bool - try_gather_evaluate_inplace(EvaluatorType &, - const VectorType &, - const EvaluationFlags::EvaluationFlags) - { - return false; - } - - /** - * Implementation for vectors that have the begin() methods. - */ - template ::value && - std::is_same().begin()), - Number *>::value, - VectorType>::type * = nullptr> - bool - try_integrate_scatter_inplace( - EvaluatorType & phi, - VectorType & destination, - const EvaluationFlags::EvaluationFlags evaluation_flag) - { - const unsigned int cell = phi.get_cell_or_face_batch_id(); - const auto & dof_info = phi.get_dof_info(); - using VectorizedArrayType = - typename std::remove_reference::type; - - // If the index storage is interleaved and contiguous and the vector storage - // has the correct alignment, we can directly pass the pointer into the - // vector to the evaluate() call, without reading the vector entries into a - // separate data field. This saves some operations. - if (std::is_same::value && - dof_info.index_storage_variants - [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell][cell] == - internal::MatrixFreeFunctions::DoFInfo::IndexStorageVariants:: - interleaved_contiguous && - reinterpret_cast( - destination.begin() + + return reinterpret_cast( + vector.begin() + dof_info.dof_indices_contiguous [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell] - [cell * VectorizedArrayType::size()]) % - sizeof(VectorizedArrayType) == - 0) - { - VectorizedArrayType *vec_values = - reinterpret_cast( - destination.begin() + - dof_info.dof_indices_contiguous - [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell] - [cell * VectorizedArrayType::size()] + - dof_info.component_dof_indices_offset - [phi.get_active_fe_index()] - [phi.get_first_selected_component()] * - VectorizedArrayType::size()); - - phi.integrate(evaluation_flag, vec_values, true); - - return true; + [cell * VectorizedArrayType::size()] + + dof_info.component_dof_indices_offset + [phi.get_active_fe_index()][phi.get_first_selected_component()] * + VectorizedArrayType::size()); } - - return false; + else + return nullptr; } /** - * Implementation for all other vectors like block vectors. + * Implementation for block vectors. */ template ().begin()), Number *>::value, VectorType>::type * = nullptr> - bool - try_integrate_scatter_inplace(EvaluatorType &, - VectorType &, - const EvaluationFlags::EvaluationFlags) + VectorizedArrayType * + check_vector_access_inplace(const EvaluatorType &, VectorType &) { - return false; + return nullptr; } } // namespace internal @@ -7317,9 +7234,12 @@ FEEvaluation(*this, - input_vector, - evaluation_flag) == false) + const VectorizedArrayType *src_ptr = + internal::check_vector_access_inplace( + *this, input_vector); + if (src_ptr != nullptr) + evaluate(src_ptr, evaluation_flag); + else { this->read_dof_values(input_vector); evaluate(this->begin_dof_values(), evaluation_flag); @@ -7524,8 +7444,12 @@ FEEvaluation( - *this, destination, integration_flag) == false) + VectorizedArrayType *dst_ptr = + internal::check_vector_access_inplace( + *this, destination); + if (dst_ptr != nullptr) + integrate(integration_flag, dst_ptr, true); + else { integrate(integration_flag, this->begin_dof_values()); this->distribute_local_to_global(destination); @@ -8121,7 +8045,9 @@ FEFaceEvaluationactive_fe_index, this->dof_info); - if (fe_degree > 0 && + if (this->data->data.front().fe_degree > 0 && + fast_evaluation_supported(this->data->data.front().fe_degree, + this->data->data.front().n_q_points_1d) && internal::FEFaceEvaluationImplGatherEvaluateSelector< dim, Number, @@ -8229,7 +8155,9 @@ FEFaceEvaluationactive_fe_index, this->dof_info); - if (fe_degree > 0 && + if (this->data->data.front().fe_degree > 0 && + fast_evaluation_supported(this->data->data.front().fe_degree, + this->data->data.front().n_q_points_1d) && internal::FEFaceEvaluationImplGatherEvaluateSelector< dim, Number, @@ -8362,7 +8290,7 @@ FEFaceEvaluation:: + internal::FEEvaluationFactory:: fast_evaluation_supported(given_degree, give_n_q_points_1d) : true; } -- 2.39.5