From f79173487af7bd55d8de6d8ff9fc54a9838fa194 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 24 Apr 2025 18:25:01 -0400 Subject: [PATCH 1/3] 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 2/3] 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 3/3] 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")