]> https://gitweb.dealii.org/ - dealii.git/commitdiff
refactor Portable::MatrixFree by combining Shared and GPU data 18267/head
authorTimo Heister <timo.heister@gmail.com>
Thu, 20 Mar 2025 19:43:23 +0000 (15:43 -0400)
committerTimo Heister <timo.heister@gmail.com>
Fri, 21 Mar 2025 02:39:26 +0000 (22:39 -0400)
examples/step-64/step-64.cc
include/deal.II/matrix_free/portable_evaluation_kernels.h
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 e51fbbfaa78252d1bf2f4842e49bf03b1c045963..c8bd2cdf3e053cd97f0fac016c345931e5dca139 100644 (file)
@@ -155,11 +155,11 @@ namespace Step64
     const int q_point) const
   {
     const int cell_index = fe_eval->get_current_cell_index();
-    const typename Portable::MatrixFree<dim, double>::Data *gpu_data =
+    const typename Portable::MatrixFree<dim, double>::Data *data =
       fe_eval->get_matrix_free_data();
 
     const unsigned int position =
-      gpu_data->local_q_point_id(cell_index, n_q_points, q_point);
+      data->local_q_point_id(cell_index, n_q_points, q_point);
     auto coeff = coef[position];
 
     auto value = fe_eval->get_value(q_point);
@@ -191,11 +191,9 @@ namespace Step64
     {}
 
     DEAL_II_HOST_DEVICE void
-    operator()(const unsigned int                                      cell,
-               const typename Portable::MatrixFree<dim, double>::Data *gpu_data,
-               Portable::SharedData<dim, double> *shared_data,
-               const double                      *src,
-               double                            *dst) const;
+    operator()(const typename Portable::MatrixFree<dim, double>::Data *data,
+               const double                                           *src,
+               double *dst) const;
 
   private:
     double *coef;
@@ -209,14 +207,12 @@ namespace Step64
   // vector.
   template <int dim, int fe_degree>
   DEAL_II_HOST_DEVICE void LocalHelmholtzOperator<dim, fe_degree>::operator()(
-    const unsigned int /*cell*/,
-    const typename Portable::MatrixFree<dim, double>::Data *gpu_data,
-    Portable::SharedData<dim, double>                      *shared_data,
+    const typename Portable::MatrixFree<dim, double>::Data *data,
     const double                                           *src,
     double                                                 *dst) const
   {
     Portable::FEEvaluation<dim, fe_degree, fe_degree + 1, 1, double> fe_eval(
-      gpu_data, shared_data);
+      data);
     fe_eval.read_dof_values(src);
     fe_eval.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
     fe_eval.apply_for_each_quad_point(
index 9879a934ed39726721253f852c012e821657eadf..28c27d61d58688a1d82db384600d042fd303b910 100644 (file)
@@ -74,31 +74,30 @@ namespace Portable
       DEAL_II_HOST_DEVICE static void
       evaluate(const unsigned int                            n_components,
                const EvaluationFlags::EvaluationFlags        evaluation_flag,
-               const typename MatrixFree<dim, Number>::Data *data,
-               SharedData<dim, Number>                      *shared_data)
+               const typename MatrixFree<dim, Number>::Data *data)
       {
         if (evaluation_flag == EvaluationFlags::nothing)
           return;
 
         // the evaluator does not need temporary storage since no in-place
         // operation takes place in this function
-        auto scratch_for_eval =
-          Kokkos::subview(shared_data->scratch_pad, Kokkos::make_pair(0, 0));
+        auto scratch_for_eval = Kokkos::subview(data->shared_data->scratch_pad,
+                                                Kokkos::make_pair(0, 0));
         EvaluatorTensorProduct<EvaluatorVariant::evaluate_general,
                                dim,
                                fe_degree + 1,
                                n_q_points_1d,
                                Number>
-          eval(shared_data->team_member,
-               data->shape_values,
-               data->shape_gradients,
-               data->co_shape_gradients,
+          eval(data->team_member,
+               data->precomputed_data->shape_values,
+               data->precomputed_data->shape_gradients,
+               data->precomputed_data->co_shape_gradients,
                scratch_for_eval);
 
         for (unsigned int c = 0; c < n_components; ++c)
           {
-            auto u      = Kokkos::subview(shared_data->values, Kokkos::ALL, c);
-            auto grad_u = Kokkos::subview(shared_data->gradients,
+            auto u = Kokkos::subview(data->shared_data->values, Kokkos::ALL, c);
+            auto grad_u = Kokkos::subview(data->shared_data->gradients,
                                           Kokkos::ALL,
                                           Kokkos::ALL,
                                           c);
@@ -106,7 +105,7 @@ namespace Portable
             if constexpr (dim == 1)
               {
                 auto temp =
-                  Kokkos::subview(shared_data->scratch_pad,
+                  Kokkos::subview(data->shared_data->scratch_pad,
                                   Kokkos::make_pair(0, n_q_points_1d));
 
                 if (evaluation_flag & EvaluationFlags::gradients)
@@ -115,7 +114,7 @@ namespace Portable
                 if (evaluation_flag & EvaluationFlags::values)
                   {
                     eval.template values<0, true, false, false>(u, temp);
-                    populate_view<false>(shared_data->team_member,
+                    populate_view<false>(data->team_member,
                                          u,
                                          temp,
                                          n_q_points_1d);
@@ -124,7 +123,7 @@ namespace Portable
             else if constexpr (dim == 2)
               {
                 constexpr int temp_size = (fe_degree + 1) * n_q_points_1d;
-                auto          temp = Kokkos::subview(shared_data->scratch_pad,
+                auto temp = Kokkos::subview(data->shared_data->scratch_pad,
                                             Kokkos::make_pair(0, temp_size));
 
                 // grad x
@@ -152,10 +151,10 @@ namespace Portable
                               temp2_size = Utilities::pow(n_q_points_1d, 2) *
                                            (fe_degree + 1);
 
-                auto temp1 = Kokkos::subview(shared_data->scratch_pad,
+                auto temp1 = Kokkos::subview(data->shared_data->scratch_pad,
                                              Kokkos::make_pair(0, temp1_size));
                 auto temp2 =
-                  Kokkos::subview(shared_data->scratch_pad,
+                  Kokkos::subview(data->shared_data->scratch_pad,
                                   Kokkos::make_pair(temp1_size,
                                                     temp1_size + temp2_size));
 
@@ -200,31 +199,30 @@ namespace Portable
       DEAL_II_HOST_DEVICE static void
       integrate(const unsigned int                            n_components,
                 const EvaluationFlags::EvaluationFlags        integration_flag,
-                const typename MatrixFree<dim, Number>::Data *data,
-                SharedData<dim, Number>                      *shared_data)
+                const typename MatrixFree<dim, Number>::Data *data)
       {
         if (integration_flag == EvaluationFlags::nothing)
           return;
 
         // the evaluator does not need temporary storage since no in-place
         // operation takes place in this function
-        auto scratch_for_eval =
-          Kokkos::subview(shared_data->scratch_pad, Kokkos::make_pair(0, 0));
+        auto scratch_for_eval = Kokkos::subview(data->shared_data->scratch_pad,
+                                                Kokkos::make_pair(0, 0));
         EvaluatorTensorProduct<EvaluatorVariant::evaluate_general,
                                dim,
                                fe_degree + 1,
                                n_q_points_1d,
                                Number>
-          eval(shared_data->team_member,
-               data->shape_values,
-               data->shape_gradients,
-               data->co_shape_gradients,
+          eval(data->team_member,
+               data->precomputed_data->shape_values,
+               data->precomputed_data->shape_gradients,
+               data->precomputed_data->co_shape_gradients,
                scratch_for_eval);
 
         for (unsigned int c = 0; c < n_components; ++c)
           {
-            auto u      = Kokkos::subview(shared_data->values, Kokkos::ALL, c);
-            auto grad_u = Kokkos::subview(shared_data->gradients,
+            auto u = Kokkos::subview(data->shared_data->values, Kokkos::ALL, c);
+            auto grad_u = Kokkos::subview(data->shared_data->gradients,
                                           Kokkos::ALL,
                                           Kokkos::ALL,
                                           c);
@@ -232,14 +230,14 @@ namespace Portable
             if constexpr (dim == 1)
               {
                 auto temp =
-                  Kokkos::subview(shared_data->scratch_pad,
+                  Kokkos::subview(data->shared_data->scratch_pad,
                                   Kokkos::make_pair(0, fe_degree + 1));
 
                 if ((integration_flag & EvaluationFlags::values) &&
                     !(integration_flag & EvaluationFlags::gradients))
                   {
                     eval.template values<0, false, false, false>(u, temp);
-                    populate_view<false>(shared_data->team_member,
+                    populate_view<false>(data->team_member,
                                          u,
                                          temp,
                                          fe_degree + 1);
@@ -251,7 +249,7 @@ namespace Portable
                         eval.template values<0, false, false, false>(u, temp);
                         eval.template gradients<0, false, true, false>(
                           Kokkos::subview(grad_u, Kokkos::ALL, 0), temp);
-                        populate_view<false>(shared_data->team_member,
+                        populate_view<false>(data->team_member,
                                              u,
                                              temp,
                                              fe_degree + 1);
@@ -264,7 +262,7 @@ namespace Portable
             else if constexpr (dim == 2)
               {
                 constexpr int temp_size = (fe_degree + 1) * n_q_points_1d;
-                auto          temp = Kokkos::subview(shared_data->scratch_pad,
+                auto temp = Kokkos::subview(data->shared_data->scratch_pad,
                                             Kokkos::make_pair(0, temp_size));
 
                 if ((integration_flag & EvaluationFlags::values) &&
@@ -292,10 +290,10 @@ namespace Portable
                               temp2_size = Utilities::pow(fe_degree + 1, 2) *
                                            n_q_points_1d;
 
-                auto temp1 = Kokkos::subview(shared_data->scratch_pad,
+                auto temp1 = Kokkos::subview(data->shared_data->scratch_pad,
                                              Kokkos::make_pair(0, temp1_size));
                 auto temp2 =
-                  Kokkos::subview(shared_data->scratch_pad,
+                  Kokkos::subview(data->shared_data->scratch_pad,
                                   Kokkos::make_pair(temp1_size,
                                                     temp1_size + temp2_size));
 
@@ -350,8 +348,7 @@ namespace Portable
       DEAL_II_HOST_DEVICE static void
       evaluate(const unsigned int                            n_components,
                const EvaluationFlags::EvaluationFlags        evaluation_flag,
-               const typename MatrixFree<dim, Number>::Data *data,
-               SharedData<dim, Number>                      *shared_data)
+               const typename MatrixFree<dim, Number>::Data *data)
       {
         // since the dof values have already been stored in
         // shared_data->values, there is nothing to do if the gradients are
@@ -360,7 +357,7 @@ namespace Portable
           return;
 
         constexpr int n_points = Utilities::pow(fe_degree + 1, dim);
-        auto scratch_for_eval  = Kokkos::subview(shared_data->scratch_pad,
+        auto scratch_for_eval  = Kokkos::subview(data->shared_data->scratch_pad,
                                                 Kokkos::make_pair(0, n_points));
 
         EvaluatorTensorProduct<EvaluatorVariant::evaluate_general,
@@ -368,16 +365,16 @@ namespace Portable
                                fe_degree + 1,
                                fe_degree + 1,
                                Number>
-          eval(shared_data->team_member,
-               data->shape_values,
-               data->shape_gradients,
-               data->co_shape_gradients,
+          eval(data->team_member,
+               data->precomputed_data->shape_values,
+               data->precomputed_data->shape_gradients,
+               data->precomputed_data->co_shape_gradients,
                scratch_for_eval);
 
         for (unsigned int c = 0; c < n_components; ++c)
           {
-            auto u      = Kokkos::subview(shared_data->values, Kokkos::ALL, c);
-            auto grad_u = Kokkos::subview(shared_data->gradients,
+            auto u = Kokkos::subview(data->shared_data->values, Kokkos::ALL, c);
+            auto grad_u = Kokkos::subview(data->shared_data->gradients,
                                           Kokkos::ALL,
                                           Kokkos::ALL,
                                           c);
@@ -397,8 +394,7 @@ namespace Portable
       DEAL_II_HOST_DEVICE static void
       integrate(const unsigned int                            n_components,
                 const EvaluationFlags::EvaluationFlags        integration_flag,
-                const typename MatrixFree<dim, Number>::Data *data,
-                SharedData<dim, Number>                      *shared_data)
+                const typename MatrixFree<dim, Number>::Data *data)
       {
         // since the quad values have already been stored in
         // shared_data->values, there is nothing to do if the gradients are
@@ -407,7 +403,7 @@ namespace Portable
           return;
 
         constexpr int n_points = Utilities::pow(fe_degree + 1, dim);
-        auto scratch_for_eval  = Kokkos::subview(shared_data->scratch_pad,
+        auto scratch_for_eval  = Kokkos::subview(data->shared_data->scratch_pad,
                                                 Kokkos::make_pair(0, n_points));
 
         EvaluatorTensorProduct<EvaluatorVariant::evaluate_general,
@@ -415,16 +411,16 @@ namespace Portable
                                fe_degree + 1,
                                fe_degree + 1,
                                Number>
-          eval(shared_data->team_member,
-               data->shape_values,
-               data->shape_gradients,
-               data->co_shape_gradients,
+          eval(data->team_member,
+               data->precomputed_data->shape_values,
+               data->precomputed_data->shape_gradients,
+               data->precomputed_data->co_shape_gradients,
                scratch_for_eval);
 
         for (unsigned int c = 0; c < n_components; ++c)
           {
-            auto u      = Kokkos::subview(shared_data->values, Kokkos::ALL, c);
-            auto grad_u = Kokkos::subview(shared_data->gradients,
+            auto u = Kokkos::subview(data->shared_data->values, Kokkos::ALL, c);
+            auto grad_u = Kokkos::subview(data->shared_data->gradients,
                                           Kokkos::ALL,
                                           Kokkos::ALL,
                                           c);
@@ -486,12 +482,11 @@ namespace Portable
       DEAL_II_HOST_DEVICE static void
       evaluate(const unsigned int                            n_components,
                const EvaluationFlags::EvaluationFlags        evaluation_flag,
-               const typename MatrixFree<dim, Number>::Data *data,
-               SharedData<dim, Number>                      *shared_data)
+               const typename MatrixFree<dim, Number>::Data *data)
       {
         constexpr int scratch_size = Utilities::pow(n_q_points_1d, dim);
         auto          scratch_for_eval =
-          Kokkos::subview(shared_data->scratch_pad,
+          Kokkos::subview(data->shared_data->scratch_pad,
                           Kokkos::make_pair(0, scratch_size));
 
         EvaluatorTensorProduct<EvaluatorVariant::evaluate_general,
@@ -499,16 +494,16 @@ namespace Portable
                                fe_degree + 1,
                                n_q_points_1d,
                                Number>
-          eval(shared_data->team_member,
-               data->shape_values,
-               data->shape_gradients,
-               data->co_shape_gradients,
+          eval(data->team_member,
+               data->precomputed_data->shape_values,
+               data->precomputed_data->shape_gradients,
+               data->precomputed_data->co_shape_gradients,
                scratch_for_eval);
 
         for (unsigned int c = 0; c < n_components; ++c)
           {
-            auto u      = Kokkos::subview(shared_data->values, Kokkos::ALL, c);
-            auto grad_u = Kokkos::subview(shared_data->gradients,
+            auto u = Kokkos::subview(data->shared_data->values, Kokkos::ALL, c);
+            auto grad_u = Kokkos::subview(data->shared_data->gradients,
                                           Kokkos::ALL,
                                           Kokkos::ALL,
                                           c);
@@ -537,12 +532,11 @@ namespace Portable
       DEAL_II_HOST_DEVICE static void
       integrate(const unsigned int                            n_components,
                 const EvaluationFlags::EvaluationFlags        integration_flag,
-                const typename MatrixFree<dim, Number>::Data *data,
-                const SharedData<dim, Number>                *shared_data)
+                const typename MatrixFree<dim, Number>::Data *data)
       {
         constexpr int scratch_size = Utilities::pow(n_q_points_1d, dim);
         auto          scratch_for_eval =
-          Kokkos::subview(shared_data->scratch_pad,
+          Kokkos::subview(data->shared_data->scratch_pad,
                           Kokkos::make_pair(0, scratch_size));
 
         EvaluatorTensorProduct<EvaluatorVariant::evaluate_general,
@@ -550,16 +544,16 @@ namespace Portable
                                fe_degree + 1,
                                n_q_points_1d,
                                Number>
-          eval(shared_data->team_member,
-               data->shape_values,
-               data->shape_gradients,
-               data->co_shape_gradients,
+          eval(data->team_member,
+               data->precomputed_data->shape_values,
+               data->precomputed_data->shape_gradients,
+               data->precomputed_data->co_shape_gradients,
                scratch_for_eval);
 
         for (unsigned int c = 0; c < n_components; ++c)
           {
-            auto u      = Kokkos::subview(shared_data->values, Kokkos::ALL, c);
-            auto grad_u = Kokkos::subview(shared_data->gradients,
+            auto u = Kokkos::subview(data->shared_data->values, Kokkos::ALL, c);
+            auto grad_u = Kokkos::subview(data->shared_data->gradients,
                                           Kokkos::ALL,
                                           Kokkos::ALL,
                                           c);
index 21de3dccb04fa5cbff27779f5119c7530e51dcce..ce75f15eb50455076b344aa2828621a83e3aa098 100644 (file)
@@ -124,7 +124,7 @@ namespace Portable
      * Constructor.
      */
     DEAL_II_HOST_DEVICE
-    FEEvaluation(const data_type *data, SharedData<dim, Number> *shdata);
+    explicit FEEvaluation(const data_type *data);
 
     /**
      * Return the index of the current cell.
@@ -264,9 +264,10 @@ namespace Portable
     apply_for_each_quad_point(const Functor &func);
 
   private:
-    const data_type         *data;
-    SharedData<dim, Number> *shared_data;
-    int                      cell_id;
+    const typename MatrixFree<dim, Number>::Data            *data;
+    const typename MatrixFree<dim, Number>::PrecomputedData *precomputed_data;
+    SharedData<dim, Number>                                 *shared_data;
+    int                                                      cell_id;
   };
 
 
@@ -278,10 +279,11 @@ namespace Portable
             typename Number>
   DEAL_II_HOST_DEVICE
   FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
-    FEEvaluation(const data_type *data, SharedData<dim, Number> *shdata)
+    FEEvaluation(const data_type *data)
     : data(data)
-    , shared_data(shdata)
-    , cell_id(shared_data->team_member.league_rank())
+    , precomputed_data(data->precomputed_data)
+    , shared_data(data->shared_data)
+    , cell_id(data->team_member.league_rank())
   {}
 
 
@@ -328,22 +330,22 @@ namespace Portable
     read_dof_values(const Number *src)
   {
     // Populate the scratch memory
-    Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
+    Kokkos::parallel_for(Kokkos::TeamThreadRange(data->team_member,
                                                  tensor_dofs_per_component),
                          [&](const int &i) {
                            for (unsigned int c = 0; c < n_components_; ++c)
                              shared_data->values(i, c) =
-                               src[data->local_to_global(
+                               src[precomputed_data->local_to_global(
                                  cell_id, i + tensor_dofs_per_component * c)];
                          });
-    shared_data->team_member.team_barrier();
+    data->team_member.team_barrier();
 
     for (unsigned int c = 0; c < n_components_; ++c)
       {
         internal::resolve_hanging_nodes<dim, fe_degree, false, Number>(
-          shared_data->team_member,
-          data->constraint_weights,
-          data->constraint_mask(cell_id * n_components + c),
+          data->team_member,
+          precomputed_data->constraint_weights,
+          precomputed_data->constraint_mask(cell_id * n_components + c),
           Kokkos::subview(shared_data->values, Kokkos::ALL, c));
       }
   }
@@ -362,32 +364,30 @@ namespace Portable
     for (unsigned int c = 0; c < n_components_; ++c)
       {
         internal::resolve_hanging_nodes<dim, fe_degree, true, Number>(
-          shared_data->team_member,
-          data->constraint_weights,
-          data->constraint_mask(cell_id * n_components + c),
+          data->team_member,
+          precomputed_data->constraint_weights,
+          precomputed_data->constraint_mask(cell_id * n_components + c),
           Kokkos::subview(shared_data->values, Kokkos::ALL, c));
       }
 
-    if (data->use_coloring)
+    if (precomputed_data->use_coloring)
       {
         Kokkos::parallel_for(
-          Kokkos::TeamThreadRange(shared_data->team_member,
-                                  tensor_dofs_per_component),
+          Kokkos::TeamThreadRange(data->team_member, tensor_dofs_per_component),
           [&](const int &i) {
             for (unsigned int c = 0; c < n_components_; ++c)
-              dst[data->local_to_global(cell_id,
-                                        i + tensor_dofs_per_component * c)] +=
+              dst[precomputed_data->local_to_global(
+                cell_id, i + tensor_dofs_per_component * c)] +=
                 shared_data->values(i, c);
           });
       }
     else
       {
         Kokkos::parallel_for(
-          Kokkos::TeamThreadRange(shared_data->team_member,
-                                  tensor_dofs_per_component),
+          Kokkos::TeamThreadRange(data->team_member, tensor_dofs_per_component),
           [&](const int &i) {
             for (unsigned int c = 0; c < n_components_; ++c)
-              Kokkos::atomic_add(&dst[data->local_to_global(
+              Kokkos::atomic_add(&dst[precomputed_data->local_to_global(
                                    cell_id, i + (tensor_dofs_per_component)*c)],
                                  shared_data->values(i, c));
           });
@@ -408,28 +408,29 @@ namespace Portable
     using ElementType = ::dealii::internal::MatrixFreeFunctions::ElementType;
 
     if (fe_degree >= 0 && fe_degree + 1 == n_q_points_1d &&
-        data->element_type == ElementType::tensor_symmetric_collocation)
+        precomputed_data->element_type ==
+          ElementType::tensor_symmetric_collocation)
       {
         internal::FEEvaluationImplCollocation<dim, fe_degree, Number>::evaluate(
-          n_components, evaluation_flag, data, shared_data);
+          n_components, evaluation_flag, data);
       }
     // '<=' on type means tensor_symmetric or tensor_symmetric_hermite, see
     // shape_info.h for more details
     else if (fe_degree >= 0 &&
              internal::use_collocation_evaluation(fe_degree, n_q_points_1d) &&
-             data->element_type <= ElementType::tensor_symmetric)
+             precomputed_data->element_type <= ElementType::tensor_symmetric)
       {
         internal::FEEvaluationImplTransformToCollocation<
           dim,
           fe_degree,
           n_q_points_1d,
-          Number>::evaluate(n_components, evaluation_flag, data, shared_data);
+          Number>::evaluate(n_components, evaluation_flag, data);
       }
-    else if (fe_degree >= 0 &&
-             data->element_type <= ElementType::tensor_symmetric_no_collocation)
+    else if (fe_degree >= 0 && precomputed_data->element_type <=
+                                 ElementType::tensor_symmetric_no_collocation)
       {
         internal::FEEvaluationImpl<dim, fe_degree, n_q_points_1d, Number>::
-          evaluate(n_components, evaluation_flag, data, shared_data);
+          evaluate(n_components, evaluation_flag, data);
       }
     else
       {
@@ -469,28 +470,29 @@ namespace Portable
     using ElementType = ::dealii::internal::MatrixFreeFunctions::ElementType;
 
     if (fe_degree >= 0 && fe_degree + 1 == n_q_points_1d &&
-        data->element_type == ElementType::tensor_symmetric_collocation)
+        precomputed_data->element_type ==
+          ElementType::tensor_symmetric_collocation)
       {
         internal::FEEvaluationImplCollocation<dim, fe_degree, Number>::
-          integrate(n_components, integration_flag, data, shared_data);
+          integrate(n_components, integration_flag, data);
       }
     // '<=' on type means tensor_symmetric or tensor_symmetric_hermite, see
     // shape_info.h for more details
     else if (fe_degree >= 0 &&
              internal::use_collocation_evaluation(fe_degree, n_q_points_1d) &&
-             data->element_type <= ElementType::tensor_symmetric)
+             precomputed_data->element_type <= ElementType::tensor_symmetric)
       {
         internal::FEEvaluationImplTransformToCollocation<
           dim,
           fe_degree,
           n_q_points_1d,
-          Number>::integrate(n_components, integration_flag, data, shared_data);
+          Number>::integrate(n_components, integration_flag, data);
       }
-    else if (fe_degree >= 0 &&
-             data->element_type <= ElementType::tensor_symmetric_no_collocation)
+    else if (fe_degree >= 0 && precomputed_data->element_type <=
+                                 ElementType::tensor_symmetric_no_collocation)
       {
         internal::FEEvaluationImpl<dim, fe_degree, n_q_points_1d, Number>::
-          integrate(n_components, integration_flag, data, shared_data);
+          integrate(n_components, integration_flag, data);
       }
     else
       {
@@ -588,13 +590,14 @@ namespace Portable
     Assert(q_point >= 0 && q_point < n_q_points, ExcInternalError());
     if constexpr (n_components_ == 1)
       {
-        shared_data->values(q_point, 0) = val_in * data->JxW(cell_id, q_point);
+        shared_data->values(q_point, 0) =
+          val_in * precomputed_data->JxW(cell_id, q_point);
       }
     else
       {
         for (unsigned int c = 0; c < n_components; ++c)
           shared_data->values(q_point, c) =
-            val_in[c] * data->JxW(cell_id, q_point);
+            val_in[c] * precomputed_data->JxW(cell_id, q_point);
       }
   }
 
@@ -645,8 +648,9 @@ namespace Portable
           {
             Number tmp = 0.;
             for (unsigned int d_2 = 0; d_2 < dim; ++d_2)
-              tmp += data->inv_jacobian(cell_id, q_point, d_2, d_1) *
-                     shared_data->gradients(q_point, d_2, 0);
+              tmp +=
+                precomputed_data->inv_jacobian(cell_id, q_point, d_2, d_1) *
+                shared_data->gradients(q_point, d_2, 0);
             grad[d_1] = tmp;
           }
       }
@@ -657,8 +661,9 @@ namespace Portable
             {
               Number tmp = 0.;
               for (unsigned int d_2 = 0; d_2 < dim; ++d_2)
-                tmp += data->inv_jacobian(cell_id, q_point, d_2, d_1) *
-                       shared_data->gradients(q_point, d_2, c);
+                tmp +=
+                  precomputed_data->inv_jacobian(cell_id, q_point, d_2, d_1) *
+                  shared_data->gradients(q_point, d_2, c);
               grad[c][d_1] = tmp;
             }
       }
@@ -685,9 +690,10 @@ namespace Portable
             Number tmp = 0.;
             for (unsigned int d_2 = 0; d_2 < dim; ++d_2)
               tmp +=
-                data->inv_jacobian(cell_id, q_point, d_1, d_2) * grad_in[d_2];
+                precomputed_data->inv_jacobian(cell_id, q_point, d_1, d_2) *
+                grad_in[d_2];
             shared_data->gradients(q_point, d_1, 0) =
-              tmp * data->JxW(cell_id, q_point);
+              tmp * precomputed_data->JxW(cell_id, q_point);
           }
       }
     else
@@ -697,10 +703,11 @@ namespace Portable
             {
               Number tmp = 0.;
               for (unsigned int d_2 = 0; d_2 < dim; ++d_2)
-                tmp += data->inv_jacobian(cell_id, q_point, d_1, d_2) *
-                       grad_in[c][d_2];
+                tmp +=
+                  precomputed_data->inv_jacobian(cell_id, q_point, d_1, d_2) *
+                  grad_in[c][d_2];
               shared_data->gradients(q_point, d_1, c) =
-                tmp * data->JxW(cell_id, q_point);
+                tmp * precomputed_data->JxW(cell_id, q_point);
             }
       }
   }
@@ -717,10 +724,9 @@ namespace Portable
   FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
     apply_for_each_quad_point(const Functor &func)
   {
-    Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
-                                                 n_q_points),
+    Kokkos::parallel_for(Kokkos::TeamThreadRange(data->team_member, n_q_points),
                          [&](const int &i) { func(this, i); });
-    shared_data->team_member.team_barrier();
+    data->team_member.team_barrier();
   }
 
 
index dee00888322b7999d1eddd40a5f77bf06b12ff51..bb621aa5b4096d41349662f61d65de00cd684801 100644 (file)
@@ -61,6 +61,8 @@ namespace Portable
     template <int dim, typename Number>
     class ReinitHelper;
   }
+  template <int dim, typename Number>
+  struct SharedData;
 #endif
 
   /**
@@ -152,9 +154,10 @@ namespace Portable
 
     /**
      * Structure which is passed to the kernel. It is used to pass all the
-     * necessary information from the CPU to the GPU.
+     * necessary information from the CPU to the GPU and is precomputed
+     * on the CPU. This data is read-only once we run on the GPU.
      */
-    struct Data
+    struct PrecomputedData
     {
       /**
        * Kokkos::View of the quadrature points.
@@ -247,6 +250,29 @@ namespace Portable
        * Size of the scratch pad for temporary storage in shared memory.
        */
       unsigned int scratch_pad_size;
+    };
+
+
+    /**
+     * A pointer to this data structure is passed to the user code in
+     * device code that computes operations in quadrature points. It
+     * can be used to construct Portable::FEEvaluation objects and contains
+     * necessary precomputed data and scratch memory space to perform
+     * matrix-free evaluations.
+     */
+    struct Data
+    {
+      using TeamHandle = Kokkos::TeamPolicy<
+        MemorySpace::Default::kokkos_space::execution_space>::member_type;
+
+      /**
+       * TeamPolicy handle.
+       */
+      TeamHandle team_member;
+
+      const int                cell_index;
+      const PrecomputedData   *precomputed_data;
+      SharedData<dim, Number> *shared_data;
 
       /**
        * Return the quadrature point index local. The index is
@@ -257,7 +283,10 @@ namespace Portable
                        const unsigned int n_q_points,
                        const unsigned int q_point) const
       {
-        return (row_start / padding_length + cell) * n_q_points + q_point;
+        return (precomputed_data->row_start / precomputed_data->padding_length +
+                cell) *
+                 n_q_points +
+               q_point;
       }
 
 
@@ -269,7 +298,7 @@ namespace Portable
       get_quadrature_point(const unsigned int cell,
                            const unsigned int q_point) const
       {
-        return q_points(cell, q_point);
+        return precomputed_data->q_points(cell, q_point);
       }
     };
 
@@ -325,7 +354,7 @@ namespace Portable
     /**
      * Return the Data structure associated with @p color.
      */
-    Data
+    PrecomputedData
     get_data(unsigned int color) const;
 
     // clang-format off
@@ -336,11 +365,9 @@ namespace Portable
      * @p func needs to define
      * \code
      * DEAL_II_HOST_DEVICE void operator()(
-     *   const unsigned int                                          cell,
-     *   const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-     *   Portable::SharedData<dim, Number> *                     shared_data,
-     *   const Number *                                              src,
-     *   Number *                                                    dst) const;
+     *   const typename Portable::MatrixFree<dim, Number>::Data *data,
+     *   const Number *                                          src,
+     *   Number *                                                dst) const;
      *   static const unsigned int n_local_dofs;
      *   static const unsigned int n_q_points;
      * \endcode
@@ -360,8 +387,7 @@ namespace Portable
      * @p func needs to define
      * \code
      *  DEAL_II_HOST_DEVICE void operator()(
-     *    const unsigned int                                          cell,
-     *    const typename Portable::MatrixFree<dim, Number>::Data *gpu_data);
+     *    const typename Portable::MatrixFree<dim, Number>::Data *data);
      * static const unsigned int n_local_dofs;
      * static const unsigned int n_q_points;
      * \endcode
@@ -642,9 +668,6 @@ namespace Portable
   template <int dim, typename Number>
   struct SharedData
   {
-    using TeamHandle = Kokkos::TeamPolicy<
-      MemorySpace::Default::kokkos_space::execution_space>::member_type;
-
     using SharedViewValues = Kokkos::View<
       Number **,
       MemorySpace::Default::kokkos_space::execution_space::scratch_memory_space,
@@ -659,21 +682,14 @@ namespace Portable
       Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
 
     DEAL_II_HOST_DEVICE
-    SharedData(const TeamHandle           &team_member,
-               const SharedViewValues     &values,
+    SharedData(const SharedViewValues     &values,
                const SharedViewGradients  &gradients,
                const SharedViewScratchPad &scratch_pad)
-      : team_member(team_member)
-      , values(values)
+      : values(values)
       , gradients(gradients)
       , scratch_pad(scratch_pad)
     {}
 
-    /**
-     * TeamPolicy handle.
-     */
-    TeamHandle team_member;
-
     /**
      * Memory for dof and quad values.
      */
@@ -792,7 +808,8 @@ namespace Portable
   template <int dim, typename Number>
   DataHost<dim, Number>
   copy_mf_data_to_host(
-    const typename dealii::Portable::MatrixFree<dim, Number>::Data &data,
+    const typename dealii::Portable::MatrixFree<dim, Number>::PrecomputedData
+                      &data,
     const UpdateFlags &update_flags)
   {
     DataHost<dim, Number> data_host;
index 1a85a456de14c958a79a33ee3acfdb63a72568de..829bf9d6745eab3f7671a3b30f0834d1cb3a2d77 100644 (file)
@@ -351,20 +351,21 @@ namespace Portable
                        scratch_memory_space,
                      Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
 
-      ApplyKernel(Functor                                      func,
-                  const typename MatrixFree<dim, Number>::Data gpu_data,
-                  Number *const                                src,
-                  Number                                      *dst)
+      ApplyKernel(
+        Functor                                                 func,
+        const typename MatrixFree<dim, Number>::PrecomputedData gpu_data,
+        Number *const                                           src,
+        Number                                                 *dst)
         : func(func)
         , gpu_data(gpu_data)
         , src(src)
         , dst(dst)
       {}
 
-      Functor                                      func;
-      const typename MatrixFree<dim, Number>::Data gpu_data;
-      Number *const                                src;
-      Number                                      *dst;
+      Functor                                                 func;
+      const typename MatrixFree<dim, Number>::PrecomputedData gpu_data;
+      Number *const                                           src;
+      Number                                                 *dst;
 
 
       // Provide the shared memory capacity. This function takes the team_size
@@ -396,11 +397,14 @@ namespace Portable
         SharedViewScratchPad scratch_pad(team_member.team_shmem(),
                                          gpu_data.scratch_pad_size);
 
-        SharedData<dim, Number> shared_data(team_member,
-                                            values,
-                                            gradients,
-                                            scratch_pad);
-        func(team_member.league_rank(), &gpu_data, &shared_data, src, dst);
+        SharedData<dim, Number> shared_data(values, gradients, scratch_pad);
+
+        typename MatrixFree<dim, Number>::Data data{team_member,
+                                                    team_member.league_rank(),
+                                                    &gpu_data,
+                                                    &shared_data};
+
+        func(&data, src, dst);
       }
     };
   } // namespace internal
@@ -487,10 +491,10 @@ namespace Portable
 
 
   template <int dim, typename Number>
-  typename MatrixFree<dim, Number>::Data
+  typename MatrixFree<dim, Number>::PrecomputedData
   MatrixFree<dim, Number>::get_data(unsigned int color) const
   {
-    Data data_copy;
+    PrecomputedData data_copy;
     if (q_points.size() > 0)
       data_copy.q_points = q_points[color];
     if (inv_jacobian.size() > 0)
@@ -622,23 +626,45 @@ namespace Portable
   void
   MatrixFree<dim, Number>::evaluate_coefficients(Functor func) const
   {
+    const unsigned int n_q_points = Functor::n_q_points;
+
     for (unsigned int i = 0; i < n_colors; ++i)
       if (n_cells[i] > 0)
         {
-          MemorySpace::Default::kokkos_space::execution_space exec;
           auto color_data = get_data(i);
-          Kokkos::parallel_for(
-            "dealii::MatrixFree::evaluate_coeff",
-            Kokkos::MDRangePolicy<
-              MemorySpace::Default::kokkos_space::execution_space,
-              Kokkos::Rank<2>>(
+
+          MemorySpace::Default::kokkos_space::execution_space exec;
+          Kokkos::TeamPolicy<
+            MemorySpace::Default::kokkos_space::execution_space>
+            team_policy(
 #if KOKKOS_VERSION >= 20900
               exec,
 #endif
-              {0, 0},
-              {n_cells[i], Functor::n_q_points}),
-            KOKKOS_LAMBDA(const int cell, const int q) {
-              func(&color_data, cell, q);
+              n_cells[i],
+              Kokkos::AUTO);
+
+          Kokkos::parallel_for(
+            "dealii::MatrixFree::evaluate_coeff_cell_loop",
+            team_policy,
+            KOKKOS_LAMBDA(const Kokkos::TeamPolicy<
+                          MemorySpace::Default::kokkos_space::execution_space>::
+                            member_type &team_member) {
+              Kokkos::parallel_for(
+#if KOKKOS_VERSION >= 20900
+                Kokkos::TeamVectorRange(team_member, n_q_points),
+#else
+                Kokkos::TeamThreadRange(team_member, n_q_points),
+#endif
+                [&](const int q_point) {
+                  const int cell = team_member.league_rank();
+
+                  Data data{team_member,
+                            cell,
+                            &color_data,
+                            /*shared_data*/ nullptr};
+
+                  func(&data, cell, q_point);
+                });
             });
         }
   }
index 00545cbfb53319143ba1c88dd9fdb516e7a5800b..3a30aa41dbebf375d4399af8aed7dfbbead6b426 100644 (file)
@@ -1397,16 +1397,16 @@ namespace MatrixFreeTools
       {}
 
       KOKKOS_FUNCTION void
-      operator()(
-        const unsigned int                                      cell,
-        const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-        Portable::SharedData<dim, Number>                      *shared_data,
-        const Number *,
-        Number *dst) const
+      operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
+                 const Number *,
+                 Number *dst) const
       {
         Portable::
           FEEvaluation<dim, fe_degree, n_q_points_1d, n_components, Number>
-            fe_eval(gpu_data, shared_data);
+            fe_eval(data);
+        const typename Portable::MatrixFree<dim, Number>::PrecomputedData
+                 *gpu_data = data->precomputed_data;
+        const int cell     = data->cell_index;
 
         constexpr int dofs_per_cell = decltype(fe_eval)::tensor_dofs_per_cell;
         typename decltype(fe_eval)::value_type
@@ -1416,7 +1416,7 @@ namespace MatrixFreeTools
             const auto c = i % n_components;
 
             Kokkos::parallel_for(
-              Kokkos::TeamThreadRange(shared_data->team_member,
+              Kokkos::TeamThreadRange(data->team_member,
                                       dofs_per_cell / n_components),
               [&](int j) {
                 typename decltype(fe_eval)::value_type val = {};
@@ -1433,14 +1433,14 @@ namespace MatrixFreeTools
                 fe_eval.submit_dof_value(val, j);
               });
 
-            shared_data->team_member.team_barrier();
+            data->team_member.team_barrier();
 
             Portable::internal::
               resolve_hanging_nodes<dim, fe_degree, false, Number>(
-                shared_data->team_member,
+                data->team_member,
                 gpu_data->constraint_weights,
                 gpu_data->constraint_mask(cell * n_components + c),
-                Kokkos::subview(shared_data->values, Kokkos::ALL, c));
+                Kokkos::subview(data->shared_data->values, Kokkos::ALL, c));
 
             fe_eval.evaluate(m_evaluation_flags);
             fe_eval.apply_for_each_quad_point(m_quad_operation);
@@ -1448,12 +1448,12 @@ namespace MatrixFreeTools
 
             Portable::internal::
               resolve_hanging_nodes<dim, fe_degree, true, Number>(
-                shared_data->team_member,
+                data->team_member,
                 gpu_data->constraint_weights,
                 gpu_data->constraint_mask(cell * n_components + c),
-                Kokkos::subview(shared_data->values, Kokkos::ALL, c));
+                Kokkos::subview(data->shared_data->values, Kokkos::ALL, c));
 
-            Kokkos::single(Kokkos::PerTeam(shared_data->team_member), [&] {
+            Kokkos::single(Kokkos::PerTeam(data->team_member), [&] {
               if constexpr (n_components == 1)
                 diagonal[i] = fe_eval.get_dof_value(i);
               else
@@ -1461,37 +1461,37 @@ namespace MatrixFreeTools
                   fe_eval.get_dof_value(i / n_components)[i % n_components];
             });
 
-            shared_data->team_member.team_barrier();
+            data->team_member.team_barrier();
           }
 
-        Kokkos::single(Kokkos::PerTeam(shared_data->team_member), [&] {
+        Kokkos::single(Kokkos::PerTeam(data->team_member), [&] {
           for (unsigned int i = 0; i < dofs_per_cell / n_components; ++i)
             fe_eval.submit_dof_value(diagonal[i], i);
         });
 
-        shared_data->team_member.team_barrier();
+        data->team_member.team_barrier();
 
         // We need to do the same as distribute_local_to_global but without
         // constraints since we have already taken care of them earlier
         if (gpu_data->use_coloring)
           {
             Kokkos::parallel_for(
-              Kokkos::TeamThreadRange(shared_data->team_member, dofs_per_cell),
+              Kokkos::TeamThreadRange(data->team_member, dofs_per_cell),
               [&](const int &i) {
                 dst[gpu_data->local_to_global(cell, i)] +=
-                  shared_data->values(i % (dofs_per_cell / n_components),
-                                      i / (dofs_per_cell / n_components));
+                  data->shared_data->values(i % (dofs_per_cell / n_components),
+                                            i / (dofs_per_cell / n_components));
               });
           }
         else
           {
             Kokkos::parallel_for(
-              Kokkos::TeamThreadRange(shared_data->team_member, dofs_per_cell),
+              Kokkos::TeamThreadRange(data->team_member, dofs_per_cell),
               [&](const int &i) {
-                Kokkos::atomic_add(
-                  &dst[gpu_data->local_to_global(cell, i)],
-                  shared_data->values(i % (dofs_per_cell / n_components),
-                                      i / (dofs_per_cell / n_components)));
+                Kokkos::atomic_add(&dst[gpu_data->local_to_global(cell, i)],
+                                   data->shared_data->values(
+                                     i % (dofs_per_cell / n_components),
+                                     i / (dofs_per_cell / n_components)));
               });
           }
       };
index a7fc778f0934c63060888e57a95c67f3d0636ff7..86450f8f22488f6ba2672a9528d1a5b7f39406e6 100644 (file)
@@ -36,11 +36,9 @@ public:
   DummyOperator() = default;
 
   DEAL_II_HOST_DEVICE void
-  operator()(const unsigned int                                      cell,
-             const typename Portable::MatrixFree<dim, double>::Data *gpu_data,
-             Portable::SharedData<dim, double> *shared_data,
-             const double                      *src,
-             double                            *dst) const;
+  operator()(const typename Portable::MatrixFree<dim, double>::Data *data,
+             const double                                           *src,
+             double                                                 *dst) const;
 
   static const unsigned int n_local_dofs =
     dealii::Utilities::pow(fe_degree + 1, dim);
@@ -53,19 +51,17 @@ public:
 template <int dim, int fe_degree>
 DEAL_II_HOST_DEVICE void
 DummyOperator<dim, fe_degree>::operator()(
-  const unsigned int                                      cell,
-  const typename Portable::MatrixFree<dim, double>::Data *gpu_data,
-  Portable::SharedData<dim, double>                      *shared_data,
+  const typename Portable::MatrixFree<dim, double>::Data *data,
   const double *,
   double *dst) const
 {
   Kokkos::parallel_for(
-    Kokkos::TeamThreadRange(shared_data->team_member, n_q_points),
+    Kokkos::TeamThreadRange(data->team_member, n_q_points),
     [&](const int q_point) {
       const unsigned int pos =
-        gpu_data->local_q_point_id(cell, n_q_points, q_point);
+        data->local_q_point_id(data->cell_index, n_q_points, q_point);
 
-      auto point = gpu_data->get_quadrature_point(cell, q_point);
+      auto point = data->get_quadrature_point(data->cell_index, q_point);
       dst[pos] =
         dim == 2 ? point[0] + point[1] : point[0] + point[1] + point[2];
     });
@@ -149,7 +145,7 @@ test()
   const unsigned int n_colors = graph.size();
   for (unsigned int color = 0; color < n_colors; ++color)
     {
-      typename Portable::MatrixFree<dim, double>::Data gpu_data =
+      typename Portable::MatrixFree<dim, double>::PrecomputedData gpu_data =
         mf_data.get_data(color);
       const unsigned int n_cells = gpu_data.n_cells;
       auto gpu_data_host         = Portable::copy_mf_data_to_host<dim, double>(
index 5bd9c5943fec11e34fef2de811ec81c2ea5db3de..908b48ffa94887b7a07389ecc576035204c17b1d 100644 (file)
@@ -39,25 +39,23 @@ public:
     : data(data_in){};
 
   DEAL_II_HOST_DEVICE void
-  operator()(const unsigned int                                      cell,
-             const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-             Portable::SharedData<dim, Number> *shared_data,
-             const Number                      *src,
-             Number                            *dst) const
+  operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
+             const Number                                           *src,
+             Number                                                 *dst) const
   {
     Portable::FEEvaluation<dim, fe_degree, n_q_points_1d, 1, Number> fe_eval(
-      gpu_data, shared_data);
+      data);
 
     // set to unit vector
     auto fe_eval_ptr = &fe_eval;
-    Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
+    Kokkos::parallel_for(Kokkos::TeamThreadRange(data->team_member,
                                                  n_local_dofs),
                          [&](int i) { fe_eval_ptr->submit_dof_value(1., i); });
-    shared_data->team_member.team_barrier();
+    data->team_member.team_barrier();
     fe_eval.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
 
 #ifndef __APPLE__
-    Kokkos::parallel_for(Kokkos::TeamThreadRange(shared_data->team_member,
+    Kokkos::parallel_for(Kokkos::TeamThreadRange(data->team_member,
                                                  n_local_dofs),
                          [&](int i) {
                            // values should evaluate to one, derivatives to zero
@@ -69,7 +67,7 @@ public:
     fe_eval.integrate(EvaluationFlags::values | EvaluationFlags::gradients);
 
     Kokkos::parallel_for(
-      Kokkos::TeamThreadRange(shared_data->team_member, n_local_dofs),
+      Kokkos::TeamThreadRange(data->team_member, n_local_dofs),
       KOKKOS_LAMBDA(int i) { assert(fe_eval_ptr->get_dof_value(i) == 1.); });
 #endif
   }
index 4e5c8524ead71061ecfc8950c64487eaf83faf68..556e56f69b9a6688b64f82ab5b894bc68a0d0d1d 100644 (file)
@@ -29,13 +29,8 @@ class HelmholtzOperatorQuad
 {
 public:
   DEAL_II_HOST_DEVICE
-  HelmholtzOperatorQuad(
-    const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-    Number                                                 *coef,
-    int                                                     cell)
-    : gpu_data(gpu_data)
-    , coef(coef)
-    , cell(cell)
+  HelmholtzOperatorQuad(Number *coef)
+    : coef(coef)
   {}
 
   DEAL_II_HOST_DEVICE void
@@ -60,7 +55,8 @@ HelmholtzOperatorQuad<dim, fe_degree, Number, n_q_points_1d>::operator()(
   Portable::FEEvaluation<dim, fe_degree, n_q_points_1d, 1, Number> *fe_eval,
   int q_point) const
 {
-  unsigned int pos = gpu_data->local_q_point_id(cell, n_q_points, q_point);
+  unsigned int pos = fe_eval->get_matrix_free_data()->local_q_point_id(
+    fe_eval->get_current_cell_index(), n_q_points, q_point);
   fe_eval->submit_value(coef[pos] * fe_eval->get_value(q_point), q_point);
   fe_eval->submit_gradient(fe_eval->get_gradient(q_point), q_point);
 }
@@ -81,11 +77,9 @@ public:
   {}
 
   DEAL_II_HOST_DEVICE void
-  operator()(const unsigned int                                      cell,
-             const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-             Portable::SharedData<dim, Number> *shared_data,
-             const Number                      *src,
-             Number                            *dst) const;
+  operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
+             const Number                                           *src,
+             Number                                                 *dst) const;
 
   Number *coef;
 };
@@ -95,20 +89,16 @@ public:
 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 unsigned int                                      cell,
-  const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-  Portable::SharedData<dim, Number>                      *shared_data,
+  const typename Portable::MatrixFree<dim, Number>::Data *data,
   const Number                                           *src,
   Number                                                 *dst) const
 {
   Portable::FEEvaluation<dim, fe_degree, n_q_points_1d, 1, Number> fe_eval(
-    gpu_data, shared_data);
+    data);
   fe_eval.read_dof_values(src);
   fe_eval.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
   fe_eval.apply_for_each_quad_point(
-    HelmholtzOperatorQuad<dim, fe_degree, Number, n_q_points_1d>(gpu_data,
-                                                                 coef,
-                                                                 cell));
+    HelmholtzOperatorQuad<dim, fe_degree, Number, n_q_points_1d>(coef));
   fe_eval.integrate(EvaluationFlags::values | EvaluationFlags::gradients);
   fe_eval.distribute_local_to_global(dst);
 }
index 3dc2b7226d6d67cb9acf5a2d58a5ac9eec32cef5..a566a4aab20a71936710684eaf97d214c262144a 100644 (file)
@@ -166,17 +166,12 @@ public:
   {}
 
   DEAL_II_HOST_DEVICE void
-  operator()(const unsigned int                                      cell,
-             const typename Portable::MatrixFree<dim, Number>::Data *gpu_data,
-             Portable::SharedData<dim, Number> *shared_data,
-             const Number                      *src,
-             Number                            *dst) const
+  operator()(const typename Portable::MatrixFree<dim, Number>::Data *data,
+             const Number                                           *src,
+             Number                                                 *dst) const
   {
-    (void)cell; // TODO?
-
     Portable::FEEvaluation<dim, fe_degree, fe_degree + 1, n_components, Number>
-      phi(
-        /*cell,*/ gpu_data, shared_data);
+      phi(data);
     phi.read_dof_values(src);
 
     if (version == 0)

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.