]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Portable::MatrixFree: introduce DeviceVector
authorTimo Heister <timo.heister@gmail.com>
Fri, 21 Mar 2025 15:16:41 +0000 (11:16 -0400)
committerTimo Heister <timo.heister@gmail.com>
Fri, 21 Mar 2025 15:53:48 +0000 (11:53 -0400)
Instead of passing double* to the user code, introduce a new type that is a type alias to a Kokkos::View.

examples/step-64/step-64.cc
include/deal.II/matrix_free/portable_fe_evaluation.h
include/deal.II/matrix_free/portable_matrix_free.h
include/deal.II/matrix_free/portable_matrix_free.templates.h
include/deal.II/matrix_free/tools.h
tests/matrix_free_kokkos/coefficient_eval_device.cc
tests/matrix_free_kokkos/matrix_free_device_no_index_initialize.cc
tests/matrix_free_kokkos/matrix_vector_device_mf.h
tests/matrix_free_kokkos/matrix_vector_host_device_multi_component.cc

index a594805b9a1ff72583f6d97503f8e041c791a071..f494ab9a9d5c856a32f69ba4a64ad355c0321cfc 100644 (file)
@@ -192,8 +192,8 @@ namespace Step64
 
     DEAL_II_HOST_DEVICE void
     operator()(const typename Portable::MatrixFree<dim, double>::Data *data,
-               const double                                           *src,
-               double *dst) const;
+               const Portable::DeviceVector<double>                   &src,
+               Portable::DeviceVector<double> &dst) const;
 
   private:
     double *coef;
