User-defined (nonlinear) objective with vectorized variables in JuMP

Question

Is it possible to use vectorized variables with user-defined objective functions in JuMP for Julia? Like so,

model = Model(GLPK.Optimizer)

A = [
    1 1 9 5
    3 5 0 8
    2 0 6 13
]

b = [7; 3; 5]

c = [1; 3; 5; 2]

@variable(model, x[1:4] >= 0)
@constraint(model, A * x .== b)

# dummy functions, could be nonlinear hypothetically
identity(x) = x
C(x, c) = c' * x

register(model, :identity, 1, identity; autodiff = true)
register(model, :C, 2, C; autodiff = true)
@NLobjective(model, Min, C(identity(x), c))

This throws the error,

ERROR: Unexpected array VariableRef[x[1], x[2], x[3], x[4]] in nonlinear expression. Nonlinear expressions may contain only scalar expression.

Which sounds like no. Is there a workaround to this? I believe scipy.optimize.minimize is capable of optimizing user-defined objectives with vectorized variables?

score 4 · Answer 1 · answered Jan 22 '22 at 23:34

No, you cannot pass vector arguments to user-defined functions.

Documentation: https://jump.dev/JuMP.jl/stable/manual/nlp/#User-defined-functions-with-vector-inputs
Issue you opened: https://github.com/jump-dev/JuMP.jl/issues/2854

The following is preferable to Prezemyslaw's answer. His suggestion to wrap things in an @expression won't work if the functions are more complicated.

using JuMP, Ipopt
model = Model(Ipopt.Optimizer)

A = [
    1 1 9 5
    3 5 0 8
    2 0 6 13
]

b = [7; 3; 5]

c = [1; 3; 5; 2]

@variable(model, x[1:4] >= 0)
@constraint(model, A * x .== b)

# dummy functions, could be nonlinear hypothetically
identity(x) = x
C(x, c) = c' * x

my_objective(x...) = C(identitiy(collect(x)), c)
register(model, :my_objective, length(x), my_objective; autodiff = true)
@NLobjective(model, Min, my_objective(x...))

what is the rule of thumb when to use `@expression` and when to switch to `register` for nonlinear models? — Przemyslaw Szufel, Jan 23 '22 at 00:27
Your method only works for linear or quadratic expressions. But a fix for that is on the roadmap. — Oscar Dowson, Jan 23 '22 at 01:14

score 3 · Answer 2 · answered Jan 22 '22 at 22:35

Firstly, use optimizer that supports nonlinear models. GLPK does not. Try Ipopt:

using Ipopt
model = Model(Ipopt.Optimizer)

Secondly, JuMP documentation reads (see https://jump.dev/JuMP.jl/stable/manual/nlp/#Syntax-notes):

The syntax accepted in nonlinear macros is more restricted than the syntax for linear and quadratic macros. (...) all expressions must be simple scalar operations. You cannot use dot, matrix-vector products, vector slices, etc.

you need wrap the goal function

@expression(model, expr,  C(identity(x), c))

Now you can do:

@NLobjective(model, Min, expr)

To show that it works I solve the model:

julia> optimize!(model)
This is Ipopt version 3.14.4, running with linear solver MUMPS 5.4.1.
...
Total seconds in IPOPT                               = 0.165

EXIT: Optimal Solution Found.

julia> value.(x)
4-element Vector{Float64}:
  0.42307697548737005
  0.3461538282496562
  0.6923076931757742
 -8.46379887234798e-9

User-defined (nonlinear) objective with vectorized variables in JuMP

2 Answers2