Matrix inverse in constraints in `@NLconstraint`

Specific Domains Optimization (Mathematical)
, ,
1

Hi, I have a question related to nonlinear optimization.

I am trying to impose a nonlinear constraint that involves a matrix inversion. The matrix contains a variable of the optimization model. When I try to do this via @NLconstraint, I get an error. An example is as follows. May I know is there another way to properly impose these types of constraints, or should I use another optimizer? Many thanks!

using JuMP, NLopt

model = Model(NLopt.Optimizer)
set_optimizer_attribute(model, "algorithm", :LD_MMA)
@variable(model, x, lower_bound = 0)
@variable(model, y, lower_bound = 0)
M1 = [1 x; 0 1]

@NLconstraint(model, inv(M1) * [y; x] == [5; 1])

The error message that I get is:

LoadError: Unexpected array AffExpr[1 x; 0 1] in nonlinear expression. Nonlinear expressions may contain only scalar expressions.
2

The nonlinear syntax is limited: Nonlinear Modeling · JuMP. In particular, you can’t use vector-valued expressions.

If it’s just the linear solve you need, why do you need the explicit inverse? If I understand correctly, your model is equivalent to:

using JuMP
model = Model()
@variable(model, x >= 0)
@variable(model, y >= 0)
@constraint(model, [y, x] .== [1 x; 0 1] * [5, 1])

This is linear, with no need for an explicit matrix inverse.

3

Thanks @odow.

The example here can be written as a linear form just for simplicity. What I was trying to do is related to a constraint that looks like A(x)B(x)^{-1}x = b, where A(x) and B(x) are two conformable matrices that contain the variable x.

Do you happen to know if there is another way to use vector-valued expressions, if it is not possible to be done in JuMP?

Thanks again!

4

Do you happen to know if there is another way to use vector-valued expressions, if it is not possible to be done in JuMP ?

You can use vector-valued linear or quadratic expressions in the @constraint macro. You just can’t use vector-valued expressions in the @NLxxx macros.

I would split out the inner inv(B(x)) * x term:

using JuMP
model = Model()
@variable(model, x[1:2] >= 0)
A(x) = [1 x; 0 1]
B(x) = [1 x; 0 1]
b = [1, 2]
@variable(model, y[1:2])
@constraint(model, B(x) * y .== x)
@constraint(model, A(x) * y .== b)

Now you have two non-convex quadratic constraints.