I am trying to solve a problem that requires a sparse matrix sparse matrix product in CUDA Fortran code.
I am trying to use the cusparse library, cusparseSpGEMM, by referring to the sample code on github (https://github.com/NVIDIA/CUDALibrarySamples/blob/master/cuSPARSE/spgemm/spgemm_example.c), but a problem has arisen.
The first problem is that the first time I do cusparseSpGEMM_workEstimation, the status becomes 7 (CUSPARSE_STATUS_INTERNAL_ERROR).
Here is my code.
The computational environment is A100 80GB with CUDA 11.0.
I would appreciate it if you could point out any problems.
Thanks.
=========================================================================
program SpGEMM
use cudafor
use cusparse
Implicit none
!!Define Matrix----------------------------------
Integer,parameter :: A_rows=4
Integer,parameter :: A_cols=4
Integer,parameter :: A_nnz=9
Integer :: Arow(A_rows+1)
Integer :: Acol(A_nnz)
Real(8) :: Aval(A_nnz)
Integer,device :: Arow_d(A_rows+1)
Integer,device :: Acol_d(A_nnz)
Real(8),device :: Aval_d(A_nnz)
Integer,parameter :: B_rows=4
Integer,parameter :: B_cols=4
Integer,parameter :: B_nnz=8
Integer :: Brow(B_rows+1)
Integer :: Bcol(B_nnz)
Real(8) :: Bval(B_nnz)
Integer,device :: Brow_d(B_rows+1)
Integer,device :: Bcol_d(B_nnz)
Real(8),device :: Bval_d(B_nnz)
Integer,allocatable :: Crow(:)
Integer,allocatable :: Ccol(:)
Integer,allocatable,device :: Crow_d(:)
Integer,allocatable,device :: Ccol_d(:)
!!Define Matrix----------------------------------
Real(8) :: alpha=1d0,beta=0d0
Integer :: status
type(cusparseHandle) :: handle
type(cusparseSpMatDescr) :: matA,matB,matC
type(cusparseSpGEMMDescr) :: SpGEMMDesc
Integer(8) :: bufferSize1
! Integer(1),pointer,device :: buffer1(:)
Integer(1),device,allocatable :: buffer1(:)
!!Define Matrix----------------------------------
Arow=(/1,4,5,8,10/)
Acol=(/1,3,4,2,1,3,4,2,4/)
Aval=(/1d0,2d0,3d0,4d0,5d0,6d0,7d0,8d0,9d0/)
Brow=(/1,3,5,8,9/)
Bcol=(/1,4,2,4,1,2,3,2/)
Bval=(/1d0,2d0,3d0,4d0,5d0,6d0,7d0,8d0/)
status=cudaDeviceSynchronize
Arow_d=Arow
Acol_d=Acol
Aval_d=Aval
Brow_d=Brow
Bcol_d=Bcol
Bval_d=Bval
status=cudaDeviceSynchronize
!!Define Matrix----------------------------------
allocate(Crow_d(A_rows+1))
! initalize CUSPARSE and matrix descriptor
status=cusparseCreate(handle)
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseCreate error: ', status
status=cusparseCreateCsr(matA,A_rows,A_cols,A_nnz, &
ARow_d,ACol_d,Aval_d, &
CUSPARSE_INDEX_32I,CUSPARSE_INDEX_32I, &
CUSPARSE_INDEX_BASE_ONE,CUDA_R_64F)
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseCreateCsr error: ', status
status=cusparseCreateCsr(matB,B_rows,B_cols,B_nnz, &
BRow_d,BCol_d,Bval_d, &
CUSPARSE_INDEX_32I,CUSPARSE_INDEX_32I, &
CUSPARSE_INDEX_BASE_ONE,CUDA_R_64F)
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseCreateCsr error: ', status
status=cusparseCreateCsr(matC,A_rows,B_cols,0, &
null(),null(),null(), &
CUSPARSE_INDEX_32I,CUSPARSE_INDEX_32I, &
CUSPARSE_INDEX_BASE_ONE,CUDA_R_64F)
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseCreateCsr error: ', status
status=cudaDeviceSynchronize
!!----------------------------------------------------------------------------------------------------
!!SpGEMM computation
status=cusparseSpGEMM_createDescr(SpGEMMDesc)
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseSpGEMM_CreateDescr error: ', status
!! ask bufferSize1 bytes for external memory
status=cusparseSpGEMM_workEstimation(handle,&
CUSPARSE_OPERATION_NON_TRANSPOSE,CUSPARSE_OPERATION_NON_TRANSPOSE,&
alpha,matA,matB,beta,matC,&
CUDA_R_64F,CUSPARSE_SPGEMM_DEFAULT,&
SpGEMMDesc,bufferSize1,null())
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseSpGEMM_workEstimation error: ', status
if(allocated(buffer1)) deallocate(buffer1)
if(bufferSize1 /= 0) allocate(buffer1(bufferSize1))
!! inspect the A and B to understand the memory requirement for the next stop
status=cusparseSpGEMM_workEstimation(handle,&
CUSPARSE_OPERATION_NON_TRANSPOSE,CUSPARSE_OPERATION_NON_TRANSPOSE,&
alpha,matA,matB,beta,matC,&
CUDA_R_64F,CUSPARSE_SPGEMM_DEFAULT,&
SpGEMMDesc,bufferSize1,buffer1)
if(status/=CUSPARSE_STATUS_SUCCESS) print *, 'cusparseSpGEMM_workEstimation error: ', status
deallocate(buffer1)
status=cusparseSpGEMM_destroyDescr(SpGEMMDesc)
status=cusparseDestroySpMat(matA)
status=cusparseDestroySpMat(matB)
status=cusparseDestroySpMat(matC)
status=cusparseDestroy(handle)
return
end program SpGEMM