]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement FEFaceEvaluation::evaluate() for ECL/unstructured meshes 13063/head
authorPeter Munch <peterrmuench@gmail.com>
Sun, 19 Dec 2021 13:02:01 +0000 (14:02 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Sun, 24 Apr 2022 07:18:42 +0000 (09:18 +0200)
include/deal.II/matrix_free/evaluation_kernels.h
include/deal.II/matrix_free/mapping_info.templates.h
tests/matrix_free/ecl.h
tests/matrix_free/ecl_01.cc
tests/matrix_free/ecl_01.with_p4est=true.mpirun=1.output
tests/matrix_free/ecl_03.cc
tests/matrix_free/ecl_04.cc

index e2dc5a51472648fcbaf53083570c8e7a4297d4c4..252c47075b960e4d8760f173b96ea10571eb9afc 100644 (file)
@@ -2861,15 +2861,20 @@ namespace internal
         const Number *                         values_dofs,
         FEEvaluationData<dim, Number, true> &  fe_eval)
     {
-      const auto &       shape_info       = fe_eval.get_shape_info();
-      const auto &       shape_data       = shape_info.data.front();
-      const unsigned int face_no          = fe_eval.get_face_no();
-      const unsigned int face_orientation = fe_eval.get_face_orientation();
+      const auto &shape_info = fe_eval.get_shape_info();
+      const auto &shape_data = shape_info.data.front();
 
       if (shape_info.element_type == MatrixFreeFunctions::tensor_none)
         {
-          const std::size_t n_dofs     = shape_info.dofs_per_component_on_cell;
-          const std::size_t n_q_points = shape_info.n_q_points_faces[face_no];
+          Assert((fe_eval.get_dof_access_index() ==
+                    MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+                  fe_eval.is_interior_face() == false) == false,
+                 ExcNotImplemented());
+
+          const unsigned int face_no          = fe_eval.get_face_no();
+          const unsigned int face_orientation = fe_eval.get_face_orientation();
+          const std::size_t  n_dofs     = shape_info.dofs_per_component_on_cell;
+          const std::size_t  n_q_points = shape_info.n_q_points_faces[face_no];
 
           using Eval =
             EvaluatorTensorProduct<evaluate_general, 1, 0, 0, Number, Number>;
@@ -2937,13 +2942,56 @@ namespace internal
       Number *temp         = fe_eval.get_scratch_data().begin();
       Number *scratch_data = temp + 3 * n_components * dofs_per_face;
 
-      FEFaceNormalEvaluationImpl<dim, fe_degree, Number>::
-        template interpolate<true, false>(n_components,
-                                          evaluation_flag,
-                                          shape_info,
-                                          values_dofs,
-                                          temp,
-                                          face_no);
+      bool use_vectorization = true;
+
+      if (fe_eval.get_dof_access_index() ==
+            MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+          fe_eval.is_interior_face() == false) // exterior faces in the ECL loop
+        use_vectorization =
+          fe_eval.get_cell_ids()[0] != numbers::invalid_unsigned_int &&
+          std::all_of(fe_eval.get_cell_ids().begin() + 1,
+                      fe_eval.get_cell_ids().end(),
+                      [&](const auto &v) {
+                        return v == fe_eval.get_cell_ids()[0] ||
+                               v == numbers::invalid_unsigned_int;
+                      });
+
+      if (use_vectorization == false)
+        {
+          for (unsigned int v = 0; v < Number::size(); ++v)
+            {
+              // the loop breaks once an invalid_unsigned_int is hit for
+              // all cases except the exterior faces in the ECL loop (where
+              // some faces might be at the boundaries but others not)
+              if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
+                {
+                  for (unsigned int i = 0; i < 3 * n_components * dofs_per_face;
+                       ++i)
+                    temp[i][v] = 0;
+                  continue;
+                }
+
+              FEFaceNormalEvaluationImpl<dim, fe_degree, Number>::
+                template interpolate<true, false>(n_components,
+                                                  evaluation_flag,
+                                                  shape_info,
+                                                  values_dofs,
+                                                  scratch_data,
+                                                  fe_eval.get_face_no(v));
+
+              for (unsigned int i = 0; i < 3 * n_components * dofs_per_face;
+                   ++i)
+                temp[i][v] = scratch_data[i][v];
+            }
+        }
+      else
+        FEFaceNormalEvaluationImpl<dim, fe_degree, Number>::
+          template interpolate<true, false>(n_components,
+                                            evaluation_flag,
+                                            shape_info,
+                                            values_dofs,
+                                            temp,
+                                            fe_eval.get_face_no());
 
       const unsigned int     subface_index = fe_eval.get_subface_index();
       constexpr unsigned int n_q_points_1d_actual =
@@ -2982,12 +3030,39 @@ namespace internal
                                                        scratch_data,
                                                        subface_index);
 
-      if (face_orientation)
+      if (use_vectorization == false)
+        {
+          for (unsigned int v = 0; v < Number::size(); ++v)
+            {
+              // the loop breaks once an invalid_unsigned_int is hit for
+              // all cases except the exterior faces in the ECL loop (where
+              // some faces might be at the boundaries but others not)
+              if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
+                continue;
+
+              if (fe_eval.get_face_orientation(v) != 0)
+                adjust_for_face_orientation_per_lane(
+                  dim,
+                  n_components,
+                  v,
+                  evaluation_flag,
+                  &fe_eval.get_shape_info().face_orientations_quad(
+                    fe_eval.get_face_orientation(v), 0),
+                  false,
+                  shape_info.n_q_points_face,
+                  &temp[0][0],
+                  fe_eval.begin_values(),
+                  fe_eval.begin_gradients(),
+                  fe_eval.begin_hessians());
+            }
+        }
+      else if (fe_eval.get_face_orientation() != 0)
         adjust_for_face_orientation(
           dim,
           n_components,
           evaluation_flag,
-          &fe_eval.get_shape_info().face_orientations_quad(face_orientation, 0),
+          &fe_eval.get_shape_info().face_orientations_quad(
+            fe_eval.get_face_orientation(), 0),
           false,
           shape_info.n_q_points_face,
           temp,
@@ -3011,15 +3086,20 @@ namespace internal
         Number *                               values_dofs,
         FEEvaluationData<dim, Number, true> &  fe_eval)
     {
-      const auto &       shape_info       = fe_eval.get_shape_info();
-      const auto &       shape_data       = shape_info.data.front();
-      const unsigned int face_no          = fe_eval.get_face_no();
-      const unsigned int face_orientation = fe_eval.get_face_orientation();
+      const auto &shape_info = fe_eval.get_shape_info();
+      const auto &shape_data = shape_info.data.front();
 
       if (shape_info.element_type == MatrixFreeFunctions::tensor_none)
         {
-          const std::size_t n_dofs     = shape_info.dofs_per_component_on_cell;
-          const std::size_t n_q_points = shape_info.n_q_points_faces[face_no];
+          Assert((fe_eval.get_dof_access_index() ==
+                    MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+                  fe_eval.is_interior_face() == false) == false,
+                 ExcNotImplemented());
+
+          const unsigned int face_no          = fe_eval.get_face_no();
+          const unsigned int face_orientation = fe_eval.get_face_orientation();
+          const std::size_t  n_dofs     = shape_info.dofs_per_component_on_cell;
+          const std::size_t  n_q_points = shape_info.n_q_points_faces[face_no];
 
           using Eval =
             EvaluatorTensorProduct<evaluate_general, 1, 0, 0, Number, Number>;
@@ -3090,12 +3170,53 @@ namespace internal
       Number *temp         = fe_eval.get_scratch_data().begin();
       Number *scratch_data = temp + 3 * n_components * dofs_per_face;
 
-      if (face_orientation)
+      bool use_vectorization = true;
+
+      if (fe_eval.get_dof_access_index() ==
+            MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+          fe_eval.is_interior_face() == false) // exterior faces in the ECL loop
+        use_vectorization =
+          fe_eval.get_cell_ids()[0] != numbers::invalid_unsigned_int &&
+          std::all_of(fe_eval.get_cell_ids().begin() + 1,
+                      fe_eval.get_cell_ids().end(),
+                      [&](const auto &v) {
+                        return v == fe_eval.get_cell_ids()[0] ||
+                               v == numbers::invalid_unsigned_int;
+                      });
+
+      if (use_vectorization == false)
+        {
+          for (unsigned int v = 0; v < Number::size(); ++v)
+            {
+              // the loop breaks once an invalid_unsigned_int is hit for
+              // all cases except the exterior faces in the ECL loop (where
+              // some faces might be at the boundaries but others not)
+              if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
+                continue;
+
+              if (fe_eval.get_face_orientation(v) != 0)
+                adjust_for_face_orientation_per_lane(
+                  dim,
+                  n_components,
+                  v,
+                  integration_flag,
+                  &fe_eval.get_shape_info().face_orientations_quad(
+                    fe_eval.get_face_orientation(v), 0),
+                  true,
+                  shape_info.n_q_points_face,
+                  &temp[0][0],
+                  fe_eval.begin_values(),
+                  fe_eval.begin_gradients(),
+                  fe_eval.begin_hessians());
+            }
+        }
+      else if (fe_eval.get_face_orientation() != 0)
         adjust_for_face_orientation(
           dim,
           n_components,
           integration_flag,
-          &fe_eval.get_shape_info().face_orientations_quad(face_orientation, 0),
+          &fe_eval.get_shape_info().face_orientations_quad(
+            fe_eval.get_face_orientation(), 0),
           true,
           shape_info.n_q_points_face,
           temp,
@@ -3141,13 +3262,37 @@ namespace internal
                                      scratch_data,
                                      subface_index);
 
-      FEFaceNormalEvaluationImpl<dim, fe_degree, Number>::
-        template interpolate<false, false>(n_components,
-                                           integration_flag,
-                                           shape_info,
-                                           temp,
-                                           values_dofs,
-                                           face_no);
+      if (use_vectorization == false)
+        {
+          for (unsigned int v = 0; v < Number::size(); ++v)
+            {
+              // the loop breaks once an invalid_unsigned_int is hit for
+              // all cases except the exterior faces in the ECL loop (where
+              // some faces might be at the boundaries but others not)
+              if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
+                continue;
+
+              FEFaceNormalEvaluationImpl<dim, fe_degree, Number>::
+                template interpolate<false, false>(n_components,
+                                                   integration_flag,
+                                                   shape_info,
+                                                   values_dofs,
+                                                   scratch_data,
+                                                   fe_eval.get_face_no(v));
+
+              for (unsigned int i = 0; i < 3 * n_components * dofs_per_face;
+                   ++i)
+                temp[i][v] = scratch_data[i][v];
+            }
+        }
+      else
+        FEFaceNormalEvaluationImpl<dim, fe_degree, Number>::
+          template interpolate<false, false>(n_components,
+                                             integration_flag,
+                                             shape_info,
+                                             temp,
+                                             values_dofs,
+                                             fe_eval.get_face_no());
       return false;
     }
   };
index 9d0cf83e770802ab7212e410553a47ecf7abf49f..209d4e65ac53fcff0e8de415faaec71929aab9fc 100644 (file)
@@ -3039,8 +3039,8 @@ namespace internal
                               for (unsigned int e = 0; e < dim; ++e)
                                 {
                                   const unsigned int ee = ExtractFaceHelper::
-                                    reorder_face_derivative_indices<dim>(face,
-                                                                         e);
+                                    reorder_face_derivative_indices<dim>(
+                                      fe_val_neigh.get_face_number(), e);
                                   face_data_by_cells[my_q]
                                     .jacobians[1][offset][d][e][v] =
                                     inv_jac[d][ee];
@@ -3081,6 +3081,23 @@ namespace internal
                                     inv_jac[d][ee];
                                 }
                           }
