From f02759eec15ba0ff49020912dba65cc04c166c89 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Mon, 3 Aug 2020 14:52:19 +0000 Subject: [PATCH] Fix compilation with CUDA 11 --- include/deal.II/lac/cuda_sparse_matrix.h | 13 +- source/lac/cuda_sparse_matrix.cu | 288 ++++++++++++++++------- tests/cuda/sparse_matrix_01.cu | 2 +- 3 files changed, 218 insertions(+), 85 deletions(-) diff --git a/include/deal.II/lac/cuda_sparse_matrix.h b/include/deal.II/lac/cuda_sparse_matrix.h index 96222a735e..45a6b6f6fe 100644 --- a/include/deal.II/lac/cuda_sparse_matrix.h +++ b/include/deal.II/lac/cuda_sparse_matrix.h @@ -320,10 +320,10 @@ namespace CUDAWrappers //@{ /** * Return a tuple containing the pointer to the values of matrix, the - * pointer to the columns indices, the pointer to the rows pointer, and - * the cuSPARSE matrix description. + * pointer to the columns indices, the pointer to the rows pointer, + * the cuSPARSE matrix description, and the cuSPARSE SP matrix description. */ - std::tuple + std::tuple get_cusparse_matrix() const; //@} @@ -364,9 +364,14 @@ namespace CUDAWrappers std::unique_ptr row_ptr_dev; /** - * cuSPARSE description of the sparse matrix. + * cuSPARSE description of the matrix. */ cusparseMatDescr_t descr; + + /** + * cuSPARSE description of the sparse matrix. + */ + cusparseSpMatDescr_t sp_descr; }; diff --git a/source/lac/cuda_sparse_matrix.cu b/source/lac/cuda_sparse_matrix.cu index 27a08e59f4..9917ae8fb0 100644 --- a/source/lac/cuda_sparse_matrix.cu +++ b/source/lac/cuda_sparse_matrix.cu @@ -44,18 +44,66 @@ namespace CUDAWrappers void - csrmv(cusparseHandle_t handle, - bool transpose, - int m, - int n, - int nnz, - const cusparseMatDescr_t descr, - const float * A_val_dev, - const int * A_row_ptr_dev, - const int * A_column_index_dev, - const float * x, - bool add, - float * y) + create_sp_mat_descr(int m, + int n, + int nnz, + const float * A_val_dev, + const int * A_row_ptr_dev, + const int * A_column_index_dev, + cusparseSpMatDescr_t &sp_descr) + { + cusparseStatus_t error_code = cusparseCreateCsr( + &sp_descr, + m, + n, + nnz, + reinterpret_cast(const_cast(A_row_ptr_dev)), + reinterpret_cast(const_cast(A_column_index_dev)), + reinterpret_cast(const_cast(A_val_dev)), + CUSPARSE_INDEX_32I, + CUSPARSE_INDEX_32I, + CUSPARSE_INDEX_BASE_ZERO, + CUDA_R_32F); + AssertCusparse(error_code); + } + + + + void + create_sp_mat_descr(int m, + int n, + int nnz, + const double * A_val_dev, + const int * A_row_ptr_dev, + const int * A_column_index_dev, + cusparseSpMatDescr_t &sp_descr) + { + cusparseStatus_t error_code = cusparseCreateCsr( + &sp_descr, + m, + n, + nnz, + reinterpret_cast(const_cast(A_row_ptr_dev)), + reinterpret_cast(const_cast(A_column_index_dev)), + reinterpret_cast(const_cast(A_val_dev)), + CUSPARSE_INDEX_32I, + CUSPARSE_INDEX_32I, + CUSPARSE_INDEX_BASE_ZERO, + CUDA_R_64F); + AssertCusparse(error_code); + } + + + + void + csrmv(cusparseHandle_t handle, + bool transpose, + int m, + int n, + const cusparseSpMatDescr_t sp_descr, + const float * x, + bool add, + float * y) { float alpha = 1.; float beta = add ? 1. : 0.; @@ -63,38 +111,72 @@ namespace CUDAWrappers transpose ? CUSPARSE_OPERATION_TRANSPOSE : CUSPARSE_OPERATION_NON_TRANSPOSE; + // Move the data to cuSPARSE data type + cusparseDnVecDescr_t x_cuvec; + cusparseStatus_t error_code = + cusparseCreateDnVec(&x_cuvec, + n, + reinterpret_cast(const_cast(x)), + CUDA_R_32F); + AssertCusparse(error_code); + + cusparseDnVecDescr_t y_cuvec; + error_code = + cusparseCreateDnVec(&y_cuvec, + m, + reinterpret_cast(const_cast(y)), + CUDA_R_32F); + AssertCusparse(error_code); + // This function performs y = alpha*op(A)*x + beta*y - cusparseStatus_t error_code = cusparseScsrmv(handle, - cusparse_operation, - m, - n, - nnz, - &alpha, - descr, - A_val_dev, - A_row_ptr_dev, - A_column_index_dev, - x, - &beta, - y); + size_t buffer_size = 0; + error_code = cusparseSpMV_bufferSize(handle, + cusparse_operation, + &alpha, + sp_descr, + x_cuvec, + &beta, + y_cuvec, + CUDA_R_32F, + CUSPARSE_MV_ALG_DEFAULT, + &buffer_size); + + void * buffer = nullptr; + cudaError_t cuda_error_code = cudaMalloc(&buffer, buffer_size); + AssertCuda(cuda_error_code); + + // execute SpMV + error_code = cusparseSpMV(handle, + cusparse_operation, + &alpha, + sp_descr, + x_cuvec, + &beta, + y_cuvec, + CUDA_R_32F, + CUSPARSE_MV_ALG_DEFAULT, + buffer); + AssertCusparse(error_code); + + cuda_error_code = cudaFree(buffer); + AssertCuda(cuda_error_code); + error_code = cusparseDestroyDnVec(x_cuvec); + AssertCusparse(error_code); + error_code = cusparseDestroyDnVec(y_cuvec); AssertCusparse(error_code); } void - csrmv(cusparseHandle_t handle, - bool transpose, - int m, - int n, - int nnz, - const cusparseMatDescr_t descr, - const double * A_val_dev, - const int * A_row_ptr_dev, - const int * A_column_index_dev, - const double * x, - bool add, - double * y) + csrmv(cusparseHandle_t handle, + bool transpose, + int m, + int n, + const cusparseSpMatDescr_t sp_descr, + const double * x, + bool add, + double * y) { double alpha = 1.; double beta = add ? 1. : 0.; @@ -102,20 +184,58 @@ namespace CUDAWrappers transpose ? CUSPARSE_OPERATION_TRANSPOSE : CUSPARSE_OPERATION_NON_TRANSPOSE; + // Move the data to cuSPARSE data type + cusparseDnVecDescr_t x_cuvec; + cusparseStatus_t error_code = + cusparseCreateDnVec(&x_cuvec, + n, + reinterpret_cast(const_cast(x)), + CUDA_R_64F); + AssertCusparse(error_code); + + cusparseDnVecDescr_t y_cuvec; + error_code = + cusparseCreateDnVec(&y_cuvec, + m, + reinterpret_cast(const_cast(y)), + CUDA_R_64F); + AssertCusparse(error_code); + // This function performs y = alpha*op(A)*x + beta*y - cusparseStatus_t error_code = cusparseDcsrmv(handle, - cusparse_operation, - m, - n, - nnz, - &alpha, - descr, - A_val_dev, - A_row_ptr_dev, - A_column_index_dev, - x, - &beta, - y); + size_t buffer_size = 0; + error_code = cusparseSpMV_bufferSize(handle, + cusparse_operation, + &alpha, + sp_descr, + x_cuvec, + &beta, + y_cuvec, + CUDA_R_64F, + CUSPARSE_MV_ALG_DEFAULT, + &buffer_size); + + void * buffer = nullptr; + cudaError_t cuda_error_code = cudaMalloc(&buffer, buffer_size); + AssertCuda(cuda_error_code); + + // execute SpMV + error_code = cusparseSpMV(handle, + cusparse_operation, + &alpha, + sp_descr, + x_cuvec, + &beta, + y_cuvec, + CUDA_R_64F, + CUSPARSE_MV_ALG_DEFAULT, + buffer); + AssertCusparse(error_code); + + cuda_error_code = cudaFree(buffer); + AssertCuda(cuda_error_code); + error_code = cusparseDestroyDnVec(x_cuvec); + AssertCusparse(error_code); + error_code = cusparseDestroyDnVec(y_cuvec); AssertCusparse(error_code); } @@ -171,6 +291,7 @@ namespace CUDAWrappers , column_index_dev(nullptr, Utilities::CUDA::delete_device_data) , row_ptr_dev(nullptr, Utilities::CUDA::delete_device_data) , descr(nullptr) + , sp_descr(nullptr) {} @@ -183,6 +304,7 @@ namespace CUDAWrappers , column_index_dev(nullptr, Utilities::CUDA::delete_device_data) , row_ptr_dev(nullptr, Utilities::CUDA::delete_device_data) , descr(nullptr) + , sp_descr(nullptr) { reinit(handle, sparse_matrix_host); } @@ -199,11 +321,13 @@ namespace CUDAWrappers , column_index_dev(std::move(other.column_index_dev)) , row_ptr_dev(std::move(other.row_ptr_dev)) , descr(other.descr) + , sp_descr(other.sp_descr) { - other.nnz = 0; - other.n_rows = 0; - other.n_cols = 0; - other.descr = nullptr; + other.nnz = 0; + other.n_rows = 0; + other.n_cols = 0; + other.descr = nullptr; + other.sp_descr = nullptr; } @@ -219,6 +343,14 @@ namespace CUDAWrappers descr = nullptr; } + if (sp_descr != nullptr) + { + const cusparseStatus_t cusparse_error_code = + cusparseDestroySpMat(sp_descr); + AssertNothrowCusparse(cusparse_error_code); + sp_descr = nullptr; + } + nnz = 0; n_rows = 0; } @@ -237,11 +369,13 @@ namespace CUDAWrappers column_index_dev = std::move(other.column_index_dev); row_ptr_dev = std::move(other.row_ptr_dev); descr = other.descr; + sp_descr = other.sp_descr; - other.nnz = 0; - other.n_rows = 0; - other.n_cols = 0; - other.descr = nullptr; + other.nnz = 0; + other.n_rows = 0; + other.n_cols = 0; + other.descr = nullptr; + other.sp_descr = nullptr; return *this; } @@ -329,6 +463,15 @@ namespace CUDAWrappers cusparse_error_code = cusparseSetMatIndexBase(descr, CUSPARSE_INDEX_BASE_ZERO); AssertCusparse(cusparse_error_code); + + // Create the sparse matrix descriptor + internal::create_sp_mat_descr(n_rows, + n_cols, + nnz, + val_dev.get(), + row_ptr_dev.get(), + column_index_dev.get(), + sp_descr); } @@ -374,11 +517,7 @@ namespace CUDAWrappers false, n_rows, n_cols, - nnz, - descr, - val_dev.get(), - row_ptr_dev.get(), - column_index_dev.get(), + sp_descr, src.get_values(), false, dst.get_values()); @@ -396,11 +535,7 @@ namespace CUDAWrappers true, n_rows, n_cols, - nnz, - descr, - val_dev.get(), - row_ptr_dev.get(), - column_index_dev.get(), + sp_descr, src.get_values(), false, dst.get_values()); @@ -418,11 +553,7 @@ namespace CUDAWrappers false, n_rows, n_cols, - nnz, - descr, - val_dev.get(), - row_ptr_dev.get(), - column_index_dev.get(), + sp_descr, src.get_values(), true, dst.get_values()); @@ -440,11 +571,7 @@ namespace CUDAWrappers true, n_rows, n_cols, - nnz, - descr, - val_dev.get(), - row_ptr_dev.get(), - column_index_dev.get(), + sp_descr, src.get_values(), true, dst.get_values()); @@ -548,13 +675,14 @@ namespace CUDAWrappers template - std::tuple + std::tuple SparseMatrix::get_cusparse_matrix() const { return std::make_tuple(val_dev.get(), column_index_dev.get(), row_ptr_dev.get(), - descr); + descr, + sp_descr); } diff --git a/tests/cuda/sparse_matrix_01.cu b/tests/cuda/sparse_matrix_01.cu index 339b688b4f..2327f9f539 100644 --- a/tests/cuda/sparse_matrix_01.cu +++ b/tests/cuda/sparse_matrix_01.cu @@ -35,7 +35,7 @@ check_matrix(SparseMatrix const & A, double * val_dev = nullptr; int * column_index_dev = nullptr; int * row_ptr_dev = nullptr; - std::tie(val_dev, column_index_dev, row_ptr_dev, std::ignore) = + std::tie(val_dev, column_index_dev, row_ptr_dev, std::ignore, std::ignore) = A_dev.get_cusparse_matrix(); int nnz = A_dev.n_nonzero_elements(); -- 2.39.5