Testing and Validating the Strided Array Interface

The strided array interface is documented in Interfaces · The Julia Language

New AbstractArray subtypes can implement this interface to enable fast BLAS, LAPACK, IO, and Python interop methods. In addition, we would like some wrappers and views of strided arrays to be strided arrays as well.

However, due to bugs and unclear documentation, this composability often fails in practice. For example,
https://github.com/JuliaPy/PythonCall.jl/issues/579 and https://github.com/JuliaLang/julia/pull/60533

I propose adding a function to Test, an exported syntax sugar function to Base, and two optional functions to the strided array interface, in order to help test if an array has a valid implementation of the strided array interface. I’m looking for feedback on these new functions.

Test.test_strided_interface(a::AbstractArray)

If isbitstype(eltype(a)) and has_strided_get(a), iterate through the array checking that elements obtained through the strided interface match elements obtained using getindex.

If isbitstype(eltype(a)) and has_strided_set(a), iterate through the array setting elements with the strided interface, and checking that elements match elements obtained using getindex.

Base.has_strided_get(a::AbstractArray)::Bool = false and Base.has_strided_set(a::AbstractArray)::Bool = false

These two functions have fallback definitions that return false and would be specialized to return true for arrays that implement the strided array interface, for example:
has_strided_get(a::Array) = true and has_strided_set(a::Array) = true

The fallback definition ensures backwards compatibility with existing array types.

strided_ptr(f, a::AbstractArray)

This is a simple syntax sugar for the cconvert, unsafe_convert, and GC.@preserve block.
The goal is to have a nice-to-use alternative to the easy-to-misuse pointer function.

function strided_ptr(f, a::AbstractArray{T})
    a_cconv = Base.cconvert(Ptr{T}, a)
    GC.@preserve a_cconv begin
        f(Base.unsafe_convert(Ptr{T}, a_cconv))
    end
end

Should these be defined on instances of AbstractArray or subtypes of AbstractArray?

I’m not sure. Using instances instead of types is more flexible. For example, in PythonCall.jl it seems like a PyArray may or may not be strided based on the value of the .strided field.

I have a draft PR New functions `has_strided_get` and `has_strided_set` by nhz2 · Pull Request #60894 · JuliaLang/julia · GitHub
The functions need to be defined on instances for the ReshapedArray and ReinterpretArray methods to work.

I also ended up putting the test functions in test/testhelpers/StridedArrays.jl.

Looking for help reviewing NFC: Add StridedArrays.jl testhelper by nhz2 · Pull Request #61250 · JuliaLang/julia · GitHub

This PR doesn’t change any functionality, it just consolidates the tests into a test helper file.

I also found a previous attempt to improve the strided array interface that seems to have stalled due to lack of testing. RFC: Make `strides` into a generic trait by andyferris · Pull Request #30432 · JuliaLang/julia · GitHub

In Strided array traits - Pull Request #60964 - JuliaLang/julia - GitHub I am trying a system of type trait functions. Multiple trait functions are needed to capture whether an array is strided at the type level because of the way the Reshape and Reinterpret array wrapper types work. The new specializable functions are:

  1. Base.is_strided(::Type{<:A}): Return true to declare that the array type follows the strided array interface.
  2. Base.is_vec_strided(::Type{<:A}): Return true to declare that the array additionally has evenly spaced elements in column-major order. Implies Base.is_strided.
  3. Base.is_contiguous(::Type{<:A}): Return true to declare that the array additionally has the same memory layout as an Array. Implies Base.is_vec_strided and provides default strides and Base.elsize definitions.
  4. is_ptr_loadable(::Type{<:A}): Return true to indicate the element pointer can be used to load elements.
  5. is_ptr_storable(::Type{<:A}): Return true to indicate the element pointer can be used to store elements.

There are are two additional functions (has_contiguous_layout and has_vec_strided_layout) that should not be specialized because they only exist to handle the fact that all strided zero dim arrays are trivially contiguous, and all strided vectors are trivially vec_strided.

Here is an example of the performance improvements we can get with a working strided array interface.

using ZipArchives, FixedSizeArrays, PythonCall, BenchmarkTools, Random
data = rand(Xoshiro(1234), UInt8, 300000000);
f_data = FixedSizeArray(copy(data));
py_data = PyArray(copy(data));
@btime zip_crc32($data)
@btime zip_crc32($f_data)
@btime zip_crc32($py_data)

Julia 1.12.6:

70.287 ms (0 allocations: 0 bytes)
82.036 ms (3 allocations: 24.07 KiB)
147.411 ms (3 allocations: 24.07 KiB)

With the strided array interface:

Implement strided array traits - Pull Request #777 - JuliaPy/PythonCall.jl - GitHub Implement strided array traits - Pull Request #191 - JuliaArrays/FixedSizeArrays.jl - GitHub and Use try_strides if available - Pull Request #106 - JuliaIO/ZipArchives.jl - GitHub

70.127 ms (0 allocations: 0 bytes)
71.035 ms (0 allocations: 0 bytes)
71.421 ms (0 allocations: 0 bytes)