Unexpected behavior of objective_value with quadratic objective and semidefinite constraint

Specific Domains Optimization (Mathematical)
1

When solving a convex quadratic optimization problem with semidefinite constraints in Mosek, the objective_value function appears to be returning an objective value that is off by a constant.

Minimal working example: when I run the following block of code in Julia 1.9.1/JuMP 1.16.0/Mosek.jl 10.1.3

using JuMP, Mosek, MosekTools, Random
Random.seed!(12345)
a=rand(10)
@show a
m=Model(Mosek.Optimizer)
@variable(m, Y[1:10, 1:10], PSD)
@objective(m, Min, sum((Y[i,i]-a[i])^2 for i=1:10))
@constraint(m, sum(Y[i,i] for i=1:10)==1.0)
optimize!(m)
@show objective_value(m)
@show objective_value(m)-sum(a[i]^2 for i=1:10)
@show sum((value.(Y)[i,i]-a[i])^2 for i=1:10)

I get the following output:

a = [0.7918054038647908, 0.1595789996994108, 0.33419142747606, 0.8113922657011057, 0.7966292033537833, 0.9178138221831842, 0.31132682723725014, 0.7529056805990061, 0.6338479453532103, 0.8999510969166794]
Problem
  Name                   :                 
  Objective sense        : minimize        
  Type                   : CONIC (conic optimization problem)
  Constraints            : 1               
  Affine conic cons.     : 1 (12 rows)
  Disjunctive cons.      : 0               
  Cones                  : 0               
  Scalar variables       : 1               
  Matrix variables       : 1 (scalarized: 55)
  Integer variables      : 0               

Optimizer started.
Presolve started.
Linear dependency checker started.
Linear dependency checker terminated.
Eliminator started.
Freed constraints in eliminator : 0
Eliminator terminated.
Eliminator - tries                  : 1                 time                   : 0.00            
Lin. dep.  - tries                  : 1                 time                   : 0.00            
Lin. dep.  - primal attempts        : 1                 successes              : 1               
Lin. dep.  - dual attempts          : 0                 successes              : 0               
Lin. dep.  - primal deps.           : 0                 dual deps.             : 0               
Presolve terminated. Time: 0.00    
Optimizer  - threads                : 12              
Optimizer  - solved problem         : the primal      
Optimizer  - Constraints            : 12              
Optimizer  - Cones                  : 2               
Optimizer  - Scalar variables       : 14                conic                  : 14              
Optimizer  - Semi-definite variables: 1                 scalarized             : 55              
Factor     - setup time             : 0.00            
Factor     - dense det. time        : 0.00              GP order time          : 0.00            
Factor     - nonzeros before factor : 78                after factor           : 78              
Factor     - dense dim.             : 0                 flops                  : 2.03e+03        
ITE PFEAS    DFEAS    GFEAS    PRSTATUS   POBJ              DOBJ              MU       TIME  
0   9.0e+00  1.0e+00  1.7e+00  0.00e+00   0.000000000e+00   -7.071067812e-01  1.0e+00  0.00  
1   3.3e+00  3.6e-01  6.9e-01  -4.16e-01  1.056923953e+00   1.164658325e+00   3.6e-01  0.00  
2   1.4e+00  1.6e-01  5.5e-02  1.32e+00   2.582687886e+00   2.298362008e+00   1.6e-01  0.00  
3   2.4e-01  2.6e-02  3.1e-03  1.57e+00   3.201035910e+00   3.169139781e+00   2.6e-02  0.00  
4   6.0e-02  6.6e-03  5.0e-04  1.16e+00   3.247189070e+00   3.240395108e+00   6.6e-03  0.00  
5   1.5e-02  1.6e-03  6.5e-05  1.06e+00   3.260054465e+00   3.258530273e+00   1.6e-03  0.00  
6   2.3e-03  2.5e-04  4.4e-06  1.01e+00   3.262749573e+00   3.262535332e+00   2.5e-04  0.00  
7   3.0e-04  3.3e-05  2.2e-07  1.00e+00   3.263244710e+00   3.263217304e+00   3.3e-05  0.00  
8   2.8e-06  3.1e-07  2.0e-10  1.00e+00   3.263311335e+00   3.263311078e+00   3.1e-07  0.00  
9   1.1e-07  1.2e-08  1.6e-12  1.00e+00   3.263311979e+00   3.263311969e+00   1.2e-08  0.00  
10  1.1e-08  1.2e-09  5.0e-14  1.00e+00   3.263312005e+00   3.263312004e+00   1.2e-09  0.00  
Optimizer terminated. Time: 0.01    

