Hi everybody!
I’m trying to find a way to obtain the dependent variables of a nonlinear system as GenericAffExpr, expressed as a function only of the independent variables of an optimization problem.

To be more precise, I have a nonlinear problem with independent variables x. The objective function is a function of both x and other dependent variables y. In order to compute the OF, is necessary to obtain the expressions of y in function of x (in type GenericAffExpr, I suppose).
These y can be found by solving a nonlinear system of equations and to do that, I wrote a second optimization model that (in my intentions) should execute this task.

To be clearer, below is reported a minimal example of what I would like to do.

using JuMP, Ipopt
# define optimization problem with independent variables (x)
m = Model(Ipopt.Optimizer)
@variable(m, x[1:2])
# solve nonlinear system to find dependent variables (y)
p = Model(Ipopt.Optimizer)
@variable(p, y[1:2])
@NLconstraint(p, x[1] * y[1] * y[2] == 2)
@NLconstraint(p, x[2] * y[1] * y[2] == 5)
@objective(p, 0.0)
JuMP.optimize!(p)
# solve original problem
@NLobjective(m, x[1] + x[2]^2 + sum(y))
JuMP.optimize!(m)

However, trying to define the ‘inner’ optimization p, JuMP provides the error Variable in nonlinear expression does not belong to the corresponding model in correspondence of the first @NLconstraint, since the variables x are inserted in the problem p but belong to the problem m.

My question is: does anyone know if exist a suitable way to solve this issue? Or, does exist another strategy to reach my purpose?

In a certain sense yes, but, since the variables x are present in the problem p and the problem p only optimize the variables y, the values of y should contain the terms of x. In other words, I would like that the problem p has the following outputs:

Ah, ok, I took another look and you’re right that wouldn’t work. You can’t mix-and-match variables and models this way. The thing to do is use a single model that contains all the necessary variables (you can add/delete constraints/variables and change the objective as necessary).

m = Model(Ipopt.Optimizer)
@variable(m, x[1:2])
@variable(m, y[1:2])
c1 = @NLconstraint(m, x[1] * y[1] * y[2] == 2)
c2 = @NLconstraint(m, x[2] * y[1] * y[2] == 5)
@objective(m, 0.0) # (note this isn't necessary)
JuMP.optimize!(m)
# should you check termination_status(m) before proceeding?
# solve original problem:
# get rid of the constraints if they don't matter any more
delete!.(m, (c1, c2))
@NLobjective(m, x[1] + x[2]^2 + sum(value.(y)))
JuMP.optimize!(m)

Alternatively, you could create a new model with only variables called x and set its objective based on the other model’s value for y.