]> https://gitweb.dealii.org/ - dealii.git/commitdiff
ECL: enable gradients for structured meshes and FE_DGQ 11747/head
authorPeter Munch <peterrmuench@gmail.com>
Sun, 14 Feb 2021 12:59:22 +0000 (13:59 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Tue, 2 Mar 2021 15:35:20 +0000 (16:35 +0100)
include/deal.II/matrix_free/evaluation_kernels.h
include/deal.II/matrix_free/fe_evaluation.h
tests/matrix_free/ecl_01.cc [new file with mode: 0644]
tests/matrix_free/ecl_01.output [new file with mode: 0644]
tests/matrix_free/ecl_02.cc
tests/matrix_free/ecl_02.output

index 2fe6ce3dee6c7fb34350f91dbcf3c5ab557c4226..12abbb4ccd0e13d0875b2109554741ad8cfaab15 100644 (file)
@@ -3523,8 +3523,6 @@ namespace internal
             else
               {
                 // case 5: default vector access
-                AssertDimension(n_face_orientations, 1);
-
                 // for the integrate_scatter path (integrate == true), we
                 // need to only prepare the data in this function for all
                 // components to later call distribute_local_to_global();
@@ -3541,8 +3539,6 @@ namespace internal
         else
           {
             // case 5: default vector access
-            AssertDimension(n_face_orientations, 1);
-
             proc.default_operation(temp1, comp);
             if (integrate)
               accesses_global_vector = false;
index 113096fedbf5d6945783aa745d80c2ce7050cb5c..dd4acbb03448c38140fdd3e898f7a7f938f6e3ff 100644 (file)
@@ -4696,10 +4696,14 @@ FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
 
   // Simple case: We have contiguous storage, so we can simply copy out the
   // data
-  if (this->dof_info->index_storage_variants[ind][this->cell] ==
-        internal::MatrixFreeFunctions::DoFInfo::IndexStorageVariants::
-          interleaved_contiguous &&
-      n_lanes == VectorizedArrayType::size())
+  if ((this->dof_info->index_storage_variants[ind][this->cell] ==
+         internal::MatrixFreeFunctions::DoFInfo::IndexStorageVariants::
+           interleaved_contiguous &&
+       n_lanes == VectorizedArrayType::size()) &&
+      !(is_face &&
+        this->dof_access_index ==
+          internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+        this->is_interior_face == false))
     {
       const unsigned int dof_index =
         dof_indices_cont[this->cell * VectorizedArrayType::size()] +
@@ -4724,23 +4728,18 @@ FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
       return;
     }
 
+  std::array<unsigned int, VectorizedArrayType::size()> cells =
+    this->get_cell_or_face_ids();
+
   // More general case: Must go through the components one by one and apply
   // some transformations
   const unsigned int n_filled_lanes =
     this->dof_info->n_vectorization_lanes_filled[ind][this->cell];
 
-  unsigned int dof_indices[VectorizedArrayType::size()];
-  for (unsigned int v = 0; v < n_filled_lanes; ++v)
-    dof_indices[v] =
-      dof_indices_cont[this->cell * VectorizedArrayType::size() + v] +
-      this->dof_info
-          ->component_dof_indices_offset[this->active_fe_index]
-                                        [this->first_selected_component] *
-        this->dof_info->dof_indices_interleave_strides
-          [ind][this->cell * VectorizedArrayType::size() + v];
-
-  for (unsigned int v = n_filled_lanes; v < VectorizedArrayType::size(); ++v)
-    dof_indices[v] = numbers::invalid_unsigned_int;
+  const bool is_ecl =
+    this->dof_access_index ==
+      internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+    this->is_interior_face == false;
 
   if (vectors_sm[0] != nullptr)
     {
@@ -4751,19 +4750,20 @@ FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
 
         for (unsigned int v = 0; v < n_filled_lanes; ++v)
           {
+            Assert(cells[v] != numbers::invalid_unsigned_int,
+                   ExcNotImplemented());
             Assert(ind < this->dof_info->dof_indices_contiguous_sm.size(),
                    ExcIndexRange(
                      ind, 0, this->dof_info->dof_indices_contiguous_sm.size()));
-            Assert(this->cell * VectorizedArrayType::size() + v <
+            Assert(cells[v] <
                      this->dof_info->dof_indices_contiguous_sm[ind].size(),
                    ExcIndexRange(
-                     this->cell * VectorizedArrayType::size() + v,
+                     cells[v],
                      0,
                      this->dof_info->dof_indices_contiguous_sm[ind].size()));
 
             const auto &temp =
-              this->dof_info->dof_indices_contiguous_sm
-                [ind][this->cell * VectorizedArrayType::size() + v];
+              this->dof_info->dof_indices_contiguous_sm[ind][cells[v]];
 
             if (temp.first != numbers::invalid_unsigned_int)
               vector_ptrs[v] = const_cast<typename VectorType::value_type *>(
@@ -4781,7 +4781,7 @@ FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
       };
 
       if (n_filled_lanes == VectorizedArrayType::size() &&
-          n_lanes == VectorizedArrayType::size())
+          n_lanes == VectorizedArrayType::size() && !is_ecl)
         {
           if (n_components == 1 || n_fe_components == 1)
             {
@@ -4842,10 +4842,26 @@ FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
       return;
     }
 
+  unsigned int dof_indices[VectorizedArrayType::size()];
+
+  for (unsigned int v = 0; v < n_filled_lanes; ++v)
+    {
+      Assert(cells[v] != numbers::invalid_unsigned_int, ExcNotImplemented());
+      dof_indices[v] =
+        dof_indices_cont[cells[v]] +
+        this->dof_info
+            ->component_dof_indices_offset[this->active_fe_index]
+                                          [this->first_selected_component] *
+          this->dof_info->dof_indices_interleave_strides[ind][cells[v]];
+    }
+
+  for (unsigned int v = n_filled_lanes; v < VectorizedArrayType::size(); ++v)
+    dof_indices[v] = numbers::invalid_unsigned_int;
+
   // In the case with contiguous cell indices, we know that there are no
   // constraints and that the indices within each element are contiguous
   if (n_filled_lanes == VectorizedArrayType::size() &&
-      n_lanes == VectorizedArrayType::size())
+      n_lanes == VectorizedArrayType::size() && !is_ecl)
     {
       if (this->dof_info->index_storage_variants[ind][this->cell] ==
           internal::MatrixFreeFunctions::DoFInfo::IndexStorageVariants::
@@ -8630,35 +8646,73 @@ FEFaceEvaluation<dim,
       !(evaluation_flag & EvaluationFlags::gradients))
     return;
 
-  if (fe_degree > -1)
-    internal::FEFaceEvaluationImplEvaluateSelector<dim, VectorizedArrayType>::
-      template run<fe_degree, n_q_points_1d>(
-        n_components,
-        *this->data,
-        values_array,
-        this->begin_values(),
-        this->begin_gradients(),
-        this->scratch_data,
-        evaluation_flag & EvaluationFlags::values,
-        evaluation_flag & EvaluationFlags::gradients,
-        this->face_no,
-        this->subface_index,
-        this->face_orientation,
-        this->descriptor->face_orientations);
+  if (this->dof_access_index ==
+        internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+      this->is_interior_face == false)
+    {
+      const auto face_nos          = this->compute_face_no_data();
+      const auto face_orientations = this->compute_face_orientations();
+
+#  ifdef DEBUG
+      // currently on structured meshes are supported -> face numers and
+      // orientations have to be the same for all filled lanes
+      for (unsigned int v = 1; v < VectorizedArrayType::size(); ++v)
+        {
+          if (face_nos[v] != numbers::invalid_unsigned_int)
+            AssertDimension(face_nos[0], face_nos[v]);
+          if (face_orientations[v] != numbers::invalid_unsigned_int)
+            AssertDimension(face_orientations[0], face_orientations[v]);
+        }
+#  endif
+
+      internal::FEFaceEvaluationImplEvaluateSelector<dim, VectorizedArrayType>::
+        template run<fe_degree, n_q_points_1d>(
+          n_components,
+          *this->data,
+          values_array,
+          this->begin_values(),
+          this->begin_gradients(),
+          this->scratch_data,
+          evaluation_flag & EvaluationFlags::values,
+          evaluation_flag & EvaluationFlags::gradients,
+          face_nos[0],
+          this->subface_index,
+          face_orientations[0],
+          this->descriptor->face_orientations);
+    }
   else
-    internal::FEFaceEvaluationFactory<dim, Number, VectorizedArrayType>::
-      evaluate(n_components,
-               *this->data,
-               values_array,
-               this->begin_values(),
-               this->begin_gradients(),
-               this->scratch_data,
-               evaluation_flag & EvaluationFlags::values,
-               evaluation_flag & EvaluationFlags::gradients,
-               this->face_no,
-               this->subface_index,
-               this->face_orientation,
-               this->descriptor->face_orientations);
+    {
+      if (fe_degree > -1)
+        internal::FEFaceEvaluationImplEvaluateSelector<dim,
+                                                       VectorizedArrayType>::
+          template run<fe_degree, n_q_points_1d>(
+            n_components,
+            *this->data,
+            values_array,
+            this->begin_values(),
+            this->begin_gradients(),
+            this->scratch_data,
+            evaluation_flag & EvaluationFlags::values,
+            evaluation_flag & EvaluationFlags::gradients,
+            this->face_no,
+            this->subface_index,
+            this->face_orientation,
+            this->descriptor->face_orientations);
+      else
+        internal::FEFaceEvaluationFactory<dim, Number, VectorizedArrayType>::
+          evaluate(n_components,
+                   *this->data,
+                   values_array,
+                   this->begin_values(),
+                   this->begin_gradients(),
+                   this->scratch_data,
+                   evaluation_flag & EvaluationFlags::values,
+                   evaluation_flag & EvaluationFlags::gradients,
+                   this->face_no,
+                   this->subface_index,
+                   this->face_orientation,
+                   this->descriptor->face_orientations);
+    }
 
 #  ifdef DEBUG
   if (evaluation_flag & EvaluationFlags::values)
@@ -9258,7 +9312,9 @@ FEFaceEvaluation<dim,
     }
   else
     {
-      std::fill(face_no_data.begin(), face_no_data.end(), 0);
+      std::fill(face_no_data.begin(),
+                face_no_data.end(),
+                numbers::invalid_unsigned_int);
 
       if (dim == 3)
         {
@@ -9308,6 +9364,16 @@ FEFaceEvaluation<dim,
               face_no_data[i] = face_orientation;
             }
         }
+      else
+        {
+          std::fill(
+            face_no_data.begin(),
+            face_no_data.begin() +
+              this->dof_info->n_vectorization_lanes_filled
+                [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell]
+                [this->cell],
+            0);
+        }
     }
 
   return face_no_data;
diff --git a/tests/matrix_free/ecl_01.cc b/tests/matrix_free/ecl_01.cc
new file mode 100644 (file)
index 0000000..ccbdc0c
--- /dev/null
@@ -0,0 +1,324 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/matrix_free/fe_evaluation.h>
+#include <deal.II/matrix_free/matrix_free.h>
+#include <deal.II/matrix_free/tools.h>
+
+#include <deal.II/multigrid/mg_constrained_dofs.h>
+
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+
+// Compare the evaluation of a SIP Laplace operator with face- and element-
+// centric loops on a structured grid.
+
+template <int dim, typename Number>
+class CosineFunction : public Function<dim, Number>
+{
+public:
+  Number
+  value(const Point<dim> &p, const unsigned int component = 0) const
+  {
+    (void)component;
+    return std::cos(2 * numbers::PI * p[0]);
+  }
+};
+
+template <int dim,
+          int fe_degree,
+          int n_points                 = fe_degree + 1,
+          typename Number              = double,
+          typename VectorizedArrayType = VectorizedArray<Number>>
+void
+test(const unsigned int n_refinements = 1)
+{
+  using VectorType = LinearAlgebra::distributed::Vector<Number>;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, 0.0, 1.0, true);
+
+  if (false)
+    {
+      std::vector<dealii::GridTools::PeriodicFacePair<
+        typename dealii::Triangulation<dim>::cell_iterator>>
+        periodic_faces;
+
+      if (dim >= 1)
+        dealii::GridTools::collect_periodic_faces(
+          tria, 0, 1, 0, periodic_faces);
+
+      if (dim >= 2)
+        dealii::GridTools::collect_periodic_faces(
+          tria, 2, 3, 1, periodic_faces);
+
+      if (dim >= 3)
+        dealii::GridTools::collect_periodic_faces(
+          tria, 4, 5, 2, periodic_faces);
+
+      tria.add_periodicity(periodic_faces);
+    }
+
+  tria.refine_global(n_refinements);
+
+  FE_DGQ<dim>     fe(fe_degree);
+  DoFHandler<dim> dof_handler(tria);
+  dof_handler.distribute_dofs(fe);
+
+  MappingQ<dim> mapping(1);
+  QGauss<1>     quad(n_points);
+
+  AffineConstraints<Number> constraint;
+
+  using MF = MatrixFree<dim, Number, VectorizedArrayType>;
+
+  typename MF::AdditionalData additional_data;
+  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.tasks_parallel_scheme =
+    MF::AdditionalData::TasksParallelScheme::none;
+  // additional_data.mapping_update_flags_faces_by_cells = update_values;
+  // additional_data.hold_all_faces_to_owned_cells       = true;
+
+  MatrixFreeTools::categorize_by_boundary_ids(tria, additional_data);
+
+  MF matrix_free;
+  matrix_free.reinit(mapping, dof_handler, constraint, quad, additional_data);
+
+  VectorType src, dst;
+
+  matrix_free.initialize_dof_vector(src);
+  matrix_free.initialize_dof_vector(dst);
+
+  FEEvaluation<dim, fe_degree, n_points, 1, Number, VectorizedArrayType> phi(
+    matrix_free);
+  FEFaceEvaluation<dim, fe_degree, n_points, 1, Number, VectorizedArrayType>
+    phi_m(matrix_free, true);
+  FEFaceEvaluation<dim, fe_degree, n_points, 1, Number, VectorizedArrayType>
+    phi_p(matrix_free, false);
+
+
+  CosineFunction<dim, Number> function;
+  VectorTools::interpolate(dof_handler, function, src);
+
+  dst = 0.0;
+
+  /**
+   * Face-centric loop
+   */
+  matrix_free.template loop<VectorType, VectorType>(
+    [&](const auto &, auto &dst, const auto &src, const auto range) {
+      for (unsigned int cell = range.first; cell < range.second; ++cell)
+        {
+          phi.reinit(cell);
+          phi.read_dof_values(src);
+          phi.evaluate(false, true, false);
+          for (unsigned int q = 0; q < phi.n_q_points; ++q)
+            phi.submit_gradient(phi.get_gradient(q), q);
+          phi.integrate(false, true);
+          phi.set_dof_values(dst);
+        }
+    },
+    [&](const auto &, auto &dst, const auto &src, const auto range) {
+      for (unsigned int face = range.first; face < range.second; face++)
+        {
+          phi_m.reinit(face);
+          phi_p.reinit(face);
+
+          phi_m.read_dof_values(src);
+          phi_m.evaluate(true, true);
+          phi_p.read_dof_values(src);
+          phi_p.evaluate(true, true);
+          VectorizedArrayType sigmaF =
+            (std::abs((phi_m.get_normal_vector(0) *
+                       phi_m.inverse_jacobian(0))[dim - 1]) +
+             std::abs((phi_m.get_normal_vector(0) *
+                       phi_p.inverse_jacobian(0))[dim - 1])) *
+            (Number)(std::max(fe_degree, 1) * (fe_degree + 1.0));
+
+          for (unsigned int q = 0; q < phi_m.n_q_points; ++q)
+            {
+              VectorizedArrayType average_value =
+                (phi_m.get_value(q) - phi_p.get_value(q)) * 0.5;
+              VectorizedArrayType average_valgrad =
+                phi_m.get_normal_derivative(q) + phi_p.get_normal_derivative(q);
+              average_valgrad =
+                average_value * 2. * sigmaF - average_valgrad * 0.5;
+              phi_m.submit_normal_derivative(-average_value, q);
+              phi_p.submit_normal_derivative(-average_value, q);
+              phi_m.submit_value(average_valgrad, q);
+              phi_p.submit_value(-average_valgrad, q);
+            }
+          phi_m.integrate(true, true);
+          phi_m.distribute_local_to_global(dst);
+          phi_p.integrate(true, true);
+          phi_p.distribute_local_to_global(dst);
+        }
+    },
+    [&](const auto &, auto &dst, const auto &src, const auto face_range) {
+      for (unsigned int face = face_range.first; face < face_range.second;
+           face++)
+        {
+          phi_m.reinit(face);
+          phi_m.read_dof_values(src);
+          phi_m.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
+          VectorizedArrayType sigmaF =
+            std::abs((phi_m.get_normal_vector(0) *
+                      phi_m.inverse_jacobian(0))[dim - 1]) *
+            Number(std::max(fe_degree, 1) * (fe_degree + 1.0)) * 2.0;
+
+          for (unsigned int q = 0; q < phi_m.n_q_points; ++q)
+            {
+              VectorizedArrayType average_value = phi_m.get_value(q);
+              VectorizedArrayType average_valgrad =
+                -phi_m.get_normal_derivative(q);
+              average_valgrad += average_value * sigmaF * 2.0;
+              phi_m.submit_normal_derivative(-average_value, q);
+              phi_m.submit_value(average_valgrad, q);
+            }
+
+          phi_m.integrate(EvaluationFlags::values | EvaluationFlags::gradients);
+          phi_m.distribute_local_to_global(dst);
+        }
+    },
+    dst,
+    src,
+    false,
+    MF::DataAccessOnFaces::gradients,
+    MF::DataAccessOnFaces::gradients);
+
+  deallog << dst.l2_norm() << std::endl;
+  dst.print(deallog.get_file_stream());
+
+  dst = 0.0;
+
+  /**
+   * Element-centric loop
+   */
+  matrix_free.template loop_cell_centric<VectorType, VectorType>(
+    [&](const auto &, auto &dst, const auto &src, const auto range) {
+      for (unsigned int cell = range.first; cell < range.second; ++cell)
+        {
+          phi.reinit(cell);
+          phi.read_dof_values(src);
+          phi.evaluate(false, true, false);
+          for (unsigned int q = 0; q < phi.n_q_points; ++q)
+            phi.submit_gradient(phi.get_gradient(q), q);
+          phi.integrate(false, true);
+
+          for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
+               ++face)
+            {
+              auto bids =
+                matrix_free.get_faces_by_cells_boundary_id(cell, face);
+
+              if (bids[0] != numbers::internal_face_boundary_id)
+                {
+                  phi_m.reinit(cell, face);
+                  phi_m.read_dof_values(src);
+                  phi_m.evaluate(EvaluationFlags::values |
+                                 EvaluationFlags::gradients);
+                  VectorizedArrayType sigmaF =
+                    std::abs((phi_m.get_normal_vector(0) *
+                              phi_m.inverse_jacobian(0))[dim - 1]) *
+                    Number(std::max(fe_degree, 1) * (fe_degree + 1.0)) * 2.0;
+
+                  for (unsigned int q = 0; q < phi_m.n_q_points; ++q)
+                    {
+                      VectorizedArrayType average_value = phi_m.get_value(q);
+                      VectorizedArrayType average_valgrad =
+                        -phi_m.get_normal_derivative(q);
+                      average_valgrad += average_value * sigmaF * 2.0;
+                      phi_m.submit_normal_derivative(-average_value, q);
+                      phi_m.submit_value(average_valgrad, q);
+                    }
+
+                  phi_m.integrate(EvaluationFlags::values |
+                                  EvaluationFlags::gradients);
+
+                  for (unsigned int q = 0; q < phi.dofs_per_cell; ++q)
+                    phi.begin_dof_values()[q] += phi_m.begin_dof_values()[q];
+                }
+              else
+                {
+                  phi_m.reinit(cell, face);
+                  phi_p.reinit(cell, face);
+
+                  phi_m.read_dof_values(src);
+                  phi_m.evaluate(EvaluationFlags::values |
+                                 EvaluationFlags::gradients);
+                  phi_p.read_dof_values(src);
+                  phi_p.evaluate(EvaluationFlags::values |
+                                 EvaluationFlags::gradients);
+
+                  VectorizedArrayType sigmaF =
+                    (std::abs((phi_m.get_normal_vector(0) *
+                               phi_m.inverse_jacobian(0))[dim - 1]) +
+                     std::abs((phi_m.get_normal_vector(0) *
+                               phi_p.inverse_jacobian(0))[dim - 1])) *
+                    (Number)(std::max(fe_degree, 1) * (fe_degree + 1.0));
+
+                  for (unsigned int q = 0; q < phi_m.n_q_points; ++q)
+                    {
+                      VectorizedArrayType average_value =
+                        (phi_m.get_value(q) - phi_p.get_value(q)) * 0.5;
+                      VectorizedArrayType average_valgrad =
+                        phi_m.get_normal_derivative(q) +
+                        phi_p.get_normal_derivative(q);
+                      average_valgrad =
+                        average_value * 2. * sigmaF - average_valgrad * 0.5;
+                      phi_m.submit_normal_derivative(-average_value, q);
+                      phi_m.submit_value(average_valgrad, q);
+                    }
+                  phi_m.integrate(EvaluationFlags::values |
+                                  EvaluationFlags::gradients);
+
+                  for (unsigned int q = 0; q < phi.dofs_per_cell; ++q)
+                    phi.begin_dof_values()[q] += phi_m.begin_dof_values()[q];
+                }
+            }
+
+          phi.set_dof_values(dst);
+        }
+    },
+    dst,
+    src,
+    false,
+    MF::DataAccessOnFaces::gradients);
+
+  deallog << dst.l2_norm() << std::endl;
+  dst.print(deallog.get_file_stream());
+}
+
+int
+main()
+{
+  initlog();
+  test<2, 2, 3, double, VectorizedArray<double>>();
+}
diff --git a/tests/matrix_free/ecl_01.output b/tests/matrix_free/ecl_01.output
new file mode 100644 (file)
index 0000000..8cd7af3
--- /dev/null
@@ -0,0 +1,11 @@
+
+DEAL::35.0111
+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::35.0111
+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 
index 8a0df9a52b31169bbec9d2d5be5e2a5eb1bdb8ec..a46b38af4aef7a81551290b0aa0fcc6a983de91f 100644 (file)
@@ -121,7 +121,7 @@ test(const unsigned int n_refinements = 1)
                 deallog << static_cast<int>(phi_m.begin_dof_values()[i][0])
                         << " ";
               deallog << std::endl;
-              phi_m.gather_evaluate(src, true, false);
+              phi_m.gather_evaluate(src, EvaluationFlags::values);
               for (unsigned int i = 0; i < phi_m.static_n_q_points; i++)
                 deallog << static_cast<int>(phi_m.begin_values()[i][0]) << " ";
               deallog << std::endl;
@@ -132,7 +132,7 @@ test(const unsigned int n_refinements = 1)
                         << " ";
               deallog << std::endl;
 
-              phi_p.gather_evaluate(src, true, false);
+              phi_p.gather_evaluate(src, EvaluationFlags::values);
               for (unsigned int i = 0; i < phi_p.static_n_q_points; i++)
                 deallog << static_cast<int>(phi_p.begin_values()[i][0]) << " ";
               deallog << std::endl << std::endl;
index d584d53cc83a6979ac6534b3d7a8562a3e5ea014..833521ad1899579bc0b3d8753b4338cb61ac8c76 100644 (file)
@@ -1,12 +1,12 @@
 
 DEAL::1 1 1 1 
 DEAL::1 1 
-DEAL::1 1 1 1 
+DEAL::2 2 2 2 
 DEAL::2 2 
 DEAL::
 DEAL::1 1 1 1 
 DEAL::1 1 
-DEAL::1 1 1 1 
+DEAL::2 2 2 2 
 DEAL::2 2 
 DEAL::
 DEAL::1 1 1 1 

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.