Containers do not support view

[odow]

From Pierre Pasquet on slack:

julia> x = Containers.@container([i=[:a, :b]], 0)
1-dimensional DenseAxisArray{Int64,1,...} with index sets:
    Dimension 1, [:a, :b]
And data, a 2-element Vector{Int64}:
 0
 0

julia> view(x, :)
ERROR: MethodError: no method matching Base.Slice(::Vector{Symbol})
Closest candidates are:
  Base.Slice(::Base.Slice) at indices.jl:353
  Base.Slice(::T) where T<:AbstractUnitRange at indices.jl:351
Stacktrace:
 [1] uncolon(inds::Tuple{Vector{Symbol}}, I::Tuple{Colon})
   @ Base ./multidimensional.jl:827
 [2] to_indices
   @ ./multidimensional.jl:822 [inlined]
 [3] to_indices
   @ ./indices.jl:325 [inlined]
 [4] view(A::JuMP.Containers.DenseAxisArray{Int64, 1, Tuple{Vector{Symbol}}, Tuple{JuMP.Containers._AxisLookup{Dict{Symbol, Int64}}}}, I::Function)
   @ Base ./subarray.jl:176
 [5] top-level scope
   @ REPL[7]:1

We should throw a nicer error message, since this might be something that people do more often in future.

[odow]

Annoyingly, it works if the indices are AbstractUnitRange:

julia> x = Containers.@container([i=1:2], 0)
2-element Vector{Int64}:
 0
 0

julia> view(x, :)
2-element view(::Vector{Int64}, :) with eltype Int64:
 0
 0

julia> x = Containers.@container([i=2:3], 0)
1-dimensional DenseAxisArray{Int64,1,...} with index sets:
    Dimension 1, 2:3
And data, a 2-element Vector{Int64}:
 0
 0

julia> v = view(x, :)
2-element view(::JuMP.Containers.DenseAxisArray{Int64, 1, Tuple{UnitRange{Int64}}, Tuple{JuMP.Containers._AxisLookup{Tuple{Int64, Int64}}}}, :) with eltype Int64 with indices 2:3:
 0
 0

julia> v[1] = 1
ERROR: BoundsError: attempt to access 2-element view(::JuMP.Containers.DenseAxisArray{Int64, 1, Tuple{UnitRange{Int64}}, Tuple{JuMP.Containers._AxisLookup{Tuple{Int64, Int64}}}}, :) with eltype Int64 with indices 2:3 at index [1]
Stacktrace:
 [1] throw_boundserror(A::SubArray{Int64, 1, JuMP.Containers.DenseAxisArray{Int64, 1, Tuple{UnitRange{Int64}}, Tuple{JuMP.Containers._AxisLookup{Tuple{Int64, Int64}}}}, Tuple{Base.Slice{UnitRange{Int64}}}, false}, I::Tuple{Int64})
   @ Base ./abstractarray.jl:651
 [2] checkbounds
   @ ./abstractarray.jl:616 [inlined]
 [3] setindex!(V::SubArray{Int64, 1, JuMP.Containers.DenseAxisArray{Int64, 1, Tuple{UnitRange{Int64}}, Tuple{JuMP.Containers._AxisLookup{Tuple{Int64, Int64}}}}, Tuple{Base.Slice{UnitRange{Int64}}}, false}, x::Int64, I::Int64)
   @ Base ./subarray.jl:316
 [4] top-level scope
   @ REPL[44]:1

julia> v[2] = 1
1

julia> v
2-element view(::JuMP.Containers.DenseAxisArray{Int64, 1, Tuple{UnitRange{Int64}}, Tuple{JuMP.Containers._AxisLookup{Tuple{Int64, Int64}}}}, :) with eltype Int64 with indices 2:3:
 1
 0

julia> x
1-dimensional DenseAxisArray{Int64,1,...} with index sets:
    Dimension 1, 2:3
And data, a 2-element Vector{Int64}:
 1
 0

[odow]

Seems like over-riding view for all instances would break existing code that currently works.

But I don't know how to intercept a method where the method signature doesn't match.

I guess the fundamental problem is that view assumes:

julia> x
1-dimensional DenseAxisArray{Float64,1,...} with index sets:
    Dimension 1, [:a, :b]
And data, a 2-element Vector{Float64}:
 1.0
 2.0

