]> https://gitweb.dealii.org/ - dealii.git/commitdiff
RPE: allow to evaluate cell-data vector 14099/head
authorPeter Munch <peterrmuench@gmail.com>
Wed, 6 Jul 2022 09:48:18 +0000 (11:48 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Fri, 8 Jul 2022 07:48:26 +0000 (09:48 +0200)
doc/news/changes/minor/20220706Munch [new file with mode: 0644]
include/deal.II/numerics/vector_tools_evaluate.h
tests/remote_point_evaluation/vector_tools_evaluate_at_points_06.cc [new file with mode: 0644]
tests/remote_point_evaluation/vector_tools_evaluate_at_points_06.with_p4est=true.mpirun=4.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20220706Munch b/doc/news/changes/minor/20220706Munch
new file mode 100644 (file)
index 0000000..76eb69b
--- /dev/null
@@ -0,0 +1,4 @@
+Improved: The function VectorTools::point_values() can now handle
+cell-data vectors.
+<br>
+(Peter Munch, 2022/07/06)
index 15ee805097a52d281ce9678b72014c553cd639b1..8660498930691bc5f7cf8acbe48543a92b0c7613 100644 (file)
@@ -21,6 +21,8 @@
 
 #include <deal.II/base/mpi_remote_point_evaluation.h>
 
+#include <deal.II/distributed/tria_base.h>
+
 #include <deal.II/dofs/dof_handler.h>
 
 #include <deal.II/matrix_free/fe_point_evaluation.h>
@@ -108,6 +110,12 @@ namespace VectorTools
    *   cache, dof_handler_2, vector_2);
    * @endcode
    *
+   * The function can also be used to evaluate cell-data vectors. For this
+   * purpose, one passes in a Triangulation instead of a DoFHandler and a
+   * vector of size Trinagulation::n_active_cells() or a vector, which
+   * has been initialized with the partitioner returned by
+   * parallel::TriangulationBase::global_active_cell_index_partitioner().
+   *
    * @note If a point cannot be found, the result for these points will
    *   be undefined (most probably 0). If you want to be sure that all
    *   points received a valid result, call `cache.all_points_found()`
