ITensor / ITensor/BlockSparseArrays.jl
[BUG] Block sparse matvec is broken
Open
Nobody has claimed this yet.
bug
- Dominant language
- Julia
- Stars
- 3
- Forks
- 3
- Avg merge
- 14m
- Merged PRs (30d)
- 3
Description
For example:
julia> using BlockSparseArrays
julia> r = [2, 2];
julia> a1 = BlockSparseArray{Float64}(undef, r, r);
julia> a2 = BlockSparseArray{Float64}(undef, r);
julia> a1 * a2
#...
julia> err
1-element ExceptionStack:
Not implemented.
Stacktrace:
[1] error(s::String)
@ Base ./error.jl:44
[2] materialize!(m::ArrayLayouts.MulAdd{BlockArrays.BlockLayout{SparseArraysBase.SparseLayout, ArrayLayouts.DenseColumnMajor}, BlockArrays.BlockLayout{SparseArraysBase.SparseLayout, ArrayLayouts.DenseColumnMajor}, BlockArrays.BlockLayout{SparseArraysBase.SparseLayout, ArrayLayouts.DenseColumnMajor}, Float64, BlockSparseMatrix{Float64, Matrix{Float64}, SparseArraysBase.SparseMatrixDOK{Matrix{Float64}, BlockSparseArrays.ZeroBlocks{2, Matrix{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}, BlockSparseVector{Float64, Vector{Float64}, SparseArraysBase.SparseVectorDOK{Vector{Float64}, BlockSparseArrays.ZeroBlocks{1, Vector{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}, BlockSparseVector{Float64, Vector{Float64}, SparseArraysBase.SparseVectorDOK{Vector{Float64}, BlockSparseArrays.ZeroBlocks{1, Vector{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}})
@ BlockSparseArrays ~/.julia/packages/BlockSparseArrays/gEM83/src/blocksparsearrayinterface/arraylayouts.jl:31
[3] muladd!(α::Float64, A::BlockSparseMatrix{Float64, Matrix{Float64}, SparseArraysBase.SparseMatrixDOK{Matrix{Float64}, BlockSparseArrays.ZeroBlocks{2, Matrix{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}, B::BlockSparseVector{Float64, Vector{Float64}, SparseArraysBase.SparseVectorDOK{Vector{Float64}, BlockSparseArrays.ZeroBlocks{1, Vector{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}, β::Float64, C::BlockSparseVector{Float64, Vector{Float64}, SparseArraysBase.SparseVectorDOK{Vector{Float64}, BlockSparseArrays.ZeroBlocks{1, Vector{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}; kw::@Kwargs{Czero::Bool})
@ ArrayLayouts ~/.julia/packages/ArrayLayouts/RCHQ8/src/muladd.jl:75
[4] copyto!
@ ~/.julia/packages/ArrayLayouts/RCHQ8/src/muladd.jl:82 [inlined]
[5] copy(M::ArrayLayouts.MulAdd{BlockArrays.BlockLayout{SparseArraysBase.SparseLayout, ArrayLayouts.DenseColumnMajor}, BlockArrays.BlockLayout{SparseArraysBase.SparseLayout, ArrayLayouts.DenseColumnMajor}, ArrayLayouts.ZerosLayout, Float64, BlockSparseMatrix{Float64, Matrix{Float64}, SparseArraysBase.SparseMatrixDOK{Matrix{Float64}, BlockSparseArrays.ZeroBlocks{2, Matrix{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}, BlockSparseVector{Float64, Vector{Float64}, SparseArraysBase.SparseVectorDOK{Vector{Float64}, BlockSparseArrays.ZeroBlocks{1, Vector{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}, FillArrays.Zeros{Float64, 1, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}})
@ ArrayLayouts ~/.julia/packages/ArrayLayouts/RCHQ8/src/muladd.jl:77
[6] copy
@ ~/.julia/packages/ArrayLayouts/RCHQ8/src/mul.jl:140 [inlined]
[7] materialize
@ ~/.julia/packages/ArrayLayouts/RCHQ8/src/mul.jl:137 [inlined]
[8] mul
@ ~/.julia/packages/ArrayLayouts/RCHQ8/src/mul.jl:138 [inlined]
[9] *(A::BlockSparseMatrix{Float64, Matrix{Float64}, SparseArraysBase.SparseMatrixDOK{Matrix{Float64}, BlockSparseArrays.ZeroBlocks{2, Matrix{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}, BlockedOneTo{Int64, Vector{Int64}}}}, B::BlockSparseVector{Float64, Vector{Float64}, SparseArraysBase.SparseVectorDOK{Vector{Float64}, BlockSparseArrays.ZeroBlocks{1, Vector{Float64}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}}}, Tuple{BlockedOneTo{Int64, Vector{Int64}}}})
@ ArrayLayouts ~/.julia/packages/ArrayLayouts/RCHQ8/src/mul.jl:226
[10] top-level scope
@ REPL[49]:1
[11] top-level scope
@ REPL:1
A workaround is adding a trivial dimension:
julia> a2 = BlockSparseArray{Float64}(undef, r, [1]);
julia> a1 * a2
2×1-blocked 4×1 BlockSparseMatrix{Float64, Matrix{Float64}, …, …}:
⋅
⋅
───
⋅
⋅
See #197.
Contributor guide
No contributing guide indexed for this repository
First steps
- Read the whole issue, then the project's contributing guide.
- Comment on the issue to say you are picking it up — it saves two people doing the same work.
- Fork the repository and make your change on a branch.
- Open a pull request that references the issue number.
Research direction
Start with src/blocksparsearrayinterface/arraylayouts.jl at the materialize! implementation named in the stack trace, then reproduce the matrix-vector multiplication example from the issue. Compare it with the workaround using a trivial dimension; done means the direct BlockSparseMatrix × BlockSparseVector operation works without raising Not implemented, with coverage for this case.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- julia
- Domain
- data
- Issue type
- Bug
- Difficulty
- 3/5
- Estimated time
- 1-2 days
- Activity status
- Stale
- Clarity
- Clearly specified
- Newbie friendliness
- 52/100