Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Implement (limited) broadcast of sparse arrays #1367

Merged
merged 4 commits into from
Feb 11, 2022
Merged

Conversation

maleadt
Copy link
Member

@maleadt maleadt commented Feb 11, 2022

Demo:

julia> cu(sparse([1f0 0f0; 1f0 1f0; 0f0 1f0]))
3×2 CUDA.CUSPARSE.CuSparseMatrixCSC{Float32, Int32} with 4 stored entries:
 1.0    
 1.0  1.0
     1.0

# best case: retains sparse output with zero-preserving function and scalar inputs
julia> cu(sparse([1f0 0f0; 1f0 1f0; 0f0 1f0])) .* 1
3×2 CUDA.CUSPARSE.CuSparseMatrixCSC{Float32, Int32} with 4 stored entries:
 1.0    
 1.0  1.0
     1.0

# with non zero-preserving function, output is dense
julia> cu(sparse([1f0 0f0; 1f0 1f0; 0f0 1f0])) .+ 1
3×2 CuArray{Float32, 2, CUDA.Mem.DeviceBuffer}:
 2.0  1.0
 2.0  2.0
 1.0  2.0

# or when a dense input is involved
julia> cu(sparse([1f0 0f0; 1f0 1f0; 0f0 1f0])) .* CUDA.ones(3,2)
3×2 CuArray{Float32, 2, CUDA.Mem.DeviceBuffer}:
 1.0  0.0
 1.0  1.0
 0.0  1.0

Limitations:

  • not implemented for BSR
  • only supports one sparse input

cc @Roger-luo @ChrisRackauckas

@maleadt maleadt added enhancement New feature or request cuda array Stuff about CuArray. labels Feb 11, 2022
@codecov
Copy link

codecov bot commented Feb 11, 2022

Codecov Report

Merging #1367 (396e6c4) into master (243c4ae) will decrease coverage by 0.78%.
The diff coverage is 62.50%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master    #1367      +/-   ##
==========================================
- Coverage   78.84%   78.06%   -0.79%     
==========================================
  Files         119      121       +2     
  Lines        8699     8929     +230     
==========================================
+ Hits         6859     6970     +111     
- Misses       1840     1959     +119     
Impacted Files Coverage Δ
lib/cusparse/CUSPARSE.jl 77.27% <ø> (ø)
lib/cusparse/device.jl 8.33% <0.00%> (ø)
lib/cusparse/generic.jl 98.26% <ø> (ø)
src/CUDA.jl 100.00% <ø> (ø)
src/random.jl 75.00% <ø> (ø)
lib/cusparse/broadcast.jl 56.81% <56.81%> (ø)
lib/cusparse/array.jl 67.07% <67.56%> (-1.70%) ⬇️
lib/cusparse/level3.jl 75.00% <83.33%> (ø)
lib/cusparse/conversions.jl 81.55% <85.71%> (ø)
lib/cusparse/level1.jl 100.00% <100.00%> (ø)
... and 9 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 243c4ae...396e6c4. Read the comment docs.

@maleadt maleadt merged commit 0aa7750 into master Feb 11, 2022
@maleadt maleadt deleted the tb/sparse_broadcast branch February 11, 2022 10:16
@@ -17,12 +18,12 @@ Base.convert(T::Type{<:AbstractCuSparseArray}, m::AbstractArray) = m isa T ? m :
mutable struct CuSparseVector{Tv, Ti} <: AbstractCuSparseVector{Tv, Ti}
iPtr::CuVector{Ti}
nzVal::CuVector{Tv}
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It would be nice if this was nzval like Base sparse has. That would make general code work.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's unlikely any generic code is actually going to work, but sure.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I just checked and generic code (SparseDiffTools.jl) almost does work, sans #101

@ChrisRackauckas
Copy link
Member

Could this be used to implement fill!(cuA,true)? That would solve the rest of most of SciML/OrdinaryDiffEq.jl#1566 .

@maleadt
Copy link
Member Author

maleadt commented Feb 11, 2022

Could this be used to implement fill!(cuA,true)?

Could you elaborate? What is cuA? fill! only makes sense for dense arrays?

@ChrisRackauckas
Copy link
Member

fill!(A,true) on a sparse matrix fills in only the non-structurally zero values, i.e. it preserves the sparsity pattern and makes all of the non-zero values be 1.

@ChrisRackauckas
Copy link
Member

cuA is just a CUDA.CUSPARSE.CuSparseMatrixCSC

@ChrisRackauckas
Copy link
Member

ChrisRackauckas commented Feb 11, 2022

using SparseArrays, CUDA
A = sprand(100,100,0.1)
fill!(A,true)
cuA = cu(A)
fill!(cuA,true) # fails

@maleadt
Copy link
Member Author

maleadt commented Feb 11, 2022

fill!(A,true) on a sparse matrix fills in only the non-structurally zero values, i.e. it preserves the sparsity pattern and makes all of the non-zero values be 1.

That doesn't look to be the case

julia> x = sprand(2,2,0.5)
2×2 SparseMatrixCSC{Float64, Int64} with 1 stored entry:
 0.507603    
            

julia> fill!(x, true)
2×2 SparseMatrixCSC{Float64, Int64} with 4 stored entries:
 1.0  1.0
 1.0  1.0

@ChrisRackauckas
Copy link
Member

peculiar...

@maleadt
Copy link
Member Author

maleadt commented Feb 11, 2022

That said, yes, such an operation could be implemented with the support I have here (although it isn't available right now, as broadcast assignment also fills zeros). But if there's an interface, I can implement it.

@ChrisRackauckas
Copy link
Member

looks like it should be fill!(nonzeros(A),true), I assume nonzeros is extended?

@ChrisRackauckas
Copy link
Member

Okay cool, I checked I think about everything is supported. But I just need to figure out how to check if something is a sparse matrix without knowing about CUDA. Is there some issparsematrix function or something that was extended somewhere?

@ChrisRackauckas
Copy link
Member

JuliaLang/julia#17670

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
cuda array Stuff about CuArray. enhancement New feature or request
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants