Add dim component test and fix strides

diff --git a/include/deal.II/matrix_free/fe_point_evaluation.h b/include/deal.II/matrix_free/fe_point_evaluation.h
index 57161ae20bff91d817af68d0700482ea0f6f9422..da845ef1bd97b2746c941bd277a74a47d3317923 100644 (file)
--- 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<n_components_, dim, spacedim, Number>::
         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<n_components_, dim, spacedim, Number>::
                 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<n_components_, dim, spacedim, Number>::
         }
       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 a0bc335459e4f8f0a4826fa4f479737b765e5adc..8833e9ca664c631f557da389efdee58c290d6dcf 100644 (file)
--- a/tests/matrix_free/point_evaluation_29.cc
+++ b/tests/matrix_free/point_evaluation_29.cc
@@ -34,7 +34,7 @@
 
 #include <deal.II/numerics/vector_tools.h>
 
-#include <../tests.h>
+#include "../tests.h"
 
 using namespace dealii;
 
@@ -42,9 +42,10 @@ template <int dim, typename Number>
 void
 do_test(unsigned int degree)
 {
-  using VectorType = LinearAlgebra::distributed::Vector<Number>;
+  using VectorType          = LinearAlgebra::distributed::Vector<Number>;
+  using VectorizedArrayType = VectorizedArray<Number>;
 
-  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  Triangulation<dim> tria;
   GridGenerator::hyper_cube(tria);
   tria.refine_global(2);
 
@@ -54,14 +55,14 @@ do_test(unsigned int degree)
   AffineConstraints<Number> constraints;
   constraints.close();
 
-  typename MatrixFree<dim, Number>::AdditionalData data;
+  typename MatrixFree<dim, Number, VectorizedArrayType>::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<dim> mapping;
 
-  MatrixFree<dim, Number> matrix_free;
+  MatrixFree<dim, Number, VectorizedArrayType> matrix_free;
   matrix_free.reinit(
     mapping, dof_handler, constraints, QGauss<dim>(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<dim - 1>(degree + 1));
+          quadrature_vector.push_back(QGauss<dim - 1>(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<dim, dim, dim, Number> 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<dim, -1, 0, dim, Number> phi(matrix_free, true, 0, 0, 1);
+      FEFaceEvaluation<dim, -1, 0, dim, Number, VectorizedArrayType> phi(
+        matrix_free, true, 0, 0, 1);
       FEFacePointEvaluation<dim, dim, dim, Number> 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<VectorizedArrayType::size()>(
+                &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<VectorizedArrayType::size()>(
+                &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<VectorType, VectorType>(
     {}, {}, 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 (file)
index 0000000..0fd8fc1
--- /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 (file)
index be8d055..0000000
--- a/tests/matrix_free/point_evaluation_29.with_p4est=true.mpirun=1.output
+++ /dev/null
@@ -1,2 +0,0 @@
-
-DEAL:0::OK