How to find the index of a Dense Axis Array?

Specific Domains Optimization (Mathematical)
1

Hello! I know I can use something like findall(x->x==1, myArray) to find the indices of elements of value 1 inside myArray.

However, my variables in my optimization problem are DenseAxisArrays, as they are indexed with numbers and strings. So when I do, for example:

findall(x->x==1,JuMP.value.(Z))

I get:

ERROR: MethodError: no method matching CartesianIndices(::Tuple{Array{Any,1},Array{Any,1},Base.OneTo{Int64}})
Closest candidates are:
  CartesianIndices(::JuMP.Containers.DenseAxisArray) at C:\Users\Alfredo\.julia\packages\JuMP\qhoVb\src\Containers\DenseAxisArray.jl:111
  CartesianIndices(::AbstractArray) at multidimensional.jl:264
  CartesianIndices(::Tuple{}) at multidimensional.jl:255

How could I get the indices?

Edit: here’s a minimal working example

using JuMP, Cbc

m = Model(Cbc.Optimizer)

TASKS         = ["TASK1"]
EQUIPMENTS          = ["EQUIPMENT1"]
@variable(m, Z[E in EQUIPMENTS, T in TASKS, 1:32], Bin)
@constraint(m, sum(Z)==1)
@objective(m, Min, sum(Z))
optimize!(m)

#output
println(JuMP.value.(Z))
println(JuMP.objective_value(m))
println(JuMP.termination_status(m))
println("Finished!")
4

You want:

[Z[i] for i in eachindex(Z) if value(Z[i]) > 0.5]

1-element Array{VariableRef,1}:
Z[EQUIPMENT1,TASK1,1]

or perhaps:

[k.I for k in keys(Z) if value(Z[k]) > 0.5]

1-element Array{Tuple{String,String,Int64},1}:
("EQUIPMENT1", "TASK1", 1)

The documentation needs to be improved on this front: https://github.com/jump-dev/JuMP.jl/issues/1663

3 Likes
5

Interesting. Apparently findall(f, A) calls pairs(A) which in turn calls CartesianIndices(axes(A)). And that errors since JuMP has defined axes(A) for a `DenseAxisArray´ to return the axes you originally specified instead of integer axes:

julia> axes(Z)
(["EQUIPMENT1"], ["TASK1"], Base.OneTo(32))

julia> CartesianIndices(axes(Z))
ERROR: MethodError: no method matching CartesianIndices(::Tuple{Vector{String}, Vector{String}, Base.OneTo{Int64}})

So JuMP probably needs to rethink this approach. A hack that would make findall work is to define pairs for a DenseAxisArray:

julia> Base.pairs(x::JuMP.Containers.DenseAxisArray) = pairs(x.data)

julia> findall(x->x==1, JuMP.value.(Z))
1-element Vector{CartesianIndex{3}}:
 CartesianIndex(1, 1, 2)

But to avoid committing type piracy, a better workaround for you might be to just do:

julia> findall(x->x==1, JuMP.value.(Z).data)
1-element Vector{CartesianIndex{3}}:
 CartesianIndex(1, 1, 2)

In any case, please open an issue for this in the JuMP repo.

EDIT: Got ninja’d by @odow there. He suggests a different approach to get the index - but wouldn’t it be better for interoperability if a DenseAxisArray worked like any other AbstractArray? Have JuMP devs now settled long term on their own container types or is the hope still that some other generic SomeAxisArray package will take care of the details?

3 Likes
6

@odow Thanks! It worked! And I agree, the documentation should be improved. Thanks for commenting the issue!

@NiclasMattsson Thank you too, your solution also worked! But I’m having trouble understanding the rest of your post, as I’m really new to Julia; I’ve only been using it from the last week.

What I’m getting is that you’re suggesting to redefine de function pairs for when the parameter is a DenseAxisArray. But what is type piracy? And why would redefining pairs would be commiting type piracy?

And what is interoperability?

Also, I will mark @odow 's reply as solution, because I realized that a cartesian index is not what I want, as it won’t give me the string of the corresponding index.

Thank you both!