From: Peter Munch Date: Sat, 21 Nov 2020 15:09:04 +0000 (+0100) Subject: Use EvaluatorTensorProduct in Simplex code path X-Git-Tag: v9.3.0-rc1~866^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ddd832d14477a2254feb5b8809d4296482076c22;p=dealii.git Use EvaluatorTensorProduct in Simplex code path --- diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index 036b628f59..bfa6319cda 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -696,35 +696,50 @@ namespace internal const unsigned int n_dofs = shape_info.dofs_per_component_on_cell; const unsigned int n_q_points = shape_info.n_q_points; + using Eval = + EvaluatorTensorProduct; + if (evaluation_flag & EvaluationFlags::values) { const auto shape_values = shape_info.data.front().shape_values.data(); + auto values_quad_ptr = values_quad; + auto values_dofs_actual_ptr = values_dofs_actual; + + Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); for (unsigned int c = 0; c < n_components; ++c) - for (unsigned int q = 0; q < n_q_points; ++q) - { - values_quad[n_q_points * c + q] = 0.0; - for (unsigned int i = 0; i < n_dofs; ++i) - values_quad[n_q_points * c + q] += - shape_values[i * n_q_points + q] * - values_dofs_actual[n_dofs * c + i]; - } + { + 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 (evaluation_flag & EvaluationFlags::gradients) { const auto shape_gradients = shape_info.data.front().shape_gradients.data(); + auto gradients_quad_ptr = gradients_quad; + auto values_dofs_actual_ptr = values_dofs_actual; + for (unsigned int c = 0; c < n_components; ++c) - for (unsigned int d = 0; d < dim; ++d) - for (unsigned int q = 0; q < n_q_points; ++q) + { + for (unsigned int d = 0; d < dim; ++d) { - gradients_quad[n_q_points * dim * c + n_q_points * d + q] = 0.0; - for (unsigned int i = 0; i < n_dofs; ++i) - gradients_quad[n_q_points * dim * c + n_q_points * d + q] += - shape_gradients[d * n_dofs * n_q_points + i * n_q_points + - q] * - values_dofs_actual[n_dofs * dim * c + i]; + Eval eval(nullptr, + shape_gradients + n_q_points * n_dofs * 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; + } } if (evaluation_flag & EvaluationFlags::hessians) @@ -757,33 +772,60 @@ namespace internal const unsigned int n_dofs = shape_info.dofs_per_component_on_cell; const unsigned int n_q_points = shape_info.n_q_points; - if (add_into_values_array == false) - for (unsigned int i = 0; i < n_components * n_dofs; ++i) - values_dofs_actual[i] = 0.0; + using Eval = + EvaluatorTensorProduct; if (integration_flag & EvaluationFlags::values) { const auto shape_values = shape_info.data.front().shape_values.data(); + auto values_quad_ptr = values_quad; + auto values_dofs_actual_ptr = values_dofs_actual; + + Eval eval(shape_values, nullptr, nullptr, n_dofs, n_q_points); for (unsigned int c = 0; c < n_components; ++c) - for (unsigned int q = 0; q < n_q_points; ++q) - for (unsigned int i = 0; i < n_dofs; ++i) - values_dofs_actual[n_dofs * c + i] += - shape_values[i * n_q_points + q] * - values_quad[n_q_points * c + q]; + { + if (add_into_values_array == false) + eval.template values<0, false, false>(values_quad_ptr, + values_dofs_actual_ptr); + else + eval.template values<0, false, true>(values_quad_ptr, + values_dofs_actual_ptr); + + values_quad_ptr += n_q_points; + values_dofs_actual_ptr += n_dofs; + } } if (integration_flag & EvaluationFlags::gradients) { const auto shape_gradients = shape_info.data.front().shape_gradients.data(); + auto gradients_quad_ptr = gradients_quad; + auto values_dofs_actual_ptr = values_dofs_actual; + for (unsigned int c = 0; c < n_components; ++c) - for (unsigned int d = 0; d < dim; ++d) - for (unsigned int q = 0; q < n_q_points; ++q) - for (unsigned int i = 0; i < n_dofs; ++i) - values_dofs_actual[n_dofs * dim * c + i] += - shape_gradients[d * n_dofs * n_q_points + i * n_q_points + - q] * - gradients_quad[n_q_points * dim * c + n_q_points * d + q]; + { + for (unsigned int d = 0; d < dim; ++d) + { + Eval eval(nullptr, + shape_gradients + n_q_points * n_dofs * d, + nullptr, + n_dofs, + n_q_points); + + if ((add_into_values_array == false && + (integration_flag & EvaluationFlags::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; + } } } diff --git a/include/deal.II/matrix_free/tensor_product_kernels.h b/include/deal.II/matrix_free/tensor_product_kernels.h index 166d7285b7..21ed078d10 100644 --- a/include/deal.II/matrix_free/tensor_product_kernels.h +++ b/include/deal.II/matrix_free/tensor_product_kernels.h @@ -625,6 +625,21 @@ namespace internal ExcDimensionMismatch(shape_hessians.size(), n_rows * n_columns)); } + /** + * Constructor, taking the data from ShapeInfo + */ + EvaluatorTensorProduct(const Number2 * shape_values, + const Number2 * shape_gradients, + const Number2 * shape_hessians, + const unsigned int n_rows, + const unsigned int n_columns) + : shape_values(shape_values) + , shape_gradients(shape_gradients) + , shape_hessians(shape_hessians) + , n_rows(n_rows) + , n_columns(n_columns) + {} + template void values(const Number *in, Number *out) const