@@ -116,7 +124,12 @@ namespace VectorTools
    * @warning This is a collective call that needs to be executed by all
    *   processors in the communicator.
    */
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
@@ -124,7 +137,7 @@ namespace VectorTools
                                typename VectorType::value_type>::value_type>
   point_values(
     const Mapping<dim> &                                  mapping,
-    const DoFHandler<dim, spacedim> &                     dof_handler,
+    const MeshType<dim, spacedim> &                       mesh,
     const VectorType &                                    vector,
     const std::vector<Point<spacedim>> &                  evaluation_points,
     Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
@@ -142,7 +155,12 @@ namespace VectorTools
    * @warning This is a collective call that needs to be executed by all
    *   processors in the communicator.
    */
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
@@ -150,7 +168,7 @@ namespace VectorTools
                                typename VectorType::value_type>::value_type>
   point_values(
     const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
-    const DoFHandler<dim, spacedim> &                           dof_handler,
+    const MeshType<dim, spacedim> &                             mesh,
     const VectorType &                                          vector,
     const EvaluationFlags::EvaluationFlags flags = EvaluationFlags::avg);
 
@@ -164,7 +182,12 @@ namespace VectorTools
    * @warning This is a collective call that needs to be executed by all
    *   processors in the communicator.
    */
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
@@ -172,7 +195,7 @@ namespace VectorTools
                                typename VectorType::value_type>::gradient_type>
   point_gradients(
     const Mapping<dim> &                                  mapping,
-    const DoFHandler<dim, spacedim> &                     dof_handler,
+    const MeshType<dim, spacedim> &                       mesh,
     const VectorType &                                    vector,
     const std::vector<Point<spacedim>> &                  evaluation_points,
     Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
@@ -189,7 +212,12 @@ namespace VectorTools
    * @warning This is a collective call that needs to be executed by all
    *   processors in the communicator.
    */
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
@@ -197,7 +225,7 @@ namespace VectorTools
                                typename VectorType::value_type>::gradient_type>
   point_gradients(
     const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
-    const DoFHandler<dim, spacedim> &                           dof_handler,
+    const MeshType<dim, spacedim> &                             mesh,
     const VectorType &                                          vector,
     const EvaluationFlags::EvaluationFlags flags = EvaluationFlags::avg);
 
@@ -207,42 +235,52 @@ namespace VectorTools
 
 
 #ifndef DOXYGEN
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   inline std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
                                spacedim,
                                typename VectorType::value_type>::value_type>
   point_values(const Mapping<dim> &                mapping,
-               const DoFHandler<dim, spacedim> &   dof_handler,
+               const MeshType<dim, spacedim> &     mesh,
                const VectorType &                  vector,
                const std::vector<Point<spacedim>> &evaluation_points,
                Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
                const EvaluationFlags::EvaluationFlags                flags)
   {
-    cache.reinit(evaluation_points, dof_handler.get_triangulation(), mapping);
+    cache.reinit(evaluation_points, mesh.get_triangulation(), mapping);
 
-    return point_values<n_components>(cache, dof_handler, vector, flags);
+    return point_values<n_components>(cache, mesh, vector, flags);
   }
 
 
 
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   inline std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
                                spacedim,
                                typename VectorType::value_type>::gradient_type>
   point_gradients(const Mapping<dim> &                mapping,
-                  const DoFHandler<dim, spacedim> &   dof_handler,
+                  const MeshType<dim, spacedim> &     mesh,
                   const VectorType &                  vector,
                   const std::vector<Point<spacedim>> &evaluation_points,
                   Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
                   const EvaluationFlags::EvaluationFlags                flags)
   {
-    cache.reinit(evaluation_points, dof_handler.get_triangulation(), mapping);
+    cache.reinit(evaluation_points, mesh.get_triangulation(), mapping);
 
-    return point_gradients<n_components>(cache, dof_handler, vector, flags);
+    return point_gradients<n_components>(cache, mesh, vector, flags);
   }
 
 
@@ -310,10 +348,179 @@ namespace VectorTools
               int spacedim,
               typename VectorType,
               typename value_type>
+    void
+    process_cell(
+      const unsigned int i,
+      const typename Utilities::MPI::RemotePointEvaluation<dim,
+                                                           spacedim>::CellData
+        &                                                         cell_data,
+      const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
+      const DoFHandler<dim, spacedim> &                           dof_handler,
+      const VectorType &                                          vector,
+      const UpdateFlags                                           update_flags,
+      const dealii::EvaluationFlags::EvaluationFlags evaluation_flags,
+      const std::function<
+        value_type(const FEPointEvaluation<n_components,
+                                           dim,
+                                           spacedim,
+                                           typename VectorType::value_type> &,
+                   const unsigned int &)>           process_quadrature_point,
+      const ArrayView<value_type> &                 values,
+      std::vector<typename VectorType::value_type> &solution_values,
+      std::vector<
+        std::unique_ptr<FEPointEvaluation<n_components,
+                                          dim,
+                                          spacedim,
+                                          typename VectorType::value_type>>>
+        &evaluators)
+    {
+      if (evaluators.size() == 0)
+        evaluators.resize(dof_handler.get_fe_collection().size());
+
+      typename DoFHandler<dim>::active_cell_iterator cell = {
+        &cache.get_triangulation(),
+        cell_data.cells[i].first,
+        cell_data.cells[i].second,
+        &dof_handler};
+
+      const ArrayView<const Point<dim>> unit_points(
+        cell_data.reference_point_values.data() +
+          cell_data.reference_point_ptrs[i],
+        cell_data.reference_point_ptrs[i + 1] -
+          cell_data.reference_point_ptrs[i]);
+
+      solution_values.resize(
+        dof_handler.get_fe(cell->active_fe_index()).n_dofs_per_cell());
+      cell->get_dof_values(vector,
+                           solution_values.begin(),
+                           solution_values.end());
+
+      if (evaluators[cell->active_fe_index()] == nullptr)
+        evaluators[cell->active_fe_index()] =
+          std::make_unique<FEPointEvaluation<n_components,
+                                             dim,
+                                             spacedim,
+                                             typename VectorType::value_type>>(
+            cache.get_mapping(), cell->get_fe(), update_flags);
+      auto &evaluator = *evaluators[cell->active_fe_index()];
+
+      evaluator.reinit(cell, unit_points);
+      evaluator.evaluate(solution_values, evaluation_flags);
+
+      for (unsigned int q = 0; q < unit_points.size(); ++q)
+        values[q + cell_data.reference_point_ptrs[i]] =
+          process_quadrature_point(evaluator, q);
+    }
+
+
+
+    template <int dim, int spacedim, typename Number>
+    Number
+    get_value(
+      const Triangulation<dim, spacedim> &                               tria,
+      const Vector<Number> &                                             vector,
+      const typename Triangulation<dim, spacedim>::active_cell_iterator &cell)
+    {
+      (void)tria;
+      AssertDimension(tria.n_active_cells(), vector.size());
+      return vector[cell->active_cell_index()];
+    }
+
+
+
+    template <int dim, int spacedim, typename Number>
+    Number
+    get_value(
+      const Triangulation<dim, spacedim> &                               tria,
+      const LinearAlgebra::distributed::Vector<Number> &                 vector,
+      const typename Triangulation<dim, spacedim>::active_cell_iterator &cell)
+    {
+      const auto distributed_tria =
+        dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(&tria);
+
+      const bool use_distributed_path =
+        (distributed_tria == nullptr) ?
+          false :
+          (vector.get_partitioner().get() ==
+           distributed_tria->global_active_cell_index_partitioner()
+             .lock()
+             .get());
+
+      if (use_distributed_path)
+        {
+          return vector[cell->global_active_cell_index()];
+        }
+      else
+        {
+          AssertDimension(tria.n_active_cells(), vector.locally_owned_size());
+          return vector[cell->active_cell_index()];
+        }
+    }
+
+
+
+    template <int n_components,
+              int dim,
+              int spacedim,
+              typename VectorType,
+              typename value_type>
+    void
+    process_cell(
+      const unsigned int i,
+      const typename Utilities::MPI::RemotePointEvaluation<dim,
+                                                           spacedim>::CellData
+        &cell_data,
+      const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &,
+      const Triangulation<dim, spacedim> &triangulation,
+      const VectorType &                  vector,
+      const UpdateFlags,
+      const dealii::EvaluationFlags::EvaluationFlags evaluation_flags,
+      const std::function<
+        value_type(const FEPointEvaluation<n_components,
+                                           dim,
+                                           spacedim,
+                                           typename VectorType::value_type> &,
+                   const unsigned int &)>,
+      const ArrayView<value_type> &values,
+      std::vector<typename VectorType::value_type> &,
+      std::vector<
+        std::unique_ptr<FEPointEvaluation<n_components,
+                                          dim,
+                                          spacedim,
+                                          typename VectorType::value_type>>> &)
+    {
+      (void)evaluation_flags;
+
+      static_assert(n_components == 1,
+                    "A cell-data vector can only have a single component.");
+
+      Assert(evaluation_flags ==
+               dealii::EvaluationFlags::EvaluationFlags::values,
+             ExcMessage("For cell-data vectors, only values can be queried."));
+
+      typename Triangulation<dim>::active_cell_iterator cell = {
+        &triangulation, cell_data.cells[i].first, cell_data.cells[i].second};
+
+      const auto value = get_value(triangulation, vector, cell);
+
+      for (unsigned int q = cell_data.reference_point_ptrs[i];
+           q < cell_data.reference_point_ptrs[i + 1];
+           ++q)
+        values[q] = value;
+    }
+
+
+
+    template <int n_components,
+              int dim,
+              int spacedim,
+              typename MeshType,
+              typename VectorType,
+              typename value_type>
     inline std::vector<value_type>
     evaluate_at_points(
       const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
-      const DoFHandler<dim, spacedim> &                           dof_handler,
+      const MeshType &                                            mesh,
       const VectorType &                                          vector,
       const EvaluationFlags::EvaluationFlags                      flags,
       const UpdateFlags                                           update_flags,
@@ -332,7 +539,7 @@ namespace VectorTools
                "yourself or another function that does this for you."));
 
       Assert(
-        &dof_handler.get_triangulation() == &cache.get_triangulation(),
+        &mesh.get_triangulation() == &cache.get_triangulation(),
         ExcMessage(
           "The provided Utilities::MPI::RemotePointEvaluation and DoFHandler "
           "object have been set up with different Triangulation objects, "
@@ -351,44 +558,21 @@ namespace VectorTools
                                               dim,
                                               spacedim,
                                               typename VectorType::value_type>>>
-            evaluators(dof_handler.get_fe_collection().size());
+            evaluators;
 
           for (unsigned int i = 0; i < cell_data.cells.size(); ++i)
-            {
-              typename DoFHandler<dim>::active_cell_iterator cell = {
-                &cache.get_triangulation(),
-                cell_data.cells[i].first,
-                cell_data.cells[i].second,
-                &dof_handler};
-
-              const ArrayView<const Point<dim>> unit_points(
-                cell_data.reference_point_values.data() +
-                  cell_data.reference_point_ptrs[i],
-                cell_data.reference_point_ptrs[i + 1] -
-                  cell_data.reference_point_ptrs[i]);
-
-              solution_values.resize(
-                dof_handler.get_fe(cell->active_fe_index()).n_dofs_per_cell());
-              cell->get_dof_values(vector,
-                                   solution_values.begin(),
-                                   solution_values.end());
-
-              if (evaluators[cell->active_fe_index()] == nullptr)
-                evaluators[cell->active_fe_index()] = std::make_unique<
-                  FEPointEvaluation<n_components,
-                                    dim,
-                                    spacedim,
-                                    typename VectorType::value_type>>(
-                  cache.get_mapping(), cell->get_fe(), update_flags);
-              auto &evaluator = *evaluators[cell->active_fe_index()];
-
-              evaluator.reinit(cell, unit_points);
-              evaluator.evaluate(solution_values, evaluation_flags);
-
-              for (unsigned int q = 0; q < unit_points.size(); ++q)
-                values[q + cell_data.reference_point_ptrs[i]] =
-                  process_quadrature_point(evaluator, q);
-            }
+            process_cell<n_components, dim, spacedim, VectorType, value_type>(
+              i,
+              cell_data,
+              cache,
+              mesh,
+              vector,
+              update_flags,
+              evaluation_flags,
+              process_quadrature_point,
+              values,
+              solution_values,
+              evaluators);
         };
 
         std::vector<value_type> evaluation_point_results;
@@ -431,7 +615,12 @@ namespace VectorTools
     }
   } // namespace internal
 
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   inline std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
@@ -439,7 +628,7 @@ namespace VectorTools
                                typename VectorType::value_type>::value_type>
   point_values(
     const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
-    const DoFHandler<dim, spacedim> &                           dof_handler,
+    const MeshType<dim, spacedim> &                             mesh,
     const VectorType &                                          vector,
     const EvaluationFlags::EvaluationFlags                      flags)
   {
@@ -447,13 +636,14 @@ namespace VectorTools
       n_components,
       dim,
       spacedim,
+      MeshType<dim, spacedim>,
       VectorType,
       typename FEPointEvaluation<n_components,
                                  dim,
                                  spacedim,
                                  typename VectorType::value_type>::value_type>(
       cache,
-      dof_handler,
+      mesh,
       vector,
       flags,
       update_values,
@@ -463,7 +653,12 @@ namespace VectorTools
       });
   }
 
