Solve minimization problem where constraint is the system of linear inequation with MathOptInterface efficiently

Specific Domains Optimization (Mathematical)
1

I want to solve

min cx
s.t. Ax + b < 0

with MathOptInterface, where A is a matrix and b and c are vectors.
However, I couldn’t understand the documentation of MathOptInterface.

Currently, my code is as follows:
(Note that A is m * n matrix, b and c are vector with n and m elements respectively.)

using MathOptInterface
const MOI = MathOptInterface
using GLPK

optimizer = GLPK.Optimizer()
x = MOI.add_variables(optimizer, m)
MOI.set(optimizer,
        MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(),
        MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(c, x), 0.0))
MOI.set(optimizer, MOI.ObjectiveSense(), MOI.MIN_SENSE)
for i in 1:n
        MOI.add_constraint(optimizer,
        MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(A[i,:],x), b[i]),
                                 MOI.LessThan(0.0))
end
MOI.optimize!(optimizer)

That is, I read the matrix with each rows and create scalar affine constraint.
Actually, this code works, but I think there is more efficient way to solve this problem, in particular using VectorAffineFunction and VectorAffineTerm. I couldn’t get it when I read the docs.

Does anyone knows the solution of this?

1 Like
2

This code looks fine to me. Why do you think that it’s inefficient?

3

Sorry, I said “inefficient” because it takes long time when the matrix A is large. (Is it true?)
Also, I think the code looks sophisticated when I use VectorAffineFunction.
I just want to know how to use VectorAffineFunction.

1 Like
4

I just want to know how to use VectorAffineFunction .

You can do it in a way similar to how you build ScalarAffineFunction with broadcasting in your first comment.
Let’s take these A and b for example (any value would work)

julia> A = reshape(1:6, 2, 3)
2×3 reshape(::UnitRange{Int64}, 2, 3) with eltype Int64:
 1  3  5
 2  4  6

julia> b = collect(1:3)
3-element Array{Int64,1}:
 1
 2
 3

You can do as follows

julia> using MathOptInterface

julia> const MOI = MathOptInterface
MathOptInterface

julia> x = MOI.VariableIndex.(1:3) # should normally comes from `MOI.add_variables`
3-element Array{MathOptInterface.VariableIndex,1}:
 MathOptInterface.VariableIndex(1)
 MathOptInterface.VariableIndex(2)
 MathOptInterface.VariableIndex(3)

julia> terms = MOI.VectorAffineTerm.(1:2, MOI.ScalarAffineTerm.(A, reshape(x, 1, 3)))
2×3 Array{MathOptInterface.VectorAffineTerm{Int64},2}:
 VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(1, VariableIndex(1)))  …  VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(5, VariableIndex(3)))
 VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(2, VariableIndex(1)))     VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(6, VariableIndex(3)))

julia> f = MOI.VectorAffineFunction(vec(terms), b)
MathOptInterface.VectorAffineFunction{Int64}(MathOptInterface.VectorAffineTerm{Int64}[VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(1, VariableIndex(1))), VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(2, VariableIndex(1))), VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(3, VariableIndex(2))), VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(4, VariableIndex(2))), VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(5, VariableIndex(3))), VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(6, VariableIndex(3)))], [1, 2, 3])
5

I got the following error

julia> using MathOptInterface

julia> const MOI = MathOptInterface
MathOptInterface

julia> using GLPK

julia> A = reshape(1:6, 2, 3)
2×3 reshape(::UnitRange{Int64}, 2, 3) with eltype Int64:
 1  3  5
 2  4  6

julia> b = collect(1:2)
2-element Array{Int64,1}:
 1
 2

julia> c = [1.0, 2.0, 3.0]
3-element Array{Float64,1}:
 1.0
 2.0
 3.0

julia> optimizer = GLPK.Optimizer()
A LinQuadOptInterface model with backend:
Prob(Ptr{Nothing} @0x00007fa352a88810)