@@ -208,8 +208,8 @@ namespace Step64
   template <int dim, int fe_degree>
   DEAL_II_HOST_DEVICE void LocalHelmholtzOperator<dim, fe_degree>::operator()(
     const typename Portable::MatrixFree<dim, double>::Data *data,
-    const double                                           *src,
-    double                                                 *dst) const
+    const Portable::DeviceVector<double>                   &src,
+    Portable::DeviceVector<double>                         &dst) const
   {
     Portable::FEEvaluation<dim, fe_degree, fe_degree + 1, 1, double> fe_eval(
       data);
index ce75f15eb50455076b344aa2828621a83e3aa098..8c47adec9e74ce503d07cd0aeffffd3c2d8f0115 100644 (file)
@@ -151,7 +151,7 @@ namespace Portable
      * AffineConstraints::read_dof_values() as well.
      */
     DEAL_II_HOST_DEVICE void
-    read_dof_values(const Number *src);
+    read_dof_values(const DeviceVector<Number> &src);
 
     /**
      * Take the value stored internally on dof values of the current cell and
@@ -160,7 +160,7 @@ namespace Portable
      * function AffineConstraints::distribute_local_to_global.
      */
     DEAL_II_HOST_DEVICE void
-    distribute_local_to_global(Number *dst) const;
+    distribute_local_to_global(DeviceVector<Number> &dst) const;
 
     /**
      * Evaluate the function values and the gradients of the FE function given
@@ -327,7 +327,7 @@ namespace Portable
             typename Number>
   DEAL_II_HOST_DEVICE void
   FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
-    read_dof_values(const Number *src)
+    read_dof_values(const DeviceVector<Number> &src)
   {
     // Populate the scratch memory
     Kokkos::parallel_for(Kokkos::TeamThreadRange(data->team_member,
@@ -359,7 +359,7 @@ namespace Portable
             typename Number>
   DEAL_II_HOST_DEVICE void
   FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
-    distribute_local_to_global(Number *dst) const
+    distribute_local_to_global(DeviceVector<Number> &dst) const
   {
     for (unsigned int c = 0; c < n_components_; ++c)
       {
index 47c1b8a9e82ceba39b55d5991887f71eed68d167..d508552604a597b43d4e2ec1af12a277b0dc8347 100644 (file)
@@ -65,6 +65,19 @@ namespace Portable
   struct SharedData;
 #endif
 
+  /**
+   * Type for source and destination vectors in device functions like
+   * MatrixFree::cell_loop().
+   *
+   * This is a type alias to a Kokkos::View to a chunk of memory, typically
+   * pointing to the local elements in the LinearAlgebra::distributed::Vector.
+   */
+  template <typename Number>
+  using DeviceVector =
+    Kokkos::View<Number *, MemorySpace::Default::kokkos_space>;
+
+
+
   /**
    * This class collects all the data that is stored for the matrix free
    * implementation. The storage scheme is tailored towards several loops
@@ -366,8 +379,8 @@ namespace Portable
      * \code
      * DEAL_II_HOST_DEVICE void operator()(
      *   const typename Portable::MatrixFree<dim, Number>::Data *data,
-     *   const Number *                                          src,
-     *   Number *                                                dst) const;
+     *   const DeviceVector<Number> &src,
+     *   DeviceVector<Number> & dst) const;
      *   static const unsigned int n_local_dofs;
      *   static const unsigned int n_q_points;
      * \endcode
index 829bf9d6745eab3f7671a3b30f0834d1cb3a2d77..1d4135d748ed44a34f10dc6ce4ecc2ca1f090468 100644 (file)
@@ -354,18 +354,19 @@ namespace Portable
       ApplyKernel(
         Functor                                                 func,
         const typename MatrixFree<dim, Number>::PrecomputedData gpu_data,
-        Number *const                                           src,
-        Number                                                 *dst)
+        const LinearAlgebra::distributed::Vector<Number, MemorySpace::Default>
+                                                                         &src,
+        LinearAlgebra::distributed::Vector<Number, MemorySpace::Default> &dst)
         : func(func)
         , gpu_data(gpu_data)
-        , src(src)
-        , dst(dst)
+        , src(src.get_values(), src.locally_owned_size())
+        , dst(dst.get_values(), dst.locally_owned_size())
       {}
 
       Functor                                                 func;
       const typename MatrixFree<dim, Number>::PrecomputedData gpu_data;
-      Number *const                                           src;
-      Number                                                 *dst;
+      const DeviceVector<Number>                              src;
+      DeviceVector<Number>                                    dst;
 
 
       // Provide the shared memory capacity. This function takes the team_size
@@ -404,7 +405,8 @@ namespace Portable
                                                     &gpu_data,
                                                     &shared_data};
 
-        func(&data, src, dst);
+        DeviceVector<Number> nonconstdst = dst;
+        func(&data, src, nonconstdst);
       }
     };
   } // namespace internal
@@ -1017,7 +1019,7 @@ namespace Portable
               Kokkos::AUTO);
 
           internal::ApplyKernel<dim, Number, Functor> apply_kernel(
-            func, get_data(color), src.get_values(), dst.get_values());
+            func, get_data(color), src, dst);
 
           Kokkos::parallel_for("dealii::MatrixFree::serial_cell_loop",
                                team_policy,
@@ -1062,7 +1064,7 @@ namespace Portable
                     Kokkos::AUTO);
 
                 internal::ApplyKernel<dim, Number, Functor> apply_kernel(
-                  func, get_data(0), src.get_values(), dst.get_values());
+                  func, get_data(0), src, dst);
 
                 Kokkos::parallel_for(
                   "dealii::MatrixFree::distributed_cell_loop_0",
@@ -1085,7 +1087,7 @@ namespace Portable
                     Kokkos::AUTO);
 
                 internal::ApplyKernel<dim, Number, Functor> apply_kernel(
-                  func, get_data(1), src.get_values(), dst.get_values());
+                  func, get_data(1), src, dst);
 
                 Kokkos::parallel_for(
                   "dealii::MatrixFree::distributed_cell_loop_1",
@@ -1113,7 +1115,7 @@ namespace Portable
                     Kokkos::AUTO);
 
                 internal::ApplyKernel<dim, Number, Functor> apply_kernel(
-                  func, get_data(2), src.get_values(), dst.get_values());
+                  func, get_data(2), src, dst);
 
                 Kokkos::parallel_for(
                   "dealii::MatrixFree::distributed_cell_loop_2",
@@ -1146,7 +1148,7 @@ namespace Portable
                       Kokkos::AUTO);
 
                   internal::ApplyKernel<dim, Number, Functor> apply_kernel(
-                    func, get_data(i), src.get_values(), dst.get_values());
+                    func, get_data(i), src, dst);
 
                   Kokkos::parallel_for(
                     "dealii::MatrixFree::distributed_cell_loop_" +
@@ -1183,10 +1185,7 @@ namespace Portable
                   Kokkos::AUTO);
 
               internal::ApplyKernel<dim, Number, Functor> apply_kernel(
-                func,
-                get_data(i),
-                ghosted_src.get_values(),
-                ghosted_dst.get_values());
+                func, get_data(i), ghosted_src, ghosted_dst);
 
               Kokkos::parallel_for(
                 "dealii::MatrixFree::distributed_cell_loop_" +
index 3a30aa41dbebf375d4399af8aed7dfbbead6b426..a675381953eae6ca25f10c8f143adcafca706845 100644 (file)
@@ -1398,8 +1398,8 @@ namespace MatrixFreeTools
 
       KOKKOS_FUNCTION void
       operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
-                 const Number *,
-                 Number *dst) const
+                 const Portable::DeviceVector<Number> &,
+                 Portable::DeviceVector<Number> &dst) const
       {
         Portable::
           FEEvaluation<dim, fe_degree, n_q_points_1d, n_components, Number>
index 86450f8f22488f6ba2672a9528d1a5b7f39406e6..a31f03a2308be4bd63132dfd2bdbb9878ed8f3a2 100644 (file)
@@ -37,8 +37,8 @@ public:
 
   DEAL_II_HOST_DEVICE void
   operator()(const typename Portable::MatrixFree<dim, double>::Data *data,
-             const double                                           *src,
-             double                                                 *dst) const;
+             const Portable::DeviceVector<double>                   &src,
+             Portable::DeviceVector<double>                         &dst) const;
 
   static const unsigned int n_local_dofs =
     dealii::Utilities::pow(fe_degree + 1, dim);
@@ -52,8 +52,8 @@ template <int dim, int fe_degree>
 DEAL_II_HOST_DEVICE void
 DummyOperator<dim, fe_degree>::operator()(
   const typename Portable::MatrixFree<dim, double>::Data *data,
-  const double *,
-  double *dst) const
+  const Portable::DeviceVector<double>                   &src,
+  Portable::DeviceVector<double>                         &dst) const
 {
   Kokkos::parallel_for(
     Kokkos::TeamThreadRange(data->team_member, n_q_points),
index 908b48ffa94887b7a07389ecc576035204c17b1d..ddc4443dd2f33f962ab81e6f0e42e26171569629 100644 (file)
@@ -40,8 +40,8 @@ public:
 
   DEAL_II_HOST_DEVICE void
   operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
-             const Number                                           *src,
-             Number                                                 *dst) const
+             const Portable::DeviceVector<Number>                   &src,
+             Portable::DeviceVector<Number>                         &dst) const
   {
     Portable::FEEvaluation<dim, fe_degree, n_q_points_1d, 1, Number> fe_eval(
       data);
index 556e56f69b9a6688b64f82ab5b894bc68a0d0d1d..36b0a6aa01ef5a587ec3ce395296ecf69e3c1e7a 100644 (file)
@@ -78,8 +78,8 @@ public:
 
   DEAL_II_HOST_DEVICE void
   operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
-             const Number                                           *src,
-             Number                                                 *dst) const;
+             const Portable::DeviceVector<Number>                   &src,
+             Portable::DeviceVector<Number>                         &dst) const;
 
   Number *coef;
 };
@@ -90,8 +90,8 @@ template <int dim, int fe_degree, typename Number, int n_q_points_1d>
 DEAL_II_HOST_DEVICE void
 HelmholtzOperator<dim, fe_degree, Number, n_q_points_1d>::operator()(
   const typename Portable::MatrixFree<dim, Number>::Data *data,
-  const Number                                           *src,
-  Number                                                 *dst) const
+  const Portable::DeviceVector<Number>                   &src,
+  Portable::DeviceVector<Number>                         &dst) const
 {
   Portable::FEEvaluation<dim, fe_degree, n_q_points_1d, 1, Number> fe_eval(
     data);
index a566a4aab20a71936710684eaf97d214c262144a..8439cc4b478e45f6d29a2c39b015c82e1ac8031b 100644 (file)
@@ -167,8 +167,8 @@ public:
 
   DEAL_II_HOST_DEVICE void
   operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
-             const Number                                           *src,
-             Number                                                 *dst) const
+             const Portable::DeviceVector<Number>                   &src,
+             Portable::DeviceVector<Number>                         &dst) const
   {
     Portable::FEEvaluation<dim, fe_degree, fe_degree + 1, n_components, Number>
       phi(data);

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.