objective_value(m) = 4.774929425505085
objective_value(m) - sum((a[i] ^ 2 for i = 1:10)) = 0.0
sum(((value.(Y)[i, i] - a[i]) ^ 2 for i = 1:10)) = 3.2633120002611253

This is quite strange because I would expect the third to last and last lines to match.

Does anyone know what’s going on? Thanks!

1 Like
2

I don’t have a Mosek license, so I can’t reproduce this, but it looks like a bug.

SCS works:

using JuMP, SCS, Random
Random.seed!(12345)
a = rand(10)
model = Model(SCS.Optimizer)
@variable(model, Y[1:10, 1:10], PSD)
@objective(model, Min, sum((Y[i,i] - a[i])^2 for i in 1:10))
@constraint(model, sum(Y[i,i] for i in 1:10) == 1)
optimize!(model)
@show objective_value(model)
@show objective_value(model) - sum(a.^2)
@show sum((value(Y[i,i]) - a[i])^2 for i in 1:10)

What are the outputs of solution_summary(m) and print_active_bridges(m)?

3

Great, thanks for the speedy reply!

Output from solution_summary(m)

* Solver : Mosek

* Status
  Result count       : 1
  Termination status : OPTIMAL
  Message from the solver:
  "Mosek.MSK_SOL_STA_OPTIMAL"

* Candidate solution (result #1)
  Primal status      : FEASIBLE_POINT
  Dual status        : FEASIBLE_POINT
  Objective value    : 3.26331e+00
  Objective bound    : 3.26331e+00
  Relative gap       : 0.00000e+00
  Dual objective value : 3.26331e+00

* Work counters
  Solve time (sec)   : 9.06205e-03
  Simplex iterations : 0
  Barrier iterations : 7
  Node count         : 0

Output from print_active_bridges(m) (this throws an error so I’m guessing it’s the issue)

* Unsupported objective: MOI.ScalarQuadraticFunction{Float64}
 |  bridged by:
 |   MOIB.Objective.SlackBridge{Float64, MOI.ScalarQuadraticFunction{Float64}, MOI.ScalarQuadraticFunction{Float64}}
 |  may introduce:
 |   * Unsupported objective: MOI.VariableIndex
 |   |  bridged by:
 |   |   MOIB.Objective.FunctionizeBridge{Float64}
 |   |  may introduce:
 |   |   * Supported objective: MOI.ScalarAffineFunction{Float64}
 |   * Unsupported constraint: MOI.ScalarQuadraticFunction{Float64}-in-MOI.GreaterThan{Float64}
 |   |  bridged by:
 |   |   MOIB.Constraint.QuadtoSOCBridge{Float64}
 |   |  may introduce:
 |   |   * Supported constraint: MOI.VectorAffineFunction{Float64}-in-MOI.RotatedSecondOrderCone
 |   * Unsupported constraint: MOI.ScalarQuadraticFunction{Float64}-in-MOI.LessThan{Float64}
 |   |  bridged by:
 |   |   MOIB.Constraint.QuadtoSOCBridge{Float64}
 |   |  may introduce:
 |   |   * Supported constraint: MOI.VectorAffineFunction{Float64}-in-MOI.RotatedSecondOrderCone
 |   * Supported variable: MOI.Reals
 * Supported constraint: MOI.ScalarAffineFunction{Float64}-in-MOI.EqualTo{Float64}
 * Supported constraint: MOI.VariableIndex-in-MOI.GreaterThan{Float64}
4

What is the error? This seems like a bug in how we’re reformulating the objective. But @blegat is the person to track it down.

For now, I’d use a different solver, or use the epigraph reformulation:

using JuMP, SCS, Random
Random.seed!(12345)
a = rand(10)
model = Model(SCS.Optimizer)
@variable(model, Y[1:10, 1:10], PSD)
# @objective(model, Min, sum((Y[i,i] - a[i])^2 for i in 1:10))
@variable(model, t)
@objective(model, Min, t)
@constraint(model, vcat(t, 0.5, [Y[i,i] - a[i] for i in 1:10]) in RotatedSecondOrderCone())
@constraint(model, sum(Y[i,i] for i in 1:10) == 1)
optimize!(model)
@show objective_value(model)
@show sum((value(Y[i,i]) - a[i])^2 for i in 1:10)
5

I don’t have time to investigate at the moment but I’ve opened an issue as a reminder:

Can you describe the error message you’re mentioning in your previous post ?

6

I think Ryan just interpreted the print output from print_active_bridges as an error.