From: Daniel Arndt Date: Mon, 22 Apr 2024 23:54:07 +0000 (-0400) Subject: Portable tensor product kernels: Simplify loop generic X-Git-Tag: v9.6.0-rc1~346^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6f40a28d09529c50f02e645960ad926627ba945e;p=dealii.git Portable tensor product kernels: Simplify loop generic --- diff --git a/include/deal.II/matrix_free/portable_tensor_product_kernels.h b/include/deal.II/matrix_free/portable_tensor_product_kernels.h index 7d94f515fb..f34bc2ff91 100644 --- a/include/deal.II/matrix_free/portable_tensor_product_kernels.h +++ b/include/deal.II/matrix_free/portable_tensor_product_kernels.h @@ -288,20 +288,22 @@ namespace Portable // This loop simply multiplies the shape function at the quadrature // point by the value finite element coefficient. - t[q_point] = 0; - for (int k = 0; k < n_q_points_1d; ++k) + const int stride_shape = dof_to_quad ? n_q_points_1d : 1; + const int stride = Utilities::pow(n_q_points_1d, direction); + const int base_shape = dof_to_quad ? q : (q * n_q_points_1d); + const int base = + (direction == 0) ? (n_q_points_1d * (i + n_q_points_1d * j)) : + (direction == 1) ? (i + n_q_points_1d * (n_q_points_1d * j)) : + (i + n_q_points_1d * j); + Number sum = + shape_data[base_shape] * (in_place ? out(base) : in(base)); + for (int k = 1; k < n_q_points_1d; ++k) { - const unsigned int shape_idx = - dof_to_quad ? (q + k * n_q_points_1d) : (k + q * n_q_points_1d); - const unsigned int source_idx = - (direction == 0) ? - (k + n_q_points_1d * (i + n_q_points_1d * j)) : - (direction == 1) ? - (i + n_q_points_1d * (k + n_q_points_1d * j)) : - (i + n_q_points_1d * (j + n_q_points_1d * k)); - t[q_point] += shape_data[shape_idx] * - (in_place ? out(source_idx) : in(source_idx)); + sum += + shape_data[base_shape + k * stride_shape] * + (in_place ? out(base + k * stride) : in(base + k * stride)); } + t[q_point] = sum; }); if (in_place)