]> https://gitweb.dealii.org/ - dealii.git/commitdiff
RemotePointEvaluation: documentation and template arg rename 18561/head
authorTimo Heister <timo.heister@gmail.com>
Fri, 13 Jun 2025 17:22:51 +0000 (13:22 -0400)
committerTimo Heister <timo.heister@gmail.com>
Sat, 14 Jun 2025 04:01:24 +0000 (00:01 -0400)
- add to the class documentation of RemotePointEvaluation
- change the template argument T to DataType

include/deal.II/base/mpi_remote_point_evaluation.h

index c23f1bc3ada94a203d345b578b851ec7d53cd2f6..c934e29d4d6069db76bc846bacb8d96c5d0d1383 100644 (file)
@@ -52,9 +52,22 @@ namespace Utilities
      * finite element mesh can be distributed among processes, which makes
      * the operations more involved due to communication.
      *
-     * The following helper functions for this class exist:
-     * VectorTools::point_values, VectorTools::point_gradients (the
-     * overloads taking an instance of RemotePointEvaluation).     *
+     * Once initialized, the data exchange is handled by the functions
+     * process_and_evaluate(), that transfers data from the Triangulation
+     * to the set of points, and evaluate_and_process(), that performs
+     * the inverse operation (taking data defined at each point and making
+     * it available in the corresponding cell of the Triangulation).
+     *
+     * The helper functions VectorTools::point_values() and
+     * VectorTools::point_gradients (the overloads taking an instance
+     * of RemotePointEvaluation) can be used in place of
+     * process_and_evaluate() and perform the evaluation of a finite
+     * element function.
+     *
+     * The template argument `DataType` of the process_and_evaluate() and
+     * evaluate_and_process() functions can be a scalar (`float`, `double`),
+     * `Tensor<rank,dim>`, or `Point<dim>`. Additionally, one can specify
+     * the number of components (in case `DataType` is a scalar).
      *
      * @note Usage and implementation details of this class and the
      *   helper functions above are explained in detail in step-87.
@@ -255,10 +268,10 @@ namespace Utilities
         /**
          * Return local view of the processed cell @p cell for the vector @p values.
          */
-        template <typename T>
-        ArrayView<T>
-        get_data_view(const unsigned int  cell,
-                      const ArrayView<T> &values) const;
+        template <typename DataType>
+        ArrayView<DataType>
+        get_data_view(const unsigned int         cell,
+                      const ArrayView<DataType> &values) const;
 
         /**
          * Level and index of processed cells.
@@ -292,9 +305,12 @@ namespace Utilities
       get_cell_data() const;
 
       /**
-       * Evaluate function @p evaluation_function in the given  points and
+       * Evaluate function @p evaluation_function in the given points and
        * triangulation. The result is stored in @p output.
        *
+       * @p buffer is a temporary buffer than can be used to avoid extra
+       * allocations.
+       *
        * @note If the map of points to cells is not a
        *   one-to-one relation (is_map_unique()==false), the result needs to be
        *   processed with the help of get_point_ptrs(). This
@@ -311,12 +327,12 @@ namespace Utilities
        *   provide unrolled tensors, which is useful, e.g., if its dimension
        *   and its rank is not known at compile time.
        */
-      template <typename T, unsigned int n_components = 1>
+      template <typename DataType, unsigned int n_components = 1>
       void
       evaluate_and_process(
-        std::vector<T> &output,
-        std::vector<T> &buffer,
-        const std::function<void(const ArrayView<T> &, const CellData &)>
+        std::vector<DataType> &output,
+        std::vector<DataType> &buffer,
+        const std::function<void(const ArrayView<DataType> &, const CellData &)>
                   &evaluation_function,
         const bool sort_data = true) const;
 
@@ -324,10 +340,10 @@ namespace Utilities
        * Same as above but with the result provided as return value and
        * without external allocation of a user-provided buffer.
        */
-      template <typename T, unsigned int n_components = 1>
-      std::vector<T>
+      template <typename DataType, unsigned int n_components = 1>
+      std::vector<DataType>
       evaluate_and_process(
-        const std::function<void(const ArrayView<T> &, const CellData &)>
+        const std::function<void(const ArrayView<DataType> &, const CellData &)>
                   &evaluation_function,
         const bool sort_data = true) const;
 