julia> x = MOI.add_variables(optimizer, 3)
3-element Array{MathOptInterface.VariableIndex,1}:
 MathOptInterface.VariableIndex(1)
 MathOptInterface.VariableIndex(2)
 MathOptInterface.VariableIndex(3)

julia> MOI.set(optimizer,
               MOI.ObjectiveFunction{MOI.ScalarAffineFunction{Float64}}(),
               MOI.ScalarAffineFunction(MOI.ScalarAffineTerm.(c, x), 0.0))

julia> MOI.set(optimizer, MOI.ObjectiveSense(), MOI.MIN_SENSE)
MIN_SENSE::OptimizationSense = 0

julia> terms = MOI.VectorAffineTerm.(1:2, MOI.ScalarAffineTerm.(A, reshape(x, 1, 3)))
2×3 Array{MathOptInterface.VectorAffineTerm{Int64},2}:
 VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(1, VariableIndex(1)))  …  VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(5, VariableIndex(3)))
 VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(2, VariableIndex(1)))     VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(6, VariableIndex(3)))

julia> f = MOI.VectorAffineFunction(vec(terms), b)
MathOptInterface.VectorAffineFunction{Int64}(MathOptInterface.VectorAffineTerm{Int64}[VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(1, VariableIndex(1))), VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(2, VariableIndex(1))), VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(3, VariableIndex(2))), VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(4, VariableIndex(2))), VectorAffineTerm{Int64}(1, ScalarAffineTerm{Int64}(5, VariableIndex(3))), VectorAffineTerm{Int64}(2, ScalarAffineTerm{Int64}(6, VariableIndex(3)))], [1, 2])

julia> MOI.add_constraint(optimizer, f, MOI.Nonpositives(2))
ERROR: MathOptInterface.UnsupportedConstraint{MathOptInterface.VectorAffineFunction{Int64},MathOptInterface.Nonpositives}: `MathOptInterface.VectorAffineFunction{Int64}`-in-`MathOptInterface.Nonpositives` constraints is not supported by the model.
Stacktrace:
 [1] #correct_throw_add_constraint_error_fallback#30(::MathOptInterface.AddConstraintNotAllowed{MathOptInterface.VectorAffineFunction{Int64},MathOptInterface.Nonpositives}, ::Function, ::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Int64}, ::MathOptInterface.Nonpositives) at /Users/myname/.julia/packages/MathOptInterface/C3lip/src/constraints.jl:122
 [2] correct_throw_add_constraint_error_fallback(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Int64}, ::MathOptInterface.Nonpositives) at /Users/myname/.julia/packages/MathOptInterface/C3lip/src/constraints.jl:119
 [3] #throw_add_constraint_error_fallback#27(::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::Function, ::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Int64}, ::MathOptInterface.Nonpositives) at /Users/myname/.julia/packages/MathOptInterface/C3lip/src/constraints.jl:99
 [4] throw_add_constraint_error_fallback(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Int64}, ::MathOptInterface.Nonpositives) at /Users/myname/.julia/packages/MathOptInterface/C3lip/src/constraints.jl:99
 [5] add_constraint(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Int64}, ::MathOptInterface.Nonpositives) at /Users/myname/.julia/packages/MathOptInterface/C3lip/src/constraints.jl:83
 [6] top-level scope at none:0

How can I fix this error? Why VectorAffineFunction doesn’t support Nonpositives? Linear constraints says it supports Nonpositives.

1 Like
6

Use Float64 instead of Int coefficients in A and b. This is slight usability bug from the MathOptInterface perspective, but isn’t visible from JuMP since it automatically converts values to Float64s.

7

Thanks, my code works well with using Float64 !

8

The example script is really helpful to understand how to use VectorAffineTerm and Function.

One thing that I cannot understand is: why do we need to apply reshape(x, 1, n) in ScalarAffineTerm. ? I can hardly couple the codes and their corresponding constraints.

Also, when matrices and vectors are given for an optimization problem (e.g. A, c, lb, ub, etc.), is there a simple way to solve it with MathOptInterface as we do so with MathProgBase?

9

