From e6e8feac9e8fe40bf3e75f94b9c26e1f9821e072 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Tue, 21 Mar 2023 09:49:59 -0400 Subject: [PATCH] Move local_to_global to Kokkos::View --- .../deal.II/matrix_free/cuda_fe_evaluation.h | 25 ++++++++++---- .../deal.II/matrix_free/cuda_matrix_free.h | 17 ++++++---- .../matrix_free/cuda_matrix_free.templates.h | 33 ++++++++++--------- 3 files changed, 47 insertions(+), 28 deletions(-) diff --git a/include/deal.II/matrix_free/cuda_fe_evaluation.h b/include/deal.II/matrix_free/cuda_fe_evaluation.h index db1b3e0638..aa1bd4e1fc 100644 --- a/include/deal.II/matrix_free/cuda_fe_evaluation.h +++ b/include/deal.II/matrix_free/cuda_fe_evaluation.h @@ -240,10 +240,19 @@ namespace CUDAWrappers apply_for_each_quad_point(const Functor &func); private: - types::global_dof_index *local_to_global; - unsigned int n_cells; - unsigned int padding_length; - const unsigned int mf_object_id; + // FIXME We would like to use + // Kokkos::Subview, int, decltype(Kokkos::ALL)> + // but we get error: incomplete type is not allowed. I cannot reproduce + // outside of deal.II. Need to investigate more. + Kokkos::Subview, + int, + Kokkos::pair> + local_to_global; + unsigned int n_cells; + unsigned int padding_length; + const unsigned int mf_object_id; const dealii::internal::MatrixFreeFunctions::ConstraintKinds constraint_mask; @@ -290,7 +299,11 @@ namespace CUDAWrappers FEEvaluation(const unsigned int cell_id, const data_type * data, SharedData *shdata) - : n_cells(data->n_cells) + : local_to_global(Kokkos::subview( + data->local_to_global, + cell_id, + Kokkos::pair(0, Utilities::pow(n_q_points_1d, dim)))) + , n_cells(data->n_cells) , padding_length(data->padding_length) , mf_object_id(data->id) , constraint_mask(data->constraint_mask[cell_id]) @@ -307,8 +320,6 @@ namespace CUDAWrappers Kokkos::pair(0, Utilities::pow(n_q_points_1d, dim)))) , values(shdata->values) { - local_to_global = data->local_to_global + padding_length * cell_id; - for (unsigned int i = 0; i < dim; ++i) gradients[i] = shdata->gradients[i]; } diff --git a/include/deal.II/matrix_free/cuda_matrix_free.h b/include/deal.II/matrix_free/cuda_matrix_free.h index d6588813f6..108274852b 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.h @@ -183,7 +183,9 @@ namespace CUDAWrappers * Map the position in the local vector to the position in the global * vector. */ - types::global_dof_index *local_to_global; + Kokkos::View + local_to_global; /** * Kokkos::View of the inverse Jacobian. @@ -514,7 +516,9 @@ namespace CUDAWrappers * Map the position in the local vector to the position in the global * vector. */ - std::vector local_to_global; + std::vector> + local_to_global; /** * Vector of Kokkos::View of the inverse Jacobian associated to the cells of @@ -737,7 +741,9 @@ namespace CUDAWrappers * Map the position in the local vector to the position in the global * vector. */ - std::vector local_to_global; + typename Kokkos::View::HostMirror + local_to_global; /** * Kokkos::View of inverse Jacobians on the host. @@ -815,9 +821,8 @@ namespace CUDAWrappers Kokkos::deep_copy(data_host.q_points, data.q_points); } - data_host.local_to_global.resize(n_elements); - Utilities::CUDA::copy_to_host(data.local_to_global, - data_host.local_to_global); + data_host.local_to_global = Kokkos::create_mirror(data.local_to_global); + Kokkos::deep_copy(data_host.local_to_global, data.local_to_global); if (update_flags & update_gradients) { diff --git a/include/deal.II/matrix_free/cuda_matrix_free.templates.h b/include/deal.II/matrix_free/cuda_matrix_free.templates.h index d23c1b82bc..7ebf79bd91 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.templates.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.templates.h @@ -232,8 +232,9 @@ namespace CUDAWrappers private: MatrixFree *data; - // Host data - std::vector local_to_global_host; + Kokkos::View + local_to_global; Kokkos::View **, MemorySpace::Default::kokkos_space> q_points; Kokkos::View JxW; @@ -348,7 +349,13 @@ namespace CUDAWrappers data->block_dim[color] = dim3(n_dofs_1d * cells_per_block, n_dofs_1d, n_dofs_1d); - local_to_global_host.resize(n_cells * padding_length); + + local_to_global = Kokkos::View( + Kokkos::view_alloc("local_to_global_" + std::to_string(color), + Kokkos::WithoutInitializing), + n_cells, + dofs_per_cell); if (update_flags & update_quadrature_points) q_points = Kokkos::View **, @@ -406,9 +413,13 @@ namespace CUDAWrappers lexicographic_dof_indices, cell_id_view); - memcpy(&local_to_global_host[cell_id * padding_length], - lexicographic_dof_indices.data(), - dofs_per_cell * sizeof(types::global_dof_index)); + // FIXME too many deep_copy + auto local_to_global_host = + Kokkos::create_mirror_view_and_copy(MemorySpace::Host::kokkos_space{}, + local_to_global); + for (unsigned int i = 0; i < dofs_per_cell; ++i) + local_to_global_host(cell_id, i) = lexicographic_dof_indices[i]; + Kokkos::deep_copy(local_to_global, local_to_global_host); fe_values.reinit(cell); @@ -460,11 +471,7 @@ namespace CUDAWrappers const unsigned int n_cells = data->n_cells[color]; // Local-to-global mapping - alloc_and_copy( - &data->local_to_global[color], - ArrayView(local_to_global_host.data(), - local_to_global_host.size()), - n_cells * padding_length); + data->local_to_global[color] = local_to_global; // Quadrature points if (update_flags & update_quadrature_points) @@ -712,10 +719,6 @@ namespace CUDAWrappers void MatrixFree::free() { - for (auto &local_to_global_color_ptr : local_to_global) - Utilities::CUDA::free(local_to_global_color_ptr); - local_to_global.clear(); - for (auto &constraint_mask_color_ptr : constraint_mask) Utilities::CUDA::free(constraint_mask_color_ptr); constraint_mask.clear(); -- 2.39.5