From 87c3d52b15513e3a8d039b84e1e1fe001e0c6c49 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Fri, 28 Aug 2020 15:50:29 +0000 Subject: [PATCH] Add functions to perform CUDA matrix-free loop on the host --- .../deal.II/matrix_free/cuda_matrix_free.h | 219 ++++++++++++++++-- .../matrix_free/cuda_matrix_free.templates.h | 6 +- 2 files changed, 198 insertions(+), 27 deletions(-) diff --git a/include/deal.II/matrix_free/cuda_matrix_free.h b/include/deal.II/matrix_free/cuda_matrix_free.h index 71e0aad607..ca8ec5895f 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.h @@ -19,22 +19,22 @@ #include -#ifdef DEAL_II_COMPILER_CUDA_AWARE +#include +#include +#include +#include -# include -# include -# include -# include +#include -# include +#include +#include +#include -# include -# include -# include +#include -# include -# include -# include +#include +#include +#include DEAL_II_NAMESPACE_OPEN @@ -42,13 +42,13 @@ DEAL_II_NAMESPACE_OPEN namespace CUDAWrappers { // forward declaration -# ifndef DOXYGEN +#ifndef DOXYGEN namespace internal { template class ReinitHelper; } -# endif +#endif /** * This class collects all the data that is stored for the matrix free @@ -82,6 +82,8 @@ namespace CUDAWrappers public: using jacobian_type = Tensor<2, dim, Tensor<1, dim, Number>>; using point_type = Point; + using CellFilter = + FilteredIterator::active_cell_iterator>; /** * Parallelization scheme used: parallel_in_elem (parallelism at the level @@ -114,12 +116,12 @@ namespace CUDAWrappers , use_coloring(use_coloring) , overlap_communication_computation(overlap_communication_computation) { -# ifndef DEAL_II_MPI_WITH_CUDA_SUPPORT +#ifndef DEAL_II_MPI_WITH_CUDA_SUPPORT AssertThrow( overlap_communication_computation == false, ExcMessage( "Overlapping communication and computation requires CUDA-aware MPI.")); -# endif +#endif if (overlap_communication_computation == true) AssertThrow( use_coloring == false || overlap_communication_computation == false, @@ -339,6 +341,12 @@ namespace CUDAWrappers initialize_dof_vector( LinearAlgebra::distributed::Vector &vec) const; + /** + * Return the colored graph of locally owned active cells. + */ + const std::vector> & + get_colored_graph() const; + /** * Return the partitioner that represents the locally owned data and the * ghost indices where access is needed to for the cell loop. The @@ -616,6 +624,11 @@ namespace CUDAWrappers */ const DoFHandler *dof_handler; + /** + * Colored graphed of locally owned active cells. + */ + std::vector> graph; + friend class internal::ReinitHelper; }; @@ -655,6 +668,8 @@ namespace CUDAWrappers +#ifdef DEAL_II_COMPILER_CUDA_AWARE + // This function determines the number of cells per block, possibly at compile // time (by virtue of being 'constexpr') // TODO this function should be rewritten using meta-programming @@ -734,20 +749,182 @@ namespace CUDAWrappers q_point_id_in_cell(n_q_points_1d)); } + /** + * Structure which is passed to the kernel. It is used to pass all the + * necessary information from the CPU to the GPU. + */ + template + struct DataHost + { + /** + * Vector of quadrature points. + */ + std::vector> q_points; + + /** + * Map the position in the local vector to the position in the global + * vector. + */ + std::vector local_to_global; + + /** + * Vector of inverse Jacobians. + */ + std::vector inv_jacobian; + + /** + * Vector of Jacobian times the weights. + */ + std::vector JxW; + + /** + * ID of the associated MatrixFree object. + */ + unsigned int id; + + /** + * Number of cells. + */ + unsigned int n_cells; + + /** + * Length of the padding. + */ + unsigned int padding_length; + + /** + * Row start (including padding). + */ + unsigned int row_start; + + /** + * Mask deciding where constraints are set on a given cell. + */ + std::vector constraint_mask; + + /** + * If true, use graph coloring has been used and we can simply add into + * the destingation vector. Otherwise, use atomic operations. + */ + bool use_coloring; + }; + + + + /** + * Copy @p data from the device to the device. @p update_flags should be + * identical to the one used in MatrixFree::AdditionalData. + * + * @relates CUDAWrappers::MatrixFree + */ + template + DataHost + copy_mf_data_to_host( + const typename dealii::CUDAWrappers::MatrixFree::Data &data, + const UpdateFlags &update_flags) + { + DataHost data_host; + + data_host.id = data.id; + data_host.n_cells = data.n_cells; + data_host.padding_length = data.padding_length; + data_host.row_start = data.row_start; + data_host.use_coloring = data.use_coloring; + + const unsigned int n_elements = + data_host.n_cells * data_host.padding_length; + if (update_flags & update_quadrature_points) + { + data_host.q_points.resize(n_elements); + Utilities::CUDA::copy_to_host(data.q_points, data_host.q_points); + } + + data_host.local_to_global.resize(n_elements); + Utilities::CUDA::copy_to_host(data.local_to_global, + data_host.local_to_global); + + if (update_flags & update_gradients) + { + data_host.inv_jacobian.resize(n_elements * dim * dim); + Utilities::CUDA::copy_to_host(data.inv_jacobian, + data_host.inv_jacobian); + } + + if (update_flags & update_JxW_values) + { + data_host.JxW.resize(n_elements); + Utilities::CUDA::copy_to_host(data.JxW, data_host.JxW); + } + + data_host.constraint_mask.resize(data_host.n_cells); + Utilities::CUDA::copy_to_host(data.constraint_mask, + data_host.constraint_mask); + + return data_host; + } + + + + /** + * This function is the host version of local_q_point_id(). + * + * @relates CUDAWrappers::MatrixFree + */ + template + inline unsigned int + local_q_point_id_host(const unsigned int cell, + const DataHost &data, + const unsigned int n_q_points, + const unsigned int i) + { + return (data.row_start / data.padding_length + cell) * n_q_points + i; + } + + + + /** + * This function is the host version of get_quadrature_point(). It assumes + * that the data in MatrixFree::Data has been copied to the host + * using copy_mf_data_to_host(). + * + * @relates CUDAWrappers::MatrixFree + */ + template + inline Point + get_quadrature_point_host(const unsigned int cell, + const DataHost &data, + const unsigned int i) + { + return data.q_points[data.padding_length * cell + i]; + } +#endif + /*----------------------- Inline functions ---------------------------------*/ -# ifndef DOXYGEN +#ifndef DOXYGEN template - const std::shared_ptr & + inline const std::vector::active_cell_iterator>>> & + MatrixFree::get_colored_graph() const + { + return graph; + } + + + + template + inline const std::shared_ptr & MatrixFree::get_vector_partitioner() const { return partitioner; } + + template - const DoFHandler & + inline const DoFHandler & MatrixFree::get_dof_handler() const { Assert(dof_handler != nullptr, ExcNotInitialized()); @@ -755,12 +932,10 @@ namespace CUDAWrappers return *dof_handler; } -# endif +#endif } // namespace CUDAWrappers DEAL_II_NAMESPACE_CLOSE #endif - -#endif 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 693016f66a..af8904b1d1 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.templates.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.templates.h @@ -31,8 +31,6 @@ # include # include -# include - # include # include @@ -952,12 +950,9 @@ namespace CUDAWrappers this, mapping, fe, quad, shape_info, *dof_handler, update_flags); // Create a graph coloring - using CellFilter = - FilteredIterator::active_cell_iterator>; CellFilter begin(IteratorFilters::LocallyOwnedCell(), dof_handler->begin_active()); CellFilter end(IteratorFilters::LocallyOwnedCell(), dof_handler->end()); - std::vector> graph; if (begin != end) { @@ -971,6 +966,7 @@ namespace CUDAWrappers } else { + graph.clear(); if (additional_data.overlap_communication_computation) { // We create one color (1) with the cells on the boundary of the -- 2.39.5