From 1dcf6de8955a569d8790e844940d89cf019749d5 Mon Sep 17 00:00:00 2001 From: =?utf8?q?Bu=C4=9Frahan=20Tem=C3=BCr?= Date: Wed, 15 Mar 2023 10:34:22 +0100 Subject: [PATCH] Apply requests from the review - Remove explicit unrolling of tensors - Use ArrayView in internal functions instead of passing around AlignedVector - Make vmult private and inline - Documentation --- .../deal.II/matrix_free/evaluation_kernels.h | 7 +- .../matrix_free/evaluation_template_factory.h | 2 +- .../evaluation_template_factory.templates.h | 19 +++-- include/deal.II/matrix_free/operators.h | 80 ++++++++++--------- 4 files changed, 58 insertions(+), 50 deletions(-) diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index f999a7eaa6..b0145fa8e4 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -5707,7 +5707,7 @@ namespace internal static bool run(const unsigned int n_desired_components, const FEEvaluationData &fe_eval, - const AlignedVector & inverse_coefficients, + const ArrayView & inverse_coefficients, const bool dyadic_coefficients, const Number * in_array, Number * out_array) @@ -5728,7 +5728,7 @@ namespace internal { if (inverse_coefficients.size() != dofs_per_component) AssertDimension(n_desired_components * dofs_per_component, - inverse_coefficients.size()) + inverse_coefficients.size()); } else { @@ -5821,8 +5821,9 @@ namespace internal return false; } + private: template - static void + static inline void vmult(const Number * inverse_coefficients, const Number * src, Number * dst, diff --git a/include/deal.II/matrix_free/evaluation_template_factory.h b/include/deal.II/matrix_free/evaluation_template_factory.h index 616386cb8b..a179402e54 100644 --- a/include/deal.II/matrix_free/evaluation_template_factory.h +++ b/include/deal.II/matrix_free/evaluation_template_factory.h @@ -107,7 +107,7 @@ namespace internal static void apply(const unsigned int n_components, const FEEvaluationData &fe_eval, - const AlignedVector & inverse_coefficients, + const ArrayView & inverse_coefficients, const bool dyadic_coefficients, const Number * in_array, Number * out_array); 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 95b93c6f29..67eba41712 100644 --- a/include/deal.II/matrix_free/evaluation_template_factory.templates.h +++ b/include/deal.II/matrix_free/evaluation_template_factory.templates.h @@ -114,7 +114,7 @@ namespace internal CellwiseInverseMassFactory::apply( const unsigned int n_components, const FEEvaluationData &fe_eval, - const AlignedVector & inverse_coefficients, + const ArrayView & inverse_coefficients, const bool dyadic_coefficients, const Number * in_array, Number * out_array) @@ -122,15 +122,14 @@ namespace internal const unsigned int fe_degree = fe_eval.get_shape_info().data[0].fe_degree; instantiation_helper_run< 1, - CellwiseInverseMassMatrixImplFlexible>( - fe_degree, - fe_degree + 1, - n_components, - fe_eval, - inverse_coefficients, - dyadic_coefficients, - in_array, - out_array); + CellwiseInverseMassMatrixImplFlexible>(fe_degree, + fe_degree + 1, + n_components, + fe_eval, + inverse_coefficients, + dyadic_coefficients, + in_array, + out_array); } diff --git a/include/deal.II/matrix_free/operators.h b/include/deal.II/matrix_free/operators.h index 480aafb288..e2c5334900 100644 --- a/include/deal.II/matrix_free/operators.h +++ b/include/deal.II/matrix_free/operators.h @@ -611,7 +611,8 @@ namespace MatrixFreeOperators * The equation may contain variable coefficients, so the user is required * to provide an array for the inverse of the local coefficient (this class * provide a helper method 'fill_inverse_JxW_values' to get the inverse of a - * constant-coefficient operator). + * constant-coefficient operator). The local coefficient can either be scalar + * in each component, or dyadic, i.e. couple between components. */ template &inverse_coefficient, const unsigned int n_actual_components, const VectorizedArrayType * in_array, - VectorizedArrayType * out_array, - const bool dyadic_coefficients = false) const; + VectorizedArrayType * out_array) const; /** * Applies the inverse @ref GlossMassMatrix "mass matrix" operation on an input array, using the @@ -672,11 +673,11 @@ namespace MatrixFreeOperators * The second-rank tensor at each quadrature point defines a linear operator * on a vector holding the dof components. It is assumed that the passed * input and output arrays are of correct size, namely - * FEEvaluation::dofs_per_cell long. + * FEEvaluation::dofs_per_cell long. The `in_array` and `out_array` + * arguments may point to the same memory position. * `inverse_dyadic_coefficients` must be dofs_per_component long, and every * element must be a second-rank tensor of dimension `n_components`. All - * entries should also contain the inverse JxW values. The `in_array` and - * `out_array` arguments may point to the same memory position. + * entries should also contain the inverse JxW values. */ void apply(const AlignedVector> @@ -1112,24 +1113,26 @@ namespace MatrixFreeOperators apply(const AlignedVector &inverse_coefficients, const unsigned int n_actual_components, const VectorizedArrayType * in_array, - VectorizedArrayType * out_array, - const bool dyadic_coefficients) const + VectorizedArrayType * out_array) const { if (fe_degree > -1) - internal::CellwiseInverseMassMatrixImplFlexible< - dim, - VectorizedArrayType>::template run(n_actual_components, - fe_eval, - inverse_coefficients, - dyadic_coefficients, - in_array, - out_array); + internal::CellwiseInverseMassMatrixImplFlexible:: + template run( + n_actual_components, + fe_eval, + ArrayView(inverse_coefficients.data(), + inverse_coefficients.size()), + false, + in_array, + out_array); else internal::CellwiseInverseMassFactory::apply( n_actual_components, fe_eval, - inverse_coefficients, - dyadic_coefficients, + ArrayView(inverse_coefficients.data(), + inverse_coefficients.size()), + false, in_array, out_array); } @@ -1150,24 +1153,29 @@ namespace MatrixFreeOperators const VectorizedArrayType *in_array, VectorizedArrayType * out_array) const { - const unsigned int dofs_per_component = inverse_dyadic_coefficients.size(); - constexpr unsigned int n_tensor_components = n_components * n_components; - - AlignedVector inverse_coefficients( - dofs_per_component * n_tensor_components); - - // Flatten the inverse dyadic coefficients into `inverse_coefficients` - { - auto begin = inverse_coefficients.begin(); - for (unsigned int q = 0; q < dofs_per_component; ++q) - { - const auto end = std::next(begin, n_tensor_components); - inverse_dyadic_coefficients[q].unroll(begin, end); - begin = end; - } - } + const unsigned int unrolled_size = + inverse_dyadic_coefficients.size() * (n_components * n_components); - apply(inverse_coefficients, n_components, in_array, out_array, true); + if (fe_degree > -1) + internal::CellwiseInverseMassMatrixImplFlexible:: + template run(n_components, + fe_eval, + ArrayView( + &inverse_dyadic_coefficients[0][0][0], + unrolled_size), + true, + in_array, + out_array); + else + internal::CellwiseInverseMassFactory::apply( + n_components, + fe_eval, + ArrayView( + &inverse_dyadic_coefficients[0][0][0], unrolled_size), + true, + in_array, + out_array); } -- 2.39.5