One thing that I cannot understand is: why do we need to apply reshape(x, 1, n) in ScalarAffineTerm. ? I can hardly couple the codes and their corresponding constraints.

You can to create a matrix of terms. The term (i, j) corresponds to the entry (i, j) of the matrix so it should have row i and variable x[j].
In the broadcast, 1:2 and x are one dimensional and A is 2-dimensional so they will be broadcasted constant along the missing dimension.
As 1:2 is constant in the second dimension, it is broadcasted correctly but x should be broadcasted constant in the first dimension so we reshape it into a row vector for that.

Also, when matrices and vectors are given for an optimization problem (e.g. A , c , lb , ub , etc.), is there a simple way to solve it with MathOptInterface as we do so with MathProgBase ?

The simple way is to do it in JuMP with

@variable(model, lb[i] <= x[i=1:n] <= ub[i])
@constraint(model, l .<= A * x .<= u)

For a way similar to MathProgBase, see a relevant discussion here:

10

Hi, This is still giving me errors. I’m running the following code:

using MathOptInterface
const MOI = MathOptInterface
using GLPK

optimizer = GLPK.Optimizer()

A = reshape(1:6, 2, 3)
A = convert(AbstractArray{Float64}, A)
b = collect(1:3)
b = convert(AbstractArray{Float64}, b)

x = MOI.VariableIndex.(1:3) # should normally comes from `MOI.add_variables`

terms = MOI.VectorAffineTerm.(1:2, MOI.ScalarAffineTerm.(A, reshape(x, 1, 3)))

f = MOI.VectorAffineFunction(vec(terms), b)

MOI.add_constraint(optimizer, f, MOI.Nonpositives(2))

and getting the following error:

ERROR: MathOptInterface.UnsupportedConstraint{MathOptInterface.VectorAffineFunction{Float64},MathOptInterface.Nonpositives}: `MathOptInterface.VectorAffineFunction{Float64}`-in-`MathOptInterface.Nonpositives` constraint is not supported by the model.
Stacktrace:
 [1] correct_throw_add_constraint_error_fallback(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Float64}, ::MathOptInterface.Nonpositives; error_if_supported::MathOptInterface.AddConstraintNotAllowed{MathOptInterface.VectorAffineFunction{Float64},MathOptInterface.Nonpositives}) at C:\Users\pnvol\.julia\packages\MathOptInterface\VjkSQ\src\constraints.jl:158
 [2] correct_throw_add_constraint_error_fallback(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Float64}, ::MathOptInterface.Nonpositives) at C:\Users\pnvol\.julia\packages\MathOptInterface\VjkSQ\src\constraints.jl:155
 [3] throw_add_constraint_error_fallback(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Float64}, ::MathOptInterface.Nonpositives; kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at C:\Users\pnvol\.julia\packages\MathOptInterface\VjkSQ\src\constraints.jl:135
 [4] throw_add_constraint_error_fallback(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Float64}, ::MathOptInterface.Nonpositives) at C:\Users\pnvol\.julia\packages\MathOptInterface\VjkSQ\src\constraints.jl:135
 [5] add_constraint(::GLPK.Optimizer, ::MathOptInterface.VectorAffineFunction{Float64}, ::MathOptInterface.Nonpositives) at C:\Users\pnvol\.julia\packages\MathOptInterface\VjkSQ\src\constraints.jl:119
 [6] top-level scope at C:\Users\pnvol\Dropbox\JYC\Code\J15_2017_chart\MWE_MOI.jl:90

Any ideas? I just really miss linprog.

1 Like
11

GLPK doesn’t support VectorAffineFunction

Either use

optimizer = MOI.Bridges.full_bridge_optimizer(GLPK.Optimizer(), Float64)

or pass ScalarAffineFunctions.

I just really miss linprog.

Use JuMP:

model = Model(GLPK.Optimizer)
@variable(model, l[i] <= x[i = 1:n] <= u[i])
@objective(model, c' * x)
@constraint(model, A * x .>= b)
optimize!(model)
12

Thanks. Both those solutions work. Much appreciated.