From: Bruno Turcksin Date: Fri, 13 Apr 2018 20:44:58 +0000 (-0400) Subject: Encapsulate all the cuda handles in a struct X-Git-Tag: v9.0.0-rc1~149^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=28ca9e96e5d34d8d9f7997baf96b156742b37ae1;p=dealii.git Encapsulate all the cuda handles in a struct --- diff --git a/include/deal.II/base/cuda.h b/include/deal.II/base/cuda.h new file mode 100644 index 0000000000..c2594dbdba --- /dev/null +++ b/include/deal.II/base/cuda.h @@ -0,0 +1,67 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#ifndef dealii_cuda_h +#define dealii_cuda_h + +#include + +#ifdef DEAL_II_WITH_CUDA + +#include + +DEAL_II_NAMESPACE_OPEN +namespace Utilities +{ + /** + * A namespace for utility structures for CUDA. + */ + namespace CUDA + { + /** + * This structure creates, stores, and destroys the handles of the different + * CUDA libraries used inside deal.II. + */ + struct Handle + { + /** + * Constructor. Create the handles for the different libraries. + */ + Handle(); + + /** + * Copy constructor is deleted. + */ + Handle(Handle const &) = delete; + + /** + * Destructor. Destroy the handles and free all the memory allocated by + * GrowingVectorMemory. + */ + ~Handle(); + + /** + * Handle to the cuSPARSE library. + */ + cusparseHandle_t cusparse_handle; + }; + } +} + +DEAL_II_NAMESPACE_CLOSE + +#endif + +#endif diff --git a/include/deal.II/lac/cuda_sparse_matrix.h b/include/deal.II/lac/cuda_sparse_matrix.h index 721363c82f..2f1a6a99fd 100644 --- a/include/deal.II/lac/cuda_sparse_matrix.h +++ b/include/deal.II/lac/cuda_sparse_matrix.h @@ -20,6 +20,7 @@ #include #ifdef DEAL_II_WITH_CUDA +#include #include #include #include @@ -74,11 +75,11 @@ namespace CUDAWrappers SparseMatrix(); /** - * Constructor. Takes a cuSPARSE handle and a sparse matrix on the host. + * Constructor. Takes a Utilities::CUDA::Handle and a sparse matrix on the host. * The sparse matrix on the host is copied on the device and the elements * are reordered according to the format supported by cuSPARSE. */ - SparseMatrix(cusparseHandle_t handle, + SparseMatrix(Utilities::CUDA::Handle &handle, const ::dealii::SparseMatrix &sparse_matrix_host); /** @@ -102,7 +103,7 @@ namespace CUDAWrappers * to the device and the elementes are reordered according to the format * supported by cuSPARSE. */ - void reinit(cusparseHandle_t handle, + void reinit(Utilities::CUDA::Handle &handle, const ::dealii::SparseMatrix &sparse_matrix_host); //@} @@ -251,8 +252,8 @@ namespace CUDAWrappers private: /** - * cuSPARSE used to call cuSPARSE function. The cuSPARSE handle needs to - * be mutable to be called in a const function. + * cuSPARSE handle used to call cuSPARSE functions. The cuSPARSE handle needs + * to be mutable to be called in a const function. */ mutable cusparseHandle_t cusparse_handle; diff --git a/source/base/CMakeLists.txt b/source/base/CMakeLists.txt index 5e80792a1a..a95b341c0d 100644 --- a/source/base/CMakeLists.txt +++ b/source/base/CMakeLists.txt @@ -88,6 +88,14 @@ SET(_separate_src symmetric_tensor.cc ) +# Add CUDA wrapper files +IF(DEAL_II_WITH_CUDA) + SET(_separate_src + ${_separate_src} + cuda.cu + ) +ENDIF() + # determined by profiling SET(_n_includes_per_unity_file 29) diff --git a/source/base/cuda.cu b/source/base/cuda.cu new file mode 100644 index 0000000000..d89c68e141 --- /dev/null +++ b/source/base/cuda.cu @@ -0,0 +1,49 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include + +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +namespace Utilities +{ + namespace CUDA + { + Handle::Handle() + { + cusparseStatus_t cusparse_error_code = cusparseCreate(&cusparse_handle); + AssertCusparse(cusparse_error_code); + } + + + + Handle::~Handle() + { + dealii::GrowingVectorMemory> + ::release_unused_memory(); + dealii::GrowingVectorMemory> + ::release_unused_memory(); + + cusparseStatus_t cusparse_error_code = cusparseDestroy(cusparse_handle); + AssertCusparse(cusparse_error_code); + } + } +} + +DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/cuda_sparse_matrix.cu b/source/lac/cuda_sparse_matrix.cu index 467d5a51f2..453330f76e 100644 --- a/source/lac/cuda_sparse_matrix.cu +++ b/source/lac/cuda_sparse_matrix.cu @@ -134,7 +134,7 @@ namespace CUDAWrappers template - SparseMatrix::SparseMatrix(cusparseHandle_t handle, + SparseMatrix::SparseMatrix(Utilities::CUDA::Handle &handle, const ::dealii::SparseMatrix &sparse_matrix_host) : val_dev(nullptr), @@ -208,10 +208,10 @@ namespace CUDAWrappers template - void SparseMatrix::reinit(cusparseHandle_t handle, + void SparseMatrix::reinit(Utilities::CUDA::Handle &handle, const ::dealii::SparseMatrix &sparse_matrix_host) { - cusparse_handle = handle; + cusparse_handle = handle.cusparse_handle; nnz = sparse_matrix_host.n_nonzero_elements(); n_rows = sparse_matrix_host.m(); n_cols = sparse_matrix_host.n(); diff --git a/tests/cuda/solver_01.cu b/tests/cuda/solver_01.cu index 48940f0a62..647aa438ec 100644 --- a/tests/cuda/solver_01.cu +++ b/tests/cuda/solver_01.cu @@ -18,6 +18,7 @@ #include "../tests.h" #include "../testmatrix.h" +#include #include #include #include @@ -26,7 +27,7 @@ #include #include -void test(cusparseHandle_t cusparse_handle) +void test(Utilities::CUDA::Handle &cuda_handle) { // Build the sparse matrix on the host const unsigned int problem_size = 10; @@ -50,7 +51,7 @@ void test(cusparseHandle_t cusparse_handle) cg_host.solve(A, sol_host, rhs_host, prec_no); // Solve on the device - CUDAWrappers::SparseMatrix A_dev(cusparse_handle, A); + CUDAWrappers::SparseMatrix A_dev(cuda_handle, A); LinearAlgebra::CUDAWrappers::Vector sol_dev(size); LinearAlgebra::CUDAWrappers::Vector rhs_dev(size); LinearAlgebra::ReadWriteVector rw_vector(size); @@ -71,17 +72,11 @@ int main() initlog(); deallog.depth_console(0); - cusparseHandle_t cusparse_handle; - cusparseStatus_t cusparse_error_code = cusparseCreate(&cusparse_handle); - AssertCusparse(cusparse_error_code); - - test(cusparse_handle); + Utilities::CUDA::Handle cuda_handle; + test(cuda_handle); GrowingVectorMemory>::release_unused_memory(); - cusparse_error_code = cusparseDestroy(cusparse_handle); - AssertCusparse(cusparse_error_code); - deallog << "OK" < #include #include #include @@ -65,7 +66,7 @@ void check_vector(Vector const &a, AssertThrow(std::abs(a[i] - b[i]) < 1e-15, ExcInternalError()); } -void test(cusparseHandle_t cusparse_handle) +void test(Utilities::CUDA::Handle &cuda_handle) { // Build the sparse matrix on the host const unsigned int size = 10; @@ -79,7 +80,7 @@ void test(cusparseHandle_t cusparse_handle) testproblem.five_point(A, true); // Create the sparse matrix on the device - CUDAWrappers::SparseMatrix A_dev(cusparse_handle, A); + CUDAWrappers::SparseMatrix A_dev(cuda_handle, A); check_matrix(A, A_dev); AssertThrow(A.m() == A_dev.m(), ExcInternalError()); @@ -190,7 +191,7 @@ void test(cusparseHandle_t cusparse_handle) if (i B_dev(cusparse_handle, B); + CUDAWrappers::SparseMatrix B_dev(cuda_handle, B); value = B.l1_norm(); value_host = B_dev.l1_norm(); AssertThrow(std::abs(value-value_host) < 1e-15, ExcInternalError()); @@ -201,14 +202,9 @@ int main() initlog(); deallog.depth_console(0); - cusparseHandle_t cusparse_handle; - cusparseStatus_t cusparse_error_code = cusparseCreate(&cusparse_handle); - AssertCusparse(cusparse_error_code); + Utilities::CUDA::Handle cuda_handle; - test(cusparse_handle); - - cusparse_error_code = cusparseDestroy(cusparse_handle); - AssertCusparse(cusparse_error_code); + test(cuda_handle); deallog << "OK" <