-  template <int n_components, int dim, int spacedim, typename VectorType>
+  template <int n_components,
+            template <int, int>
+            class MeshType,
+            int dim,
+            int spacedim,
+            typename VectorType>
   inline std::vector<
     typename FEPointEvaluation<n_components,
                                dim,
@@ -471,7 +666,7 @@ namespace VectorTools
                                typename VectorType::value_type>::gradient_type>
   point_gradients(
     const Utilities::MPI::RemotePointEvaluation<dim, spacedim> &cache,
-    const DoFHandler<dim, spacedim> &                           dof_handler,
+    const MeshType<dim, spacedim> &                             mesh,
     const VectorType &                                          vector,
     const EvaluationFlags::EvaluationFlags                      flags)
   {
@@ -479,6 +674,7 @@ namespace VectorTools
       n_components,
       dim,
       spacedim,
+      MeshType<dim, spacedim>,
       VectorType,
       typename FEPointEvaluation<
         n_components,
@@ -486,12 +682,12 @@ namespace VectorTools
         spacedim,
         typename VectorType::value_type>::gradient_type>(
       cache,
-      dof_handler,
+      mesh,
       vector,
       flags,
       update_gradients,
       dealii::EvaluationFlags::gradients,
-      [](const auto &evaluator, const auto &q) {
+      [](const auto &evaluator, const unsigned &q) {
         return evaluator.get_gradient(q);
       });
   }
diff --git a/tests/remote_point_evaluation/vector_tools_evaluate_at_points_06.cc b/tests/remote_point_evaluation/vector_tools_evaluate_at_points_06.cc
new file mode 100644 (file)
index 0000000..81563f4
--- /dev/null
@@ -0,0 +1,145 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 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.
+//
+// ---------------------------------------------------------------------
+
+// Test VectorTools::point_values() for cell-data vectors.
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/mapping_fe.h>
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/lac/la_parallel_vector.h>
+#include <deal.II/lac/precondition.h>
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/lac/solver_control.h>
+
+#include <deal.II/matrix_free/fe_point_evaluation.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/numerics/vector_tools_evaluate.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim, int spacedim>
+std::shared_ptr<const Utilities::MPI::Partitioner>
+create_partitioner(const DoFHandler<dim, spacedim> &dof_handler)
+{
+  IndexSet locally_relevant_dofs;
+
+  DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs);
+
+  return std::make_shared<const Utilities::MPI::Partitioner>(
+    dof_handler.locally_owned_dofs(),
+    locally_relevant_dofs,
+    dof_handler.get_communicator());
+}
+
+void
+test()
+{
+  const unsigned int dim             = 2;
+  const unsigned int n_refinements_1 = 3;
+
+  parallel::distributed::Triangulation<dim> tria_1(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria_1, -1.0, +1.0);
+  tria_1.refine_global(n_refinements_1);
+  DoFHandler<dim> dof_handler_1(tria_1);
+
+  FE_DGQ<dim> fe_1(0);
+  dof_handler_1.distribute_dofs(fe_1);
+
+  // approach 1:
+  LinearAlgebra::distributed::Vector<double> vector_1(
+    create_partitioner(dof_handler_1));
+
+  // approach 2:
+  LinearAlgebra::distributed::Vector<double> vector_2(
+    tria_1.global_active_cell_index_partitioner().lock());
+
+  // approach 3:
+  Vector<double> vector_3(tria_1.n_active_cells());
+
+  for (const auto &cell : dof_handler_1.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        const auto value = cell->center()[0];
+
+        // approach 1:
+        Vector<double> temp(cell->get_fe().n_dofs_per_cell());
+        temp = value;
+        cell->set_dof_values(temp, vector_1);
+
+        // approach 2:
+        vector_2[cell->global_active_cell_index()] = value;
+
+        // approach 3:
+        vector_3[cell->active_cell_index()] = value;
+      }
+
+  const MappingQ1<dim> mapping_1;
+
+  std::vector<Point<dim>> evaluation_points;
+
+  const unsigned int n_intervals = 100;
+
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+    for (unsigned int i = 0; i < n_intervals; ++i)
+      {
+        const double rad = 2.0 * numbers::PI / n_intervals * i;
+        evaluation_points.emplace_back(0.75 * std::cos(rad),
+                                       0.75 * std::sin(rad));
+      }
+
+  Utilities::MPI::RemotePointEvaluation<dim> evaluation_cache;
+  const auto evaluation_point_results_1 = VectorTools::point_values<1>(
+    mapping_1, dof_handler_1, vector_1, evaluation_points, evaluation_cache);
+  const auto evaluation_point_results_2 =
+    VectorTools::point_values<1>(evaluation_cache, tria_1, vector_2);
+  const auto evaluation_point_results_3 =
+    VectorTools::point_values<1>(evaluation_cache, tria_1, vector_3);
+
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+    for (unsigned int i = 0; i < n_intervals; ++i)
+      {
+        if (std::abs(evaluation_point_results_1[i] -
+                     evaluation_point_results_2[i]) > 1e-10 ||
+            std::abs(evaluation_point_results_1[i] -
+                     evaluation_point_results_3[i]) > 1e-10)
+          Assert(false, ExcNotImplemented());
+      }
+
+  deallog << "OK!" << std::endl;
+}
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+  initlog();
+
+  test();
+}
diff --git a/tests/remote_point_evaluation/vector_tools_evaluate_at_points_06.with_p4est=true.mpirun=4.output b/tests/remote_point_evaluation/vector_tools_evaluate_at_points_06.with_p4est=true.mpirun=4.output
new file mode 100644 (file)
index 0000000..5cfb783
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK!

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.