From 1cf310c4611dd9d91ad761f5bb9513428b6201c5 Mon Sep 17 00:00:00 2001 From: Maximilian Bergbauer Date: Fri, 5 Jan 2024 18:05:12 +0100 Subject: [PATCH] Add dim component test and fix strides --- .../deal.II/matrix_free/fe_point_evaluation.h | 49 ++++++++++--------- tests/matrix_free/point_evaluation_29.cc | 32 ++++++------ tests/matrix_free/point_evaluation_29.output | 2 + ...luation_29.with_p4est=true.mpirun=1.output | 2 - 4 files changed, 46 insertions(+), 39 deletions(-) create mode 100644 tests/matrix_free/point_evaluation_29.output delete mode 100644 tests/matrix_free/point_evaluation_29.with_p4est=true.mpirun=1.output diff --git a/include/deal.II/matrix_free/fe_point_evaluation.h b/include/deal.II/matrix_free/fe_point_evaluation.h index 57161ae20b..da845ef1bd 100644 --- a/include/deal.II/matrix_free/fe_point_evaluation.h +++ b/include/deal.II/matrix_free/fe_point_evaluation.h @@ -3420,13 +3420,17 @@ FEFacePointEvaluation:: is_linear ? Utilities::pow(2, dim - 1) : this->dofs_per_component_face; for (unsigned int comp = 0; comp < n_components; ++comp) for (unsigned int i = 0; i < 2 * dofs_per_comp_face; ++i) - ETT::read_value(face_dof_values[i + comp * 3 * dofs_per_comp_face], + ETT::read_value(face_dof_values[(i + comp * 3 * dofs_per_comp_face) * + stride_face_dof], comp, this->solution_renumbered[i]); face_dof_values_ptr = this->solution_renumbered.data(); } + constexpr int stride_face_dof_actual = + n_components == 1 ? stride_face_dof : 1; + // loop over quadrature batches qb const unsigned int n_shapes = is_linear ? 2 : this->poly.size(); @@ -3444,17 +3448,18 @@ FEFacePointEvaluation:: scalar_value_type, VectorizedArrayType, 2, - stride_face_dof>(face_dof_values_ptr, - this->unit_point_faces_ptr[qb]) : + stride_face_dof_actual>(face_dof_values_ptr, + this->unit_point_faces_ptr[qb]) : internal::evaluate_tensor_product_value_and_gradient_shapes< dim - 1, scalar_value_type, VectorizedArrayType, 2, false, - stride_face_dof>(this->shapes_faces.data() + qb * n_shapes, - n_shapes, - face_dof_values_ptr); + stride_face_dof_actual>(this->shapes_faces.data() + + qb * n_shapes, + n_shapes, + face_dof_values_ptr); value = interpolated_value[dim - 1]; // reorder derivative from tangential/normal derivatives into tensor @@ -3497,22 +3502,22 @@ FEFacePointEvaluation:: } else { - value = - is_linear ? - internal::evaluate_tensor_product_value_linear< - dim - 1, - scalar_value_type, - VectorizedArrayType, - stride_face_dof>(face_dof_values_ptr, - this->unit_point_faces_ptr[qb]) : - internal::evaluate_tensor_product_value_shapes< - dim - 1, - scalar_value_type, - VectorizedArrayType, - false, - stride_face_dof>(this->shapes_faces.data() + qb * n_shapes, - n_shapes, - face_dof_values_ptr); + value = is_linear ? + internal::evaluate_tensor_product_value_linear< + dim - 1, + scalar_value_type, + VectorizedArrayType, + stride_face_dof_actual>(face_dof_values_ptr, + this->unit_point_faces_ptr[qb]) : + internal::evaluate_tensor_product_value_shapes< + dim - 1, + scalar_value_type, + VectorizedArrayType, + false, + stride_face_dof_actual>(this->shapes_faces.data() + + qb * n_shapes, + n_shapes, + face_dof_values_ptr); } if (evaluation_flags & EvaluationFlags::values) diff --git a/tests/matrix_free/point_evaluation_29.cc b/tests/matrix_free/point_evaluation_29.cc index a0bc335459..8833e9ca66 100644 --- a/tests/matrix_free/point_evaluation_29.cc +++ b/tests/matrix_free/point_evaluation_29.cc @@ -34,7 +34,7 @@ #include -#include <../tests.h> +#include "../tests.h" using namespace dealii; @@ -42,9 +42,10 @@ template void do_test(unsigned int degree) { - using VectorType = LinearAlgebra::distributed::Vector; + using VectorType = LinearAlgebra::distributed::Vector; + using VectorizedArrayType = VectorizedArray; - parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + Triangulation tria; GridGenerator::hyper_cube(tria); tria.refine_global(2); @@ -54,14 +55,14 @@ do_test(unsigned int degree) AffineConstraints constraints; constraints.close(); - typename MatrixFree::AdditionalData data; + typename MatrixFree::AdditionalData data; data.mapping_update_flags = update_values; data.mapping_update_flags_inner_faces = update_values; data.mapping_update_flags_boundary_faces = update_values; MappingQ1 mapping; - MatrixFree matrix_free; + MatrixFree matrix_free; matrix_free.reinit( mapping, dof_handler, constraints, QGauss(degree + 1), data); @@ -96,7 +97,7 @@ do_test(unsigned int degree) vector_face_accessors.push_back( matrix_free.get_face_iterator(face_batch, 0)); - quadrature_vector.emplace_back(Quadrature(degree + 1)); + quadrature_vector.push_back(QGauss(degree + 1)); } } @@ -135,7 +136,7 @@ do_test(unsigned int degree) src); const auto boundary_function_1 = [&](const auto &data, - auto & dst, + auto &dst, const auto &src, const auto face_range) { FEFacePointEvaluation phi_pnt( @@ -168,7 +169,8 @@ do_test(unsigned int degree) const auto boundary_function_2 = [&](const auto &data, auto &dst, const auto &src, const auto face_range) { - FEFaceEvaluation phi(matrix_free, true, 0, 0, 1); + FEFaceEvaluation phi( + matrix_free, true, 0, 0, 1); FEFacePointEvaluation phi_pnt( nm_mapping_info, data.get_dof_handler().get_fe(), true, 1); @@ -183,14 +185,14 @@ do_test(unsigned int degree) ++v) { phi_pnt.reinit(face * n_lanes + v); - phi_pnt.evaluate_in_face(&phi.get_scratch_data().begin()[0][v], - EvaluationFlags::values); + phi_pnt.template evaluate_in_face( + &phi.get_scratch_data().begin()[0][v], EvaluationFlags::values); for (unsigned int q : phi_pnt.quadrature_point_indices()) phi_pnt.submit_value(phi_pnt.get_value(q), q); - phi_pnt.integrate_in_face(&phi.get_scratch_data().begin()[0][v], - EvaluationFlags::values); + phi_pnt.template integrate_in_face( + &phi.get_scratch_data().begin()[0][v], EvaluationFlags::values); } phi.collect_from_face(EvaluationFlags::values); phi.distribute_local_to_global(dst); @@ -204,14 +206,14 @@ do_test(unsigned int degree) matrix_free.template loop( {}, {}, boundary_function_2, dst_2, src, true); - Assert(std::abs(dst_1.l2_norm() - dst_2.l2_norm()) < 1.0e-12); + Assert(std::abs(dst_1.l2_norm() - dst_2.l2_norm()) < 1.0e-12, + ExcInternalError()); } int main(int argc, char *argv[]) { - Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1); - MPILogInitAll all; + initlog(); do_test<2, double>(2); diff --git a/tests/matrix_free/point_evaluation_29.output b/tests/matrix_free/point_evaluation_29.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/matrix_free/point_evaluation_29.output @@ -0,0 +1,2 @@ + +DEAL::OK diff --git a/tests/matrix_free/point_evaluation_29.with_p4est=true.mpirun=1.output b/tests/matrix_free/point_evaluation_29.with_p4est=true.mpirun=1.output deleted file mode 100644 index be8d055f86..0000000000 --- a/tests/matrix_free/point_evaluation_29.with_p4est=true.mpirun=1.output +++ /dev/null @@ -1,2 +0,0 @@ - -DEAL:0::OK -- 2.39.5