From d04693f18dea018e8ac2d91e3d05ac7913a8cd7c Mon Sep 17 00:00:00 2001 From: mtfishman Date: Fri, 7 Feb 2025 10:29:10 -0500 Subject: [PATCH 1/7] Fix a bug in blocktype --- Project.toml | 2 +- .../blocksparsearrayinterface.jl | 35 +++++++++---------- test/test_basics.jl | 3 ++ 3 files changed, 20 insertions(+), 20 deletions(-) diff --git a/Project.toml b/Project.toml index 50581f2b..b611987f 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "BlockSparseArrays" uuid = "2c9a651f-6452-4ace-a6ac-809f4280fbb4" authors = ["ITensor developers and contributors"] -version = "0.2.13" +version = "0.2.14" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" diff --git a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl index b907a6ad..2f79c1e7 100644 --- a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl +++ b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl @@ -41,29 +41,26 @@ function eachstoredblock(a::AbstractArray) return storedvalues(blocks(a)) end -# TODO: Generalize this, this catches simple cases -# where the more general definition isn't specific enough. -blocktype(a::Array) = typeof(a) -# TODO: Maybe unwrap SubArrays? function blocktype(a::AbstractArray) - # TODO: Unfortunately, this doesn't always give - # a concrete type, even when it could be concrete, i.e. - #= - ```julia - julia> eltype(blocks(BlockArray(randn(2, 2), [1, 1], [1, 1]))) - Matrix{Float64} (alias for Array{Float64, 2}) - - julia> eltype(blocks(BlockedArray(randn(2, 2), [1, 1], [1, 1]))) - AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) - - julia> eltype(blocks(randn(2, 2))) - AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) - ``` - =# if isempty(blocks(a)) + # TODO: Unfortunately, this doesn't always give + # a concrete type, even when it could be concrete, i.e. + #= + ```julia + julia> eltype(blocks(BlockArray(randn(2, 2), [1, 1], [1, 1]))) + Matrix{Float64} (alias for Array{Float64, 2}) + + julia> eltype(blocks(BlockedArray(randn(2, 2), [1, 1], [1, 1]))) + AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) + + julia> eltype(blocks(randn(2, 2))) + AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) + ``` + =# + # That is an issue in BlockArrays.jl that we should address. return eltype(blocks(a)) end - return eltype(first(blocks(a))) + return mapreduce(typeof, promote_type, blocks(a)) end abstract type AbstractBlockSparseArrayInterface <: AbstractSparseArrayInterface end diff --git a/test/test_basics.jl b/test/test_basics.jl index 17544287..e62f7cae 100644 --- a/test/test_basics.jl +++ b/test/test_basics.jl @@ -132,6 +132,9 @@ arrayts = (Array, JLArray) end end end + @testset "blocktype" begin + @test blocktype(arrayt(randn(elt, 2, 2))) <: SubArray{elt,2,arrayt{elt,2}} + end @testset "Basics" begin a = dev(BlockSparseArray{elt}([2, 3], [2, 3])) @allowscalar @test a == dev( From 833d4a7056d6132a46294b99b99bd49322358652 Mon Sep 17 00:00:00 2001 From: mtfishman Date: Fri, 7 Feb 2025 11:26:14 -0500 Subject: [PATCH 2/7] Add more tests for `blocktype` and `blockstype` --- .../blocksparsearrayinterface.jl | 9 ++++++ test/test_basics.jl | 30 +++++++++++++++++-- 2 files changed, 37 insertions(+), 2 deletions(-) diff --git a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl index 2f79c1e7..0bf5c148 100644 --- a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl +++ b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl @@ -41,6 +41,9 @@ function eachstoredblock(a::AbstractArray) return storedvalues(blocks(a)) end +function blockstype(a::AbstractArray) + return typeof(blocks(a)) +end function blocktype(a::AbstractArray) if isempty(blocks(a)) # TODO: Unfortunately, this doesn't always give @@ -63,6 +66,12 @@ function blocktype(a::AbstractArray) return mapreduce(typeof, promote_type, blocks(a)) end +using BlockArrays: BlockArray +blockstype(::Type{<:BlockArray{<:Any,<:Any,B}}) where {B} = B +blockstype(a::BlockArray) = blockstype(typeof(a)) +blocktype(arraytype::Type{<:BlockArray}) = eltype(blockstype(arraytype)) +blocktype(a::BlockArray) = eltype(blocks(a)) + abstract type AbstractBlockSparseArrayInterface <: AbstractSparseArrayInterface end # TODO: Also support specifying the `blocktype` along with the `eltype`. diff --git a/test/test_basics.jl b/test/test_basics.jl index e62f7cae..941992e6 100644 --- a/test/test_basics.jl +++ b/test/test_basics.jl @@ -132,8 +132,34 @@ arrayts = (Array, JLArray) end end end - @testset "blocktype" begin - @test blocktype(arrayt(randn(elt, 2, 2))) <: SubArray{elt,2,arrayt{elt,2}} + @testset "blockstype, blocktype" begin + a = arrayt(randn(elt, 2, 2)) + @test blockstype(a) <: BlockArrays.BlocksView{elt,2} + # TODO: This is difficult to determine just from type information. + @test_broken blockstype(typeof(a)) <: BlockArrays.BlocksView{elt,2} + @test blocktype(a) <: SubArray{elt,2,arrayt{elt,2}} + # TODO: This is difficult to determine just from type information. + @test_broken blocktype(typeof(a)) <: SubArray{elt,2,arrayt{elt,2}} + + a = BlockSparseMatrix{elt,arrayt{elt,2}}([1, 1], [1, 1]) + @test blockstype(a) <: SparseMatrixDOK{arrayt{elt,2}} + @test blockstype(typeof(a)) <: SparseMatrixDOK{arrayt{elt,2}} + @test blocktype(a) <: arrayt{elt,2} + @test blocktype(typeof(a)) <: arrayt{elt,2} + + a = BlockArray(arrayt(randn(elt, (2, 2))), [1, 1], [1, 1]) + @test blockstype(a) === Matrix{arrayt{elt,2}} + @test blockstype(typeof(a)) === Matrix{arrayt{elt,2}} + @test blocktype(a) <: arrayt{elt,2} + @test blocktype(typeof(a)) <: arrayt{elt,2} + + a = BlockedArray(arrayt(randn(elt, 2, 2)), [1, 1], [1, 1]) + @test blockstype(a) <: BlockArrays.BlocksView{elt,2} + # TODO: This is difficult to determine just from type information. + @test_broken blockstype(typeof(a)) <: BlockArrays.BlocksView{elt,2} + @test blocktype(a) <: SubArray{elt,2,arrayt{elt,2}} + # TODO: This is difficult to determine just from type information. + @test_broken blocktype(typeof(a)) <: SubArray{elt,2,arrayt{elt,2}} end @testset "Basics" begin a = dev(BlockSparseArray{elt}([2, 3], [2, 3])) From 174248f6ac3a75ca2e895e23254ed99dd0b9f837 Mon Sep 17 00:00:00 2001 From: mtfishman Date: Fri, 7 Feb 2025 11:37:36 -0500 Subject: [PATCH 3/7] Delete extraneous definition of `blockstype` --- src/blocksparsearrayinterface/blocksparsearrayinterface.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl index 0bf5c148..24da762c 100644 --- a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl +++ b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl @@ -84,8 +84,6 @@ struct BlockSparseArrayInterface <: AbstractBlockSparseArrayInterface end @interface ::AbstractBlockSparseArrayInterface BlockArrays.blocks(a::AbstractArray) = error("Not implemented") -blockstype(a::AbstractArray) = blockstype(typeof(a)) - @interface ::AbstractBlockSparseArrayInterface function Base.getindex( a::AbstractArray{<:Any,N}, I::Vararg{Int,N} ) where {N} From 0b5fcdbc55e790de23f9873a739a73eec77d76de Mon Sep 17 00:00:00 2001 From: mtfishman Date: Fri, 7 Feb 2025 11:41:16 -0500 Subject: [PATCH 4/7] Fix missing imports in test --- test/test_basics.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/test/test_basics.jl b/test/test_basics.jl index 941992e6..7c2fa0d4 100644 --- a/test/test_basics.jl +++ b/test/test_basics.jl @@ -1,13 +1,16 @@ using Adapt: adapt using ArrayLayouts: zero! using BlockArrays: + BlockArrays, Block, + BlockArray, BlockIndexRange, BlockRange, BlockSlice, BlockVector, BlockedOneTo, BlockedUnitRange, + BlockedArray, BlockedVector, blockedrange, blocklength, From 1eac901e0900dcf117d792d5c46ec1ecc66e02a5 Mon Sep 17 00:00:00 2001 From: mtfishman Date: Fri, 7 Feb 2025 12:20:32 -0500 Subject: [PATCH 5/7] Test inference --- test/Project.toml | 1 + test/test_basics.jl | 25 +++++++++++++------------ 2 files changed, 14 insertions(+), 12 deletions(-) diff --git a/test/Project.toml b/test/Project.toml index 0cc36241..e5418bd8 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -18,3 +18,4 @@ Suppressor = "fd094767-a336-5f1f-9728-57cf17d0bbfb" SymmetrySectors = "f8a8ad64-adbc-4fce-92f7-ffe2bb36a86e" TensorAlgebra = "68bd88dc-f39d-4e12-b2ca-f046b68fcc6a" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +TestExtras = "5ed8adda-3752-4e41-b88a-e8b09835ee3a" diff --git a/test/test_basics.jl b/test/test_basics.jl index 7c2fa0d4..6120bb95 100644 --- a/test/test_basics.jl +++ b/test/test_basics.jl @@ -37,6 +37,7 @@ using LinearAlgebra: Adjoint, Transpose, dot, mul!, norm using SparseArraysBase: SparseArrayDOK, SparseMatrixDOK, SparseVectorDOK, storedlength using TensorAlgebra: contract using Test: @test, @test_broken, @test_throws, @testset, @inferred +using TestExtras: @constinferred include("TestBlockSparseArraysUtils.jl") arrayts = (Array, JLArray) @@ -137,30 +138,30 @@ arrayts = (Array, JLArray) end @testset "blockstype, blocktype" begin a = arrayt(randn(elt, 2, 2)) - @test blockstype(a) <: BlockArrays.BlocksView{elt,2} + @test (@constinferred blockstype(a)) <: BlockArrays.BlocksView{elt,2} # TODO: This is difficult to determine just from type information. @test_broken blockstype(typeof(a)) <: BlockArrays.BlocksView{elt,2} - @test blocktype(a) <: SubArray{elt,2,arrayt{elt,2}} + @test (@constinferred blocktype(a)) <: SubArray{elt,2,arrayt{elt,2}} # TODO: This is difficult to determine just from type information. @test_broken blocktype(typeof(a)) <: SubArray{elt,2,arrayt{elt,2}} a = BlockSparseMatrix{elt,arrayt{elt,2}}([1, 1], [1, 1]) - @test blockstype(a) <: SparseMatrixDOK{arrayt{elt,2}} - @test blockstype(typeof(a)) <: SparseMatrixDOK{arrayt{elt,2}} - @test blocktype(a) <: arrayt{elt,2} - @test blocktype(typeof(a)) <: arrayt{elt,2} + @test (@constinferred blockstype(a)) <: SparseMatrixDOK{arrayt{elt,2}} + @test (@constinferred blockstype(typeof(a))) <: SparseMatrixDOK{arrayt{elt,2}} + @test (@constinferred blocktype(a)) <: arrayt{elt,2} + @test (@constinferred blocktype(typeof(a))) <: arrayt{elt,2} a = BlockArray(arrayt(randn(elt, (2, 2))), [1, 1], [1, 1]) - @test blockstype(a) === Matrix{arrayt{elt,2}} - @test blockstype(typeof(a)) === Matrix{arrayt{elt,2}} - @test blocktype(a) <: arrayt{elt,2} - @test blocktype(typeof(a)) <: arrayt{elt,2} + @test (@constinferred blockstype(a)) === Matrix{arrayt{elt,2}} + @test (@constinferred blockstype(typeof(a))) === Matrix{arrayt{elt,2}} + @test (@constinferred blocktype(a)) <: arrayt{elt,2} + @test (@constinferred blocktype(typeof(a))) <: arrayt{elt,2} a = BlockedArray(arrayt(randn(elt, 2, 2)), [1, 1], [1, 1]) - @test blockstype(a) <: BlockArrays.BlocksView{elt,2} + @test (@constinferred blockstype(a)) <: BlockArrays.BlocksView{elt,2} # TODO: This is difficult to determine just from type information. @test_broken blockstype(typeof(a)) <: BlockArrays.BlocksView{elt,2} - @test blocktype(a) <: SubArray{elt,2,arrayt{elt,2}} + @test (@constinferred blocktype(a)) <: SubArray{elt,2,arrayt{elt,2}} # TODO: This is difficult to determine just from type information. @test_broken blocktype(typeof(a)) <: SubArray{elt,2,arrayt{elt,2}} end From 9ecd84230a3eef470a7ea04a25bfc47cdb99a5e6 Mon Sep 17 00:00:00 2001 From: mtfishman Date: Fri, 7 Feb 2025 12:31:51 -0500 Subject: [PATCH 6/7] Improve type stability --- .../blocksparsearrayinterface.jl | 36 ++++++++++--------- 1 file changed, 20 insertions(+), 16 deletions(-) diff --git a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl index 24da762c..04fd135e 100644 --- a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl +++ b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl @@ -44,24 +44,28 @@ end function blockstype(a::AbstractArray) return typeof(blocks(a)) end + +#= +Ideally this would just be defined as `eltype(blockstype(a))`. +However, BlockArrays.jl doesn't make `eltype(blocks(a))` concrete +even when it could be +(https://github.com/JuliaArrays/BlockArrays.jl/blob/v1.4.0/src/blocks.jl#L71-L74): +```julia +julia> eltype(blocks(BlockArray(randn(2, 2), [1, 1], [1, 1]))) +Matrix{Float64} (alias for Array{Float64, 2}) + +julia> eltype(blocks(BlockedArray(randn(2, 2), [1, 1], [1, 1]))) +AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) + +julia> eltype(blocks(randn(2, 2))) +AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) +``` +Also note the current definition doesn't handle the limit +when `blocks(a)` is empty +=# function blocktype(a::AbstractArray) if isempty(blocks(a)) - # TODO: Unfortunately, this doesn't always give - # a concrete type, even when it could be concrete, i.e. - #= - ```julia - julia> eltype(blocks(BlockArray(randn(2, 2), [1, 1], [1, 1]))) - Matrix{Float64} (alias for Array{Float64, 2}) - - julia> eltype(blocks(BlockedArray(randn(2, 2), [1, 1], [1, 1]))) - AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) - - julia> eltype(blocks(randn(2, 2))) - AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) - ``` - =# - # That is an issue in BlockArrays.jl that we should address. - return eltype(blocks(a)) + error("`blocktype` can't be determined if `isempty(blocks(a))`.") end return mapreduce(typeof, promote_type, blocks(a)) end From 767d82bcf72c81e3bf14bcf6848da8e4f79e1c95 Mon Sep 17 00:00:00 2001 From: Matt Fishman Date: Fri, 7 Feb 2025 13:59:35 -0500 Subject: [PATCH 7/7] Expand comment a bit --- .../blocksparsearrayinterface.jl | 15 +++++++++++++-- 1 file changed, 13 insertions(+), 2 deletions(-) diff --git a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl index 04fd135e..6f4200aa 100644 --- a/src/blocksparsearrayinterface/blocksparsearrayinterface.jl +++ b/src/blocksparsearrayinterface/blocksparsearrayinterface.jl @@ -60,8 +60,19 @@ AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) julia> eltype(blocks(randn(2, 2))) AbstractMatrix{Float64} (alias for AbstractArray{Float64, 2}) ``` -Also note the current definition doesn't handle the limit -when `blocks(a)` is empty +Also note the current definition errors in the limit +when `blocks(a)` is empty, but even empty arrays generally +have at least one block: +```julia +julia> length(blocks(randn(0))) +1 + +julia> length(blocks(BlockVector{Float64}(randn(0)))) +1 + +julia> length(blocks(BlockedVector{Float64}(randn(0)))) +1 +``` =# function blocktype(a::AbstractArray) if isempty(blocks(a))