julia> eachindex(IndexLinear(), x)
2-element Vector{Symbol}:
 :a
 :b

returns an AbtractUnitRange.

[odow]

For example, this is clearly wrong:

julia> c = Containers.@container([i=1:4, j=2:3], i + j)
2-dimensional DenseAxisArray{Int64,2,...} with index sets:
    Dimension 1, Base.OneTo(4)
    Dimension 2, 2:3
And data, a 4×2 Matrix{Int64}:
 3  4
 4  5
 5  6
 6  7

julia> eachindex(IndexLinear(), c)
Base.OneTo(8)

[metab0t]

function Base.view(A::DenseAxisArray{T,N}, idx...) where {T,N}
    new_indices = Base.to_index(A, idx)
    new_axes = _getindex_recurse(A.axes, new_indices, _is_range)
    return DenseAxisArray(view(A.data, new_indices...), new_axes...) 
end

But the code does not work because we use data::Array{T,N} in DenseAxisArray. Maybe we should use AbstractArray to allow SubArray as data storage?

[odow]

Maybe we should use AbstractArray to allow SubArray as data storage?

I'd prefer not too, because we'd need a new type parameter in DenseAxisArray. But I'm not really sure what the right solution is.

[odow]

@metab0t what is the motivation for using a view? When is it useful?

[metab0t]

For performance reasons:

julia> using JuMP                                                                               
                                                                                                
julia> A = Containers.@container([i=1:100, j=1:200], i+j, container=Array);                     
                                                                                                
julia> DA = Containers.@container([i=1:100, j=1:200], i+j, container=Containers.DenseAxisArray);
                                                                                                
julia> size(A)                                                                                  
(100, 200)                                                                                      
                                                                                                
julia> size(DA)                                                                                 
(100, 200)                                                                                      
                                                                                                
julia> using BenchmarkTools

julia> @btime A[:, 50];
  212.571 ns (2 allocations: 912 bytes)

julia> @btime @view A[:, 50];
  94.044 ns (2 allocations: 64 bytes)

julia> @btime DA[:, 50];
  257.864 ns (5 allocations: 1008 bytes)                                                                

[odow]

Sure. But where does this show up as a problem when modeling?

…

On Tue, 2 Aug 2022, 7:18 pm Y. Yang, ***@***.***> wrote: For performance reasons: julia> using JuMP julia> A = ***@***.***([i=1:100, j=1:200], i+j, container=Array); julia> DA = ***@***.***([i=1:100, j=1:200], i+j, container=Containers.DenseAxisArray); julia> size(A) (100, 200) julia> size(DA) (100, 200) julia> using BenchmarkTools julia> @Btime A[:, 50]; 212.571 ns (2 allocations: 912 bytes) julia> @Btime @view A[:, 50]; 94.044 ns (2 allocations: 64 bytes) julia> @Btime DA[:, 50]; 257.864 ns (5 allocations: 1008 bytes) — Reply to this email directly, view it on GitHub <#2998 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AB6MQJN23ATEOYF6HEHFFQDVXDDS7ANCNFSM5YERJXQA> . You are receiving this because you authored the thread.Message ID: ***@***.***>

[metab0t]

Sometimes we want to extract a row or column of variables to treat it like a 1-d array.
For example:

for i in axes(variables, 1)
    @constraint(LinearAlgebra.dot(f(i), variables[i, :]) <= 0.0)
end

[odow]

Less convenient, but you could go something like

for (i, key) in enumerate(axes(variables, 1))
    @constraint(
        model, 
        LinearAlgebra.dot(f(key), view(variables.data, i, :)) <= 0.0,
    )
end

# or

for i in axes(x, 1)
    c = f(i)
    @constraint(model, sum(c[k] * x[i, k] for k in axes(x, 2)) <= 0)
end

[metab0t]

Yes, we can bypass DenseAxisArray to select the view of its underlying array directly as an alternative.

In fact, view of DenseAxisArray is just a wrapper built upon view of Array which will be nice to have as a shortcut.

The problem is what view of DenseAxisArray should return. Maybe a DenseAxisArray wrapping SubArray?

[odow]

It would have to be something like ViewedDenseAxisArray. But that's getting complicated. Do you have a realistic benchmark where this is a problem? (And not just an artificial one. I get the issue. The question is a trade-off in added complexity to JuMP.)