]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove temporary Kokkos::View
authorBruno Turcksin <bruno.turcksin@gmail.com>
Thu, 18 May 2023 15:51:44 +0000 (15:51 +0000)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Mon, 22 May 2023 13:01:21 +0000 (13:01 +0000)
include/deal.II/matrix_free/cuda_matrix_free.templates.h

index 421e0abbcc5cedc91976bc0faa86265efed671e4..455dedf809c2905ce1d7a4e34ad16093b1952c6c 100644 (file)
@@ -712,24 +712,22 @@ namespace CUDAWrappers
     dofs_per_cell     = fe.n_dofs_per_cell();
     q_points_per_cell = Utilities::fixed_power<dim>(n_q_points_1d);
 
-    const ::dealii::internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info(
-      quad, fe);
+    ::dealii::internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info(quad,
+                                                                          fe);
 
     unsigned int size_shape_values = n_dofs_1d * n_q_points_1d;
 
-    FE_DGQArbitraryNodes<1> fe_quad_co(quad);
-    const ::dealii::internal::MatrixFreeFunctions::ShapeInfo<Number>
-      shape_info_co(quad, fe_quad_co);
+    FE_DGQArbitraryNodes<1>                                    fe_quad_co(quad);
+    ::dealii::internal::MatrixFreeFunctions::ShapeInfo<Number> shape_info_co(
+      quad, fe_quad_co);
 
     shape_values = Kokkos::View<Number *, MemorySpace::Default::kokkos_space>(
       Kokkos::view_alloc("shape_values", Kokkos::WithoutInitializing),
       size_shape_values);
-    auto shape_values_host = Kokkos::create_mirror_view(shape_values);
-    for (unsigned int i = 0; i < size_shape_values; ++i)
-      {
-        shape_values_host[i] = shape_info.data.front().shape_values[i];
-      }
-    Kokkos::deep_copy(shape_values, shape_values_host);
+    Kokkos::deep_copy(shape_values,
+                      Kokkos::View<Number *, Kokkos::HostSpace>(
+                        shape_info.data.front().shape_values.data(),
+                        size_shape_values));
 
     if (update_flags & update_gradients)
       {
@@ -737,13 +735,10 @@ namespace CUDAWrappers
           Kokkos::View<Number *, MemorySpace::Default::kokkos_space>(
             Kokkos::view_alloc("shape_gradients", Kokkos::WithoutInitializing),
             size_shape_values);
-        auto shape_gradients_host = Kokkos::create_mirror_view(shape_gradients);
-        for (unsigned int i = 0; i < size_shape_values; ++i)
-          {
-            shape_gradients_host[i] =
-              shape_info.data.front().shape_gradients[i];
-          }
-        Kokkos::deep_copy(shape_gradients, shape_gradients_host);
+        Kokkos::deep_copy(shape_gradients,
+                          Kokkos::View<Number *, Kokkos::HostSpace>(
+                            shape_info.data.front().shape_gradients.data(),
+                            size_shape_values));
 
 
         co_shape_gradients =
@@ -751,14 +746,10 @@ namespace CUDAWrappers
             Kokkos::view_alloc("co_shape_gradients",
                                Kokkos::WithoutInitializing),
             size_shape_values);
-        auto co_shape_gradients_host =
-          Kokkos::create_mirror_view(co_shape_gradients);
-        for (unsigned int i = 0; i < size_shape_values; ++i)
-          {
-            co_shape_gradients_host[i] =
-              shape_info_co.data.front().shape_gradients[i];
-          }
-        Kokkos::deep_copy(co_shape_gradients, co_shape_gradients_host);
+        Kokkos::deep_copy(co_shape_gradients,
+                          Kokkos::View<Number *, Kokkos::HostSpace>(
+                            shape_info_co.data.front().shape_gradients.data(),
+                            size_shape_values));
       }
 
     internal::ReinitHelper<dim, Number> helper(

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.