Non-convex optimization

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
  1. How can I make use of NLparameter and NLexpression or at least NLexpression for the example above?
  2. 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.

I recently revamped the NLP documentation. Here’s a sneak preview of what will be included in the next release: Nonlinear Modeling · JuMP.

I didn’t test this, so there may be typos, but you want something like:

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)
    
    @NLparameter(model, p1[j=1:k, i=1:n] == c1[j, i])
    @NLparameter(model, p2[j=1:k, i=1:n] == c1[j, i])

    @NLexpression(model, xp, sum(x[i] * p[i] for i in 1:n))
    @NLexpression(model, c1_expr[j=1:k], sum(p1[j, i] * x[i] for i in 1:n))
    @NLexpression(model, c2_expr[j=1:k], sum(p2[j, i] * x[i] for i in 1:n))

    @NLconstraint(model, con3, xp - sum(b[j]       / (1+var1)^j for j in 1:k) == 0)
    @NLconstraint(model, con5, xp - sum(c1_expr[j] / (1+var2)^j for j in 1:k) == 0)
    @NLconstraint(model, con6, xp - sum(c2_expr[j] / (1+var3)^j for j in 1:k) == 0)

    @constraint(model, con1, c1*x .>= b)
    return model
end
1 Like

Thank you

Updates to the documentation looks really good.

One area where it will be useful to have more examples in the documentation is DenseAxisArray, e.g., how to write this non-linear expression if x is defined with non-integer indices: @NLexpression(model, c1_expr[j=1:k], sum(p1[j, i] * x[i] for i in 1:n))