+                      if (is_local && (update_flags & update_jacobians))
+                        for (unsigned int q = 0; q < fe_val.n_quadrature_points;
+                             ++q)
+                          {
+                            DerivativeForm<1, dim, dim> inv_jac =
+                              fe_val_neigh.jacobian(q).covariant_form();
+                            for (unsigned int d = 0; d < dim; ++d)
+                              for (unsigned int e = 0; e < dim; ++e)
+                                {
+                                  const unsigned int ee = ExtractFaceHelper::
+                                    reorder_face_derivative_indices<dim>(
+                                      fe_val_neigh.get_face_number(), e);
+                                  face_data_by_cells[my_q]
+                                    .jacobians[1][offset + q][d][e][v] =
+                                    inv_jac[d][ee];
+                                }
+                          }
                       if (update_flags & update_jacobian_grads)
                         {
                           Assert(false, ExcNotImplemented());
index e877335bb0a1ee7ffa87baf3597c9a393cd5b45e..99d92576e138ba6d3e7340d372ef10d3e6f33b93 100644 (file)
@@ -51,14 +51,25 @@ template <int dim,
           typename Number              = double,
           typename VectorizedArrayType = VectorizedArray<Number>>
 void
-test(const unsigned int n_refinements = 1,
+test(const unsigned int geometry      = 0,
+     const unsigned int n_refinements = 1,
      const bool         print_vector  = true,
      const MPI_Comm     comm          = MPI_COMM_SELF)
 {
   using VectorType = LinearAlgebra::distributed::Vector<Number>;
 
   parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
-  GridGenerator::hyper_cube(tria, 0.0, 1.0, true);
+
+  if (geometry == 0)
+    GridGenerator::hyper_cube(tria, 0.0, 1.0, true);
+  else if (geometry == 1)
+    GridGenerator::hyper_shell(tria, Point<dim>(), 0.5, 1.0);
+  else if (geometry == 2)
+    GridGenerator::hyper_ball(tria);
+  else
+    Assert(false, ExcNotImplemented());
+
+  tria.reset_all_manifolds();
 
   if (false)
     {
index a335c637d5654c3e91949affeb967497206b4b66..c514146ea8d7864b5eba3b0bf13910b4792b9120 100644 (file)
@@ -34,7 +34,7 @@
 
 
 // Compare the evaluation of a SIP Laplace operator with face- and element-
-// centric loops on a structured grid.
+// centric loops on a hypercube, hypershell, and hyperball.
 
 int
 main(int argc, char **argv)
@@ -42,5 +42,7 @@ main(int argc, char **argv)
   Utilities::MPI::MPI_InitFinalize mpi_init(argc, argv, 1);
 
   mpi_initlog();
-  test<2, 2, 3, double, VectorizedArray<double>>();
+  test<2, 2, 3, double, VectorizedArray<double>>(0);
+  test<2, 2, 3, double, VectorizedArray<double>>(1);
+  test<2, 2, 3, double, VectorizedArray<double>>(2);
 }
index 8cd7af3429e442a9dab82d6fdc37bb79126c5905..18b2bc574ed5cd77ab2e37906f3f7744d35266e9 100644 (file)
@@ -9,3 +9,23 @@ Process #0
 Local range: [0, 36), global size: 36
 Vector data:
 7.000e+00 6.667e-01 -4.000e+00 1.467e+01 2.667e+00 -2.667e+00 3.333e+00 6.667e-01 -3.333e-01 -4.000e+00 6.667e-01 7.000e+00 -2.667e+00 2.667e+00 1.467e+01 -3.333e-01 6.667e-01 3.333e+00 3.333e+00 6.667e-01 -3.333e-01 1.467e+01 2.667e+00 -2.667e+00 7.000e+00 6.667e-01 -4.000e+00 -3.333e-01 6.667e-01 3.333e+00 -2.667e+00 2.667e+00 1.467e+01 -4.000e+00 6.667e-01 7.000e+00 
+DEAL::72.3168
+Process #0
+Local range: [0, 360), global size: 360
+Vector data:
+4.941e+00 1.716e+01 3.811e+00 1.083e+00 3.815e+00 8.082e-01 -1.810e-01 -7.247e-01 -1.784e-01 3.781e+00 1.116e+01 1.837e+00 7.886e-01 2.139e+00 3.138e-01 -1.804e-01 -6.036e-01 -1.160e-01 5.725e-02 2.190e-01 4.532e-02 -5.636e-01 -2.347e+00 -5.921e-01 -2.040e+00 -9.318e+00 -2.273e+00 4.317e-02 1.340e-01 2.360e-02 -5.950e-01 -2.307e+00 -5.633e-01 -2.271e+00 -8.507e+00 -1.698e+00 1.413e+00 -7.166e+00 -4.468e+00 1.531e-01 -2.185e+00 -1.153e+00 -1.109e-01 8.951e-02 1.347e-01 -4.755e+00 -1.585e+01 -2.389e+00 -1.189e+00 -3.772e+00 -6.212e-01 1.747e-01 6.227e-01 7.765e-02 2.275e-02 -1.590e-02 -4.084e-02 -5.566e-01 -1.628e+00 -2.593e-01 -1.754e+00 -5.194e+00 -4.948e-01 -2.041e-02 -1.789e-01 -8.169e-02 -1.506e-01 2.347e-01 2.377e-01 -4.155e-01 1.993e+00 1.197e+00 -1.167e+00 9.612e+00 5.016e+00 -1.064e-01 2.498e+00 1.263e+00 9.472e-02 -2.707e-01 -1.725e-01 5.016e+00 9.612e+00 -1.167e+00 1.263e+00 2.498e+00 -1.064e-01 -1.725e-01 -2.707e-01 9.472e-02 -2.628e-02 -1.973e-01 -4.778e-02 4.129e-01 2.142e+00 6.641e-01 1.220e+00 8.134e+00 2.470e+00 -4.778e-02 -1.973e-01 -2.628e-02 6.641e-01 2.142e+00 4.129e-01 2.470e+00 8.134e+00 1.220e+00 -2.389e+00 -1.585e+01 -4.755e+00 -6.212e-01 -3.772e+00 -1.189e+00 7.765e-02 6.227e-01 1.747e-01 -4.468e+00 -7.166e+00 1.413e+00 -1.153e+00 -2.185e+00 1.531e-01 1.347e-01 8.951e-02 -1.109e-01 -8.169e-02 -1.789e-01 -2.041e-02 2.377e-01 2.347e-01 -1.506e-01 1.197e+00 1.993e+00 -4.155e-01 -4.084e-02 -1.590e-02 2.275e-02 -2.593e-01 -1.628e+00 -5.566e-01 -4.948e-01 -5.194e+00 -1.754e+00 1.837e+00 1.116e+01 3.781e+00 3.138e-01 2.139e+00 7.886e-01 -1.160e-01 -6.036e-01 -1.804e-01 3.811e+00 1.716e+01 4.941e+00 8.082e-01 3.815e+00 1.083e+00 -1.784e-01 -7.247e-01 -1.810e-01 2.360e-02 1.340e-01 4.317e-02 -5.633e-01 -2.307e+00 -5.950e-01 -1.698e+00 -8.507e+00 -2.271e+00 4.532e-02 2.190e-01 5.725e-02 -5.921e-01 -2.347e+00 -5.636e-01 -2.273e+00 -9.318e+00 -2.040e+00 4.941e+00 1.716e+01 3.811e+00 1.083e+00 3.815e+00 8.082e-01 -1.810e-01 -7.247e-01 -1.784e-01 3.781e+00 1.116e+01 1.837e+00 7.886e-01 2.139e+00 3.138e-01 -1.804e-01 -6.036e-01 -1.160e-01 5.725e-02 2.190e-01 4.532e-02 -5.636e-01 -2.347e+00 -5.921e-01 -2.040e+00 -9.318e+00 -2.273e+00 4.317e-02 1.340e-01 2.360e-02 -5.950e-01 -2.307e+00 -5.633e-01 -2.271e+00 -8.507e+00 -1.698e+00 1.413e+00 -7.166e+00 -4.468e+00 1.531e-01 -2.185e+00 -1.153e+00 -1.109e-01 8.951e-02 1.347e-01 -4.755e+00 -1.585e+01 -2.389e+00 -1.189e+00 -3.772e+00 -6.212e-01 1.747e-01 6.227e-01 7.765e-02 2.275e-02 -1.590e-02 -4.084e-02 -5.566e-01 -1.628e+00 -2.593e-01 -1.754e+00 -5.194e+00 -4.948e-01 -2.041e-02 -1.789e-01 -8.169e-02 -1.506e-01 2.347e-01 2.377e-01 -4.155e-01 1.993e+00 1.197e+00 -1.167e+00 9.612e+00 5.016e+00 -1.064e-01 2.498e+00 1.263e+00 9.472e-02 -2.707e-01 -1.725e-01 5.016e+00 9.612e+00 -1.167e+00 1.263e+00 2.498e+00 -1.064e-01 -1.725e-01 -2.707e-01 9.472e-02 -2.628e-02 -1.973e-01 -4.778e-02 4.129e-01 2.142e+00 6.641e-01 1.220e+00 8.134e+00 2.470e+00 -4.778e-02 -1.973e-01 -2.628e-02 6.641e-01 2.142e+00 4.129e-01 2.470e+00 8.134e+00 1.220e+00 -2.389e+00 -1.585e+01 -4.755e+00 -6.212e-01 -3.772e+00 -1.189e+00 7.765e-02 6.227e-01 1.747e-01 -4.468e+00 -7.166e+00 1.413e+00 -1.153e+00 -2.185e+00 1.531e-01 1.347e-01 8.951e-02 -1.109e-01 -8.169e-02 -1.789e-01 -2.041e-02 2.377e-01 2.347e-01 -1.506e-01 1.197e+00 1.993e+00 -4.155e-01 -4.084e-02 -1.590e-02 2.275e-02 -2.593e-01 -1.628e+00 -5.566e-01 -4.948e-01 -5.194e+00 -1.754e+00 1.837e+00 1.116e+01 3.781e+00 3.138e-01 2.139e+00 7.886e-01 -1.160e-01 -6.036e-01 -1.804e-01 3.811e+00 1.716e+01 4.941e+00 8.082e-01 3.815e+00 1.083e+00 -1.784e-01 -7.247e-01 -1.810e-01 2.360e-02 1.340e-01 4.317e-02 -5.633e-01 -2.307e+00 -5.950e-01 -1.698e+00 -8.507e+00 -2.271e+00 4.532e-02 2.190e-01 5.725e-02 -5.921e-01 -2.347e+00 -5.636e-01 -2.273e+00 -9.318e+00 -2.040e+00 
+DEAL::72.3168
+Process #0
+Local range: [0, 360), global size: 360
+Vector data:
+4.941e+00 1.716e+01 3.811e+00 1.083e+00 3.815e+00 8.082e-01 -1.810e-01 -7.247e-01 -1.784e-01 3.781e+00 1.116e+01 1.837e+00 7.886e-01 2.139e+00 3.138e-01 -1.804e-01 -6.036e-01 -1.160e-01 5.725e-02 2.190e-01 4.532e-02 -5.636e-01 -2.347e+00 -5.921e-01 -2.040e+00 -9.318e+00 -2.273e+00 4.317e-02 1.340e-01 2.360e-02 -5.950e-01 -2.307e+00 -5.633e-01 -2.271e+00 -8.507e+00 -1.698e+00 1.413e+00 -7.166e+00 -4.468e+00 1.531e-01 -2.185e+00 -1.153e+00 -1.109e-01 8.951e-02 1.347e-01 -4.755e+00 -1.585e+01 -2.389e+00 -1.189e+00 -3.772e+00 -6.212e-01 1.747e-01 6.227e-01 7.765e-02 2.275e-02 -1.590e-02 -4.084e-02 -5.566e-01 -1.628e+00 -2.593e-01 -1.754e+00 -5.194e+00 -4.948e-01 -2.041e-02 -1.789e-01 -8.169e-02 -1.506e-01 2.347e-01 2.377e-01 -4.155e-01 1.993e+00 1.197e+00 -1.167e+00 9.612e+00 5.016e+00 -1.064e-01 2.498e+00 1.263e+00 9.472e-02 -2.707e-01 -1.725e-01 5.016e+00 9.612e+00 -1.167e+00 1.263e+00 2.498e+00 -1.064e-01 -1.725e-01 -2.707e-01 9.472e-02 -2.628e-02 -1.973e-01 -4.778e-02 4.129e-01 2.142e+00 6.641e-01 1.220e+00 8.134e+00 2.470e+00 -4.778e-02 -1.973e-01 -2.628e-02 6.641e-01 2.142e+00 4.129e-01 2.470e+00 8.134e+00 1.220e+00 -2.389e+00 -1.585e+01 -4.755e+00 -6.212e-01 -3.772e+00 -1.189e+00 7.765e-02 6.227e-01 1.747e-01 -4.468e+00 -7.166e+00 1.413e+00 -1.153e+00 -2.185e+00 1.531e-01 1.347e-01 8.951e-02 -1.109e-01 -8.169e-02 -1.789e-01 -2.041e-02 2.377e-01 2.347e-01 -1.506e-01 1.197e+00 1.993e+00 -4.155e-01 -4.084e-02 -1.590e-02 2.275e-02 -2.593e-01 -1.628e+00 -5.566e-01 -4.948e-01 -5.194e+00 -1.754e+00 1.837e+00 1.116e+01 3.781e+00 3.138e-01 2.139e+00 7.886e-01 -1.160e-01 -6.036e-01 -1.804e-01 3.811e+00 1.716e+01 4.941e+00 8.082e-01 3.815e+00 1.083e+00 -1.784e-01 -7.247e-01 -1.810e-01 2.360e-02 1.340e-01 4.317e-02 -5.633e-01 -2.307e+00 -5.950e-01 -1.698e+00 -8.507e+00 -2.271e+00 4.532e-02 2.190e-01 5.725e-02 -5.921e-01 -2.347e+00 -5.636e-01 -2.273e+00 -9.318e+00 -2.040e+00 4.941e+00 1.716e+01 3.811e+00 1.083e+00 3.815e+00 8.082e-01 -1.810e-01 -7.247e-01 -1.784e-01 3.781e+00 1.116e+01 1.837e+00 7.886e-01 2.139e+00 3.138e-01 -1.804e-01 -6.036e-01 -1.160e-01 5.725e-02 2.190e-01 4.532e-02 -5.636e-01 -2.347e+00 -5.921e-01 -2.040e+00 -9.318e+00 -2.273e+00 4.317e-02 1.340e-01 2.360e-02 -5.950e-01 -2.307e+00 -5.633e-01 -2.271e+00 -8.507e+00 -1.698e+00 1.413e+00 -7.166e+00 -4.468e+00 1.531e-01 -2.185e+00 -1.153e+00 -1.109e-01 8.951e-02 1.347e-01 -4.755e+00 -1.585e+01 -2.389e+00 -1.189e+00 -3.772e+00 -6.212e-01 1.747e-01 6.227e-01 7.765e-02 2.275e-02 -1.590e-02 -4.084e-02 -5.566e-01 -1.628e+00 -2.593e-01 -1.754e+00 -5.194e+00 -4.948e-01 -2.041e-02 -1.789e-01 -8.169e-02 -1.506e-01 2.347e-01 2.377e-01 -4.155e-01 1.993e+00 1.197e+00 -1.167e+00 9.612e+00 5.016e+00 -1.064e-01 2.498e+00 1.263e+00 9.472e-02 -2.707e-01 -1.725e-01 5.016e+00 9.612e+00 -1.167e+00 1.263e+00 2.498e+00 -1.064e-01 -1.725e-01 -2.707e-01 9.472e-02 -2.628e-02 -1.973e-01 -4.778e-02 4.129e-01 2.142e+00 6.641e-01 1.220e+00 8.134e+00 2.470e+00 -4.778e-02 -1.973e-01 -2.628e-02 6.641e-01 2.142e+00 4.129e-01 2.470e+00 8.134e+00 1.220e+00 -2.389e+00 -1.585e+01 -4.755e+00 -6.212e-01 -3.772e+00 -1.189e+00 7.765e-02 6.227e-01 1.747e-01 -4.468e+00 -7.166e+00 1.413e+00 -1.153e+00 -2.185e+00 1.531e-01 1.347e-01 8.951e-02 -1.109e-01 -8.169e-02 -1.789e-01 -2.041e-02 2.377e-01 2.347e-01 -1.506e-01 1.197e+00 1.993e+00 -4.155e-01 -4.084e-02 -1.590e-02 2.275e-02 -2.593e-01 -1.628e+00 -5.566e-01 -4.948e-01 -5.194e+00 -1.754e+00 1.837e+00 1.116e+01 3.781e+00 3.138e-01 2.139e+00 7.886e-01 -1.160e-01 -6.036e-01 -1.804e-01 3.811e+00 1.716e+01 4.941e+00 8.082e-01 3.815e+00 1.083e+00 -1.784e-01 -7.247e-01 -1.810e-01 2.360e-02 1.340e-01 4.317e-02 -5.633e-01 -2.307e+00 -5.950e-01 -1.698e+00 -8.507e+00 -2.271e+00 4.532e-02 2.190e-01 5.725e-02 -5.921e-01 -2.347e+00 -5.636e-01 -2.273e+00 -9.318e+00 -2.040e+00 
+DEAL::54.4721
+Process #0
+Local range: [0, 180), global size: 180
+Vector data:
+-8.357e+00 -1.936e+01 7.502e+00 -1.714e+00 -3.825e+00 1.988e+00 2.975e-01 1.023e+00 -2.189e-01 7.502e+00 -1.936e+01 -8.357e+00 1.988e+00 -3.825e+00 -1.714e+00 -2.189e-01 1.023e+00 2.975e-01 -6.685e-02 1.331e-01 1.371e-01 -1.531e-01 6.871e-01 4.191e-01 -4.180e-02 1.877e-01 8.959e-02 1.371e-01 1.331e-01 -6.685e-02 4.191e-01 6.871e-01 -1.531e-01 8.959e-02 1.877e-01 -4.180e-02 -3.630e+00 -1.540e+00 -1.219e-01 -1.304e+01 -4.109e+00 1.228e-01 -3.261e+00 -1.027e+00 3.070e-02 -2.549e-01 -6.908e-01 -1.212e-01 -4.832e-01 -1.106e+00 -2.925e-01 -1.208e-01 -2.766e-01 -7.311e-02 -3.261e+00 -1.027e+00 3.070e-02 -1.304e+01 -4.109e+00 1.228e-01 -3.630e+00 -1.540e+00 -1.219e-01 -1.208e-01 -2.766e-01 -7.311e-02 -4.832e-01 -1.106e+00 -2.925e-01 -2.549e-01 -6.908e-01 -1.212e-01 1.772e-02 2.378e-01 1.136e-01 1.242e-01 8.492e-01 4.196e-01 3.105e-02 2.123e-01 1.049e-01 1.136e-01 2.378e-01 1.772e-02 4.196e-01 8.492e-01 1.242e-01 1.049e-01 2.123e-01 3.105e-02 3.105e-02 2.123e-01 1.049e-01 1.242e-01 8.492e-01 4.196e-01 1.772e-02 2.378e-01 1.136e-01 1.049e-01 2.123e-01 3.105e-02 4.196e-01 8.492e-01 1.242e-01 1.136e-01 2.378e-01 1.772e-02 -3.630e+00 -1.304e+01 -3.261e+00 -1.540e+00 -4.109e+00 -1.027e+00 -1.219e-01 1.228e-01 3.070e-02 -3.261e+00 -1.304e+01 -3.630e+00 -1.027e+00 -4.109e+00 -1.540e+00 3.070e-02 1.228e-01 -1.219e-01 -2.549e-01 -4.832e-01 -1.208e-01 -6.908e-01 -1.106e+00 -2.766e-01 -1.212e-01 -2.925e-01 -7.311e-02 -1.208e-01 -4.832e-01 -2.549e-01 -2.766e-01 -1.106e+00 -6.908e-01 -7.311e-02 -2.925e-01 -1.212e-01 -8.357e+00 -1.714e+00 2.975e-01 -1.936e+01 -3.825e+00 1.023e+00 7.502e+00 1.988e+00 -2.189e-01 -6.685e-02 -1.531e-01 -4.180e-02 1.331e-01 6.871e-01 1.877e-01 1.371e-01 4.191e-01 8.959e-02 7.502e+00 1.988e+00 -2.189e-01 -1.936e+01 -3.825e+00 1.023e+00 -8.357e+00 -1.714e+00 2.975e-01 1.371e-01 4.191e-01 8.959e-02 1.331e-01 6.871e-01 1.877e-01 -6.685e-02 -1.531e-01 -4.180e-02 
+DEAL::54.4721
+Process #0
+Local range: [0, 180), global size: 180
+Vector data:
+-8.357e+00 -1.936e+01 7.502e+00 -1.714e+00 -3.825e+00 1.988e+00 2.975e-01 1.023e+00 -2.189e-01 7.502e+00 -1.936e+01 -8.357e+00 1.988e+00 -3.825e+00 -1.714e+00 -2.189e-01 1.023e+00 2.975e-01 -6.685e-02 1.331e-01 1.371e-01 -1.531e-01 6.871e-01 4.191e-01 -4.180e-02 1.877e-01 8.959e-02 1.371e-01 1.331e-01 -6.685e-02 4.191e-01 6.871e-01 -1.531e-01 8.959e-02 1.877e-01 -4.180e-02 -3.630e+00 -1.540e+00 -1.219e-01 -1.304e+01 -4.109e+00 1.228e-01 -3.261e+00 -1.027e+00 3.070e-02 -2.549e-01 -6.908e-01 -1.212e-01 -4.832e-01 -1.106e+00 -2.925e-01 -1.208e-01 -2.766e-01 -7.311e-02 -3.261e+00 -1.027e+00 3.070e-02 -1.304e+01 -4.109e+00 1.228e-01 -3.630e+00 -1.540e+00 -1.219e-01 -1.208e-01 -2.766e-01 -7.311e-02 -4.832e-01 -1.106e+00 -2.925e-01 -2.549e-01 -6.908e-01 -1.212e-01 1.772e-02 2.378e-01 1.136e-01 1.242e-01 8.492e-01 4.196e-01 3.105e-02 2.123e-01 1.049e-01 1.136e-01 2.378e-01 1.772e-02 4.196e-01 8.492e-01 1.242e-01 1.049e-01 2.123e-01 3.105e-02 3.105e-02 2.123e-01 1.049e-01 1.242e-01 8.492e-01 4.196e-01 1.772e-02 2.378e-01 1.136e-01 1.049e-01 2.123e-01 3.105e-02 4.196e-01 8.492e-01 1.242e-01 1.136e-01 2.378e-01 1.772e-02 -3.630e+00 -1.304e+01 -3.261e+00 -1.540e+00 -4.109e+00 -1.027e+00 -1.219e-01 1.228e-01 3.070e-02 -3.261e+00 -1.304e+01 -3.630e+00 -1.027e+00 -4.109e+00 -1.540e+00 3.070e-02 1.228e-01 -1.219e-01 -2.549e-01 -4.832e-01 -1.208e-01 -6.908e-01 -1.106e+00 -2.766e-01 -1.212e-01 -2.925e-01 -7.311e-02 -1.208e-01 -4.832e-01 -2.549e-01 -2.766e-01 -1.106e+00 -6.908e-01 -7.311e-02 -2.925e-01 -1.212e-01 -8.357e+00 -1.714e+00 2.975e-01 -1.936e+01 -3.825e+00 1.023e+00 7.502e+00 1.988e+00 -2.189e-01 -6.685e-02 -1.531e-01 -4.180e-02 1.331e-01 6.871e-01 1.877e-01 1.371e-01 4.191e-01 8.959e-02 7.502e+00 1.988e+00 -2.189e-01 -1.936e+01 -3.825e+00 1.023e+00 -8.357e+00 -1.714e+00 2.975e-01 1.371e-01 4.191e-01 8.959e-02 1.331e-01 6.871e-01 1.877e-01 -6.685e-02 -1.531e-01 -4.180e-02 
index 62ce5fbb8db059ffe08d3333868ba17dc2cada23..c55064bb871d47b6a662fc57b19be0dee34dd344 100644 (file)
@@ -44,7 +44,7 @@ main(int argc, char **argv)
 
   // test without shared memory
   {
-    test<2, 2, 3, double, VectorizedArray<double>>(6, false, MPI_COMM_SELF);
+    test<2, 2, 3, double, VectorizedArray<double>>(0, 6, false, MPI_COMM_SELF);
   }
 
   // test with shared memory
@@ -58,7 +58,10 @@ main(int argc, char **argv)
                         MPI_INFO_NULL,
                         &subcommunicator);
 
-    test<2, 2, 3, double, VectorizedArray<double>>(6, false, subcommunicator);
+    test<2, 2, 3, double, VectorizedArray<double>>(0,
+                                                   6,
+                                                   false,
+                                                   subcommunicator);
 
     MPI_Comm_free(&subcommunicator);
   }
index e746eca8349063ea9a1696e9c5df6a6dc7864e39..b53a850454831db8bd978e8e9430216358323394 100644 (file)
@@ -108,12 +108,17 @@ test(const unsigned int geometry, const MPI_Comm comm = MPI_COMM_SELF)
   using MF = MatrixFree<dim, Number, VectorizedArrayType>;
 
   typename MF::AdditionalData additional_data;
-  additional_data.mapping_update_flags                = update_values;
-  additional_data.mapping_update_flags_inner_faces    = update_values;
-  additional_data.mapping_update_flags_boundary_faces = update_values;
-  additional_data.mapping_update_flags_faces_by_cells = update_values;
-  additional_data.hold_all_faces_to_owned_cells       = true;
-  additional_data.communicator_sm                     = comm;
+  additional_data.mapping_update_flags = update_values | update_gradients;
+  additional_data.mapping_update_flags_inner_faces =
+    update_values | update_gradients;
+  additional_data.mapping_update_flags_boundary_faces =
+    update_values | update_gradients;
+  additional_data.mapping_update_flags_faces_by_cells =
+    update_values | update_gradients;
+  additional_data.hold_all_faces_to_owned_cells = true;
+  additional_data.communicator_sm               = comm;
+  additional_data.tasks_parallel_scheme =
+    MF::AdditionalData::TasksParallelScheme::none;
 
   MatrixFreeTools::categorize_by_boundary_ids(tria, additional_data);
 
@@ -168,6 +173,35 @@ test(const unsigned int geometry, const MPI_Comm comm = MPI_COMM_SELF)
                              ExcMessage("Entries do not match!"));
                     }
                 }
+
+              phi_m.gather_evaluate(src,
+                                    EvaluationFlags::values |
+                                      EvaluationFlags::gradients);
+              phi_p.gather_evaluate(src,
+                                    EvaluationFlags::values |
+                                      EvaluationFlags::gradients);
+
+              for (unsigned int q = 0; q < phi_m.n_q_points; ++q)
+                {
+                  const auto u_minus = phi_m.get_value(q);
+                  const auto u_plus  = phi_p.get_value(q);
+
+                  const auto grad_u_minus = phi_m.get_gradient(q);
+                  const auto grad_u_plus  = phi_p.get_gradient(q);
+
+                  for (unsigned int v = 0; v < VectorizedArray<double>::size();
+                       ++v)
+                    {
+                      Assert(std::abs(u_minus[v] - u_plus[v]) < 1e-10,
+                             ExcMessage("Entries do not match!"));
+
+                      if (false)
+                        for (int d = 0; d < dim; ++d)
+                          Assert(std::abs(grad_u_minus[d][v] -
+                                          grad_u_plus[d][v]) < 1e-6,
+                                 ExcMessage("Entries do not match!"));
+                    }
+                }
             }
         }
     },

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.