Julia's Extended Array Indexing Interface

The following ArrayInterface functions extend Julia's Base LinearAlgebra interface to improve the ability to write code for generic array types.

Indexing Traits

The following traits allow for one to accurately determine the type of indexing allowed on arrays in order to write optimal code for generic array types.

ArrayInterface.can_avxFunction
can_avx(f) -> Bool

Returns true if the function f is guaranteed to be compatible with LoopVectorization.@avx for supported element and array types. While a return value of false does not indicate the function isn't supported, this allows a library to conservatively apply @avx only when it is known to be safe to do so.

function mymap!(f, y, args...)
    if can_avx(f)
        @avx @. y = f(args...)
    else
        @. y = f(args...)
    end
end
source
ArrayInterface.can_change_sizeFunction
can_change_size(::Type{T}) -> Bool

Returns true if the Base.size of T can change, in which case operations such as pop! and popfirst! are available for collections of type T.

source
ArrayInterface.ismutableFunction
ismutable(::Type{T}) -> Bool

Query whether instances of type T are mutable or not, see https://github.com/JuliaDiffEq/RecursiveArrayTools.jl/issues/19.

source
ArrayInterface.ndims_indexFunction
ndims_index(::Type{I}) -> Int

Returns the number of dimensions that an instance of I indexes into. If this method is not explicitly defined, then 1 is returned.

See also ndims_shape

Examples

julia> ArrayInterface.ndims_index(Int)
1

julia> ArrayInterface.ndims_index(CartesianIndex(1, 2, 3))
3

julia> ArrayInterface.ndims_index([CartesianIndex(1, 2), CartesianIndex(1, 3)])
2
source
ArrayInterface.ndims_shapeFunction
ndims_shape(::Type{I}) -> Union{Int,Tuple{Vararg{Int}}}

Returns the number of dimension that are represented in the shape of the returned array when indexing with an instance of I.

See also ndims_index

Examples

```julia julia> ArrayInterface.ndims_shape([CartesianIndex(1, 1), CartesianIndex(1, 2)]) 1

julia> ndims(CartesianIndices((2,2))[[CartesianIndex(1, 1), CartesianIndex(1, 2)]]) 1

source
ArrayInterface.defines_stridesFunction
defines_strides(::Type{T}) -> Bool

Is strides(::T) defined? It is assumed that types returning true also return a valid pointer on pointer(::T).

source
ArrayInterface.ensures_all_uniqueFunction
ensures_all_unique(T::Type) -> Bool

Returns true if all instances of type T are composed of a unique set of elements. This does not require that T subtypes AbstractSet or implements the AbstractSet interface.

Examples

julia> ArrayInterface.ensures_all_unique(BitSet())
true

julia> ArrayInterface.ensures_all_unique([])
false

julia> ArrayInterface.ensures_all_unique(typeof(1:10))
true

julia> ArrayInterface.ensures_all_unique(LinRange(1, 1, 10))
false
source
ArrayInterface.ensures_sortedFunction
ensures_sorted(T::Type) -> Bool

Returns true if all instances of T are sorted.

Examples

julia> ArrayInterface.ensures_sorted(BitSet())
true

julia> ArrayInterface.ensures_sorted([])
false

julia> ArrayInterface.ensures_sorted(1:10)
true
source
ArrayInterface.indices_do_not_aliasFunction
indices_do_not_alias(::Type{T<:AbstractArray}) -> Bool

Is it safe to ivdep arrays of type T? That is, would it be safe to write to an array of type T in parallel? Examples where this is not true are BitArrays or view(rand(6), [1,2,3,1,2,3]). That is, it is not safe whenever different indices may alias the same memory.

source
ArrayInterface.instances_do_not_aliasFunction
instances_do_not_alias(::Type{T}) -> Bool

Is it safe to ivdep arrays containing elements of type T? That is, would it be safe to write to an array full of T in parallel? This is not true for mutable structs in general, where editing one index could edit other indices. That is, it is not safe when different instances may alias the same memory.

source
ArrayInterface.deviceFunction
device(::Type{T}) -> AbstractDevice

Indicates the most efficient way to access elements from the collection in low-level code. For GPUArrays, will return ArrayInterface.GPU(). For AbstractArray supporting a pointer method, returns ArrayInterface.CPUPointer(). For other AbstractArrays and Tuples, returns ArrayInterface.CPUIndex(). Otherwise, returns nothing.

source

Allowed Indexing Functions

These are generic functions for forced "allowed indexing". For example, with CUDA.jl's CuArrays a mode can be enabled such that allowscalar(false) forces errors to be thrown if a GPU array is scalar indexed. Instead of using the CUDA-specific CUDA.@allowscalar on an operation, these functions allow for a general generic "allowed indexing" for all array types.

Indexing Type Buffers

The following indexing types allow for generically controlling bounds checking and index translations.

ArrayInterface.ArrayIndexType
ArrayIndex{N}

Subtypes of ArrayIndex represent series of transformations for a provided index to some buffer which is typically accomplished with square brackets (e.g., buffer[index[inds...]]). The only behavior that is required of a subtype of ArrayIndex is the ability to transform individual index elements (i.e. not collections). This does not guarantee bounds checking or the ability to iterate (although additional functionality may be provided for specific types).

source
ArrayInterface.GetIndexType
GetIndex(buffer) = GetIndex{true}(buffer)
GetIndex{check}(buffer) -> g

Wraps an indexable buffer in a function type that is indexed when called, so that g(inds..) is equivalent to buffer[inds...]. If check is false, then all indexing arguments are considered in-bounds. The default value for check is true, requiring bounds checking for each index.

See also SetIndex!

Warning

Passing false as check may result in incorrect results/crashes/corruption for out-of-bounds indices, similar to inappropriate use of @inbounds. The user is responsible for ensuring this is correctly used.

Examples

julia> ArrayInterface.GetIndex(1:10)(3)
3

julia> ArrayInterface.GetIndex{false}(1:10)(11)  # shouldn't be in-bounds
11
source
ArrayInterface.SetIndex!Type
SetIndex!(buffer) = SetIndex!{true}(buffer)
SetIndex!{check}(buffer) -> g

Wraps an indexable buffer in a function type that sets a value at an index when called, so that g(val, inds..) is equivalent to setindex!(buffer, val, inds...). If check is false, then all indexing arguments are considered in-bounds. The default value for check is true, requiring bounds checking for each index.

See also GetIndex

Warning

Passing false as check may result in incorrect results/crashes/corruption for out-of-bounds indices, similar to inappropriate use of @inbounds. The user is responsible for ensuring this is correctly used.

Examples


julia> x = [1, 2, 3, 4];

julia> ArrayInterface.SetIndex!(x)(10, 2);

julia> x[2]
10
source