From 7269531bf8bfd0ddf6f478978b2135986edcf1a2 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sat, 28 Nov 2020 18:50:38 +0100 Subject: [PATCH] Use ShapeInfo::shape_values_face for evaluation/integration --- .../deal.II/matrix_free/evaluation_kernels.h | 147 ++++++++++++++++-- 1 file changed, 138 insertions(+), 9 deletions(-) diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index 7634a69293..2fe6ce3dee 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -2546,8 +2546,70 @@ namespace internal const unsigned int face_orientation, const Table<2, unsigned int> &orientation_map) { - Assert(data.element_type != MatrixFreeFunctions::tensor_none, - ExcNotImplemented()); + if (data.element_type == MatrixFreeFunctions::tensor_none) + { + const unsigned int n_dofs = data.dofs_per_component_on_cell; + const unsigned int n_q_points = data.n_q_points_face; + const auto shape_info = data.data.front(); + + using Eval = EvaluatorTensorProduct; + + if (evaluate_values) + { + const auto shape_values = + &shape_info.shape_values_face(face_no, face_orientation, 0); + + auto values_quad_ptr = values_quad; + auto values_dofs_actual_ptr = values_array; + + Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); + for (unsigned int c = 0; c < n_components; ++c) + { + eval.template values<0, true, false>(values_dofs_actual_ptr, + values_quad_ptr); + + values_quad_ptr += n_q_points; + values_dofs_actual_ptr += n_dofs; + } + } + + if (evaluate_gradients) + { + auto gradients_quad_ptr = gradients_quad; + auto values_dofs_actual_ptr = values_array; + + std::array shape_gradients; + for (unsigned int d = 0; d < dim; ++d) + shape_gradients[d] = &shape_info.shape_gradients_face( + face_no, face_orientation, d, 0); + + for (unsigned int c = 0; c < n_components; ++c) + { + for (unsigned int d = 0; d < dim; ++d) + { + Eval eval(nullptr, + shape_gradients[d], + nullptr, + n_dofs, + n_q_points); + + eval.template gradients<0, true, false>( + values_dofs_actual_ptr, gradients_quad_ptr); + + gradients_quad_ptr += n_q_points; + } + values_dofs_actual_ptr += n_dofs; + } + } + + + return true; + } constexpr unsigned int static_dofs_per_face = fe_degree > -1 ? Utilities::pow(fe_degree + 1, dim - 1) : @@ -2639,8 +2701,74 @@ namespace internal const unsigned int face_orientation, const Table<2, unsigned int> &orientation_map) { - Assert(data.element_type != MatrixFreeFunctions::tensor_none, - ExcNotImplemented()); + if (data.element_type == MatrixFreeFunctions::tensor_none) + { + const unsigned int n_dofs = data.dofs_per_component_on_cell; + const unsigned int n_q_points = data.n_q_points_face; + const auto shape_info = data.data.front(); + + using Eval = EvaluatorTensorProduct; + + if (integrate_values) + { + const auto shape_values = + &shape_info.shape_values_face(face_no, face_orientation, 0); + + auto values_quad_ptr = values_quad; + auto values_dofs_actual_ptr = values_array; + + Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); + for (unsigned int c = 0; c < n_components; ++c) + { + eval.template values<0, false, false>(values_quad_ptr, + values_dofs_actual_ptr); + + values_quad_ptr += n_q_points; + values_dofs_actual_ptr += n_dofs; + } + } + + if (integrate_gradients) + { + auto gradients_quad_ptr = gradients_quad; + auto values_dofs_actual_ptr = values_array; + + std::array shape_gradients; + for (unsigned int d = 0; d < dim; ++d) + shape_gradients[d] = &shape_info.shape_gradients_face( + face_no, face_orientation, d, 0); + + for (unsigned int c = 0; c < n_components; ++c) + { + for (unsigned int d = 0; d < dim; ++d) + { + Eval eval(nullptr, + shape_gradients[d], + nullptr, + n_dofs, + n_q_points); + + if ((integrate_values == false) && d == 0) + eval.template gradients<0, false, false>( + gradients_quad_ptr, values_dofs_actual_ptr); + else + eval.template gradients<0, false, true>( + gradients_quad_ptr, values_dofs_actual_ptr); + + gradients_quad_ptr += n_q_points; + } + values_dofs_actual_ptr += n_dofs; + } + } + + + return true; + } if (face_orientation) adjust_for_face_orientation(dim, @@ -3479,9 +3607,10 @@ namespace internal const Table<2, unsigned int> &orientation_map) { if (src_ptr == nullptr) - { - return false; - } + return false; + + if (data.element_type == MatrixFreeFunctions::tensor_none) + return false; (void)sm_ptr; @@ -3735,10 +3864,10 @@ namespace internal { (void)sm_ptr; - if (dst_ptr == nullptr) + if (dst_ptr == nullptr || + data.element_type == MatrixFreeFunctions::tensor_none) { AssertDimension(n_face_orientations, 1); - AssertDimension(n_face_orientations, 1); // for block vectors simply integrate FEFaceEvaluationImplIntegrateSelector:: -- 2.39.5