From a87efd8ce95c07199f8b1604b56d8f1825490f8e Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Tue, 11 Jun 2024 18:20:50 +0530 Subject: [PATCH 01/94] test: fix symbolic indexing adjoint test --- test/downstream/symbol_indexing.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/downstream/symbol_indexing.jl b/test/downstream/symbol_indexing.jl index d398a886..1dc45dd3 100644 --- a/test/downstream/symbol_indexing.jl +++ b/test/downstream/symbol_indexing.jl @@ -40,7 +40,7 @@ gs, = Zygote.gradient(sol) do sol end @testset "Symbolic Indexing ADjoint" begin - @test all(all.(isone, gs.u)) + @test all(all.(isone, gs)) end # Tables interface From 320b577c7459ef8a0ec42ce47122d65246d0d6d0 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 11 Jun 2024 15:22:58 -0400 Subject: [PATCH 02/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index c04a8cc5..15c2e649 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.23.0" +version = "3.23.1" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 721c2a7007511f45cf59e826d7da2790971e69ef Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Fri, 3 May 2024 20:12:45 +0530 Subject: [PATCH 03/94] feat: add parameter timeseries support to `AbstractDiffEqArray` --- src/vector_of_array.jl | 88 +++++++++++++++++++++--------------------- 1 file changed, 44 insertions(+), 44 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 2c13a573..36f95ed2 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -60,11 +60,13 @@ A[1, :] # all time periods for f(t) A.t ``` """ -mutable struct DiffEqArray{T, N, A, B, F, S} <: AbstractDiffEqArray{T, N, A} +mutable struct DiffEqArray{T, N, A, B, F, S, D <: Union{Nothing, ParameterTimeseriesCollection}} <: + AbstractDiffEqArray{T, N, A} u::A # A <: AbstractVector{<: AbstractArray{T, N - 1}} t::B p::F sys::S + discretes::D end ### Abstract Interface struct AllObserved @@ -174,29 +176,32 @@ function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p = nothing, - sys = nothing) where {T, N} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys)}(vec, + sys = nothing; discretes = nothing) where {T, N} + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, ts, p, - sys) + sys, + discretes) end # ambiguity resolution function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, ::NTuple{N, Int}) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing}(vec, + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing, Nothing}(vec, ts, nothing, + nothing, nothing) end function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, - ::NTuple{N, Int}, p) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing}(vec, + ::NTuple{N, Int}, p; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, ts, p, - nothing) + nothing, + discretes) end # Assume that the first element is representative of all other elements @@ -206,7 +211,8 @@ function DiffEqArray(vec::AbstractVector, sys = nothing; variables = nothing, parameters = nothing, - independent_variables = nothing) + independent_variables = nothing, + discretes = nothing) sys = something(sys, SymbolCache(something(variables, []), something(parameters, []), @@ -219,11 +225,13 @@ function DiffEqArray(vec::AbstractVector, typeof(vec), typeof(ts), typeof(p), - typeof(sys) + typeof(sys), + typeof(discretes) }(vec, ts, p, - sys) + sys, + discretes) end function DiffEqArray(vec::AbstractVector{VT}, @@ -232,7 +240,8 @@ function DiffEqArray(vec::AbstractVector{VT}, sys = nothing; variables = nothing, parameters = nothing, - independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} + independent_variables = nothing, + discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} sys = something(sys, SymbolCache(something(variables, []), something(parameters, []), @@ -243,18 +252,30 @@ function DiffEqArray(vec::AbstractVector{VT}, typeof(vec), typeof(ts), typeof(p), - typeof(sys) + typeof(sys), + typeof(discretes), }(vec, ts, p, - sys) + sys, + discretes) end +has_discretes(::T) where {T <: AbstractDiffEqArray} = hasfield(T, :discretes) +get_discretes(x) = getfield(x, :discretes) + SymbolicIndexingInterface.is_timeseries(::Type{<:AbstractVectorOfArray}) = Timeseries() +function SymbolicIndexingInterface.is_parameter_timeseries(::Type{DiffEqArray{T, N, A, B, + F, S, D}}) where {T, N, A, B, F, S, D <: ParameterIndexingProxy} + Timeseries() +end SymbolicIndexingInterface.state_values(A::AbstractDiffEqArray) = A.u SymbolicIndexingInterface.current_time(A::AbstractDiffEqArray) = A.t SymbolicIndexingInterface.parameter_values(A::AbstractDiffEqArray) = A.p SymbolicIndexingInterface.symbolic_container(A::AbstractDiffEqArray) = A.sys +function SymbolicIndexingInterface.get_parameter_timeseries_collection(A::AbstractDiffEqArray) + return get_discretes(A) +end Base.IndexStyle(A::AbstractVectorOfArray) = Base.IndexStyle(typeof(A)) Base.IndexStyle(::Type{<:AbstractVectorOfArray}) = IndexCartesian() @@ -363,39 +384,18 @@ end # Symbolic Indexing Methods for (symtype, elsymtype, valtype, errcheck) in [ - (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym))), - (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym))), + (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), + (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), (NotSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Union{<:Tuple, <:AbstractArray}, - :(all(x -> is_parameter(A, x), sym))), + :(all(x -> is_parameter(A, x) && !is_timeseries_parameter(A, x), sym))), ] -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym)(A) -end -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, arg) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym)(A, arg) -end -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, arg::Union{AbstractArray{Int}, AbstractArray{Bool}}) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym).((A,), arg) -end -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, ::Colon) - if $errcheck - throw(ParameterIndexingError(sym)) + @eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, + ::$elsymtype, sym::$valtype, arg...) + if $errcheck + throw(ParameterIndexingError(sym)) + end + getu(A, sym)(A, arg...) end - getu(A, sym)(A) -end end Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::ScalarSymbolic, From 33b2fe33d24f4c4dcd1d8bdaf4c726f548716fae Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Fri, 14 Jun 2024 12:41:56 +0530 Subject: [PATCH 04/94] build: bump SII compat --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 15c2e649..bb434c6e 100644 --- a/Project.toml +++ b/Project.toml @@ -59,7 +59,7 @@ StaticArrays = "1.6" StaticArraysCore = "1.4" Statistics = "1.10" StructArrays = "0.6.11" -SymbolicIndexingInterface = "0.3.20" +SymbolicIndexingInterface = "0.3.23" Tables = "1.11" Test = "1" Tracker = "0.2.15" From f41afbbd6ed7cdd16e32025f928ac29976671d97 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 14 Jun 2024 09:28:31 -0400 Subject: [PATCH 05/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index bb434c6e..dcd25c4a 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.23.1" +version = "3.24.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 7d912cb3f243dd7ae7034b6b741f94f77813a57c Mon Sep 17 00:00:00 2001 From: Bagaev Dmitry Date: Thu, 27 Jun 2024 16:35:19 +0200 Subject: [PATCH 06/94] fix any --- src/array_partition.jl | 4 ++-- test/partitions_test.jl | 20 ++++++++++++++++++++ 2 files changed, 22 insertions(+), 2 deletions(-) diff --git a/src/array_partition.jl b/src/array_partition.jl index b4047d3b..5b3a1c09 100644 --- a/src/array_partition.jl +++ b/src/array_partition.jl @@ -169,8 +169,8 @@ function Base.mapreduce(f, op, A::ArrayPartition{T}; kwargs...) where {T} mapreduce(f, op, (i for i in A); kwargs...) end Base.filter(f, A::ArrayPartition) = ArrayPartition(map(x -> filter(f, x), A.x)) -Base.any(f, A::ArrayPartition) = any(f, (any(f, x) for x in A.x)) -Base.any(f::Function, A::ArrayPartition) = any(f, (any(f, x) for x in A.x)) +Base.any(f, A::ArrayPartition) = any((any(f, x) for x in A.x)) +Base.any(f::Function, A::ArrayPartition) = any((any(f, x) for x in A.x)) Base.any(A::ArrayPartition) = any(identity, A) Base.all(f, A::ArrayPartition) = all(f, (all(f, x) for x in A.x)) Base.all(f::Function, A::ArrayPartition) = all((all(f, x) for x in A.x)) diff --git a/test/partitions_test.jl b/test/partitions_test.jl index e31e9da5..1b0cf3f8 100644 --- a/test/partitions_test.jl +++ b/test/partitions_test.jl @@ -136,6 +136,26 @@ y = ArrayPartition(ArrayPartition([1], [2.0]), ArrayPartition([3], [4.0])) @inferred mapreduce(string, *, x) @test mapreduce(i -> string(i) * "q", *, x) == "1q2q3.0q4.0q" +# any +@test !any(isnan, ArrayPartition([1, 2], [3.0, 4.0])) +@test !any(isnan, ArrayPartition([3.0, 4.0])) +@test any(isnan, ArrayPartition([NaN], [3.0, 4.0])) +@test any(isnan, ArrayPartition([NaN])) +@test any(isnan, ArrayPartition(ArrayPartition([NaN]))) +@test any(isnan, ArrayPartition([2], [NaN])) +@test any(isnan, ArrayPartition([2], ArrayPartition([NaN]))) + +# all +@test !all(isnan, ArrayPartition([1, 2], [3.0, 4.0])) +@test !all(isnan, ArrayPartition([3.0, 4.0])) +@test !all(isnan, ArrayPartition([NaN], [3.0, 4.0])) +@test all(isnan, ArrayPartition([NaN])) +@test all(isnan, ArrayPartition(ArrayPartition([NaN]))) +@test !all(isnan, ArrayPartition([2], [NaN])) +@test all(isnan, ArrayPartition([NaN], [NaN])) +@test all(isnan, ArrayPartition([NaN], ArrayPartition([NaN]))) + + # broadcasting _scalar_op(y) = y + 1 # Can't do `@inferred(_scalar_op.(x))` so we wrap that in a function: From 7517384593a3e3d933ed648fa505b255860062a7 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Thu, 27 Jun 2024 21:15:37 +0530 Subject: [PATCH 07/94] Revert "feat: add parameter timeseries support to `AbstractDiffEqArray`" --- Project.toml | 2 +- src/vector_of_array.jl | 88 +++++++++++++++++++++--------------------- 2 files changed, 45 insertions(+), 45 deletions(-) diff --git a/Project.toml b/Project.toml index dcd25c4a..b8b9127b 100644 --- a/Project.toml +++ b/Project.toml @@ -59,7 +59,7 @@ StaticArrays = "1.6" StaticArraysCore = "1.4" Statistics = "1.10" StructArrays = "0.6.11" -SymbolicIndexingInterface = "0.3.23" +SymbolicIndexingInterface = "0.3.20" Tables = "1.11" Test = "1" Tracker = "0.2.15" diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 36f95ed2..2c13a573 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -60,13 +60,11 @@ A[1, :] # all time periods for f(t) A.t ``` """ -mutable struct DiffEqArray{T, N, A, B, F, S, D <: Union{Nothing, ParameterTimeseriesCollection}} <: - AbstractDiffEqArray{T, N, A} +mutable struct DiffEqArray{T, N, A, B, F, S} <: AbstractDiffEqArray{T, N, A} u::A # A <: AbstractVector{<: AbstractArray{T, N - 1}} t::B p::F sys::S - discretes::D end ### Abstract Interface struct AllObserved @@ -176,32 +174,29 @@ function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p = nothing, - sys = nothing; discretes = nothing) where {T, N} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, + sys = nothing) where {T, N} + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys)}(vec, ts, p, - sys, - discretes) + sys) end # ambiguity resolution function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, ::NTuple{N, Int}) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing, Nothing}(vec, + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing}(vec, ts, nothing, - nothing, nothing) end function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, - ::NTuple{N, Int}, p; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, + ::NTuple{N, Int}, p) where {T, N, VT <: AbstractArray{T, N}} + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing}(vec, ts, p, - nothing, - discretes) + nothing) end # Assume that the first element is representative of all other elements @@ -211,8 +206,7 @@ function DiffEqArray(vec::AbstractVector, sys = nothing; variables = nothing, parameters = nothing, - independent_variables = nothing, - discretes = nothing) + independent_variables = nothing) sys = something(sys, SymbolCache(something(variables, []), something(parameters, []), @@ -225,13 +219,11 @@ function DiffEqArray(vec::AbstractVector, typeof(vec), typeof(ts), typeof(p), - typeof(sys), - typeof(discretes) + typeof(sys) }(vec, ts, p, - sys, - discretes) + sys) end function DiffEqArray(vec::AbstractVector{VT}, @@ -240,8 +232,7 @@ function DiffEqArray(vec::AbstractVector{VT}, sys = nothing; variables = nothing, parameters = nothing, - independent_variables = nothing, - discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} + independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} sys = something(sys, SymbolCache(something(variables, []), something(parameters, []), @@ -252,30 +243,18 @@ function DiffEqArray(vec::AbstractVector{VT}, typeof(vec), typeof(ts), typeof(p), - typeof(sys), - typeof(discretes), + typeof(sys) }(vec, ts, p, - sys, - discretes) + sys) end -has_discretes(::T) where {T <: AbstractDiffEqArray} = hasfield(T, :discretes) -get_discretes(x) = getfield(x, :discretes) - SymbolicIndexingInterface.is_timeseries(::Type{<:AbstractVectorOfArray}) = Timeseries() -function SymbolicIndexingInterface.is_parameter_timeseries(::Type{DiffEqArray{T, N, A, B, - F, S, D}}) where {T, N, A, B, F, S, D <: ParameterIndexingProxy} - Timeseries() -end SymbolicIndexingInterface.state_values(A::AbstractDiffEqArray) = A.u SymbolicIndexingInterface.current_time(A::AbstractDiffEqArray) = A.t SymbolicIndexingInterface.parameter_values(A::AbstractDiffEqArray) = A.p SymbolicIndexingInterface.symbolic_container(A::AbstractDiffEqArray) = A.sys -function SymbolicIndexingInterface.get_parameter_timeseries_collection(A::AbstractDiffEqArray) - return get_discretes(A) -end Base.IndexStyle(A::AbstractVectorOfArray) = Base.IndexStyle(typeof(A)) Base.IndexStyle(::Type{<:AbstractVectorOfArray}) = IndexCartesian() @@ -384,18 +363,39 @@ end # Symbolic Indexing Methods for (symtype, elsymtype, valtype, errcheck) in [ - (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), - (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), + (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym))), + (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym))), (NotSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Union{<:Tuple, <:AbstractArray}, - :(all(x -> is_parameter(A, x) && !is_timeseries_parameter(A, x), sym))), + :(all(x -> is_parameter(A, x), sym))), ] - @eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, arg...) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym)(A, arg...) +@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, + ::$elsymtype, sym::$valtype) + if $errcheck + throw(ParameterIndexingError(sym)) end + getu(A, sym)(A) +end +@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, + ::$elsymtype, sym::$valtype, arg) + if $errcheck + throw(ParameterIndexingError(sym)) + end + getu(A, sym)(A, arg) +end +@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, + ::$elsymtype, sym::$valtype, arg::Union{AbstractArray{Int}, AbstractArray{Bool}}) + if $errcheck + throw(ParameterIndexingError(sym)) + end + getu(A, sym).((A,), arg) +end +@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, + ::$elsymtype, sym::$valtype, ::Colon) + if $errcheck + throw(ParameterIndexingError(sym)) + end + getu(A, sym)(A) +end end Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::ScalarSymbolic, From 91e05cda4da4e7d5b82f3775f70c240243a4d768 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 29 Jun 2024 00:01:18 +0800 Subject: [PATCH 08/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index b8b9127b..53ef6bf6 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.24.0" +version = "3.25.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From eabcbf1be5fa9b5529e26e3246a1ec84acbadb0a Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Fri, 28 Jun 2024 11:23:40 +0530 Subject: [PATCH 09/94] Revert "Revert "feat: add parameter timeseries support to `AbstractDiffEqArray`"" --- Project.toml | 2 +- src/vector_of_array.jl | 88 +++++++++++++++++++++--------------------- 2 files changed, 45 insertions(+), 45 deletions(-) diff --git a/Project.toml b/Project.toml index 53ef6bf6..9cbe8cf4 100644 --- a/Project.toml +++ b/Project.toml @@ -59,7 +59,7 @@ StaticArrays = "1.6" StaticArraysCore = "1.4" Statistics = "1.10" StructArrays = "0.6.11" -SymbolicIndexingInterface = "0.3.20" +SymbolicIndexingInterface = "0.3.23" Tables = "1.11" Test = "1" Tracker = "0.2.15" diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 2c13a573..36f95ed2 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -60,11 +60,13 @@ A[1, :] # all time periods for f(t) A.t ``` """ -mutable struct DiffEqArray{T, N, A, B, F, S} <: AbstractDiffEqArray{T, N, A} +mutable struct DiffEqArray{T, N, A, B, F, S, D <: Union{Nothing, ParameterTimeseriesCollection}} <: + AbstractDiffEqArray{T, N, A} u::A # A <: AbstractVector{<: AbstractArray{T, N - 1}} t::B p::F sys::S + discretes::D end ### Abstract Interface struct AllObserved @@ -174,29 +176,32 @@ function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p = nothing, - sys = nothing) where {T, N} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys)}(vec, + sys = nothing; discretes = nothing) where {T, N} + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, ts, p, - sys) + sys, + discretes) end # ambiguity resolution function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, ::NTuple{N, Int}) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing}(vec, + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing, Nothing}(vec, ts, nothing, + nothing, nothing) end function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, - ::NTuple{N, Int}, p) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing}(vec, + ::NTuple{N, Int}, p; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} + DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, ts, p, - nothing) + nothing, + discretes) end # Assume that the first element is representative of all other elements @@ -206,7 +211,8 @@ function DiffEqArray(vec::AbstractVector, sys = nothing; variables = nothing, parameters = nothing, - independent_variables = nothing) + independent_variables = nothing, + discretes = nothing) sys = something(sys, SymbolCache(something(variables, []), something(parameters, []), @@ -219,11 +225,13 @@ function DiffEqArray(vec::AbstractVector, typeof(vec), typeof(ts), typeof(p), - typeof(sys) + typeof(sys), + typeof(discretes) }(vec, ts, p, - sys) + sys, + discretes) end function DiffEqArray(vec::AbstractVector{VT}, @@ -232,7 +240,8 @@ function DiffEqArray(vec::AbstractVector{VT}, sys = nothing; variables = nothing, parameters = nothing, - independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} + independent_variables = nothing, + discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} sys = something(sys, SymbolCache(something(variables, []), something(parameters, []), @@ -243,18 +252,30 @@ function DiffEqArray(vec::AbstractVector{VT}, typeof(vec), typeof(ts), typeof(p), - typeof(sys) + typeof(sys), + typeof(discretes), }(vec, ts, p, - sys) + sys, + discretes) end +has_discretes(::T) where {T <: AbstractDiffEqArray} = hasfield(T, :discretes) +get_discretes(x) = getfield(x, :discretes) + SymbolicIndexingInterface.is_timeseries(::Type{<:AbstractVectorOfArray}) = Timeseries() +function SymbolicIndexingInterface.is_parameter_timeseries(::Type{DiffEqArray{T, N, A, B, + F, S, D}}) where {T, N, A, B, F, S, D <: ParameterIndexingProxy} + Timeseries() +end SymbolicIndexingInterface.state_values(A::AbstractDiffEqArray) = A.u SymbolicIndexingInterface.current_time(A::AbstractDiffEqArray) = A.t SymbolicIndexingInterface.parameter_values(A::AbstractDiffEqArray) = A.p SymbolicIndexingInterface.symbolic_container(A::AbstractDiffEqArray) = A.sys +function SymbolicIndexingInterface.get_parameter_timeseries_collection(A::AbstractDiffEqArray) + return get_discretes(A) +end Base.IndexStyle(A::AbstractVectorOfArray) = Base.IndexStyle(typeof(A)) Base.IndexStyle(::Type{<:AbstractVectorOfArray}) = IndexCartesian() @@ -363,39 +384,18 @@ end # Symbolic Indexing Methods for (symtype, elsymtype, valtype, errcheck) in [ - (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym))), - (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym))), + (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), + (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), (NotSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Union{<:Tuple, <:AbstractArray}, - :(all(x -> is_parameter(A, x), sym))), + :(all(x -> is_parameter(A, x) && !is_timeseries_parameter(A, x), sym))), ] -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym)(A) -end -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, arg) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym)(A, arg) -end -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, arg::Union{AbstractArray{Int}, AbstractArray{Bool}}) - if $errcheck - throw(ParameterIndexingError(sym)) - end - getu(A, sym).((A,), arg) -end -@eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, - ::$elsymtype, sym::$valtype, ::Colon) - if $errcheck - throw(ParameterIndexingError(sym)) + @eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, + ::$elsymtype, sym::$valtype, arg...) + if $errcheck + throw(ParameterIndexingError(sym)) + end + getu(A, sym)(A, arg...) end - getu(A, sym)(A) -end end Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::ScalarSymbolic, From 4292ab56a560f969829153aa55b0f497d696a2f6 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Thu, 4 Jul 2024 14:44:03 +0530 Subject: [PATCH 10/94] build: bump SII compat --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 9cbe8cf4..b71909ad 100644 --- a/Project.toml +++ b/Project.toml @@ -59,7 +59,7 @@ StaticArrays = "1.6" StaticArraysCore = "1.4" Statistics = "1.10" StructArrays = "0.6.11" -SymbolicIndexingInterface = "0.3.23" +SymbolicIndexingInterface = "0.3.25" Tables = "1.11" Test = "1" Tracker = "0.2.15" From def0d0a087ae4daabb2920a5da90c849a23a7a85 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 5 Jul 2024 10:33:50 -0400 Subject: [PATCH 11/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index b71909ad..c6cf5fac 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.25.0" +version = "3.26.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 8acba2b76479be99592e0a582e5380b25c07c29e Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 5 Jul 2024 10:49:06 -0400 Subject: [PATCH 12/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index c6cf5fac..b71909ad 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.26.0" +version = "3.25.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From e0f1b1ddc98671ade78fb4fd89339224e145fb95 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Fri, 5 Jul 2024 14:18:09 -0400 Subject: [PATCH 13/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index b71909ad..c6cf5fac 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.25.0" +version = "3.26.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 8e66b0ae16abdf131d4328d3a51e9abb183c73e8 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 10 Aug 2024 04:57:09 -0400 Subject: [PATCH 14/94] Make sparsearrays into a weak dependency --- Project.toml | 8 +++++--- ext/RecursiveArrayToolsSparseArraysExt.jl | 20 ++++++++++++++++++++ src/RecursiveArrayTools.jl | 1 - 3 files changed, 25 insertions(+), 4 deletions(-) create mode 100644 ext/RecursiveArrayToolsSparseArraysExt.jl diff --git a/Project.toml b/Project.toml index c6cf5fac..6a4a3dc1 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.26.0" +version = "3.27.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" @@ -11,7 +11,6 @@ GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" IteratorInterfaceExtensions = "82899510-4779-5014-852e-03e436cf321d" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" -SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" StaticArraysCore = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5" @@ -22,6 +21,7 @@ FastBroadcast = "7034ab61-46d4-4ed7-9d0f-46aef9175898" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" +SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" @@ -32,6 +32,7 @@ RecursiveArrayToolsForwardDiffExt = "ForwardDiff" RecursiveArrayToolsMeasurementsExt = "Measurements" RecursiveArrayToolsMonteCarloMeasurementsExt = "MonteCarloMeasurements" RecursiveArrayToolsReverseDiffExt = ["ReverseDiff", "Zygote"] +RecursiveArrayToolsSparseArraysExt = ["SparseArrays"] RecursiveArrayToolsTrackerExt = "Tracker" RecursiveArrayToolsZygoteExt = "Zygote" @@ -78,6 +79,7 @@ OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f" +SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" StructArrays = "09ab397b-f2b6-538f-b94a-2f83cf4a842a" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" @@ -86,4 +88,4 @@ Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" [targets] -test = ["SafeTestsets", "Aqua", "FastBroadcast", "ForwardDiff", "NLsolve", "OrdinaryDiffEq", "Pkg", "Test", "Unitful", "Random", "StaticArrays", "StructArrays", "Zygote", "Measurements"] +test = ["SafeTestsets", "Aqua", "FastBroadcast", "SparseArrays", "ForwardDiff", "NLsolve", "OrdinaryDiffEq", "Pkg", "Test", "Unitful", "Random", "StaticArrays", "StructArrays", "Zygote", "Measurements"] diff --git a/ext/RecursiveArrayToolsSparseArraysExt.jl b/ext/RecursiveArrayToolsSparseArraysExt.jl new file mode 100644 index 00000000..48e196ee --- /dev/null +++ b/ext/RecursiveArrayToolsSparseArraysExt.jl @@ -0,0 +1,20 @@ +module RecursiveArrayToolsSparseArrays + +import SparseArrays +import RecursiveArrayTools + +function Base.copyto!(dest::SparseArrays.AbstractCompressedVector, A::RecursiveArrayToolsArrayPartition) + @assert length(dest) == length(A) + cur = 1 + @inbounds for i in 1:length(A.x) + if A.x[i] isa Number + dest[cur:(cur + length(A.x[i]) - 1)] .= A.x[i] + else + dest[cur:(cur + length(A.x[i]) - 1)] .= vec(A.x[i]) + end + cur += length(A.x[i]) + end + dest +end + +end \ No newline at end of file diff --git a/src/RecursiveArrayTools.jl b/src/RecursiveArrayTools.jl index 71cc01f9..9e251ac0 100644 --- a/src/RecursiveArrayTools.jl +++ b/src/RecursiveArrayTools.jl @@ -8,7 +8,6 @@ using DocStringExtensions using RecipesBase, StaticArraysCore, Statistics, ArrayInterface, LinearAlgebra using SymbolicIndexingInterface -using SparseArrays import Adapt From 1e91f0e49dcce036ea9d075c4921df2ccac00ea5 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 10 Aug 2024 05:17:12 -0400 Subject: [PATCH 15/94] don't forget to remove --- src/array_partition.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/array_partition.jl b/src/array_partition.jl index 5b3a1c09..7bd8bb52 100644 --- a/src/array_partition.jl +++ b/src/array_partition.jl @@ -176,7 +176,7 @@ Base.all(f, A::ArrayPartition) = all(f, (all(f, x) for x in A.x)) Base.all(f::Function, A::ArrayPartition) = all((all(f, x) for x in A.x)) Base.all(A::ArrayPartition) = all(identity, A) -for type in [AbstractArray, SparseArrays.AbstractCompressedVector, PermutedDimsArray] +for type in [AbstractArray, PermutedDimsArray] @eval function Base.copyto!(dest::$(type), A::ArrayPartition) @assert length(dest) == length(A) cur = 1 From 11837c0b4df8b71c435002a29eeb7b16959e4889 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 10 Aug 2024 05:58:57 -0400 Subject: [PATCH 16/94] fix ext? --- ext/RecursiveArrayToolsSparseArraysExt.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ext/RecursiveArrayToolsSparseArraysExt.jl b/ext/RecursiveArrayToolsSparseArraysExt.jl index 48e196ee..3a40a70d 100644 --- a/ext/RecursiveArrayToolsSparseArraysExt.jl +++ b/ext/RecursiveArrayToolsSparseArraysExt.jl @@ -1,4 +1,4 @@ -module RecursiveArrayToolsSparseArrays +module RecursiveArrayToolsSparseArraysExt import SparseArrays import RecursiveArrayTools From a3a9d32c598f11ba54e185242a391383479f9f45 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 10 Aug 2024 10:22:05 -0400 Subject: [PATCH 17/94] Update ext/RecursiveArrayToolsSparseArraysExt.jl --- ext/RecursiveArrayToolsSparseArraysExt.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ext/RecursiveArrayToolsSparseArraysExt.jl b/ext/RecursiveArrayToolsSparseArraysExt.jl index 3a40a70d..f9261bac 100644 --- a/ext/RecursiveArrayToolsSparseArraysExt.jl +++ b/ext/RecursiveArrayToolsSparseArraysExt.jl @@ -3,7 +3,7 @@ module RecursiveArrayToolsSparseArraysExt import SparseArrays import RecursiveArrayTools -function Base.copyto!(dest::SparseArrays.AbstractCompressedVector, A::RecursiveArrayToolsArrayPartition) +function Base.copyto!(dest::SparseArrays.AbstractCompressedVector, A::RecursiveArrayTools.ArrayPartition) @assert length(dest) == length(A) cur = 1 @inbounds for i in 1:length(A.x) From 1af01c93149b14d1c3dab21cf5c24b8bdbcf6352 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Wed, 11 Sep 2024 22:14:48 +0800 Subject: [PATCH 18/94] Fix similar for VectorOfArray Signed-off-by: ErikQQY <2283984853@qq.com> --- ext/RecursiveArrayToolsReverseDiffExt.jl | 3 +- ext/RecursiveArrayToolsSparseArraysExt.jl | 5 +- ext/RecursiveArrayToolsZygoteExt.jl | 1 - src/vector_of_array.jl | 60 ++++++++++++++--------- test/basic_indexing.jl | 4 +- test/interface_tests.jl | 4 ++ test/partitions_test.jl | 7 ++- 7 files changed, 50 insertions(+), 34 deletions(-) diff --git a/ext/RecursiveArrayToolsReverseDiffExt.jl b/ext/RecursiveArrayToolsReverseDiffExt.jl index d551992d..c82d1c79 100644 --- a/ext/RecursiveArrayToolsReverseDiffExt.jl +++ b/ext/RecursiveArrayToolsReverseDiffExt.jl @@ -24,7 +24,8 @@ end return Array(VA), Array_adjoint end -@adjoint function Base.view(A::AbstractVectorOfArray{<:ReverseDiff.TrackedReal, N}, I::Colon...) where {N} +@adjoint function Base.view( + A::AbstractVectorOfArray{<:ReverseDiff.TrackedReal, N}, I::Colon...) where {N} view_adjoint = let A = A, I = I function (y) A = recursivecopy(A) diff --git a/ext/RecursiveArrayToolsSparseArraysExt.jl b/ext/RecursiveArrayToolsSparseArraysExt.jl index f9261bac..929f63f5 100644 --- a/ext/RecursiveArrayToolsSparseArraysExt.jl +++ b/ext/RecursiveArrayToolsSparseArraysExt.jl @@ -3,7 +3,8 @@ module RecursiveArrayToolsSparseArraysExt import SparseArrays import RecursiveArrayTools -function Base.copyto!(dest::SparseArrays.AbstractCompressedVector, A::RecursiveArrayTools.ArrayPartition) +function Base.copyto!( + dest::SparseArrays.AbstractCompressedVector, A::RecursiveArrayTools.ArrayPartition) @assert length(dest) == length(A) cur = 1 @inbounds for i in 1:length(A.x) @@ -17,4 +18,4 @@ function Base.copyto!(dest::SparseArrays.AbstractCompressedVector, A::RecursiveA dest end -end \ No newline at end of file +end diff --git a/ext/RecursiveArrayToolsZygoteExt.jl b/ext/RecursiveArrayToolsZygoteExt.jl index 415c7f92..14cb9e87 100644 --- a/ext/RecursiveArrayToolsZygoteExt.jl +++ b/ext/RecursiveArrayToolsZygoteExt.jl @@ -146,7 +146,6 @@ end view(A, I...), view_adjoint end - @adjoint function Broadcast.broadcasted(::typeof(+), x::AbstractVectorOfArray, y::Union{Zygote.Numeric, AbstractVectorOfArray}) broadcast(+, x, y), ȳ -> (nothing, map(x -> Zygote.unbroadcast(x, ȳ), (x, y))...) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 36f95ed2..04c19fd0 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -29,6 +29,7 @@ returns a vector of the series for each component, that is, `A[i,:]` for each `i A plot recipe is provided, which plots the `A[i,:]` series. There is also support for `VectorOfArray` constructed from multi-dimensional arrays + ```julia VectorOfArray(u::AbstractArray{AT}) where {T, N, AT <: AbstractArray{T, N}} ``` @@ -60,8 +61,9 @@ A[1, :] # all time periods for f(t) A.t ``` """ -mutable struct DiffEqArray{T, N, A, B, F, S, D <: Union{Nothing, ParameterTimeseriesCollection}} <: - AbstractDiffEqArray{T, N, A} +mutable struct DiffEqArray{ + T, N, A, B, F, S, D <: Union{Nothing, ParameterTimeseriesCollection}} <: + AbstractDiffEqArray{T, N, A} u::A # A <: AbstractVector{<: AbstractArray{T, N - 1}} t::B p::F @@ -177,7 +179,9 @@ function DiffEqArray(vec::AbstractVector{T}, ::NTuple{N, Int}, p = nothing, sys = nothing; discretes = nothing) where {T, N} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, + DiffEqArray{ + eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}( + vec, ts, p, sys, @@ -197,7 +201,8 @@ end function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, ::NTuple{N, Int}, p; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, + DiffEqArray{ + eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, ts, p, nothing, @@ -253,7 +258,7 @@ function DiffEqArray(vec::AbstractVector{VT}, typeof(ts), typeof(p), typeof(sys), - typeof(discretes), + typeof(discretes) }(vec, ts, p, @@ -375,19 +380,23 @@ Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::NotSymboli end struct ParameterIndexingError <: Exception - sym + sym::Any end function Base.showerror(io::IO, pie::ParameterIndexingError) - print(io, "Indexing with parameters is deprecated. Use `getp(A, $(pie.sym))` for parameter indexing.") + print(io, + "Indexing with parameters is deprecated. Use `getp(A, $(pie.sym))` for parameter indexing.") end # Symbolic Indexing Methods for (symtype, elsymtype, valtype, errcheck) in [ - (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), - (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), - (NotSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Union{<:Tuple, <:AbstractArray}, - :(all(x -> is_parameter(A, x) && !is_timeseries_parameter(A, x), sym))), + (ScalarSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, + :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), + (ArraySymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, Any, + :(is_parameter(A, sym) && !is_timeseries_parameter(A, sym))), + (NotSymbolic, SymbolicIndexingInterface.SymbolicTypeTrait, + Union{<:Tuple, <:AbstractArray}, + :(all(x -> is_parameter(A, x) && !is_timeseries_parameter(A, x), sym))) ] @eval Base.@propagate_inbounds function _getindex(A::AbstractDiffEqArray, ::$symtype, ::$elsymtype, sym::$valtype, arg...) @@ -413,8 +422,9 @@ Base.@propagate_inbounds function Base.getindex(A::AbstractVectorOfArray, _arg, elsymtype = symbolic_type(eltype(_arg)) if symtype == NotSymbolic() && elsymtype == NotSymbolic() - if _arg isa Union{Tuple, AbstractArray} && any(x -> symbolic_type(x) != NotSymbolic(), _arg) - _getindex(A, symtype, elsymtype, _arg, args...) + if _arg isa Union{Tuple, AbstractArray} && + any(x -> symbolic_type(x) != NotSymbolic(), _arg) + _getindex(A, symtype, elsymtype, _arg, args...) else _getindex(A, symtype, _arg, args...) end @@ -707,30 +717,32 @@ end # Tools for creating similar objects Base.eltype(::Type{<:AbstractVectorOfArray{T}}) where {T} = T -# TODO: Is there a better way to do this? + @inline function Base.similar(VA::AbstractVectorOfArray, args...) if args[end] isa Type - return Base.similar(eltype(VA)[], args..., size(VA)) + return Base.similar(VA.u, args..., size(VA)) else - return Base.similar(eltype(VA)[], args...) + return Base.similar(VA.u, args...) end end -@inline function Base.similar(VA::VectorOfArray, ::Type{T} = eltype(VA)) where {T} - VectorOfArray([similar(VA[:, i], T) for i in eachindex(VA.u)]) -end -# for VectorOfArray with multi-dimensional parent arrays of arrays where all elements are the same type function Base.similar(vec::VectorOfArray{ T, N, AT}) where {T, N, AT <: AbstractArray{<:AbstractArray{T}}} return VectorOfArray(similar.(Base.parent(vec))) end -# special-case when the multi-dimensional parent array is just an AbstractVector (call the old method) -function Base.similar(vec::VectorOfArray{ - T, N, AT}) where {T, N, AT <: AbstractVector{<:AbstractArray{T}}} - return Base.similar(vec, eltype(vec)) +@inline function Base.similar(VA::VectorOfArray, ::Type{T}) where {T} + VectorOfArray(similar(VA.u, T)) +end + +@inline function Base.similar(VA::VectorOfArray, dims::N) where {N} + VectorOfArray(similar(VA.u, dims)) end +@inline function Base.similar(VA::VectorOfArray{T, N, AT}, + dims::Tuple) where {T, N, AT <: AbstractArray{<:AbstractArray{T}}} + VectorOfArray(similar(VA.u, dims)) +end # fill! # For DiffEqArray it ignores ts and fills only u diff --git a/test/basic_indexing.jl b/test/basic_indexing.jl index c2053159..a71100b3 100644 --- a/test/basic_indexing.jl +++ b/test/basic_indexing.jl @@ -248,13 +248,13 @@ function foo!(u) end foo!(u_matrix) foo!(u_vector) -@test all(u_matrix .== [3, 10]) +@test all(u_matrix .== [3, 10]) @test all(vec(u_matrix) .≈ vec(u_vector)) # test that, for VectorOfArray with multi-dimensional parent arrays, # broadcast and `similar` preserve the structure of the parent array @test typeof(parent(similar(u_matrix))) == typeof(parent(u_matrix)) -@test typeof(parent((x->x).(u_matrix))) == typeof(parent(u_matrix)) +@test typeof(parent((x -> x).(u_matrix))) == typeof(parent(u_matrix)) # test efficiency num_allocs = @allocations foo!(u_matrix) diff --git a/test/interface_tests.jl b/test/interface_tests.jl index 51b92a30..f370e21d 100644 --- a/test/interface_tests.jl +++ b/test/interface_tests.jl @@ -145,6 +145,10 @@ testva2 = similar(testva) @test typeof(testva2) == typeof(testva) @test size(testva2) == size(testva) +testva3 = similar(testva, 10) +@test typeof(testva3) == typeof(testva) +@test length(testva3) == 10 + # Fill AbstractVectorOfArray and check all testval = 3.0 fill!(testva2, testval) diff --git a/test/partitions_test.jl b/test/partitions_test.jl index 1b0cf3f8..2625941e 100644 --- a/test/partitions_test.jl +++ b/test/partitions_test.jl @@ -155,7 +155,6 @@ y = ArrayPartition(ArrayPartition([1], [2.0]), ArrayPartition([3], [4.0])) @test all(isnan, ArrayPartition([NaN], [NaN])) @test all(isnan, ArrayPartition([NaN], ArrayPartition([NaN]))) - # broadcasting _scalar_op(y) = y + 1 # Can't do `@inferred(_scalar_op.(x))` so we wrap that in a function: @@ -303,7 +302,7 @@ end end @testset "Scalar copyto!" begin - u = [2.0,1.0] - copyto!(u, ArrayPartition(1.0,-1.2)) - @test u == [1.0,-1.2] + u = [2.0, 1.0] + copyto!(u, ArrayPartition(1.0, -1.2)) + @test u == [1.0, -1.2] end From fdd3e66e9dda739e002aee62340918de9fd323e1 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Wed, 11 Sep 2024 23:03:16 +0800 Subject: [PATCH 19/94] Fix similar for a specific type Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 4 ++-- test/basic_indexing.jl | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 04c19fd0..1588d693 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -731,8 +731,8 @@ function Base.similar(vec::VectorOfArray{ return VectorOfArray(similar.(Base.parent(vec))) end -@inline function Base.similar(VA::VectorOfArray, ::Type{T}) where {T} - VectorOfArray(similar(VA.u, T)) +@inline function Base.similar(VA::VectorOfArray{T1, N, AT}, ::Type{T2}) where {T1, T2, N, AT <: AbstractArray{<:AbstractArray{T1}}} + eltype(VA.u) <: Vector ? VectorOfArray(similar(VA.u, Vector{T1})) : VectorOfArray(similar(VA.u, Matrix{T1})) end @inline function Base.similar(VA::VectorOfArray, dims::N) where {N} diff --git a/test/basic_indexing.jl b/test/basic_indexing.jl index a71100b3..4a1cb4a0 100644 --- a/test/basic_indexing.jl +++ b/test/basic_indexing.jl @@ -45,7 +45,7 @@ testvasim = similar(testva) @test size(testvasim) == size(testva) @test eltype(testvasim) == eltype(testva) testvasim = similar(testva, Float32) -@test size(testvasim) == size(testva) +#@test size(testvasim) == size(testva) @test eltype(testvasim) == Float32 testvb = deepcopy(testva) @test testva == testvb == recs From c87bc1cb01352041d7edb057de599cb0ef6ed04a Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Wed, 11 Sep 2024 23:12:50 +0800 Subject: [PATCH 20/94] Convert to the right type in similar Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 1588d693..12b9be6b 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -732,7 +732,7 @@ function Base.similar(vec::VectorOfArray{ end @inline function Base.similar(VA::VectorOfArray{T1, N, AT}, ::Type{T2}) where {T1, T2, N, AT <: AbstractArray{<:AbstractArray{T1}}} - eltype(VA.u) <: Vector ? VectorOfArray(similar(VA.u, Vector{T1})) : VectorOfArray(similar(VA.u, Matrix{T1})) + eltype(VA.u) <: Vector ? VectorOfArray(similar(VA.u, Vector{T2})) : VectorOfArray(similar(VA.u, Matrix{T2})) end @inline function Base.similar(VA::VectorOfArray, dims::N) where {N} From 3047711e8349839f0a4e7b8c05b62a0f9183aed7 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Thu, 12 Sep 2024 00:06:02 +0800 Subject: [PATCH 21/94] Use the original type converter Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 4 ++-- test/basic_indexing.jl | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 12b9be6b..b77d30bb 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -731,8 +731,8 @@ function Base.similar(vec::VectorOfArray{ return VectorOfArray(similar.(Base.parent(vec))) end -@inline function Base.similar(VA::VectorOfArray{T1, N, AT}, ::Type{T2}) where {T1, T2, N, AT <: AbstractArray{<:AbstractArray{T1}}} - eltype(VA.u) <: Vector ? VectorOfArray(similar(VA.u, Vector{T2})) : VectorOfArray(similar(VA.u, Matrix{T2})) +@inline function Base.similar(VA::VectorOfArray, ::Type{T} = eltype(VA)) where {T} + VectorOfArray([similar(VA[:, i], T) for i in eachindex(VA.u)]) end @inline function Base.similar(VA::VectorOfArray, dims::N) where {N} diff --git a/test/basic_indexing.jl b/test/basic_indexing.jl index 4a1cb4a0..a71100b3 100644 --- a/test/basic_indexing.jl +++ b/test/basic_indexing.jl @@ -45,7 +45,7 @@ testvasim = similar(testva) @test size(testvasim) == size(testva) @test eltype(testvasim) == eltype(testva) testvasim = similar(testva, Float32) -#@test size(testvasim) == size(testva) +@test size(testvasim) == size(testva) @test eltype(testvasim) == Float32 testvb = deepcopy(testva) @test testva == testvb == recs From 8ff79f52fea39d957cef0b3934bdca4a020b3f80 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Thu, 12 Sep 2024 01:14:10 +0800 Subject: [PATCH 22/94] Avoid undef error Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index b77d30bb..68c4fe2f 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -736,7 +736,21 @@ end end @inline function Base.similar(VA::VectorOfArray, dims::N) where {N} - VectorOfArray(similar(VA.u, dims)) + l = length(VA) + if dims <= l + VectorOfArray([similar(VA[:, i]) for i in 1:l]) + else + VectorOfArray([[similar(VA[:, i]) for i in 1:l]; [similar(VA.u[end]) for _ in (l+1):dims]]) + end +end + +@inline function Base.similar(VA::VectorOfArray, ::Type{T}, dims::N) where {N, T} + l = length(VA) + if dims <= l + VectorOfArray([similar(VA[:, i], T) for i in 1:l]) + else + VectorOfArray([[similar(VA[:, i], T) for i in 1:l]; [similar(VA.u[end], T) for _ in (l+1):dims]]) + end end @inline function Base.similar(VA::VectorOfArray{T, N, AT}, From c4e0076cdb40765add869ff957689cbb0e7bbfab Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Thu, 12 Sep 2024 01:33:41 +0800 Subject: [PATCH 23/94] Fix similar for VectorOfArray when both type and dims are specified Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 11 +---------- 1 file changed, 1 insertion(+), 10 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 68c4fe2f..a662abf8 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -735,7 +735,7 @@ end VectorOfArray([similar(VA[:, i], T) for i in eachindex(VA.u)]) end -@inline function Base.similar(VA::VectorOfArray, dims::N) where {N} +@inline function Base.similar(VA::VectorOfArray, dims::N) where {N <: Number} l = length(VA) if dims <= l VectorOfArray([similar(VA[:, i]) for i in 1:l]) @@ -744,15 +744,6 @@ end end end -@inline function Base.similar(VA::VectorOfArray, ::Type{T}, dims::N) where {N, T} - l = length(VA) - if dims <= l - VectorOfArray([similar(VA[:, i], T) for i in 1:l]) - else - VectorOfArray([[similar(VA[:, i], T) for i in 1:l]; [similar(VA.u[end], T) for _ in (l+1):dims]]) - end -end - @inline function Base.similar(VA::VectorOfArray{T, N, AT}, dims::Tuple) where {T, N, AT <: AbstractArray{<:AbstractArray{T}}} VectorOfArray(similar(VA.u, dims)) From 02501d62cdd2198e83bb736cbf528ca16269e1d5 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Thu, 12 Sep 2024 12:08:00 +0800 Subject: [PATCH 24/94] Change to broadcast Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index a662abf8..f4336df7 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -546,9 +546,7 @@ end function Base.zero(VA::AbstractVectorOfArray) val = copy(VA) - for i in eachindex(VA.u) - val.u[i] = zero(VA.u[i]) - end + val.u = zero.(VA.u) return val end @@ -732,15 +730,15 @@ function Base.similar(vec::VectorOfArray{ end @inline function Base.similar(VA::VectorOfArray, ::Type{T} = eltype(VA)) where {T} - VectorOfArray([similar(VA[:, i], T) for i in eachindex(VA.u)]) + VectorOfArray(similar.(VA.u, T)) end @inline function Base.similar(VA::VectorOfArray, dims::N) where {N <: Number} l = length(VA) if dims <= l - VectorOfArray([similar(VA[:, i]) for i in 1:l]) + VectorOfArray(similar.(VA.u[1:dims])) else - VectorOfArray([[similar(VA[:, i]) for i in 1:l]; [similar(VA.u[end]) for _ in (l+1):dims]]) + VectorOfArray([similar.(VA.u); [similar(VA.u[end]) for _ in (l + 1):dims]]) end end From e8f1fc380161688d8eb3490fcbcc7686e52e8ff4 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Thu, 12 Sep 2024 13:26:17 +0800 Subject: [PATCH 25/94] Revert rewrap changes Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index f4336df7..506c0a81 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -718,9 +718,9 @@ Base.eltype(::Type{<:AbstractVectorOfArray{T}}) where {T} = T @inline function Base.similar(VA::AbstractVectorOfArray, args...) if args[end] isa Type - return Base.similar(VA.u, args..., size(VA)) + return return Base.similar(eltype(VA)[], args..., size(VA)) else - return Base.similar(VA.u, args...) + return Base.similar(eltype(VA)[], args...) end end From da37b6f474b84ede79e50da5e0022b3480ff3255 Mon Sep 17 00:00:00 2001 From: ErikQQY <2283984853@qq.com> Date: Thu, 12 Sep 2024 15:29:13 +0800 Subject: [PATCH 26/94] Fix similar for specific size Signed-off-by: ErikQQY <2283984853@qq.com> --- src/vector_of_array.jl | 7 +------ 1 file changed, 1 insertion(+), 6 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 506c0a81..9271eead 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -718,7 +718,7 @@ Base.eltype(::Type{<:AbstractVectorOfArray{T}}) where {T} = T @inline function Base.similar(VA::AbstractVectorOfArray, args...) if args[end] isa Type - return return Base.similar(eltype(VA)[], args..., size(VA)) + return Base.similar(eltype(VA)[], args..., size(VA)) else return Base.similar(eltype(VA)[], args...) end @@ -742,11 +742,6 @@ end end end -@inline function Base.similar(VA::VectorOfArray{T, N, AT}, - dims::Tuple) where {T, N, AT <: AbstractArray{<:AbstractArray{T}}} - VectorOfArray(similar(VA.u, dims)) -end - # fill! # For DiffEqArray it ignores ts and fills only u function Base.fill!(VA::AbstractVectorOfArray, x) From 161d3995c2d21452d84ed2340efbdde80f983153 Mon Sep 17 00:00:00 2001 From: Anant Thazhemadam Date: Thu, 17 Oct 2024 10:59:20 +0530 Subject: [PATCH 27/94] ci: test with `1`, `lts` and `pre` versions of julia --- .github/workflows/Tests.yml | 2 ++ 1 file changed, 2 insertions(+) diff --git a/.github/workflows/Tests.yml b/.github/workflows/Tests.yml index a9e4259e..08d0f0a5 100644 --- a/.github/workflows/Tests.yml +++ b/.github/workflows/Tests.yml @@ -24,6 +24,8 @@ jobs: matrix: version: - "1" + - "lts" + - "pre" group: - "Core" - "Downstream" From 765f3783cf06712e28549bdafc00516f5446ad1a Mon Sep 17 00:00:00 2001 From: CompatHelper Julia Date: Thu, 17 Oct 2024 17:33:39 +0000 Subject: [PATCH 28/94] CompatHelper: bump compat for GPUArraysCore to 0.2, (keep existing compat) --- Project.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index 6a4a3dc1..d6510f85 100644 --- a/Project.toml +++ b/Project.toml @@ -21,8 +21,8 @@ FastBroadcast = "7034ab61-46d4-4ed7-9d0f-46aef9175898" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" -SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" +SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" @@ -43,7 +43,7 @@ ArrayInterface = "7.6" DocStringExtensions = "0.9" FastBroadcast = "0.2.8, 0.3" ForwardDiff = "0.10.19" -GPUArraysCore = "0.1.1" +GPUArraysCore = "0.1.1, 0.2" IteratorInterfaceExtensions = "1" LinearAlgebra = "1.10" Measurements = "2.3" From 52f86198f244d1501fe0e9140e860b717ea72fcc Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Oct 2024 10:20:05 +0000 Subject: [PATCH 29/94] Overload issingular for AbstractVectorOfArray Fixes https://github.com/SciML/DifferentialEquations.jl/issues/1050 --- src/vector_of_array.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 9271eead..352d4ae4 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -792,6 +792,9 @@ end @inline Statistics.cor(VA::AbstractVectorOfArray; kwargs...) = cor(Array(VA); kwargs...) @inline Base.adjoint(VA::AbstractVectorOfArray) = Adjoint(VA) +# linear algebra +ArrayInterface.issingular(va::AbstractVectorOfArray) = ArrayInterface.issingular(Matrix(va)) + # make it show just like its data function Base.show(io::IO, m::MIME"text/plain", x::AbstractVectorOfArray) (println(io, summary(x), ':'); show(io, m, x.u)) From 02048c30da23ca19990b3cfbee64dd3580daf5dd Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Oct 2024 14:15:51 +0000 Subject: [PATCH 30/94] add issingular test --- test/linalg.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/test/linalg.jl b/test/linalg.jl index cd778aa5..dddbcd56 100644 --- a/test/linalg.jl +++ b/test/linalg.jl @@ -55,3 +55,5 @@ mat = Array(va) @test size(va') == (size(va', 1), size(va', 2)) == (size(va, 2), size(va, 1)) @test all(va'[i] == mat'[i] for i in eachindex(mat')) @test Array(va') == mat' + +@test !ArrayInterface.issingular(VectorOfArray([rand(2),rand(2)])) From ace1e4cacb8878d2f939d1dafb57eb824fac13b6 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 28 Oct 2024 22:35:18 -0100 Subject: [PATCH 31/94] Fix issingular test --- test/linalg.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/linalg.jl b/test/linalg.jl index dddbcd56..23fedffb 100644 --- a/test/linalg.jl +++ b/test/linalg.jl @@ -1,5 +1,5 @@ using RecursiveArrayTools, Test, Random -using LinearAlgebra +using LinearAlgebra, ArrayInterface n, m = 5, 6 bb = rand(n), rand(m) From 23549bdc7c265ef593303aa24e290bf9afcfe106 Mon Sep 17 00:00:00 2001 From: huiyuxie Date: Sun, 27 Oct 2024 17:10:05 -1000 Subject: [PATCH 32/94] Fix --- src/vector_of_array.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 352d4ae4..a9b72e5f 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -906,7 +906,7 @@ for (type, N_expr) in [ unpacked = unpack_voa(bc, i) arr_type = StaticArraysCore.similar_type(dest[:, i]) dest[:, i] = if length(unpacked) == 1 - fill(copy(unpacked), arr_type) + arr_type(unpacked[1]) else arr_type(unpacked[j] for j in eachindex(unpacked)) end From f362c88f9a199cbdeeb1bd8c300b61a07d95c32a Mon Sep 17 00:00:00 2001 From: huiyuxie Date: Sun, 27 Oct 2024 17:30:36 -1000 Subject: [PATCH 33/94] Add test --- test/copy_static_array_test.jl | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/test/copy_static_array_test.jl b/test/copy_static_array_test.jl index 76b9668b..c0aea1d2 100644 --- a/test/copy_static_array_test.jl +++ b/test/copy_static_array_test.jl @@ -82,3 +82,23 @@ b = recursivecopy(a) @test a[1] == b[1] a[1] *= 2 @test a[1] != b[1] + +# Broadcasting when SVector{1} +a = [SVector(0.0) for _ in 1:2] +a_voa = VectorOfArray(a) +b_voa = copy(a_voa) +a_voa[1] = SVector(1.0) +a_voa[2] = SVector(1.0) +@. b_voa = a_voa +@test b_voa[1] == a_voa[1] +@test b_voa[2] == a_voa[2] + +# Broadcasting when SVector{N} where N > 1 +a = [SVector(0.0, 0.0) for _ in 1:2] +a_voa = VectorOfArray(a) +b_voa = copy(a_voa) +a_voa[1] = SVector(1.0, 1.0) +a_voa[2] = SVector(1.0, 1.0) +@. b_voa = a_voa +@test b_voa[1] == a_voa[1] +@test b_voa[2] == a_voa[2] From baa2af98f698ebd0959c36bc5723719a332c4dc8 Mon Sep 17 00:00:00 2001 From: huiyuxie Date: Mon, 28 Oct 2024 10:46:46 -1000 Subject: [PATCH 34/94] Fix again --- src/vector_of_array.jl | 4 +++- test/copy_static_array_test.jl | 14 +++++++++++++- 2 files changed, 16 insertions(+), 2 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index a9b72e5f..eef8625f 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -905,8 +905,10 @@ for (type, N_expr) in [ else unpacked = unpack_voa(bc, i) arr_type = StaticArraysCore.similar_type(dest[:, i]) - dest[:, i] = if length(unpacked) == 1 + dest[:, i] = if length(unpacked) == 1 && length(dest[:, i]) == 1 arr_type(unpacked[1]) + elseif length(unpacked) == 1 + fill(copy(unpacked), arr_type) else arr_type(unpacked[j] for j in eachindex(unpacked)) end diff --git a/test/copy_static_array_test.jl b/test/copy_static_array_test.jl index c0aea1d2..ffcfa52c 100644 --- a/test/copy_static_array_test.jl +++ b/test/copy_static_array_test.jl @@ -83,7 +83,7 @@ b = recursivecopy(a) a[1] *= 2 @test a[1] != b[1] -# Broadcasting when SVector{1} +# Broadcasting when SVector{N} where N = 1 a = [SVector(0.0) for _ in 1:2] a_voa = VectorOfArray(a) b_voa = copy(a_voa) @@ -93,6 +93,12 @@ a_voa[2] = SVector(1.0) @test b_voa[1] == a_voa[1] @test b_voa[2] == a_voa[2] +a = [SVector(0.0) for _ in 1:2] +a_voa = VectorOfArray(a) +a_voa .= 1.0 +@test a_voa[1] == SVector(1.0) +@test a_voa[2] == SVector(1.0) + # Broadcasting when SVector{N} where N > 1 a = [SVector(0.0, 0.0) for _ in 1:2] a_voa = VectorOfArray(a) @@ -102,3 +108,9 @@ a_voa[2] = SVector(1.0, 1.0) @. b_voa = a_voa @test b_voa[1] == a_voa[1] @test b_voa[2] == a_voa[2] + +a = [SVector(0.0, 0.0) for _ in 1:2] +a_voa = VectorOfArray(a) +a_voa .= 1.0 +@test a_voa[1] == SVector(1.0, 1.0) +@test a_voa[2] == SVector(1.0, 1.0) From b7de81ed8c3b4571703b3bde0379551e73592b6c Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 28 Oct 2024 23:23:23 -0100 Subject: [PATCH 35/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index d6510f85..32e4a821 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.27.0" +version = "3.27.1" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From c207bc0f34b5fb107a030201834ee51e3a9e2a90 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 30 Oct 2024 06:46:48 -0100 Subject: [PATCH 36/94] Simplify VectorOfArray indexing Relying on stack removes the need to `adapt`, which should make the GPUs more efficient. With `stack`, those extra dispatches were unnecessary. They were rarely hit and one had a bug too! So they can just be removed. Seems to fix downstream. --- Project.toml | 2 +- src/vector_of_array.jl | 25 ------------------------- 2 files changed, 1 insertion(+), 26 deletions(-) diff --git a/Project.toml b/Project.toml index 32e4a821..c0b801fb 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.27.1" +version = "3.27.2" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index eef8625f..95a2a930 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -321,32 +321,7 @@ end @deprecate Base.getindex(A::AbstractDiffEqArray, i::Int) Base.getindex(A, :, i) false __parameterless_type(T) = Base.typename(T).wrapper -Base.@propagate_inbounds function _getindex(A::AbstractVectorOfArray{T, N}, - ::NotSymbolic, I::Colon...) where {T, N} - @assert length(I) == ndims(A.u[1]) + 1 - vecs = if N == 1 - A.u - else - vec.(A.u) - end - return Adapt.adapt(__parameterless_type(T), - reshape(reduce(hcat, vecs), size(A.u[1])..., length(A.u))) -end -Base.@propagate_inbounds function _getindex(A::AbstractVectorOfArray{T, N}, - ::NotSymbolic, I::Colon...) where {T <: Number, N} - @assert length(I) == ndims(A.u) - return A.u[I...] -end - -Base.@propagate_inbounds function _getindex(A::AbstractVectorOfArray{T, N}, - ::NotSymbolic, I::AbstractArray{Bool}, - J::Colon...) where {T, N} - @assert length(J) == ndims(A.u[1]) + 1 - ndims(I) - @assert size(I) == size(A)[1:(ndims(A) - length(J))] - return A[ntuple(x -> Colon(), ndims(A))...][I, J...] -end -# Need two of each methods to avoid ambiguities Base.@propagate_inbounds function _getindex( A::AbstractVectorOfArray, ::NotSymbolic, ::Colon, I::Int) A.u[I] From 72ac0ce6f9b6ba9944dcc3c35f4f998673c85449 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Wed, 30 Oct 2024 16:35:45 +0530 Subject: [PATCH 37/94] fix: do not assume `AbstractVectorOfArray` is mutable in `Base.zero` --- src/vector_of_array.jl | 2 +- test/interface_tests.jl | 10 ++++++++++ 2 files changed, 11 insertions(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 95a2a930..fc2fce0b 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -521,7 +521,7 @@ end function Base.zero(VA::AbstractVectorOfArray) val = copy(VA) - val.u = zero.(VA.u) + val.u .= zero.(VA.u) return val end diff --git a/test/interface_tests.jl b/test/interface_tests.jl index f370e21d..a1266082 100644 --- a/test/interface_tests.jl +++ b/test/interface_tests.jl @@ -269,3 +269,13 @@ end f3!(z, zz) @test z == VectorOfArray([fill(4, SVector{2, Float64}), fill(2, SVector{2, Float64})]) @test (@allocated f3!(z, zz)) == 0 + +struct ImmutableVectorOfArray{T, N, A} <: AbstractVectorOfArray{T, N, A} + u::A # A <: AbstractArray{<: AbstractArray{T, N - 1}} +end + +@testset "Base.zero does not assume mutable struct" begin + voa = ImmutableVectorOfArray{Float64, 2, Vector{Vector{Float64}}}([ones(3), 2ones(3)]) + zvoa = zero(voa) + @test zvoa.u[1] == zvoa.u[2] == zeros(3) +end From 8b4f756ead697407d1bf796f57f412c85c682587 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Wed, 30 Oct 2024 16:35:49 +0530 Subject: [PATCH 38/94] refactor: format --- test/linalg.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/linalg.jl b/test/linalg.jl index 23fedffb..33a8c4ec 100644 --- a/test/linalg.jl +++ b/test/linalg.jl @@ -56,4 +56,4 @@ mat = Array(va) @test all(va'[i] == mat'[i] for i in eachindex(mat')) @test Array(va') == mat' -@test !ArrayInterface.issingular(VectorOfArray([rand(2),rand(2)])) +@test !ArrayInterface.issingular(VectorOfArray([rand(2), rand(2)])) From 439354c5e2910e9b21ebbed35e42c7f2b37ecc3d Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Wed, 30 Oct 2024 17:20:02 +0530 Subject: [PATCH 39/94] test: remove working `@test_broken` --- test/partitions_test.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/test/partitions_test.jl b/test/partitions_test.jl index 2625941e..d444a4ff 100644 --- a/test/partitions_test.jl +++ b/test/partitions_test.jl @@ -161,7 +161,6 @@ _scalar_op(y) = y + 1 _broadcast_wrapper(y) = _scalar_op.(y) # Issue #8 @inferred _broadcast_wrapper(x) -@test_broken @inferred _broadcast_wrapper(y) # Testing map @test map(x -> x^2, x) == ArrayPartition(x.x[1] .^ 2, x.x[2] .^ 2) From 64f696414ec602ee151b4fa66e438b894ada3581 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 4 Nov 2024 07:22:17 -0100 Subject: [PATCH 40/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index c0b801fb..10afd1ec 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.27.2" +version = "3.27.3" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From e29918d3c275d7e381130ac6885e4486fed0bc75 Mon Sep 17 00:00:00 2001 From: "dependabot[bot]" <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon, 18 Nov 2024 11:34:34 +0000 Subject: [PATCH 41/94] Bump codecov/codecov-action from 4 to 5 Bumps [codecov/codecov-action](https://github.com/codecov/codecov-action) from 4 to 5. - [Release notes](https://github.com/codecov/codecov-action/releases) - [Changelog](https://github.com/codecov/codecov-action/blob/main/CHANGELOG.md) - [Commits](https://github.com/codecov/codecov-action/compare/v4...v5) --- updated-dependencies: - dependency-name: codecov/codecov-action dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] --- .github/workflows/Downstream.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/Downstream.yml b/.github/workflows/Downstream.yml index 3bbfea5d..bcbe7788 100644 --- a/.github/workflows/Downstream.yml +++ b/.github/workflows/Downstream.yml @@ -62,7 +62,7 @@ jobs: exit(0) # Exit immediately, as a success end - uses: julia-actions/julia-processcoverage@v1 - - uses: codecov/codecov-action@v4 + - uses: codecov/codecov-action@v5 with: file: lcov.info token: ${{ secrets.CODECOV_TOKEN }} From 149083bbd8a15e4c04676fb02f9dbfa23ca720e4 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 18 Nov 2024 10:40:22 -0100 Subject: [PATCH 42/94] Update .github/workflows/Downstream.yml --- .github/workflows/Downstream.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/Downstream.yml b/.github/workflows/Downstream.yml index bcbe7788..4e4647f5 100644 --- a/.github/workflows/Downstream.yml +++ b/.github/workflows/Downstream.yml @@ -64,6 +64,6 @@ jobs: - uses: julia-actions/julia-processcoverage@v1 - uses: codecov/codecov-action@v5 with: - file: lcov.info + files: lcov.info token: ${{ secrets.CODECOV_TOKEN }} fail_ci_if_error: true From 83c990ac98bef49f39107dc0310342e3522794a5 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 19 Nov 2024 21:36:38 -0100 Subject: [PATCH 43/94] Fix StructArray broadcast in VectorOfArray Fixes https://github.com/SciML/RecursiveArrayTools.jl/issues/410 This specializes so that if `u.u` is not a vector, it will convert the broadcast to fix that. I couldn't find a nice generic way to use `map` so the fallback is to build the vector and convert, which seems to not be a big performance issue. For StructArrays, `convert(typeof(x), Vector(x))` fails, and so this case is specialized. --- Project.toml | 2 ++ ext/RecursiveArrayToolsStructArraysExt.jl | 6 ++++++ src/vector_of_array.jl | 17 +++++++++++------ test/copy_static_array_test.jl | 7 +++++++ 4 files changed, 26 insertions(+), 6 deletions(-) create mode 100644 ext/RecursiveArrayToolsStructArraysExt.jl diff --git a/Project.toml b/Project.toml index 10afd1ec..88ce5b46 100644 --- a/Project.toml +++ b/Project.toml @@ -23,6 +23,7 @@ Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" +StructArrays = "09ab397b-f2b6-538f-b94a-2f83cf4a842a" Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" @@ -33,6 +34,7 @@ RecursiveArrayToolsMeasurementsExt = "Measurements" RecursiveArrayToolsMonteCarloMeasurementsExt = "MonteCarloMeasurements" RecursiveArrayToolsReverseDiffExt = ["ReverseDiff", "Zygote"] RecursiveArrayToolsSparseArraysExt = ["SparseArrays"] +RecursiveArrayToolsStructArraysExt = "StructArrays" RecursiveArrayToolsTrackerExt = "Tracker" RecursiveArrayToolsZygoteExt = "Zygote" diff --git a/ext/RecursiveArrayToolsStructArraysExt.jl b/ext/RecursiveArrayToolsStructArraysExt.jl new file mode 100644 index 00000000..b4a5d07c --- /dev/null +++ b/ext/RecursiveArrayToolsStructArraysExt.jl @@ -0,0 +1,6 @@ +module RecursiveArrayToolsStructArraysExt + +import RecursiveArrayTools, StructArrays +RecursiveArrayTools.rewrap(::StructArrays.StructArray, u) = StructArrays.StructArray(u) + +end \ No newline at end of file diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index fc2fce0b..78442274 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -849,28 +849,33 @@ end @inline function Base.copy(bc::Broadcast.Broadcasted{<:VectorOfArrayStyle}) bc = Broadcast.flatten(bc) - parent = find_VoA_parent(bc.args) - if parent isa AbstractVector + u = if parent isa AbstractVector # this is the default behavior in v3.15.0 N = narrays(bc) - return VectorOfArray(map(1:N) do i + map(1:N) do i copy(unpack_voa(bc, i)) - end) + end else # if parent isa AbstractArray - return VectorOfArray(map(enumerate(Iterators.product(axes(parent)...))) do (i, _) + map(enumerate(Iterators.product(axes(parent)...))) do (i, _) copy(unpack_voa(bc, i)) - end) + end end + VectorOfArray(rewrap(parent, u)) end +rewrap(::Array,u) = u +rewrap(parent, u) = convert(typeof(parent), u) + for (type, N_expr) in [ (Broadcast.Broadcasted{<:VectorOfArrayStyle}, :(narrays(bc))), (Broadcast.Broadcasted{<:Broadcast.DefaultArrayStyle}, :(length(dest.u))) ] @eval @inline function Base.copyto!(dest::AbstractVectorOfArray, bc::$type) + @show typeof(dest) + error() bc = Broadcast.flatten(bc) N = $N_expr @inbounds for i in 1:N diff --git a/test/copy_static_array_test.jl b/test/copy_static_array_test.jl index ffcfa52c..d3346593 100644 --- a/test/copy_static_array_test.jl +++ b/test/copy_static_array_test.jl @@ -114,3 +114,10 @@ a_voa = VectorOfArray(a) a_voa .= 1.0 @test a_voa[1] == SVector(1.0, 1.0) @test a_voa[2] == SVector(1.0, 1.0) + +#Broadcast Copy of StructArray +x = StructArray{SVector{2, Float64}}((randn(2), randn(2))) +vx = VectorOfArray(x) +vx2 = copy(vx) .+ 1 +ans = vx .+ vx2 +@test ans.u isa StructArray \ No newline at end of file From dfdf6caec6454619937a11e00d2efeab620cccce Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 20 Nov 2024 09:45:54 -0500 Subject: [PATCH 44/94] Update src/vector_of_array.jl --- src/vector_of_array.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 78442274..eef23eac 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -874,8 +874,6 @@ for (type, N_expr) in [ ] @eval @inline function Base.copyto!(dest::AbstractVectorOfArray, bc::$type) - @show typeof(dest) - error() bc = Broadcast.flatten(bc) N = $N_expr @inbounds for i in 1:N From 98c9f157d655c20824b430e3350701e3819a3326 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 20 Nov 2024 12:23:08 -0500 Subject: [PATCH 45/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 88ce5b46..ce574609 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.27.3" +version = "3.27.4" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 69b86619a522af6b3076bfb71cc843625ecd30b2 Mon Sep 17 00:00:00 2001 From: CompatHelper Julia Date: Wed, 20 Nov 2024 17:42:02 +0000 Subject: [PATCH 46/94] CompatHelper: bump compat for StructArrays in [weakdeps] to 0.7, (keep existing compat) --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index ce574609..130fa233 100644 --- a/Project.toml +++ b/Project.toml @@ -61,7 +61,7 @@ SparseArrays = "1.10" StaticArrays = "1.6" StaticArraysCore = "1.4" Statistics = "1.10" -StructArrays = "0.6.11" +StructArrays = "0.6.11, 0.7" SymbolicIndexingInterface = "0.3.25" Tables = "1.11" Test = "1" From 39ec49410024f647b0b44e8d06df74785f047a42 Mon Sep 17 00:00:00 2001 From: CompatHelper Julia Date: Sat, 4 Jan 2025 08:03:57 +0000 Subject: [PATCH 47/94] CompatHelper: bump compat for Zygote in [weakdeps] to 0.7, (keep existing compat) --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 130fa233..32a5b499 100644 --- a/Project.toml +++ b/Project.toml @@ -67,7 +67,7 @@ Tables = "1.11" Test = "1" Tracker = "0.2.15" Unitful = "1" -Zygote = "0.6.67" +Zygote = "0.6.67, 0.7" julia = "1.10" [extras] From 987f8357a057299b5b942a66c3a78a56929acf22 Mon Sep 17 00:00:00 2001 From: anicusan Date: Mon, 27 Jan 2025 13:02:41 +0000 Subject: [PATCH 48/94] Added unrolled implementation of recursivefill! which works on GPUs and avoids recomputing global indices for each setindex! --- src/array_partition.jl | 7 +++++++ src/utils.jl | 5 +++++ 2 files changed, 12 insertions(+) diff --git a/src/array_partition.jl b/src/array_partition.jl index 7bd8bb52..4ee1a9ea 100644 --- a/src/array_partition.jl +++ b/src/array_partition.jl @@ -209,6 +209,13 @@ function Base.copyto!(A::ArrayPartition, src::ArrayPartition) A end +function recursivefill!(b::ArrayPartition, a::T2) where {T2 <: Union{Number, Bool}} + unrolled_foreach!(b.x) do x + fill!(x, a) + end +end + + ## indexing # Interface for the linear indexing. This is just a view of the underlying nested structure diff --git a/src/utils.jl b/src/utils.jl index ae4f1d2f..58945065 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -1,3 +1,7 @@ +unrolled_foreach!(f, t::Tuple) = (f(t[1]); unrolled_foreach!(f, Base.tail(t))) +unrolled_foreach!(f, ::Tuple{}) = nothing + + """ ```julia recursivecopy(a::Union{AbstractArray{T, N}, AbstractVectorOfArray{T, N}}) @@ -127,6 +131,7 @@ function recursivefill!(bs::AbstractVectorOfArray{T, N}, end end + for type in [AbstractArray, AbstractVectorOfArray] @eval function recursivefill!(b::$type{T, N}, a::T2) where {T <: Enum, T2 <: Enum, N} fill!(b, a) From 3a2b75fb87b5dc8fa8c90c254cf04d336a22bd4c Mon Sep 17 00:00:00 2001 From: anicusan Date: Mon, 27 Jan 2025 13:17:47 +0000 Subject: [PATCH 49/94] Added GPU tests for ArrayPartition --- test/gpu/arraypartition_gpu.jl | 16 ++++++++++++++++ test/runtests.jl | 1 + 2 files changed, 17 insertions(+) create mode 100644 test/gpu/arraypartition_gpu.jl diff --git a/test/gpu/arraypartition_gpu.jl b/test/gpu/arraypartition_gpu.jl new file mode 100644 index 00000000..c9a87dc8 --- /dev/null +++ b/test/gpu/arraypartition_gpu.jl @@ -0,0 +1,16 @@ +using RecursiveArrayTools, CUDA, Test +CUDA.allowscalar(false) + + +# Test indexing with colon +a = (CUDA.zeros(5), CUDA.zeros(5)) +pA = ArrayPartition(a) +pA[:, :] + +# Indexing with boolean masks does not work yet +mask = pA .> 0 +# pA[mask] + +# Test recursive filling is done using GPU kernels and not scalar indexing +RecursiveArrayTools.recursivefill!(pA, true) +@test all(pA .== true) diff --git a/test/runtests.jl b/test/runtests.jl index 819e40f3..4ec9d6f4 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -54,5 +54,6 @@ end if GROUP == "GPU" activate_gpu_env() @time @safetestset "VectorOfArray GPU" include("gpu/vectorofarray_gpu.jl") + @time @safetestset "ArrayPartition GPU" include("gpu/arraypartition_gpu.jl") end end From 81847679811a385ab4c91d2aff18c05758484b32 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 28 Jan 2025 05:54:20 -0500 Subject: [PATCH 50/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 32a5b499..4a4b9322 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.27.4" +version = "3.28.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 6bf7697de0c71513439b8c991c64bd0b9fc16f51 Mon Sep 17 00:00:00 2001 From: Erik Faulhaber <44124897+efaulhaber@users.noreply.github.com> Date: Tue, 28 Jan 2025 18:12:53 +0100 Subject: [PATCH 51/94] Fix `fill!(::ArrayPartition)` --- src/array_partition.jl | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/array_partition.jl b/src/array_partition.jl index 7bd8bb52..ce93d555 100644 --- a/src/array_partition.jl +++ b/src/array_partition.jl @@ -209,6 +209,13 @@ function Base.copyto!(A::ArrayPartition, src::ArrayPartition) A end +function Base.fill!(A::ArrayPartition, x) + for i in eachindex(A.x) + fill!(A.x[i], x) + end + A +end + ## indexing # Interface for the linear indexing. This is just a view of the underlying nested structure From fa3ffc2b9ad6e4ab87b03e11ff361a6a87159e03 Mon Sep 17 00:00:00 2001 From: Erik Faulhaber <44124897+efaulhaber@users.noreply.github.com> Date: Wed, 29 Jan 2025 11:30:53 +0100 Subject: [PATCH 52/94] Use for type stability --- src/array_partition.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/array_partition.jl b/src/array_partition.jl index 4ad46df2..b0325fe0 100644 --- a/src/array_partition.jl +++ b/src/array_partition.jl @@ -210,8 +210,8 @@ function Base.copyto!(A::ArrayPartition, src::ArrayPartition) end function Base.fill!(A::ArrayPartition, x) - for i in eachindex(A.x) - fill!(A.x[i], x) + unrolled_foreach!(A.x) do x_ + fill!(x_, x) end A end @@ -221,7 +221,7 @@ function recursivefill!(b::ArrayPartition, a::T2) where {T2 <: Union{Number, Boo fill!(x, a) end end - + ## indexing # Interface for the linear indexing. This is just a view of the underlying nested structure From 28aff8208e0df4efe7d13a21a6fb2147c1cceefc Mon Sep 17 00:00:00 2001 From: Erik Faulhaber <44124897+efaulhaber@users.noreply.github.com> Date: Wed, 29 Jan 2025 11:32:10 +0100 Subject: [PATCH 53/94] Add test --- test/gpu/arraypartition_gpu.jl | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/test/gpu/arraypartition_gpu.jl b/test/gpu/arraypartition_gpu.jl index c9a87dc8..3b335855 100644 --- a/test/gpu/arraypartition_gpu.jl +++ b/test/gpu/arraypartition_gpu.jl @@ -14,3 +14,7 @@ mask = pA .> 0 # Test recursive filling is done using GPU kernels and not scalar indexing RecursiveArrayTools.recursivefill!(pA, true) @test all(pA .== true) + +# Test that regular filling is done using GPU kernels and not scalar indexing +fill!(pA, false) +@test all(pA .== false) From 9454936d1b2f20046f2367e5db022978269790e5 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 29 Jan 2025 05:53:30 -0500 Subject: [PATCH 54/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 4a4b9322..8e5a9c38 100644 --- a/Project.toml +++ b/Project.toml @@ -60,7 +60,7 @@ SafeTestsets = "0.1" SparseArrays = "1.10" StaticArrays = "1.6" StaticArraysCore = "1.4" -Statistics = "1.10" +Statistics = "1.10, 1.11" StructArrays = "0.6.11, 0.7" SymbolicIndexingInterface = "0.3.25" Tables = "1.11" From b9b35d93f8d2bf6b68ea1d109c61027d2e85df21 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Mon, 3 Feb 2025 13:35:09 -0600 Subject: [PATCH 55/94] remove recursive `similar` broadcast --- src/vector_of_array.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index eef23eac..2e31e69e 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -701,7 +701,7 @@ end function Base.similar(vec::VectorOfArray{ T, N, AT}) where {T, N, AT <: AbstractArray{<:AbstractArray{T}}} - return VectorOfArray(similar.(Base.parent(vec))) + return VectorOfArray(similar(Base.parent(vec))) end @inline function Base.similar(VA::VectorOfArray, ::Type{T} = eltype(VA)) where {T} From 3f658973c320abac22f3ad8d629dad0817b3c2f0 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Mon, 3 Feb 2025 14:09:24 -0600 Subject: [PATCH 56/94] restore original Base.similar behavior, specialize on StaticArrays --- src/vector_of_array.jl | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 2e31e69e..f1947de4 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -701,6 +701,12 @@ end function Base.similar(vec::VectorOfArray{ T, N, AT}) where {T, N, AT <: AbstractArray{<:AbstractArray{T}}} + return VectorOfArray(similar.(Base.parent(vec))) +end + +function Base.similar(vec::VectorOfArray{ + T, N, AT}) where {T, N, AT <: AbstractArray{<:StaticArraysCore.StaticVecOrMat{T}}} + # this avoids behavior such as similar(SVector) returning an MVector return VectorOfArray(similar(Base.parent(vec))) end From 2ee17246ed9e17081594a06523555d35cebe3280 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Mon, 3 Feb 2025 14:09:28 -0600 Subject: [PATCH 57/94] add tests --- test/copy_static_array_test.jl | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/test/copy_static_array_test.jl b/test/copy_static_array_test.jl index d3346593..de0132ce 100644 --- a/test/copy_static_array_test.jl +++ b/test/copy_static_array_test.jl @@ -120,4 +120,14 @@ x = StructArray{SVector{2, Float64}}((randn(2), randn(2))) vx = VectorOfArray(x) vx2 = copy(vx) .+ 1 ans = vx .+ vx2 -@test ans.u isa StructArray \ No newline at end of file +@test ans.u isa StructArray + +# check that Base.similar(VectorOfArray{<:StaticArray}) returns the +# same type as the original VectorOfArray +x_staticvector = [SVector(0.0, 0.0) for _ in 1:2] +x_structarray = StructArray{SVector{2, Float64}}((randn(2), randn(2))) +x_mutablefv = [MutableFV(1.0, 2.0)] +x_immutablefv = [ImmutableFV(1.0, 2.0)] +for vec in [x_staticvector, x_structarray, x_mutablefv, x_immutablefv] + @test typeof(similar(VectorOfArray(vec))) === typeof(VectorOfArray(vec)) +end \ No newline at end of file From 7fdbfbae0636b8e82b2237b44af02e78237da9cd Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 4 Feb 2025 07:17:56 -0800 Subject: [PATCH 58/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 8e5a9c38..b1be4604 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.28.0" +version = "3.29.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 20b90c17be238368029219e7b191a76752c244d6 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Thu, 20 Feb 2025 16:35:29 +0530 Subject: [PATCH 59/94] fix: fix DiffEqArray constructors ambiguity --- src/vector_of_array.jl | 164 ++++++++++++++++++++++++++++++----------- 1 file changed, 119 insertions(+), 45 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index f1947de4..2609d9ff 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -174,54 +174,77 @@ end Base.parent(vec::VectorOfArray) = vec.u -function DiffEqArray(vec::AbstractVector{T}, - ts::AbstractVector, - ::NTuple{N, Int}, - p = nothing, - sys = nothing; discretes = nothing) where {T, N} - DiffEqArray{ - eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}( - vec, +#### 2-argument + +# first element representative +function DiffEqArray(vec::AbstractVector, ts::AbstractVector; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) + sys = SymbolCache(something(variables, []), + something(parameters, []), + something(independent_variables, [])) + _size = size(vec[1]) + T = eltype(vec[1]) + return DiffEqArray{ + T, + length(_size) + 1, + typeof(vec), + typeof(ts), + Nothing, + typeof(sys), + typeof(discretes) + }(vec, ts, - p, + nothing, + sys, + discretes) +end + +# T and N from type +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} + sys = SymbolCache(something(variables, []), + something(parameters, []), + something(independent_variables, [])) + return DiffEqArray{ + eltype(eltype(vec)), + N + 1, + typeof(vec), + typeof(ts), + Nothing, + typeof(sys), + typeof(discretes) + }(vec, + ts, + nothing, sys, discretes) end -# ambiguity resolution -function DiffEqArray(vec::AbstractVector{VT}, - ts::AbstractVector, - ::NTuple{N, Int}) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing, Nothing}(vec, +#### 3-argument + +# NTuple, T from type +function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}; discretes = nothing) where {T, N} + DiffEqArray{ + eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing, typeof(discretes)}( + vec, ts, nothing, nothing, - nothing) + discretes) end -function DiffEqArray(vec::AbstractVector{VT}, - ts::AbstractVector, - ::NTuple{N, Int}, p; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{ - eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, + +# NTuple parameter +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, Int}; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}, N2} + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, ts, p, nothing, discretes) end -# Assume that the first element is representative of all other elements -function DiffEqArray(vec::AbstractVector, - ts::AbstractVector, - p = nothing, - sys = nothing; - variables = nothing, - parameters = nothing, - independent_variables = nothing, - discretes = nothing) - sys = something(sys, - SymbolCache(something(variables, []), +# first element representative +function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) + sys = SymbolCache(something(variables, []), something(parameters, []), - something(independent_variables, []))) + something(independent_variables, [])) _size = size(vec[1]) T = eltype(vec[1]) return DiffEqArray{ @@ -239,21 +262,50 @@ function DiffEqArray(vec::AbstractVector, discretes) end -function DiffEqArray(vec::AbstractVector{VT}, - ts::AbstractVector, - p = nothing, - sys = nothing; - variables = nothing, - parameters = nothing, - independent_variables = nothing, - discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} - sys = something(sys, - SymbolCache(something(variables, []), +# T and N from type +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} + sys = SymbolCache(something(variables, []), something(parameters, []), - something(independent_variables, []))) + something(independent_variables, [])) + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, + ts, + p, + sys, + discretes) +end + +#### 4-argument + +# NTuple, T from type +function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p; discretes = nothing) where {T, N} + DiffEqArray{ + eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}( + vec, + ts, + p, + nothing, + discretes) +end + +# NTuple parameter +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, Int}, sys; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}, N2} + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, + ts, + p, + sys, + discretes) +end + +# first element representative +function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p, sys; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) + sys = SymbolCache(something(variables, []), + something(parameters, []), + something(independent_variables, [])) + _size = size(vec[1]) + T = eltype(vec[1]) return DiffEqArray{ - eltype(eltype(vec)), - N + 1, + T, + length(_size) + 1, typeof(vec), typeof(ts), typeof(p), @@ -266,6 +318,28 @@ function DiffEqArray(vec::AbstractVector{VT}, discretes) end +# T and N from type +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p, sys; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, + ts, + p, + sys, + discretes) +end + +#### 5-argument + +# NTuple, T from type +function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p, sys; discretes = nothing) where {T, N} + DiffEqArray{ + eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}( + vec, + ts, + p, + sys, + discretes) +end + has_discretes(::T) where {T <: AbstractDiffEqArray} = hasfield(T, :discretes) get_discretes(x) = getfield(x, :discretes) From 720735d3b0b0506516d43c1c38085916624f60bc Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Thu, 20 Feb 2025 16:37:12 +0530 Subject: [PATCH 60/94] test: test some ambiguous `DiffEqArray` constructors --- test/interface_tests.jl | 19 +++++++++++++++++++ 1 file changed, 19 insertions(+) diff --git a/test/interface_tests.jl b/test/interface_tests.jl index a1266082..674d933c 100644 --- a/test/interface_tests.jl +++ b/test/interface_tests.jl @@ -279,3 +279,22 @@ end zvoa = zero(voa) @test zvoa.u[1] == zvoa.u[2] == zeros(3) end + +@testset "Issue SciMLBase#889: `DiffEqArray` constructor ambiguity" begin + darr = DiffEqArray([[1.0, 1.0]], [1.0], ()) + @test darr.p == () + @test darr.sys === nothing + @test ndims(darr) == 2 + darr = DiffEqArray([[1.0, 1.0]], [1.0], (), "A") + @test darr.p == () + @test darr.sys == "A" + @test ndims(darr) == 2 + darr = DiffEqArray([ones(2, 2)], [1.0], (1, 1, 1)) + @test darr.p == (1, 1, 1) + @test darr.sys === nothing + @test ndims(darr) == 3 + darr = DiffEqArray([ones(2, 2)], [1.0], (1, 1, 1), "A") + @test darr.p == (1, 1, 1) + @test darr.sys == "A" + @test ndims(darr) == 3 +end From d6f139a95b3f0b788d99d6acb08f0f91fb7d19e4 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Thu, 20 Feb 2025 17:13:14 +0530 Subject: [PATCH 61/94] test: fix test to account for MTK's `Initial` parameters --- test/downstream/symbol_indexing.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/downstream/symbol_indexing.jl b/test/downstream/symbol_indexing.jl index 1dc45dd3..557b65be 100644 --- a/test/downstream/symbol_indexing.jl +++ b/test/downstream/symbol_indexing.jl @@ -27,7 +27,7 @@ sol_new = DiffEqArray(sol.u[1:10], @test all(isequal.(all_variable_symbols(sol), all_variable_symbols(sol_new))) @test all(isequal.(all_variable_symbols(sol), [x, RHS])) @test all(isequal.(all_symbols(sol), all_symbols(sol_new))) -@test all(isequal.(all_symbols(sol), [x, RHS, τ, t])) +@test all([any(isequal(sym), all_symbols(sol)) for sym in [x, RHS, τ, t, Initial(x), Initial(RHS)]]) @test sol[solvedvariables] == sol[[x]] @test sol_new[solvedvariables] == sol_new[[x]] @test sol[allvariables] == sol[[x, RHS]] From d2e036072ecf313b9fd2ccc31040510261449bc1 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 22 Feb 2025 00:50:55 -0500 Subject: [PATCH 62/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index b1be4604..3092e971 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.29.0" +version = "3.30.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From abc5b39bef81d402900070ad23ba1e6c56723ce8 Mon Sep 17 00:00:00 2001 From: Aayush Sabharwal Date: Mon, 24 Feb 2025 14:23:24 +0530 Subject: [PATCH 63/94] fix: fix system not being retained in 4-argument constructor --- src/vector_of_array.jl | 5 +---- test/interface_tests.jl | 5 +++++ 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 2609d9ff..cf32d89d 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -297,10 +297,7 @@ function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, end # first element representative -function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p, sys; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) - sys = SymbolCache(something(variables, []), - something(parameters, []), - something(independent_variables, [])) +function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p, sys; discretes = nothing) _size = size(vec[1]) T = eltype(vec[1]) return DiffEqArray{ diff --git a/test/interface_tests.jl b/test/interface_tests.jl index 674d933c..3384848a 100644 --- a/test/interface_tests.jl +++ b/test/interface_tests.jl @@ -298,3 +298,8 @@ end @test darr.sys == "A" @test ndims(darr) == 3 end + +@testset "System retained in 4-argument constructor" begin + darr = DiffEqArray([ones(2)], [1.0], :params, :sys) + @test darr.sys == :sys +end From 451c79de6714b853c4b019deae364bf90ec65e8a Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 25 Feb 2025 05:39:35 -0800 Subject: [PATCH 64/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 3092e971..841aa97e 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.30.0" +version = "3.31.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 90586127782fcef21a20f435a6100420a624586c Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Fri, 14 Mar 2025 15:15:53 -0500 Subject: [PATCH 65/94] add `setindex!` specialization for VectorOfArray{StructArray} --- ext/RecursiveArrayToolsStructArraysExt.jl | 23 ++++++++++++++++++++++- 1 file changed, 22 insertions(+), 1 deletion(-) diff --git a/ext/RecursiveArrayToolsStructArraysExt.jl b/ext/RecursiveArrayToolsStructArraysExt.jl index b4a5d07c..80c8b71f 100644 --- a/ext/RecursiveArrayToolsStructArraysExt.jl +++ b/ext/RecursiveArrayToolsStructArraysExt.jl @@ -3,4 +3,25 @@ module RecursiveArrayToolsStructArraysExt import RecursiveArrayTools, StructArrays RecursiveArrayTools.rewrap(::StructArrays.StructArray, u) = StructArrays.StructArray(u) -end \ No newline at end of file +using RecursiveArrayTools: VectorOfArray +using StructArrays: StructArray + +const VectorOfStructArray{T, N} = VectorOfArray{T, N, <:StructArray} + +# Since `StructArray` lazily materializes struct entries, the general `setindex!(x, val, I)` +# operation `VA.u[I[end]][Base.front(I)...]` will only update a lazily materialized struct +# entry of `u`, but will not actually mutate `x::StructArray`. See the StructArray documentation +# for more details: +# +# https://juliaarrays.github.io/StructArrays.jl/stable/counterintuitive/#Modifying-a-field-of-a-struct-element +# +# To avoid this, we can materialize a struct entry, modify it, and then use `setindex!` +# with the modified struct entry. +function Base.setindex!(VA::VectorOfStructArray{T, N}, v, + I::Int...) where {T, N} + u_I = VA.u[I[end]] + u_I[Base.front(I)...] = v + return VA.u[I[end]] = u_I +end + +end From 107fe3c419b44fd8fa102809ad4baf3c84c4f757 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Fri, 14 Mar 2025 15:15:58 -0500 Subject: [PATCH 66/94] add tests --- test/basic_indexing.jl | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/test/basic_indexing.jl b/test/basic_indexing.jl index a71100b3..37c24857 100644 --- a/test/basic_indexing.jl +++ b/test/basic_indexing.jl @@ -262,3 +262,12 @@ num_allocs = @allocations foo!(u_matrix) # issue 354 @test VectorOfArray(ones(1))[:] == ones(1) + +# check VectorOfArray indexing for a StructArray of mutable structs +using StructArrays +using StaticArrays: MVector +x = VectorOfArray(StructArray{MVector{1, Float64}}(ntuple(_ -> [1.0, 2.0], 1))) + +# check VectorOfArray assignment +x[1, 1] = 10 +@test x[1, 1] == 10 From 9d0bd3507a3c031bf377e8f6d75516384c93deea Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Fri, 14 Mar 2025 17:44:11 -0500 Subject: [PATCH 67/94] specialize broadcasting for StructArrays --- ext/RecursiveArrayToolsStructArraysExt.jl | 37 ++++++++++++++++++++++- 1 file changed, 36 insertions(+), 1 deletion(-) diff --git a/ext/RecursiveArrayToolsStructArraysExt.jl b/ext/RecursiveArrayToolsStructArraysExt.jl index 80c8b71f..a77ca20a 100644 --- a/ext/RecursiveArrayToolsStructArraysExt.jl +++ b/ext/RecursiveArrayToolsStructArraysExt.jl @@ -3,7 +3,8 @@ module RecursiveArrayToolsStructArraysExt import RecursiveArrayTools, StructArrays RecursiveArrayTools.rewrap(::StructArrays.StructArray, u) = StructArrays.StructArray(u) -using RecursiveArrayTools: VectorOfArray +using RecursiveArrayTools: VectorOfArray, VectorOfArrayStyle, ArrayInterface, unpack_voa, + narrays, StaticArraysCore using StructArrays: StructArray const VectorOfStructArray{T, N} = VectorOfArray{T, N, <:StructArray} @@ -17,6 +18,7 @@ const VectorOfStructArray{T, N} = VectorOfArray{T, N, <:StructArray} # # To avoid this, we can materialize a struct entry, modify it, and then use `setindex!` # with the modified struct entry. +# function Base.setindex!(VA::VectorOfStructArray{T, N}, v, I::Int...) where {T, N} u_I = VA.u[I[end]] @@ -24,4 +26,37 @@ function Base.setindex!(VA::VectorOfStructArray{T, N}, v, return VA.u[I[end]] = u_I end +for (type, N_expr) in [ + (Broadcast.Broadcasted{<:VectorOfArrayStyle}, :(narrays(bc))), + (Broadcast.Broadcasted{<:Broadcast.DefaultArrayStyle}, :(length(dest.u))) +] + @eval @inline function Base.copyto!(dest::VectorOfStructArray, + bc::$type) + bc = Broadcast.flatten(bc) + N = $N_expr + @inbounds for i in 1:N + dest_i = dest[:, i] + if dest_i isa AbstractArray + if ArrayInterface.ismutable(dest_i) + copyto!(dest_i, unpack_voa(bc, i)) + else + unpacked = unpack_voa(bc, i) + arr_type = StaticArraysCore.similar_type(dest_i) + dest_i = if length(unpacked) == 1 && length(dest_i) == 1 + arr_type(unpacked[1]) + elseif length(unpacked) == 1 + fill(copy(unpacked), arr_type) + else + arr_type(unpacked[j] for j in eachindex(unpacked)) + end + end + else + dest_i = copy(unpack_voa(bc, i)) + end + dest[:, i] = dest_i + end + dest + end +end + end From 96776b97f95ca75d7f2a98f93e09959718c246a5 Mon Sep 17 00:00:00 2001 From: Jesse Chan Date: Fri, 14 Mar 2025 17:44:20 -0500 Subject: [PATCH 68/94] add tests --- test/basic_indexing.jl | 13 +++++++++++-- 1 file changed, 11 insertions(+), 2 deletions(-) diff --git a/test/basic_indexing.jl b/test/basic_indexing.jl index 37c24857..442e761f 100644 --- a/test/basic_indexing.jl +++ b/test/basic_indexing.jl @@ -265,9 +265,18 @@ num_allocs = @allocations foo!(u_matrix) # check VectorOfArray indexing for a StructArray of mutable structs using StructArrays -using StaticArrays: MVector +using StaticArrays: MVector, SVector x = VectorOfArray(StructArray{MVector{1, Float64}}(ntuple(_ -> [1.0, 2.0], 1))) +y = 2 * x -# check VectorOfArray assignment +# check mutable VectorOfArray assignment and broadcast x[1, 1] = 10 @test x[1, 1] == 10 +@. x = y +@test all(all.(y .== x)) + +# check immutable VectorOfArray broadcast +x = VectorOfArray(StructArray{SVector{1, Float64}}(ntuple(_ -> [1.0, 2.0], 1))) +y = 2 * x +@. x = y +@test all(all.(y .== x)) From 6784255b26a4d083ba422b0ec868b18df89bca85 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Sat, 15 Mar 2025 04:26:16 -0100 Subject: [PATCH 69/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 841aa97e..8b008a92 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.31.0" +version = "3.31.1" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 1654ac058f25126b958e5b16a5dfa92455ab7c48 Mon Sep 17 00:00:00 2001 From: Oscar Smith Date: Thu, 27 Mar 2025 09:27:35 -0400 Subject: [PATCH 70/94] Update ForwardDiff compat --- Project.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index 8b008a92..b07b08f2 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.31.1" +version = "3.31.2" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" @@ -44,7 +44,7 @@ Aqua = "0.8" ArrayInterface = "7.6" DocStringExtensions = "0.9" FastBroadcast = "0.2.8, 0.3" -ForwardDiff = "0.10.19" +ForwardDiff = "0.10.19, 1" GPUArraysCore = "0.1.1, 0.2" IteratorInterfaceExtensions = "1" LinearAlgebra = "1.10" From 59ae7b8b7753b312b9b13f786c771cc45a51ebe3 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Fri, 28 Mar 2025 11:37:54 -0400 Subject: [PATCH 71/94] improve tests --- Project.toml | 6 +++--- test/adjoints.jl | 2 +- test/copy_static_array_test.jl | 28 ++++++++++++++-------------- 3 files changed, 18 insertions(+), 18 deletions(-) diff --git a/Project.toml b/Project.toml index b07b08f2..a4c6c955 100644 --- a/Project.toml +++ b/Project.toml @@ -51,12 +51,12 @@ LinearAlgebra = "1.10" Measurements = "2.3" MonteCarloMeasurements = "1.1" NLsolve = "4.5" -OrdinaryDiffEq = "6.62" Pkg = "1" Random = "1" RecipesBase = "1.1" ReverseDiff = "1.15" SafeTestsets = "0.1" +SciMLBase = "2" SparseArrays = "1.10" StaticArrays = "1.6" StaticArraysCore = "1.4" @@ -77,10 +77,10 @@ ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56" -OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f" +SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" StructArrays = "09ab397b-f2b6-538f-b94a-2f83cf4a842a" @@ -90,4 +90,4 @@ Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" [targets] -test = ["SafeTestsets", "Aqua", "FastBroadcast", "SparseArrays", "ForwardDiff", "NLsolve", "OrdinaryDiffEq", "Pkg", "Test", "Unitful", "Random", "StaticArrays", "StructArrays", "Zygote", "Measurements"] +test = ["SafeTestsets", "Aqua", "FastBroadcast", "SparseArrays", "ForwardDiff", "NLsolve", "Pkg", "Test", "Unitful", "Random", "SciMLBase", "StaticArrays", "StructArrays", "Zygote", "Measurements"] diff --git a/test/adjoints.jl b/test/adjoints.jl index e5a1fc50..2fd48fd1 100644 --- a/test/adjoints.jl +++ b/test/adjoints.jl @@ -1,5 +1,5 @@ using RecursiveArrayTools, Zygote, ForwardDiff, Test -using OrdinaryDiffEq +using SciMLBase function loss(x) sum(abs2, Array(VectorOfArray([x .* i for i in 1:5]))) diff --git a/test/copy_static_array_test.jl b/test/copy_static_array_test.jl index de0132ce..43539dcf 100644 --- a/test/copy_static_array_test.jl +++ b/test/copy_static_array_test.jl @@ -87,33 +87,33 @@ a[1] *= 2 a = [SVector(0.0) for _ in 1:2] a_voa = VectorOfArray(a) b_voa = copy(a_voa) -a_voa[1] = SVector(1.0) -a_voa[2] = SVector(1.0) +a_voa[:, 1] = SVector(1.0) +a_voa[:, 2] = SVector(1.0) @. b_voa = a_voa -@test b_voa[1] == a_voa[1] -@test b_voa[2] == a_voa[2] +@test b_voa[:, 1] == a_voa[1] +@test b_voa[:, 2] == a_voa[2] a = [SVector(0.0) for _ in 1:2] a_voa = VectorOfArray(a) a_voa .= 1.0 -@test a_voa[1] == SVector(1.0) -@test a_voa[2] == SVector(1.0) +@test a_voa[:, 1] == SVector(1.0) +@test a_voa[:, 2] == SVector(1.0) # Broadcasting when SVector{N} where N > 1 a = [SVector(0.0, 0.0) for _ in 1:2] a_voa = VectorOfArray(a) b_voa = copy(a_voa) -a_voa[1] = SVector(1.0, 1.0) -a_voa[2] = SVector(1.0, 1.0) +a_voa[:, 1] = SVector(1.0, 1.0) +a_voa[:, 2] = SVector(1.0, 1.0) @. b_voa = a_voa -@test b_voa[1] == a_voa[1] -@test b_voa[2] == a_voa[2] +@test b_voa[:, 1] == a_voa[:, 1] +@test b_voa[:, 2] == a_voa[:, 2] a = [SVector(0.0, 0.0) for _ in 1:2] a_voa = VectorOfArray(a) a_voa .= 1.0 -@test a_voa[1] == SVector(1.0, 1.0) -@test a_voa[2] == SVector(1.0, 1.0) +@test a_voa[:, 1] == SVector(1.0, 1.0) +@test a_voa[:, 2] == SVector(1.0, 1.0) #Broadcast Copy of StructArray x = StructArray{SVector{2, Float64}}((randn(2), randn(2))) @@ -122,7 +122,7 @@ vx2 = copy(vx) .+ 1 ans = vx .+ vx2 @test ans.u isa StructArray -# check that Base.similar(VectorOfArray{<:StaticArray}) returns the +# check that Base.similar(VectorOfArray{<:StaticArray}) returns the # same type as the original VectorOfArray x_staticvector = [SVector(0.0, 0.0) for _ in 1:2] x_structarray = StructArray{SVector{2, Float64}}((randn(2), randn(2))) @@ -130,4 +130,4 @@ x_mutablefv = [MutableFV(1.0, 2.0)] x_immutablefv = [ImmutableFV(1.0, 2.0)] for vec in [x_staticvector, x_structarray, x_mutablefv, x_immutablefv] @test typeof(similar(VectorOfArray(vec))) === typeof(VectorOfArray(vec)) -end \ No newline at end of file +end From c0468b6e8c2fb06efd2ad1b2dd6c97c286949851 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Fri, 28 Mar 2025 13:12:41 -0400 Subject: [PATCH 72/94] test that downstream packages aren't added --- Project.toml | 2 +- test/qa.jl | 10 +++++++++- 2 files changed, 10 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index a4c6c955..bd393da9 100644 --- a/Project.toml +++ b/Project.toml @@ -90,4 +90,4 @@ Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" [targets] -test = ["SafeTestsets", "Aqua", "FastBroadcast", "SparseArrays", "ForwardDiff", "NLsolve", "Pkg", "Test", "Unitful", "Random", "SciMLBase", "StaticArrays", "StructArrays", "Zygote", "Measurements"] +test = ["Aqua", "FastBroadcast", "ForwardDiff", "Measurements", "NLsolve", "Pkg", "Random", "SafeTestsets", "SciMLBase", "SparseArrays", "StaticArrays", "StructArrays", "Test", "Unitful", "Zygote"] diff --git a/test/qa.jl b/test/qa.jl index 29323ac1..e061435f 100644 --- a/test/qa.jl +++ b/test/qa.jl @@ -1,4 +1,12 @@ -using RecursiveArrayTools, Aqua +using RecursiveArrayTools, Aqua, Pkg + +# yes this is horrible, we'll fix it when Pkg or Base provides a decent API +manifest = Pkg.Types.EnvCache().manifest +# these are good sentinels to test whether someone has added a heavy SciML package to the test deps +if haskey(manifest.deps, "NonlinearSolveBase") || haskey(manifest.deps, "DiffEqBase") + error("Don't put Downstream Packages in non Downstream CI") +end + @testset "Aqua" begin Aqua.find_persistent_tasks_deps(RecursiveArrayTools) Aqua.test_ambiguities(RecursiveArrayTools; recursive = false, broken = true) From e13b491bdc404d1aff3fd44df736ffed979f48fc Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 3 Apr 2025 11:13:26 +0200 Subject: [PATCH 73/94] Use similar methods from ArrayPartion for NamedArrayPartition --- src/named_array_partition.jl | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/named_array_partition.jl b/src/named_array_partition.jl index d0737e38..3fd662e9 100644 --- a/src/named_array_partition.jl +++ b/src/named_array_partition.jl @@ -26,6 +26,11 @@ end # fields except through `getfield` and accessor functions. ArrayPartition(x::NamedArrayPartition) = getfield(x, :array_partition) +function Base.similar(x::NamedArrayPartition, args...) + NamedArrayPartition( + similar(getfield(x, :array_partition), args...), getfield(x, :names_to_indices)) +end + Base.Array(x::NamedArrayPartition) = Array(ArrayPartition(x)) function Base.zero(x::NamedArrayPartition{T, S, TN}) where {T, S, TN} From 6fdcbb2a4c13d8437343a20b28a7db78acb58887 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 3 Apr 2025 11:29:41 +0200 Subject: [PATCH 74/94] Explicitly write out methods --- src/named_array_partition.jl | 29 +++++++++++++++++++++++++++-- test/named_array_partition_tests.jl | 2 ++ 2 files changed, 29 insertions(+), 2 deletions(-) diff --git a/src/named_array_partition.jl b/src/named_array_partition.jl index 3fd662e9..2babcc86 100644 --- a/src/named_array_partition.jl +++ b/src/named_array_partition.jl @@ -26,9 +26,34 @@ end # fields except through `getfield` and accessor functions. ArrayPartition(x::NamedArrayPartition) = getfield(x, :array_partition) -function Base.similar(x::NamedArrayPartition, args...) +function Base.similar(A::NamedArrayPartition{T, S}) where {T, S} + NamedArrayPartition(ArrayPartition{T, S}(similar.(getfield(A, :array_partition))), + getfield(A, :names_to_indices)) +end + +# return ArrayPartition when possible, otherwise next best thing of the correct size +function Base.similar(A::NamedArrayPartition, dims::NTuple{N, Int}) where {N} + NamedArrayPartition( + similar(getfield(A, :array_partition), dims), getfield(A, :names_to_indices)) +end + +# similar array partition of common type +@inline function Base.similar(A::NamedArrayPartition, ::Type{T}) where {T} + NamedArrayPartition( + similar(getfield(A, :array_partition), T), getfield(A, :names_to_indices)) +end + +# return ArrayPartition when possible, otherwise next best thing of the correct size +function Base.similar(A::NamedArrayPartition, ::Type{T}, dims::NTuple{N, Int}) where {T, N} + NamedArrayPartition( + similar(getfield(A, :array_partition), T, dims), getfield(A, :names_to_indices)) +end + +# similar array partition with different types +function Base.similar( + A::NamedArrayPartition, ::Type{T}, ::Type{S}, R::DataType...) where {T, S} NamedArrayPartition( - similar(getfield(x, :array_partition), args...), getfield(x, :names_to_indices)) + similar(getfield(A, :array_partition), T, S, R), getfield(A, :names_to_indices)) end Base.Array(x::NamedArrayPartition) = Array(ArrayPartition(x)) diff --git a/test/named_array_partition_tests.jl b/test/named_array_partition_tests.jl index d8164edf..d5647bad 100644 --- a/test/named_array_partition_tests.jl +++ b/test/named_array_partition_tests.jl @@ -4,6 +4,8 @@ using RecursiveArrayTools, Test x = NamedArrayPartition(a = ones(10), b = rand(20)) @test typeof(@. sin(x * x^2 / x - 1)) <: NamedArrayPartition @test typeof(x .^ 2) <: NamedArrayPartition + @test typeof(similar(x)) <: NamedArrayPartition + @test typeof(similar(x, Int)) <: NamedArrayPartition @test x.a ≈ ones(10) @test typeof(x .+ x[1:end]) <: Vector # test broadcast precedence @test all(x .== x[1:end]) From fe41c3c93976ced3565843e9188190364f363e67 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 3 Apr 2025 11:33:54 +0200 Subject: [PATCH 75/94] Fix --- src/named_array_partition.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/named_array_partition.jl b/src/named_array_partition.jl index 2babcc86..de8fa91a 100644 --- a/src/named_array_partition.jl +++ b/src/named_array_partition.jl @@ -26,9 +26,9 @@ end # fields except through `getfield` and accessor functions. ArrayPartition(x::NamedArrayPartition) = getfield(x, :array_partition) -function Base.similar(A::NamedArrayPartition{T, S}) where {T, S} - NamedArrayPartition(ArrayPartition{T, S}(similar.(getfield(A, :array_partition))), - getfield(A, :names_to_indices)) +function Base.similar(A::NamedArrayPartition) + NamedArrayPartition( + similar(getfield(A, :array_partition)), getfield(A, :names_to_indices)) end # return ArrayPartition when possible, otherwise next best thing of the correct size From 6a66fc8f700b883a2f118f0cb97079f8a81599e5 Mon Sep 17 00:00:00 2001 From: Bart de Koning Date: Thu, 3 Apr 2025 12:17:16 +0200 Subject: [PATCH 76/94] Formatting --- src/utils.jl | 2 - src/vector_of_array.jl | 63 +++++++++++++++++++----------- test/downstream/symbol_indexing.jl | 3 +- test/gpu/arraypartition_gpu.jl | 1 - 4 files changed, 42 insertions(+), 27 deletions(-) diff --git a/src/utils.jl b/src/utils.jl index 58945065..8c34c316 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -1,7 +1,6 @@ unrolled_foreach!(f, t::Tuple) = (f(t[1]); unrolled_foreach!(f, Base.tail(t))) unrolled_foreach!(f, ::Tuple{}) = nothing - """ ```julia recursivecopy(a::Union{AbstractArray{T, N}, AbstractVectorOfArray{T, N}}) @@ -131,7 +130,6 @@ function recursivefill!(bs::AbstractVectorOfArray{T, N}, end end - for type in [AbstractArray, AbstractVectorOfArray] @eval function recursivefill!(b::$type{T, N}, a::T2) where {T <: Enum, T2 <: Enum, N} fill!(b, a) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index cf32d89d..408d5bf5 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -177,10 +177,11 @@ Base.parent(vec::VectorOfArray) = vec.u #### 2-argument # first element representative -function DiffEqArray(vec::AbstractVector, ts::AbstractVector; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) +function DiffEqArray(vec::AbstractVector, ts::AbstractVector; discretes = nothing, + variables = nothing, parameters = nothing, independent_variables = nothing) sys = SymbolCache(something(variables, []), - something(parameters, []), - something(independent_variables, [])) + something(parameters, []), + something(independent_variables, [])) _size = size(vec[1]) T = eltype(vec[1]) return DiffEqArray{ @@ -199,10 +200,12 @@ function DiffEqArray(vec::AbstractVector, ts::AbstractVector; discretes = nothin end # T and N from type -function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector; + discretes = nothing, variables = nothing, parameters = nothing, + independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} sys = SymbolCache(something(variables, []), - something(parameters, []), - something(independent_variables, [])) + something(parameters, []), + something(independent_variables, [])) return DiffEqArray{ eltype(eltype(vec)), N + 1, @@ -221,7 +224,8 @@ end #### 3-argument # NTuple, T from type -function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}; discretes = nothing) where {T, N} +function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, + ::NTuple{N, Int}; discretes = nothing) where {T, N} DiffEqArray{ eltype(T), N, typeof(vec), typeof(ts), Nothing, Nothing, typeof(discretes)}( vec, @@ -232,8 +236,11 @@ function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int end # NTuple parameter -function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, Int}; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}, N2} - DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}(vec, +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, Int}; + discretes = nothing) where {T, N, VT <: AbstractArray{T, N}, N2} + DiffEqArray{ + eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}( + vec, ts, p, nothing, @@ -241,10 +248,11 @@ function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, end # first element representative -function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) +function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p; discretes = nothing, + variables = nothing, parameters = nothing, independent_variables = nothing) sys = SymbolCache(something(variables, []), - something(parameters, []), - something(independent_variables, [])) + something(parameters, []), + something(independent_variables, [])) _size = size(vec[1]) T = eltype(vec[1]) return DiffEqArray{ @@ -263,11 +271,14 @@ function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p; discretes = not end # T and N from type -function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p; discretes = nothing, variables = nothing, parameters = nothing, independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p; + discretes = nothing, variables = nothing, parameters = nothing, + independent_variables = nothing) where {T, N, VT <: AbstractArray{T, N}} sys = SymbolCache(something(variables, []), - something(parameters, []), - something(independent_variables, [])) - DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, + something(parameters, []), + something(independent_variables, [])) + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), + typeof(p), typeof(sys), typeof(discretes)}(vec, ts, p, sys, @@ -277,7 +288,8 @@ end #### 4-argument # NTuple, T from type -function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p; discretes = nothing) where {T, N} +function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, + ::NTuple{N, Int}, p; discretes = nothing) where {T, N} DiffEqArray{ eltype(T), N, typeof(vec), typeof(ts), typeof(p), Nothing, typeof(discretes)}( vec, @@ -288,8 +300,10 @@ function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int end # NTuple parameter -function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, Int}, sys; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}, N2} - DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p::NTuple{N2, Int}, sys; + discretes = nothing) where {T, N, VT <: AbstractArray{T, N}, N2} + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), + typeof(p), typeof(sys), typeof(discretes)}(vec, ts, p, sys, @@ -316,8 +330,10 @@ function DiffEqArray(vec::AbstractVector, ts::AbstractVector, p, sys; discretes end # T and N from type -function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p, sys; discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} - DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}(vec, +function DiffEqArray(vec::AbstractVector{VT}, ts::AbstractVector, p, sys; + discretes = nothing) where {T, N, VT <: AbstractArray{T, N}} + DiffEqArray{eltype(T), N + 1, typeof(vec), typeof(ts), + typeof(p), typeof(sys), typeof(discretes)}(vec, ts, p, sys, @@ -327,7 +343,8 @@ end #### 5-argument # NTuple, T from type -function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, ::NTuple{N, Int}, p, sys; discretes = nothing) where {T, N} +function DiffEqArray(vec::AbstractVector{T}, ts::AbstractVector, + ::NTuple{N, Int}, p, sys; discretes = nothing) where {T, N} DiffEqArray{ eltype(T), N, typeof(vec), typeof(ts), typeof(p), typeof(sys), typeof(discretes)}( vec, @@ -942,7 +959,7 @@ end VectorOfArray(rewrap(parent, u)) end -rewrap(::Array,u) = u +rewrap(::Array, u) = u rewrap(parent, u) = convert(typeof(parent), u) for (type, N_expr) in [ diff --git a/test/downstream/symbol_indexing.jl b/test/downstream/symbol_indexing.jl index 557b65be..5d11fab0 100644 --- a/test/downstream/symbol_indexing.jl +++ b/test/downstream/symbol_indexing.jl @@ -27,7 +27,8 @@ sol_new = DiffEqArray(sol.u[1:10], @test all(isequal.(all_variable_symbols(sol), all_variable_symbols(sol_new))) @test all(isequal.(all_variable_symbols(sol), [x, RHS])) @test all(isequal.(all_symbols(sol), all_symbols(sol_new))) -@test all([any(isequal(sym), all_symbols(sol)) for sym in [x, RHS, τ, t, Initial(x), Initial(RHS)]]) +@test all([any(isequal(sym), all_symbols(sol)) + for sym in [x, RHS, τ, t, Initial(x), Initial(RHS)]]) @test sol[solvedvariables] == sol[[x]] @test sol_new[solvedvariables] == sol_new[[x]] @test sol[allvariables] == sol[[x, RHS]] diff --git a/test/gpu/arraypartition_gpu.jl b/test/gpu/arraypartition_gpu.jl index 3b335855..80f9a8d6 100644 --- a/test/gpu/arraypartition_gpu.jl +++ b/test/gpu/arraypartition_gpu.jl @@ -1,7 +1,6 @@ using RecursiveArrayTools, CUDA, Test CUDA.allowscalar(false) - # Test indexing with colon a = (CUDA.zeros(5), CUDA.zeros(5)) pA = ArrayPartition(a) From 78820ed9f320712dfbfb71fac858d1b882feb0f1 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 3 Apr 2025 09:18:06 -0400 Subject: [PATCH 77/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index bd393da9..9971174f 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.31.2" +version = "3.32.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 2dae2f134cf5f1af3578e90c90a115477dd9443e Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 21 Apr 2025 04:27:03 -0400 Subject: [PATCH 78/94] Recurse adapt Needs a test case to really know if this is right. Thinking of GPUArrays, you can't have a GPUArray of GPUArrays so it should be an Array of GPUArrays --- src/vector_of_array.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 408d5bf5..eb6e095e 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -143,7 +143,7 @@ end Base.convert(::Type{AbstractArray}, VA::AbstractVectorOfArray) = stack(VA.u) function Adapt.adapt_structure(to, VA::AbstractVectorOfArray) - Adapt.adapt(to, Array(VA)) + Adapt.adapt.(to, VA.u) end function VectorOfArray(vec::AbstractVector{T}, ::NTuple{N}) where {T, N} From b34fd8ff3756a0235e34022c3444036e7b5c369a Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 21 Apr 2025 22:33:36 -0400 Subject: [PATCH 79/94] Update vector_of_array.jl --- src/vector_of_array.jl | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index eb6e095e..3fa52087 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -143,7 +143,11 @@ end Base.convert(::Type{AbstractArray}, VA::AbstractVectorOfArray) = stack(VA.u) function Adapt.adapt_structure(to, VA::AbstractVectorOfArray) - Adapt.adapt.(to, VA.u) + VectorOfArray(Adapt.adapt.(to, VA.u)) +end + +function Adapt.adapt_structure(to, VA::AbstractDiffEqArray) + DiffEqArray(Adapt.adapt.(to, VA.u), Adapt.adapt(to, VA.t)) end function VectorOfArray(vec::AbstractVector{T}, ::NTuple{N}) where {T, N} From e812893306e63e8aa7965f9f01cea1b657efeac1 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 21 Apr 2025 22:37:35 -0400 Subject: [PATCH 80/94] Update gpu.jl --- test/gpu.jl | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/test/gpu.jl b/test/gpu.jl index 2ed5211c..88097a98 100644 --- a/test/gpu.jl +++ b/test/gpu.jl @@ -1,4 +1,15 @@ -using RecursiveArrayTools, CuArrays +using RecursiveArrayTools, Adapt, CuArrays, Test a = ArrayPartition(([1.0f0] |> cu, [2.0f0] |> cu, [3.0f0] |> cu)) b = ArrayPartition(([0.0f0] |> cu, [0.0f0] |> cu, [0.0f0] |> cu)) @. a + b + +a = VectorOfArray([ones(2) for i in 1:3]) +_a = Adapt.adapt(CuArray,a) +@test _a isa VectorOfArray +@test _a.u isa Vector{<:CuArray} + +b = DiffEqArray([ones(2) for i in 1:3],ones(2)) +_b = Adapt.adapt(CuArray,b) +@test _b isa DiffEqArray +@test _b.u isa Vector{<:CuArray} +@test _b.t isa CuArray From 2419648cf928003d682822f2fd265c6baf98f85d Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 21 Apr 2025 23:33:43 -0400 Subject: [PATCH 81/94] Update vectorofarray_gpu.jl --- test/gpu/vectorofarray_gpu.jl | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/test/gpu/vectorofarray_gpu.jl b/test/gpu/vectorofarray_gpu.jl index c97beca2..7422a6c1 100644 --- a/test/gpu/vectorofarray_gpu.jl +++ b/test/gpu/vectorofarray_gpu.jl @@ -1,4 +1,4 @@ -using RecursiveArrayTools, CUDA, Test, Zygote +using RecursiveArrayTools, CUDA, Test, Zygote, Adapt CUDA.allowscalar(false) # Test indexing with colon @@ -37,3 +37,14 @@ va_cu = convert(AbstractArray, va) @test va_cu isa CuArray @test size(va_cu) == size(x) + +a = VectorOfArray([ones(2) for i in 1:3]) +_a = Adapt.adapt(CuArray,a) +@test _a isa VectorOfArray +@test _a.u isa Vector{<:CuArray} + +b = DiffEqArray([ones(2) for i in 1:3],ones(2)) +_b = Adapt.adapt(CuArray,b) +@test _b isa DiffEqArray +@test _b.u isa Vector{<:CuArray} +@test _b.t isa CuArray From 72fbbf499d1a8508b91b9813dbae316ee4db82e5 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 21 Apr 2025 23:34:03 -0400 Subject: [PATCH 82/94] Update arraypartition_gpu.jl --- test/gpu/arraypartition_gpu.jl | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/test/gpu/arraypartition_gpu.jl b/test/gpu/arraypartition_gpu.jl index 80f9a8d6..08fdb69a 100644 --- a/test/gpu/arraypartition_gpu.jl +++ b/test/gpu/arraypartition_gpu.jl @@ -17,3 +17,7 @@ RecursiveArrayTools.recursivefill!(pA, true) # Test that regular filling is done using GPU kernels and not scalar indexing fill!(pA, false) @test all(pA .== false) + +a = ArrayPartition(([1.0f0] |> cu, [2.0f0] |> cu, [3.0f0] |> cu)) +b = ArrayPartition(([0.0f0] |> cu, [0.0f0] |> cu, [0.0f0] |> cu)) +@. a + b From 63be2a17f8c325cecd801785a34cc860615e1e21 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 21 Apr 2025 23:34:11 -0400 Subject: [PATCH 83/94] Delete test/gpu.jl --- test/gpu.jl | 15 --------------- 1 file changed, 15 deletions(-) delete mode 100644 test/gpu.jl diff --git a/test/gpu.jl b/test/gpu.jl deleted file mode 100644 index 88097a98..00000000 --- a/test/gpu.jl +++ /dev/null @@ -1,15 +0,0 @@ -using RecursiveArrayTools, Adapt, CuArrays, Test -a = ArrayPartition(([1.0f0] |> cu, [2.0f0] |> cu, [3.0f0] |> cu)) -b = ArrayPartition(([0.0f0] |> cu, [0.0f0] |> cu, [0.0f0] |> cu)) -@. a + b - -a = VectorOfArray([ones(2) for i in 1:3]) -_a = Adapt.adapt(CuArray,a) -@test _a isa VectorOfArray -@test _a.u isa Vector{<:CuArray} - -b = DiffEqArray([ones(2) for i in 1:3],ones(2)) -_b = Adapt.adapt(CuArray,b) -@test _b isa DiffEqArray -@test _b.u isa Vector{<:CuArray} -@test _b.t isa CuArray From 71a6428da90b5c2ecff578711ef7131ef74654eb Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 22 Apr 2025 00:03:18 -0400 Subject: [PATCH 84/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 9971174f..c3dbc132 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.32.0" +version = "3.33.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 65a6d617f198deea4100dfbe71ce0d27bbd83520 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 22 Apr 2025 00:04:43 -0400 Subject: [PATCH 85/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index c3dbc132..9971174f 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.33.0" +version = "3.32.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 1146cbeda2bc58d37bf9cb6001f2ddcbe227da25 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 23 Apr 2025 08:06:26 -0400 Subject: [PATCH 86/94] Make VectorOfArray adapt not broadcast the `to` --- src/vector_of_array.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 3fa52087..29bc3143 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -143,11 +143,11 @@ end Base.convert(::Type{AbstractArray}, VA::AbstractVectorOfArray) = stack(VA.u) function Adapt.adapt_structure(to, VA::AbstractVectorOfArray) - VectorOfArray(Adapt.adapt.(to, VA.u)) + VectorOfArray(Adapt.adapt.((to,), VA.u)) end function Adapt.adapt_structure(to, VA::AbstractDiffEqArray) - DiffEqArray(Adapt.adapt.(to, VA.u), Adapt.adapt(to, VA.t)) + DiffEqArray(Adapt.adapt.((to,), VA.u), Adapt.adapt((to,), VA.t)) end function VectorOfArray(vec::AbstractVector{T}, ::NTuple{N}) where {T, N} From a2a0e41f8f5909dd265dc5d475324c7c0318c384 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 23 Apr 2025 08:18:41 -0400 Subject: [PATCH 87/94] Update src/vector_of_array.jl --- src/vector_of_array.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 29bc3143..7e38fb89 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -147,7 +147,7 @@ function Adapt.adapt_structure(to, VA::AbstractVectorOfArray) end function Adapt.adapt_structure(to, VA::AbstractDiffEqArray) - DiffEqArray(Adapt.adapt.((to,), VA.u), Adapt.adapt((to,), VA.t)) + DiffEqArray(Adapt.adapt.((to,), VA.u), Adapt.adapt(to, VA.t)) end function VectorOfArray(vec::AbstractVector{T}, ::NTuple{N}) where {T, N} From 7ba9b51b718c750f7dd1bc9dae3f5a3e474f6087 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 23 Apr 2025 08:29:01 -0400 Subject: [PATCH 88/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 9971174f..11132fdf 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.32.0" +version = "3.32.1" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From f79173487af7bd55d8de6d8ff9fc54a9838fa194 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Apr 2025 18:25:01 -0400 Subject: [PATCH 89/94] More precise deprecations of linear indexing Noticed in https://github.com/SciML/Catalyst.jl/pull/1248, the linear indexing is actually fine when done on scalar arrays because that matches the linear indexing that you get from an AbstractArray interface. So instead the deprecations are made more precise to be just about the array of array case. --- src/vector_of_array.jl | 34 ++++++++++++++++++++-------------- 1 file changed, 20 insertions(+), 14 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 7e38fb89..1a7d06c6 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -388,29 +388,32 @@ end @inline Base.IteratorSize(::Type{<:AbstractVectorOfArray}) = Base.HasLength() @inline Base.first(VA::AbstractVectorOfArray) = first(VA.u) @inline Base.last(VA::AbstractVectorOfArray) = last(VA.u) -function Base.firstindex(VA::AbstractVectorOfArray) - Base.depwarn( +function Base.firstindex(VA::AbstractVectorOfArray{T,N,A}) where {T,N,A} + N > 1 && Base.depwarn( "Linear indexing of `AbstractVectorOfArray` is deprecated. Change `A[i]` to `A.u[i]` ", :firstindex) return firstindex(VA.u) end -function Base.lastindex(VA::AbstractVectorOfArray) - Base.depwarn( +function Base.lastindex(VA::AbstractVectorOfArray{T,N,A}) where {T,N,A} + N > 1 && Base.depwarn( "Linear indexing of `AbstractVectorOfArray` is deprecated. Change `A[i]` to `A.u[i]` ", :lastindex) return lastindex(VA.u) end -@deprecate Base.getindex(A::AbstractVectorOfArray, I::Int) Base.getindex(A, :, I) false +Base.getindex(A::AbstractVectorOfArray, I::Int) = A.u[I] +Base.getindex(A::AbstractVectorOfArray, I::AbstractArray{Int}) = A.u[I] +Base.getindex(A::AbstractDiffEqArray, I::Int) = A.u[I] +Base.getindex(A::AbstractDiffEqArray, I::AbstractArray{Int}) = A.u[I] -@deprecate Base.getindex(A::AbstractVectorOfArray, I::AbstractArray{Int}) Base.getindex( - A, :, I) false +@deprecate Base.getindex(A::AbstractVectorOfArray{T,N,A}, I::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[I] false -@deprecate Base.getindex(A::AbstractDiffEqArray, I::AbstractArray{Int}) Base.getindex( - A, :, I) false +@deprecate Base.getindex(A::AbstractVectorOfArray{T,N,A}, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[I] false -@deprecate Base.getindex(A::AbstractDiffEqArray, i::Int) Base.getindex(A, :, i) false +@deprecate Base.getindex(A::AbstractDiffEqArray{T,N,A}, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[I] false + +@deprecate Base.getindex(A::AbstractDiffEqArray{T,N,A}, i::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[i] false __parameterless_type(T) = Base.typename(T).wrapper @@ -520,22 +523,25 @@ Base.@propagate_inbounds function Base.setindex!(VA::AbstractVectorOfArray{T, N} VA.u[I] = v end -@deprecate Base.setindex!(VA::AbstractVectorOfArray, v, I::Int) Base.setindex!(VA, v, :, I) false +Base.@propagate_inbounds Base.setindex!(VA::AbstractVectorOfArray, v, I::Int) = Base.setindex!(VA.u, v, I) +@deprecate Base.setindex!(VA::AbstractVectorOfArray{T,N,A}, v, I::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} Base.setindex!(VA.u, v, I) false Base.@propagate_inbounds function Base.setindex!(VA::AbstractVectorOfArray{T, N}, v, ::Colon, I::Colon) where {T, N} VA.u[I] = v end -@deprecate Base.setindex!(VA::AbstractVectorOfArray, v, I::Colon) Base.setindex!( - VA, v, :, I) false +Base.@propagate_inbounds Base.setindex!(VA::AbstractVectorOfArray, v, I::Colon) = Base.setindex!(VA.u, v, I) +@deprecate Base.setindex!(VA::AbstractVectorOfArray{T,N,A}, v, I::Colon) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} Base.setindex!( + VA.u, v, I) false Base.@propagate_inbounds function Base.setindex!(VA::AbstractVectorOfArray{T, N}, v, ::Colon, I::AbstractArray{Int}) where {T, N} VA.u[I] = v end -@deprecate Base.setindex!(VA::AbstractVectorOfArray, v, I::AbstractArray{Int}) Base.setindex!( +Base.@propagate_inbounds Base.setindex!(VA::AbstractVectorOfArray, v, I::AbstractArray{Int}) = Base.setindex!(VA.u, v, I) +@deprecate Base.setindex!(VA::AbstractVectorOfArray{T,N,A}, v, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} Base.setindex!( VA, v, :, I) false Base.@propagate_inbounds function Base.setindex!( From a2d0d5e5d582d423ac001a1760fe23c64ac1086e Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Apr 2025 18:37:30 -0400 Subject: [PATCH 90/94] Update vector_of_array.jl --- src/vector_of_array.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index 1a7d06c6..c05a35c2 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -407,13 +407,13 @@ Base.getindex(A::AbstractVectorOfArray, I::AbstractArray{Int}) = A.u[I] Base.getindex(A::AbstractDiffEqArray, I::Int) = A.u[I] Base.getindex(A::AbstractDiffEqArray, I::AbstractArray{Int}) = A.u[I] -@deprecate Base.getindex(A::AbstractVectorOfArray{T,N,A}, I::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[I] false +@deprecate Base.getindex(VA::AbstractVectorOfArray{T,N,A}, I::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} VA.u[I] false -@deprecate Base.getindex(A::AbstractVectorOfArray{T,N,A}, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[I] false +@deprecate Base.getindex(VA::AbstractVectorOfArray{T,N,A}, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} VA.u[I] false -@deprecate Base.getindex(A::AbstractDiffEqArray{T,N,A}, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[I] false +@deprecate Base.getindex(VA::AbstractDiffEqArray{T,N,A}, I::AbstractArray{Int}) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} VA.u[I] false -@deprecate Base.getindex(A::AbstractDiffEqArray{T,N,A}, i::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} A.u[i] false +@deprecate Base.getindex(VA::AbstractDiffEqArray{T,N,A}, i::Int) where {T,N,A<:Union{AbstractArray, AbstractVectorOfArray}} VA.u[i] false __parameterless_type(T) = Base.typename(T).wrapper From 50b064f3f05dcd7079bd06353634f88007ed7b3b Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Apr 2025 18:57:29 -0400 Subject: [PATCH 91/94] add a documentation page for the interface --- docs/pages.jl | 1 + docs/src/AbstractVectorOfArrayInterface.md | 6 ++ src/RecursiveArrayTools.jl | 113 +++++++++++++++++++++ 3 files changed, 120 insertions(+) create mode 100644 docs/src/AbstractVectorOfArrayInterface.md diff --git a/docs/pages.jl b/docs/pages.jl index 5d9ec136..16338f15 100644 --- a/docs/pages.jl +++ b/docs/pages.jl @@ -2,6 +2,7 @@ pages = [ "Home" => "index.md", + "AbstractVectorOfArrayInterface.md", "array_types.md", "recursive_array_functions.md" ] diff --git a/docs/src/AbstractVectorOfArrayInterface.md b/docs/src/AbstractVectorOfArrayInterface.md new file mode 100644 index 00000000..907cc511 --- /dev/null +++ b/docs/src/AbstractVectorOfArrayInterface.md @@ -0,0 +1,6 @@ +# The AbstractVectorOfArray and AbstractDiffEqArray Interfaces + +```@doc +AbstractVectorOfArray +AbstractDiffEqArray +``` diff --git a/src/RecursiveArrayTools.jl b/src/RecursiveArrayTools.jl index 9e251ac0..a616c592 100644 --- a/src/RecursiveArrayTools.jl +++ b/src/RecursiveArrayTools.jl @@ -13,7 +13,120 @@ import Adapt import Tables, IteratorInterfaceExtensions +""" + AbstractVectorOfArray{T, N, A} + +An AbstractVectorOfArray is an object which represents arrays of arrays, +and arbitrary recursive nesting of arrays, as a single array-like object. +Thus a canonical example of an AbstractVectorOfArray is something of the +form `VectorOfArray([[1,2],[3,4]])`, which "acts" like the matrix [1 3; 2 4] +where the data is stored and accessed in a column-ordered fashion (as is typical +in Julia), but the actual matrix is never constructed and instead lazily represented +through the type. + +An AbstractVectorOfArray subtype should match the following behaviors. + +!!! note + + In 2023 the linear indexing `A[i]`` was deprecated. It previously had the behavior that + `A[i] = A.u[i]`. However, this is incompatible with standard `AbstractArray` interfaces, + Since if `A = VectorOfArray([[1,2],[3,4]])` and `A` is supposed to act like `[1 3; 2 4]`, + then there is a difference `A[1] = [1,2]` for the VectorOfArray while `A[1] = 1` for the + matrix. This causes many issues if `AbstractVectorOfArray <: AbstractArray`. Thus we + plan in 2026 to complete the deprecation and thus have a breaking update where `A[i]` + matches the linear indexing of an `AbstractArray`, and then making + `AbstractVectorOfArray <: AbstractArray`. Until then, `AbstractVectorOfArray` due to + this interface break but manaully implements an AbstractArray-like interface for + future compatability. + +## Fields + +An AbstractVectorOfArray has the following fields: + +* `u` which holds the Vector of values at each timestep + +## Array Interface + +The general operations are as follows. Use + +```julia +A.u[j] +``` + +to access the `j`th array. For multidimensional systems, this +will address first by component and lastly by time, and thus + +```julia +A[i, j] +``` + +will be the `i`th component at array `j`. Hence, `A[j][i] == A[i, j]`. This is done +because Julia is column-major, so the leading dimension should be contiguous in memory. +If the independent variables had shape (for example, was a matrix), then `i` is the +linear index. We can also access solutions with shape: + +```julia +A[i, k, j] +``` + +gives the `[i,k]` component of the system at array `j`. The colon operator is +supported, meaning that + +```julia +A[i, :] +``` + +gives the timeseries for the `i`th component. + +## Using the AbstractArray Interface + +The `AbstractArray` interface can be directly used. For example, for a vector +system of variables `A[i,j]` is a matrix with rows being the variables and +columns being the timepoints. Operations like `A'` will +transpose the solution type. Functionality written for `AbstractArray`s can +directly use this. For example, the Base `cov` function computes correlations +amongst columns, and thus: + +```julia +cov(A) +``` + +computes the correlation of the system state in time, whereas + +```julia +cov(A, 2) +``` + +computes the correlation between the variables. Similarly, `mean(A,2)` is the +mean of the variable in time, and `var(A,2)` is the variance. Other statistical +functions and packages which work on `AbstractArray` types will work on the +solution type. + +## Conversions + +At anytime, a true `Array` can be created using `Array(A)`, or more generally `stack(A)` +to make the array type match the internal array type (for example, if `A` is an array +of GPU arrays, `stack(A)` will be a GPU array). +""" abstract type AbstractVectorOfArray{T, N, A} end + +""" + AbstractDiffEqArray{T, N, A} <: AbstractVectorOfArray{T, N, A} + +An AbstractVectorOfArray object which has extra information of a time array `A.t` +in order to specify a time series. A canonical AbstractDiffEqArray is for example +the pairing `DiffEqArray([[1,2],[3,4]],[1.0,2.0])` which means that at time 1.0 +the values were `[1,2]` and at time 2.0 the values were `[3,4]`. + +An AbstractDiffEqArray has all of the same behaviors as an AbstractVectorOfArray with the +additional properties: + +## Fields + +An AbstractDiffEqArray adds the following fields: + +* `t` which holds the times of each timestep. +""" abstract type AbstractDiffEqArray{T, N, A} <: AbstractVectorOfArray{T, N, A} end include("utils.jl") From 32da6b1ce867ff804411af21cfe3daff46f102aa Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Apr 2025 19:07:24 -0400 Subject: [PATCH 92/94] Update Project.toml --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 11132fdf..c3dbc132 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "RecursiveArrayTools" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" authors = ["Chris Rackauckas "] -version = "3.32.1" +version = "3.33.0" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" From 0f2f3a3810eedb14430087abbd46d1d0a76320f5 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 15 May 2025 20:22:23 +0000 Subject: [PATCH 93/94] Fix interface doc builds --- docs/src/AbstractVectorOfArrayInterface.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/AbstractVectorOfArrayInterface.md b/docs/src/AbstractVectorOfArrayInterface.md index 907cc511..c93c6cca 100644 --- a/docs/src/AbstractVectorOfArrayInterface.md +++ b/docs/src/AbstractVectorOfArrayInterface.md @@ -1,6 +1,6 @@ # The AbstractVectorOfArray and AbstractDiffEqArray Interfaces -```@doc +```@docs AbstractVectorOfArray AbstractDiffEqArray ``` From 4b9d9d4897ddf253dd369dc3a4d6bb1029f46815 Mon Sep 17 00:00:00 2001 From: adienes <51664769+adienes@users.noreply.github.com> Date: Mon, 19 May 2025 16:51:25 -0400 Subject: [PATCH 94/94] Update vector_of_array.jl --- src/vector_of_array.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/vector_of_array.jl b/src/vector_of_array.jl index c05a35c2..98af2ddd 100644 --- a/src/vector_of_array.jl +++ b/src/vector_of_array.jl @@ -914,7 +914,7 @@ end ## broadcasting struct VectorOfArrayStyle{N} <: Broadcast.AbstractArrayStyle{N} end # N is only used when voa sees other abstract arrays -VectorOfArrayStyle(::Val{N}) where {N} = VectorOfArrayStyle{N}() +VectorOfArrayStyle{N}(::Val{N}) where {N} = VectorOfArrayStyle{N}() # The order is important here. We want to override Base.Broadcast.DefaultArrayStyle to return another Base.Broadcast.DefaultArrayStyle. Broadcast.BroadcastStyle(a::VectorOfArrayStyle, ::Base.Broadcast.DefaultArrayStyle{0}) = a