I am using Julia v1.0.2 with JuMP v0.18.5 and Ipopt v3.12.12. After solving my problem with juMP I want to estimate my objective to know how good is my solution (closer to 0 is better).
The problem is that the objective in the Ipopt solver is super different from the one given by getobjectivevalue() and the one I compute manually (given the solutions) while I expect all three to be equal.
-I suspect that the objective displayed by Ipopt has some sort of offset because it is far from zero. Would you know how to suppress it?
With the formulation 2, this offset disappear and objective of Ipopt goes to zero.
-In formulation 1, I sometime have a small but negative value with getobjectivevalue() which should not be possible because objective is a sum of square terms.
-getobjectivevalue() and the manually computed objective are both close to zero but different in general. Which objective should I believe?
-Also the two formulation should be equivalent, which one should I trust more?
I reproduced the problem I have with this reduced example. I try two formulations:
Formulation 1 & 2: I have square variables that I want to square in the objective, so I use auxilary variable
+
Formulation 1: The affine equations are directly squared in the objective
Formulation 2: The affine equations are equal to an auxiliary variable that is squared in the objective
using JuMP,Ipopt
N=5;
L=50;
xi=randn(N,L)
println("######Formulation 1########\n\n\n")
model = Model(solver = IpoptSolver())
@variable(model, x[1:N,1:L])
@variable(model, y[1:N,1:L])
#Quadratic equations
NL_const = Array{JuMP.GenericQuadExpr,2}(undef, N,L)
for n=1:N,k=1:L
NL_const[n,k] = 3x[n,k]^2
end
@variable(model, auxNL[1:N,1:L])
for n=1:N,k=1:L
@constraint(model, auxNL[n,k] == NL_const[n,k])
end
#Linear equations
L_const = sum( (2y[n,k] + 3x[n,k] + xi[n,k])^2 for n in 1:N, k in 1:L)
#objective
@objective(model, Min, sum(auxNL.^2) + L_const)
@time sol = solve(model)
println("Objective value with getobjectivevalue(): ", getobjectivevalue(model))
xsol1 = getvalue(x);
ysol1 = getvalue(y);
println("Objective value computed manually: ", sum( (2ysol1[n,k] + 3xsol1[n,k] + xi[n,k])^2 for n in 1:N, k in 1:L) + sum( (3xsol1[n,k]^2)^2 for n in 1:N, k in 1:L))
println("\n\n\n######Formulation 2########\n\n\n")
model = Model(solver = IpoptSolver())
@variable(model, x[1:N,1:L])
@variable(model, y[1:N,1:L])
#Quadratic equations
NL_const = Array{JuMP.GenericQuadExpr,2}(undef, N,L)
for n=1:N,k=1:L
NL_const[n,k] = 3x[n,k]^2
end
@variable(model, auxNL[1:N,1:L])
for n=1:N,k=1:L
@constraint(model, auxNL[n,k] == NL_const[n,k])
end
#Linear equations
L_const = Array{JuMP.GenericAffExpr,2}(undef, N,L)
for n=1:N,k=1:L
L_const[n,k] = 2y[n,k] + 3x[n,k] + xi[n,k]
end
@variable(model, auxL[1:N,1:L])
for n=1:N,k=1:L
@constraint(model, auxL[n,k] == L_const[n,k])
end
#objective
@objective(model, Min, sum(auxNL.^2) + sum(auxL.^2))
@time sol = solve(model)
println("Objective value with getobjectivevalue(): ", getobjectivevalue(model))
xsol2 = getvalue(x);
ysol2 = getvalue(y);
println("Objective value computed manually: ", sum( (2ysol2[n,k] + 3xsol2[n,k] + xi[n,k])^2 for n in 1:N, k in 1:L) + sum( (3xsol2[n,k]^2)^2 for n in 1:N, k in 1:L))
Here is a typical output (it depends on the random matrix xi)
######Formulation 1########
This is Ipopt version 3.12.12, running with linear solver ma27.
Number of nonzeros in equality constraint Jacobian...: 750
Number of nonzeros in inequality constraint Jacobian.: 0
Number of nonzeros in Lagrangian Hessian.............: 1500
Total number of variables............................: 750
variables with only lower bounds: 0
variables with lower and upper bounds: 0
variables with only upper bounds: 0
Total number of equality constraints.................: 250
Total number of inequality constraints...............: 0
inequality constraints with only lower bounds: 0
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 0
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 0.0000000e+00 0.00e+00 1.69e+01 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 -2.6797426e+02 1.27e+00 6.50e-05 -1.0 6.50e-01 -4.0 1.00e+00 1.00e+00f 1
2 -2.6797426e+02 3.17e-01 3.11e-15 -1.0 4.87e-01 - 1.00e+00 1.00e+00h 1
3 -2.6797426e+02 7.92e-02 2.72e-15 -1.7 2.44e-01 - 1.00e+00 1.00e+00h 1
4 -2.6797426e+02 1.98e-02 2.02e-15 -2.5 1.22e-01 - 1.00e+00 1.00e+00h 1
5 -2.6797426e+02 4.95e-03 1.88e-15 -3.8 6.09e-02 - 1.00e+00 1.00e+00h 1
6 -2.6797426e+02 1.24e-03 2.54e-15 -3.8 3.05e-02 - 1.00e+00 1.00e+00h 1
7 -2.6797426e+02 3.09e-04 2.53e-15 -5.7 1.52e-02 - 1.00e+00 1.00e+00h 1
8 -2.6797426e+02 7.73e-05 2.14e-15 -5.7 7.62e-03 - 1.00e+00 1.00e+00h 1
9 -2.6797426e+02 1.93e-05 2.85e-15 -5.7 3.81e-03 - 1.00e+00 1.00e+00h 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
10 -2.6797426e+02 4.83e-06 2.49e-15 -5.7 1.90e-03 - 1.00e+00 1.00e+00h 1
11 -2.6797426e+02 1.21e-06 1.92e-15 -8.6 9.52e-04 - 1.00e+00 1.00e+00h 1
12 -2.6797426e+02 1.35e-08 3.35e-09 -8.6 1.00e-04 -4.5 1.00e+00 1.00e+00h 1
13 -2.6797426e+02 2.60e-08 1.55e-09 -9.0 1.40e-04 -5.0 1.00e+00 1.00e+00h 1
14 -2.6797426e+02 3.42e-08 7.51e-10 -9.0 1.60e-04 -5.4 1.00e+00 1.00e+00h 1
15 -2.6797426e+02 2.90e-08 3.37e-10 -9.0 1.48e-04 -5.9 1.00e+00 1.00e+00h 1
16 -2.6797426e+02 1.91e-08 1.28e-10 -9.0 1.20e-04 -6.4 1.00e+00 1.00e+00h 1
17 -2.6797426e+02 1.08e-08 4.35e-11 -9.0 8.98e-05 -6.9 1.00e+00 1.00e+00h 1
18 -2.6797426e+02 5.50e-09 1.36e-11 -9.0 6.42e-05 -7.3 1.00e+00 1.00e+00h 1
Number of Iterations....: 18
(scaled) (unscaled)
Objective...............: -2.6797425556154758e+02 -2.6797425556154758e+02
Dual infeasibility......: 1.3560987184059177e-11 1.3560987184059177e-11
Constraint violation....: 5.4983582362510339e-09 5.4983582362510339e-09
Complementarity.........: 0.0000000000000000e+00 0.0000000000000000e+00
Overall NLP error.......: 5.4983582362510339e-09 5.4983582362510339e-09
Number of objective function evaluations = 19
Number of objective gradient evaluations = 19
Number of equality constraint evaluations = 19
Number of inequality constraint evaluations = 0
Number of equality constraint Jacobian evaluations = 19
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations = 18
Total CPU secs in IPOPT (w/o function evaluations) = 0.009
Total CPU secs in NLP function evaluations = 0.001
EXIT: Optimal Solution Found.
0.021784 seconds (14.67 k allocations: 634.977 KiB, 30.99% gc time)
Objective value with getobjectivevalue(): -8.526512829121202e-13
Objective value computed manually: 3.0007444607343187e-14
######Formulation 2########
This is Ipopt version 3.12.12, running with linear solver ma27.
Number of nonzeros in equality constraint Jacobian...: 1500
Number of nonzeros in inequality constraint Jacobian.: 0
Number of nonzeros in Lagrangian Hessian.............: 750
Total number of variables............................: 1000
variables with only lower bounds: 0
variables with lower and upper bounds: 0
variables with only upper bounds: 0
Total number of equality constraints.................: 500
Total number of inequality constraints...............: 0
inequality constraints with only lower bounds: 0
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 0
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 0.0000000e+00 2.82e+00 0.00e+00 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 3.9636893e-09 1.27e+00 6.50e-05 -1.0 6.50e-01 -4.0 1.00e+00 1.00e+00h 1
2 1.2965119e-40 3.17e-01 4.19e-21 -1.0 4.87e-01 - 1.00e+00 1.00e+00h 1
3 3.2412798e-41 7.92e-02 1.05e-21 -1.7 2.44e-01 - 1.00e+00 1.00e+00h 1
4 8.1031995e-42 1.98e-02 2.62e-22 -2.5 1.22e-01 - 1.00e+00 1.00e+00h 1
5 2.0257999e-42 4.95e-03 6.55e-23 -3.8 6.09e-02 - 1.00e+00 1.00e+00h 1
6 5.0644997e-43 1.24e-03 1.64e-23 -3.8 3.05e-02 - 1.00e+00 1.00e+00h 1
7 1.2661249e-43 3.09e-04 4.10e-24 -5.7 1.52e-02 - 1.00e+00 1.00e+00h 1
8 3.1653123e-44 7.73e-05 1.02e-24 -5.7 7.62e-03 - 1.00e+00 1.00e+00h 1
9 7.9132807e-45 1.93e-05 2.56e-25 -5.7 3.81e-03 - 1.00e+00 1.00e+00h 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
10 1.9783202e-45 4.83e-06 6.40e-26 -5.7 1.90e-03 - 1.00e+00 1.00e+00h 1
11 4.9458005e-46 1.21e-06 1.60e-26 -8.6 9.52e-04 - 1.00e+00 1.00e+00h 1
12 1.2364501e-46 3.02e-07 4.00e-27 -8.6 4.76e-04 - 1.00e+00 1.00e+00h 1
13 3.0911253e-47 7.55e-08 1.00e-27 -8.6 2.38e-04 - 1.00e+00 1.00e+00h 1
14 7.7278132e-48 1.89e-08 2.50e-28 -8.6 1.19e-04 - 1.00e+00 1.00e+00h 1
15 1.9319533e-48 4.72e-09 6.25e-29 -9.0 5.95e-05 - 1.00e+00 1.00e+00h 1
Number of Iterations....: 15
(scaled) (unscaled)
Objective...............: 1.9319533031007233e-48 1.9319533031007233e-48
Dual infeasibility......: 6.2490981135057184e-29 6.2490981135057184e-29
Constraint violation....: 4.7193765244310251e-09 4.7193765244310251e-09
Complementarity.........: 0.0000000000000000e+00 0.0000000000000000e+00
Overall NLP error.......: 4.7193765244310251e-09 4.7193765244310251e-09
Number of objective function evaluations = 16
Number of objective gradient evaluations = 16
Number of equality constraint evaluations = 16
Number of inequality constraint evaluations = 0
Number of equality constraint Jacobian evaluations = 16
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations = 15
Total CPU secs in IPOPT (w/o function evaluations) = 0.005
Total CPU secs in NLP function evaluations = 0.000
EXIT: Optimal Solution Found.
0.009841 seconds (11.39 k allocations: 661.414 KiB)
Objective value with getobjectivevalue(): 1.9319533031007233e-48
Objective value computed manually: 2.9111864471415683e-16
You see here
Formulation 1:
Objective Ipopt...............: -2.6797425556154758e+02 (not close to zero, must be an offset?? But why?)
Objective value with getobjectivevalue(): -8.526512829121202e-13 (Negative!! But close to zero)
Objective value computed manually: 3.0007444607343187e-14 (makes sense = close to zero and positive)
Formulation 2:
Objective Ipopt...............: 1.9319533031007233e-48 (Makes sense but super small)
Objective value with getobjectivevalue():1.9319533031007233e-48 (Makes sense and equal to Ipopt objective)
Objective value computed manually: 2.9111864471415683e-16 (Makes sense but different from the other two)