From 856007a732bef1484b2eb540714e9c78b2ad3c08 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 6 Feb 2019 12:53:59 +0100 Subject: [PATCH] Do not use transformation to collocation for large n_q_points. --- .../deal.II/matrix_free/evaluation_kernels.h | 14 +++++- .../deal.II/matrix_free/evaluation_selector.h | 44 ++++++++++++++----- 2 files changed, 46 insertions(+), 12 deletions(-) diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index 7cfa435b6a..3f47fbe06e 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -1337,6 +1337,16 @@ namespace internal typename Number> struct FEFaceEvaluationImpl { + // We enable a transformation to collocation for derivatives if it gives + // correct results (first two conditions), if it is the most efficient + // choice in terms of operation counts (third condition) and if we were + // able to initialize the fields in shape_info.templates.h from the + // polynomials (fourth condition). + static constexpr bool use_collocation = + symmetric_evaluate && + n_q_points_1d > fe_degree &&n_q_points_1d <= 3 * fe_degree / 2 + 1 && + n_q_points_1d < 200; + static void evaluate_in_face(const MatrixFreeFunctions::ShapeInfo &data, Number * values_dofs, @@ -1431,7 +1441,7 @@ namespace internal switch (dim) { case 3: - if (symmetric_evaluate && n_q_points_1d > fe_degree) + if (use_collocation) { eval1.template values<0, true, false>(values_dofs, values_quad); @@ -1590,7 +1600,7 @@ namespace internal 2 * n_q_points); eval1.template values<0, false, false>( gradients_quad + 2 * n_q_points, values_dofs + size_deg); - if (symmetric_evaluate && n_q_points_1d > fe_degree) + if (use_collocation) { internal::EvaluatorTensorProduct< internal::evaluate_evenodd, diff --git a/include/deal.II/matrix_free/evaluation_selector.h b/include/deal.II/matrix_free/evaluation_selector.h index 2fd382aa06..30115e481b 100644 --- a/include/deal.II/matrix_free/evaluation_selector.h +++ b/include/deal.II/matrix_free/evaluation_selector.h @@ -251,6 +251,17 @@ namespace internal n_q_points_1d, typename std::enable_if<(n_q_points_1d < degree + 3)>::type> { + /** + * We enable a transformation to collocation for derivatives if it gives + * correct results (first condition), if it is the most efficient choice + * in terms of operation counts (second condition) and if we were able + * to initialize the fields in shape_info.templates.h from the + * polynomials (third condition). + */ + static constexpr bool use_collocation = + n_q_points_1d > degree &&n_q_points_1d <= 3 * degree / 2 + 1 && + n_q_points_1d < 200; + static inline void evaluate( const internal::MatrixFreeFunctions::ShapeInfo &shape_info, @@ -280,7 +291,7 @@ namespace internal evaluate_values, evaluate_gradients, evaluate_hessians); - else if (degree < n_q_points_1d) + else if (use_collocation) internal::FEEvaluationImplTransformToCollocation< dim, degree, @@ -352,7 +363,7 @@ namespace internal integrate_values, integrate_gradients, sum_into_values_array); - else if (degree < n_q_points_1d) + else if (use_collocation) internal::FEEvaluationImplTransformToCollocation< dim, degree, @@ -481,6 +492,17 @@ template struct SelectEvaluator { + /** + * We enable a transformation to collocation for derivatives if it gives + * correct results (first condition), if it is the most efficient choice in + * terms of operation counts (second condition) and if we were able to + * initialize the fields in shape_info.templates.h from the polynomials + * (third condition). + */ + static constexpr bool use_collocation = + n_q_points_1d > fe_degree &&n_q_points_1d <= 3 * fe_degree / 2 + 1 && + n_q_points_1d < 200; + /** * Chooses an appropriate evaluation strategy for the evaluate function, i.e. * this calls internal::FEEvaluationImpl::evaluate(), @@ -606,9 +628,10 @@ SelectEvaluator::evaluate( evaluate_gradients, evaluate_hessians); } - else if (fe_degree < n_q_points_1d && - shape_info.element_type <= - internal::MatrixFreeFunctions::tensor_symmetric) + // '<=' on type means tensor_symmetric or tensor_symmetric_hermite, see + // shape_info.h for more details + else if (use_collocation && shape_info.element_type <= + internal::MatrixFreeFunctions::tensor_symmetric) { internal::FEEvaluationImplTransformToCollocation< dim, @@ -625,7 +648,7 @@ SelectEvaluator::evaluate( evaluate_gradients, evaluate_hessians); } - else if (shape_info.element_type == + else if (shape_info.element_type <= internal::MatrixFreeFunctions::tensor_symmetric) { internal::FEEvaluationImpl< @@ -739,9 +762,10 @@ SelectEvaluator::integrate( integrate_gradients, sum_into_values_array); } - else if (fe_degree < n_q_points_1d && - shape_info.element_type <= - internal::MatrixFreeFunctions::tensor_symmetric) + // '<=' on type means tensor_symmetric or tensor_symmetric_hermite, see + // shape_info.h for more details + else if (use_collocation && shape_info.element_type <= + internal::MatrixFreeFunctions::tensor_symmetric) { internal::FEEvaluationImplTransformToCollocation< dim, @@ -757,7 +781,7 @@ SelectEvaluator::integrate( integrate_gradients, sum_into_values_array); } - else if (shape_info.element_type == + else if (shape_info.element_type <= internal::MatrixFreeFunctions::tensor_symmetric) { internal::FEEvaluationImpl< -- 2.39.5