@@ -345,26 +361,26 @@ namespace Utilities
        *   provide unrolled tensors, which is useful, e.g., if its dimension
        *   and its rank is not known at compile time.
        */
-      template <typename T, unsigned int n_components = 1>
+      template <typename DataType, unsigned int n_components = 1>
       void
       process_and_evaluate(
-        const std::vector<T> &input,
-        std::vector<T>       &buffer,
-        const std::function<void(const ArrayView<const T> &, const CellData &)>
-                  &evaluation_function,
-        const bool sort_data = true) const;
+        const std::vector<DataType>                 &input,
+        std::vector<DataType>                       &buffer,
+        const std::function<void(const ArrayView<const DataType> &,
+                                 const CellData &)> &evaluation_function,
+        const bool                                   sort_data = true) const;
 
       /**
        * Same as above but without external allocation of a user-provided
        * buffer.
        */
-      template <typename T, unsigned int n_components = 1>
+      template <typename DataType, unsigned int n_components = 1>
       void
       process_and_evaluate(
-        const std::vector<T> &input,
-        const std::function<void(const ArrayView<const T> &, const CellData &)>
-                  &evaluation_function,
-        const bool sort_data = true) const;
+        const std::vector<DataType>                 &input,
+        const std::function<void(const ArrayView<const DataType> &,
+                                 const CellData &)> &evaluation_function,
+        const bool                                   sort_data = true) const;
 
       /**
        * Return a CRS-like data structure to determine the position of the
@@ -705,11 +721,11 @@ namespace Utilities
 
 
     template <int dim, int spacedim>
-    template <typename T>
-    ArrayView<T>
+    template <typename DataType>
+    ArrayView<DataType>
     RemotePointEvaluation<dim, spacedim>::CellData::get_data_view(
-      const unsigned int  cell,
-      const ArrayView<T> &values) const
+      const unsigned int         cell,
+      const ArrayView<DataType> &values) const
     {
       AssertIndexRange(cell, cells.size());
       return {values.data() + reference_point_ptrs[cell],
@@ -719,12 +735,12 @@ namespace Utilities
 
 
     template <int dim, int spacedim>
-    template <typename T, unsigned int n_components>
+    template <typename DataType, unsigned int n_components>
     void
     RemotePointEvaluation<dim, spacedim>::evaluate_and_process(
-      std::vector<T> &output,
-      std::vector<T> &buffer,
-      const std::function<void(const ArrayView<T> &, const CellData &)>
+      std::vector<DataType> &output,
+      std::vector<DataType> &buffer,
+      const std::function<void(const ArrayView<DataType> &, const CellData &)>
                 &evaluation_function,
       const bool sort_data) const
     {
@@ -749,15 +765,14 @@ namespace Utilities
         n_components);
 
       // ... for evaluation
-      ArrayView<T> buffer_eval(buffer.data(),
-                               send_permutation.size() * n_components);
+      ArrayView<DataType> buffer_eval(buffer.data(),
+                                      send_permutation.size() * n_components);
 
       // ... for communication (send)
-      ArrayView<T> buffer_send(sort_data ?
-                                 (buffer.data() +
-                                  send_permutation.size() * n_components) :
-                                 buffer.data(),
-                               send_permutation.size() * n_components);
+      ArrayView<DataType> buffer_send(
+        sort_data ? (buffer.data() + send_permutation.size() * n_components) :
+                    buffer.data(),
+        send_permutation.size() * n_components);
 
       // more arrays
       std::vector<MPI_Request>       send_requests;
@@ -819,7 +834,7 @@ namespace Utilities
             continue;
 
           internal::pack_and_isend(
-            ArrayView<const T>(
+            ArrayView<const DataType>(
               buffer_send.begin() + send_ptrs[i] * n_components,
               (send_ptrs[i + 1] - send_ptrs[i]) * n_components),
             send_ranks[i],
@@ -875,7 +890,7 @@ namespace Utilities
           const unsigned int j = std::distance(recv_ranks.begin(), ptr);
 
           // ... for communication (recv)
-          ArrayView<T> recv_buffer(
+          ArrayView<DataType> recv_buffer(
             sort_data ?
               (buffer.data() + send_permutation.size() * 2 * n_components) :
               (output.data() + recv_ptrs[j] * n_components),
@@ -910,20 +925,20 @@ namespace Utilities
 
 
     template <int dim, int spacedim>
-    template <typename T, unsigned int n_components>
-    std::vector<T>
+    template <typename DataType, unsigned int n_components>
+    std::vector<DataType>
     RemotePointEvaluation<dim, spacedim>::evaluate_and_process(
-      const std::function<void(const ArrayView<T> &, const CellData &)>
+      const std::function<void(const ArrayView<DataType> &, const CellData &)>
                 &evaluation_function,
       const bool sort_data) const
     {
-      std::vector<T> output;
-      std::vector<T> buffer;
+      std::vector<DataType> output;
+      std::vector<DataType> buffer;
 
-      this->evaluate_and_process<T, n_components>(output,
-                                                  buffer,
-                                                  evaluation_function,
-                                                  sort_data);
+      this->evaluate_and_process<DataType, n_components>(output,
+                                                         buffer,
+                                                         evaluation_function,
+                                                         sort_data);
 
       return output;
     }
@@ -931,14 +946,14 @@ namespace Utilities
 
 
     template <int dim, int spacedim>
-    template <typename T, unsigned int n_components>
+    template <typename DataType, unsigned int n_components>
     void
     RemotePointEvaluation<dim, spacedim>::process_and_evaluate(
-      const std::vector<T> &input,
-      std::vector<T>       &buffer,
-      const std::function<void(const ArrayView<const T> &, const CellData &)>
-                &evaluation_function,
-      const bool sort_data) const
+      const std::vector<DataType>                 &input,
+      std::vector<DataType>                       &buffer,
+      const std::function<void(const ArrayView<const DataType> &,
+                               const CellData &)> &evaluation_function,
+      const bool                                   sort_data) const
     {
 #ifndef DEAL_II_WITH_MPI
       Assert(false, ExcNeedsMPI());
@@ -964,15 +979,14 @@ namespace Utilities
         n_components);
 
       // ... for evaluation
-      ArrayView<T> buffer_eval(buffer.data(),
-                               send_permutation.size() * n_components);
+      ArrayView<DataType> buffer_eval(buffer.data(),
+                                      send_permutation.size() * n_components);
 
       // ... for communication (send)
-      ArrayView<T> buffer_send(sort_data ?
-                                 (buffer.data() +
-                                  send_permutation.size() * n_components) :
-                                 const_cast<T *>(input.data()),
-                               point_ptrs.back() * n_components);
+      ArrayView<DataType> buffer_send(
+        sort_data ? (buffer.data() + send_permutation.size() * n_components) :
+                    const_cast<DataType *>(input.data()),
+        point_ptrs.back() * n_components);
 
       // more arrays
       std::vector<MPI_Request>       send_requests;
@@ -1041,7 +1055,7 @@ namespace Utilities
             continue;
 
           internal::pack_and_isend(
-            ArrayView<const T>(
+            ArrayView<const DataType>(
               buffer_send.begin() + recv_ptrs[i] * n_components,
               (recv_ptrs[i + 1] - recv_ptrs[i]) * n_components),
             recv_ranks[i],
@@ -1096,7 +1110,7 @@ namespace Utilities
 
           const unsigned int j = std::distance(send_ranks.begin(), ptr);
 
-          ArrayView<T> recv_buffer(
+          ArrayView<DataType> recv_buffer(
             sort_data ?
               (buffer.data() +
                (point_ptrs.back() + send_permutation.size()) * n_components) :
@@ -1134,19 +1148,19 @@ namespace Utilities
 
 
     template <int dim, int spacedim>
-    template <typename T, unsigned int n_components>
+    template <typename DataType, unsigned int n_components>
     void
     RemotePointEvaluation<dim, spacedim>::process_and_evaluate(
-      const std::vector<T> &input,
-      const std::function<void(const ArrayView<const T> &, const CellData &)>
-                &evaluation_function,
-      const bool sort_data) const
+      const std::vector<DataType>                 &input,
+      const std::function<void(const ArrayView<const DataType> &,
+                               const CellData &)> &evaluation_function,
+      const bool                                   sort_data) const
     {
-      std::vector<T> buffer;
-      this->process_and_evaluate<T, n_components>(input,
-                                                  buffer,
-                                                  evaluation_function,
-                                                  sort_data);
+      std::vector<DataType> buffer;
+      this->process_and_evaluate<DataType, n_components>(input,
+                                                         buffer,
+                                                         evaluation_function,
+                                                         sort_data);
     }
 
   } // end of namespace MPI

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.