Why does the following setindex call allocate?

General Usage Performance
,
1 
julia> function f!(v, A, r)
           v[eachindex(r)] = @view A[r]
           v
       end
f! (generic function with 1 method)

julia> v = zeros(10); A = ones(4,10);

julia> @btime f!($v, $A, 4:3:13);
  56.050 ns (2 allocations: 80 bytes)

julia> @btime f!($v, $A, 4:4);
  37.753 ns (2 allocations: 80 bytes)
2 Likes
2

Seems likely that it allocates a subarray and fails to optimize it away. Cf

julia> g(v, A, r) = @view A[r]
g (generic function with 1 method)

julia> @btime g($v, $A, 4:3:13);
  11.263 ns (2 allocations: 80 bytes)
1 Like
3

This is also “strange”:

julia> g(A, r) = @view A[r]
g (generic function with 2 methods)

julia> A = reshape(1:(4*10),4,10)
4Ă—10 reshape(::UnitRange{Int64}, 4, 10) with eltype Int64:
 1  5   9  13  17  21  25  29  33  37
 2  6  10  14  18  22  26  30  34  38
 3  7  11  15  19  23  27  31  35  39
 4  8  12  16  20  24  28  32  36  40

julia> B = collect(reshape(1:(4*10),4,10))
4Ă—10 Matrix{Int64}:
 1  5   9  13  17  21  25  29  33  37
 2  6  10  14  18  22  26  30  34  38
 3  7  11  15  19  23  27  31  35  39
 4  8  12  16  20  24  28  32  36  40

julia> @btime g($A, $(4:2:13));
  10.065 ns (0 allocations: 0 bytes)

julia> @btime g($B, $(4:2:13));
  37.539 ns (2 allocations: 80 bytes)
1 Like
4

But a SubArray shouldn’t be allocated in the first place. Unless you mean another Array that is a subset, but if that were to happen, I’d expect the reported amount of memory to vary by the slice’s length, but trying 4:2:13 and 4:1:13 still reports 80 bytes.

5 
 @view A[r]

That seems like it would create a SubArray, no?

6

As far as I can trace it, the view does linear indexing into a 2D-array, which therefore needs to be reshaped into a vector, which eventually lands in julia/reshapedarray.jl at master · JuliaLang/julia · GitHub and the ccall needs to allocate memory.

3 Likes
7

Oh now I notice the view/SubArray’s parent has a different type. The linear indexing reshaped a Matrix, and it allocated a Vector instead of wrapping the Matrix in a ReshapedArray. The allocated memory doesn’t vary with slice length because the Vector and the Matrix share the same buffer.

julia> time_view(A, I) = @time @view(A[I]);
julia> A = ones(4,10); x = time_view(A, 4:3:13)
  0.000004 seconds (2 allocations: 80 bytes)
4-element view(::Vector{Float64}, 4:3:13) with eltype Float64:
 1.0
 1.0
 1.0
 1.0
julia> typeof(reshape(A, length(A)))
Vector{Float64} (alias for Array{Float64, 1})
julia> Base.mightalias(x.parent, A)
true
julia> x.parent[1] = 0; (x.parent[1], A[1])
(0.0, 0.0)

Wonder why it works this way though. A Matrix can be linearly indexed directly, so it seems fine as a parent for ReshapedArray. Why choose the overhead of heap-allocating an aliased Vector?

julia> time_oddview(A, I) = @time view(Base.ReshapedArray(A, (length(A),), ()), I);
julia> A = ones(4,10); x = time_oddview(A, 4:3:13)
  0.000002 seconds
4-element view(reshape(::Matrix{Float64}, 40), 4:3:13) with eltype Float64:
 1.0
 1.0
 1.0
 1.0