I am solving a non linear optimization using Ipopt. The problem is non-convex but Ipopt gives a good enough answer.
The documentation on use of nonlinear parameters and non-linear expressions is somewhat thin.
Suppose I have: sum(c1[j,i] *x[i] / (1+var2)^j for i in 1:n, j in 1:k)
where c1 is a j x i matrix of constants
x is a vector n x 1
that I want to write as expression. Here is a MWE without expression and NLparameters:
using Ipopt, JuMP
function nlp_model(p, c1, c2,b, k, n)
model = Model(Ipopt.Optimizer)
@variable(model, 0 <= x[i = 1:n] <= 1)
@variable(model, 0<=var1<=1)
@variable(model, 0<=var2<=1)
@variable(model, 0<=var3<=1)
@objective(model,Max,var1-var2+var3)
@NLconstraint(model,con3,sum(x[i]*p[i] for i in 1:n)-sum(b[j]/(1+var1)^j for j in 1:k) == 0)
@NLconstraint(model,con5,sum(x[i]*p[i] for i in 1:n)-sum(c1[j,i] *x[i] / (1+var2)^j for i in 1:n, j in 1:k) == 0)
@NLconstraint(model,con6,sum(x[i]*p[i] for i in 1:n)-sum(c2[j,i] *x[i] / (1+var3)^j for i in 1:n, j in 1:k) == 0)
@constraint(model, con1, c1*x .>= b)
JuMP.optimize!(model)
println("")
println("variable l = $(value(var1)) , variable 2 = $(value(var2)) , variable 3 = $(value(var3))")
end
function data()
p = rand(400:700,5,1); c1= (rand(100:200,4,5))'; c2 = 0.9 .* c1
b =rand(150:250,5,1); n= size(c1,2); k = size(c1,1)
return p, c1, c2, b, n, k
end
p, c1, c2, b, n, k = data()
nlp_model(p, c1, c2,b, k, n)
Here is my attempt at writing it using NLparameter
and NLexpression
function nlp_model_expr(p, c1, c2,b, k, n)
model = Model(Ipopt.Optimizer)
@variable(model, 0<= x[i = 1:n] <= 1)
@variable(model, 0<=var1<=1)
@variable(model, 0<=var2<=1)
@variable(model, 0<=var3<=1)
@NLparameter(model,para1[j=1:k,i=1:n]=c1[j,i])
@NLparameter(model,para2[j=1:k,i=1:n]=c2[j,i])
@NLexpression(model,my_expr1[j=1:k],para1[j,:]*x)
@NLexpression(model,my_expr2[j=1:k], para2[j,:]*x)
@NLconstraint(model,con3,sum(x[i]*p[i] for i in 1:n)-sum(b[j]/(1+var1)^j for j in 1:k) == 0)
@NLconstraint(model,con5,sum(x[i]*p[i] for i in 1:n)-sum(my_expr1[j=1:k] / (1+var2)^j for j in 1:k) == 0)
@NLconstraint(model,con6,sum(x[i]*p[i] for i in 1:n)-sum(my_expr2[j=1:k] *x[i] / (1+var3)^j for j in 1:k) == 0)
@constraint(model, con1, c1*x .>= b)
JuMP.optimize!(model)
println("")
println("variable l = $(value(var1)) , variable 2 = $(value(var2)) , variable 3 = $(value(var3))")
end
- How can I make use of
NLparameter
andNLexpression
or at leastNLexpression
for the example above? - The code above without NLexpression works but when I use DenseAxisArray then I can’t use
sum(c1[j,i] *x[i] / (1+var2)^j for i in 1:n, j in 1:k)
in constraints. I don’t know if there is workaround for it.