diff --git a/lib/cusparse/conversions.jl b/lib/cusparse/conversions.jl index fc1f2c75a3..e477e049d4 100644 --- a/lib/cusparse/conversions.jl +++ b/lib/cusparse/conversions.jl @@ -332,7 +332,7 @@ for elty in (:Float32, :Float64, :ComplexF32, :ComplexF64) @eval begin function CuSparseMatrixCSC{$elty}(csr::CuSparseMatrixCSR{$elty}; index::SparseChar='O', action::cusparseAction_t=CUSPARSE_ACTION_NUMERIC, algo::cusparseCsr2CscAlg_t=CUSPARSE_CSR2CSC_ALG1) m,n = size(csr) - colPtr = CUDA.zeros(Cint, n+1) + colPtr = (index == 'O') ? CUDA.ones(Cint, n+1) : CUDA.zeros(Cint, n+1) rowVal = CUDA.zeros(Cint, nnz(csr)) nzVal = CUDA.zeros($elty, nnz(csr)) function bufferSize() @@ -352,7 +352,7 @@ for elty in (:Float32, :Float64, :ComplexF32, :ComplexF64) function CuSparseMatrixCSR{$elty}(csc::CuSparseMatrixCSC{$elty}; index::SparseChar='O', action::cusparseAction_t=CUSPARSE_ACTION_NUMERIC, algo::cusparseCsr2CscAlg_t=CUSPARSE_CSR2CSC_ALG1) m,n = size(csc) - rowPtr = CUDA.zeros(Cint,m+1) + rowPtr = (index == 'O') ? CUDA.ones(Cint, m+1) : CUDA.zeros(Cint, m+1) colVal = CUDA.zeros(Cint,nnz(csc)) nzVal = CUDA.zeros($elty,nnz(csc)) function bufferSize() @@ -379,7 +379,7 @@ for (elty, welty) in ((:Float16, :Float32), @eval begin function CuSparseMatrixCSC{$elty}(csr::CuSparseMatrixCSR{$elty}; index::SparseChar='O', action::cusparseAction_t=CUSPARSE_ACTION_NUMERIC, algo::cusparseCsr2CscAlg_t=CUSPARSE_CSR2CSC_ALG1) m,n = size(csr) - colPtr = CUDA.zeros(Cint, n+1) + colPtr = (index == 'O') ? CUDA.ones(Cint, n+1) : CUDA.zeros(Cint, n+1) rowVal = CUDA.zeros(Cint, nnz(csr)) nzVal = CUDA.zeros($elty, nnz(csr)) if $elty == Float16 #broken for ComplexF16? @@ -405,7 +405,7 @@ for (elty, welty) in ((:Float16, :Float32), function CuSparseMatrixCSR{$elty}(csc::CuSparseMatrixCSC{$elty}; index::SparseChar='O', action::cusparseAction_t=CUSPARSE_ACTION_NUMERIC, algo::cusparseCsr2CscAlg_t=CUSPARSE_CSR2CSC_ALG1) m,n = size(csc) - rowPtr = CUDA.zeros(Cint,m+1) + rowPtr = (index == 'O') ? CUDA.ones(Cint, m+1) : CUDA.zeros(Cint, m+1) colVal = CUDA.zeros(Cint,nnz(csc)) nzVal = CUDA.zeros($elty,nnz(csc)) if $elty == Float16 #broken for ComplexF16? @@ -586,10 +586,10 @@ end ## CSR to COO and vice-versa function CuSparseMatrixCSR{Tv}(coo::CuSparseMatrixCOO{Tv}; index::SparseChar='O') where {Tv} - m,n = size(coo) - nnz(coo) == 0 && return CuSparseMatrixCSR{Tv}(CUDA.ones(Cint, m+1), coo.colInd, nonzeros(coo), size(coo)) + m, n = size(coo) + csrRowPtr = (index == 'O') ? CUDA.ones(Cint, m + 1) : CUDA.zeros(Cint, m + 1) + nnz(coo) == 0 && return CuSparseMatrixCSR{Tv}(csrRowPtr, coo.colInd, nonzeros(coo), size(coo)) coo = sort_coo(coo, 'R') - csrRowPtr = CuVector{Cint}(undef, m+1) cusparseXcoo2csr(handle(), coo.rowInd, nnz(coo), m, csrRowPtr, index) CuSparseMatrixCSR{Tv}(csrRowPtr, coo.colInd, nonzeros(coo), size(coo)) end @@ -605,10 +605,10 @@ end ### CSC to COO and viceversa function CuSparseMatrixCSC{Tv}(coo::CuSparseMatrixCOO{Tv}; index::SparseChar='O') where {Tv} - m,n = size(coo) - nnz(coo) == 0 && return CuSparseMatrixCSC{Tv}(CUDA.ones(Cint, n+1), coo.rowInd, nonzeros(coo), size(coo)) + m, n = size(coo) + cscColPtr = (index == 'O') ? CUDA.ones(Cint, n + 1) : CUDA.zeros(Cint, n + 1) + nnz(coo) == 0 && return CuSparseMatrixCSC{Tv}(cscColPtr, coo.rowInd, nonzeros(coo), size(coo)) coo = sort_coo(coo, 'C') - cscColPtr = CuVector{Cint}(undef, n+1) cusparseXcoo2csr(handle(), coo.colInd, nnz(coo), n, cscColPtr, index) CuSparseMatrixCSC{Tv}(cscColPtr, coo.rowInd, nonzeros(coo), size(coo)) end diff --git a/test/libraries/cusparse/conversions.jl b/test/libraries/cusparse/conversions.jl index 870d591fa0..f946fd1967 100644 --- a/test/libraries/cusparse/conversions.jl +++ b/test/libraries/cusparse/conversions.jl @@ -113,17 +113,19 @@ if !(v"12.0" <= CUSPARSE.version() < v"12.1") end end -for (n, bd, p) in [(100, 5, 0.02), (5, 1, 0.8), (4, 2, 0.5)] +for (n, bd, p) in [(100, 5, 0.02), (5, 1, 0.8), (4, 2, 0.5), (0, 1, 0.0)] v"12.0" <= CUSPARSE.version() < v"12.1" && n == 4 && continue @testset "conversions between CuSparseMatrices (n, bd, p) = ($n, $bd, $p)" begin A = sprand(n, n, p) blockdim = bd for CuSparseMatrixType1 in (CuSparseMatrixCSC, CuSparseMatrixCSR, CuSparseMatrixCOO, CuSparseMatrixBSR) + if CuSparseMatrixType1 == CuSparseMatrixBSR && n == 0 continue end # TODO: conversion to CuSparseMatrixBSR breaks with (0x0) matrices dA1 = CuSparseMatrixType1 == CuSparseMatrixBSR ? CuSparseMatrixType1(A, blockdim) : CuSparseMatrixType1(A) @testset "conversion $CuSparseMatrixType1 --> SparseMatrixCSC" begin @test SparseMatrixCSC(dA1) ≈ A end for CuSparseMatrixType2 in (CuSparseMatrixCSC, CuSparseMatrixCSR, CuSparseMatrixCOO, CuSparseMatrixBSR) + if CuSparseMatrixType2 == CuSparseMatrixBSR && n == 0 continue end # TODO: conversion to CuSparseMatrixBSR breaks with (0x0) matrices CuSparseMatrixType1 == CuSparseMatrixType2 && continue dA2 = CuSparseMatrixType2 == CuSparseMatrixBSR ? CuSparseMatrixType2(dA1, blockdim) : CuSparseMatrixType2(dA1) @testset "conversion $CuSparseMatrixType1 --> $CuSparseMatrixType2" begin