# JuMP : not getting the correct solution of the optimization problem

Dear everyone

I am solving a minimization problem under equality and inequality constraints. I solve this problem with COSMO and I get a solution which is the vector alpha_opt=[0 0 0]
But I found another vector which is alpha_test=[1 1 0] which also respects the constraints and such that the objective function is equal to -5.5<0

julia> sum(0.5*alpha_opt[k]*G[k,i]*y[k]*y[i]*alpha_opt[i] -alpha_opt[k] for i=1:n for k=1:n)
0.0

julia> sum(0.5*alpha_test[k]*G[k,i]*y[k]*y[i]*alpha_test[i] -alpha_test[k] for i=1:n for k=1:n)
-5.5

However, I do not understand what I did wrong when I defined my optimization problem.

Could you show what’s in alpha_opt, G? Based on your question on Slack, your constraints are all linear, is that right? Your objective function is probably nonconvex and your problem has local minima (and/or saddle points). Most solver only look for a stationary point. In order to find the point you’re looking for, you can start your solver from a different initial guess or use a global optimization solver.

1 Like

1 Like

@dpo : this is right : my constraints are linear and my objective function is quadratic.

Here is the minimal working example

using LinearAlgebra
using JuMP
using COSMO

X=[1 1; 1 2; 2 2 ];
y=[+1 -1 -1 ];
n=length(y);
p=size(X,2);
G=X*X';

model = Model(with_optimizer(COSMO.Optimizer))
@variables model begin
alpha[i=1:length(y)]
end
@NLobjective(model, Min, sum(0.5*alpha[k]*G[k,i]*y[k]*y[i]*alpha[i] -alpha[k] for i=1:n for k=1:n))
@constraint(model, sum( y[i] * alpha[i] for i=1:n ) == 0)
@constraint(model, [i=1:n], -alpha[i]<=0)

optimize!(model);
alpha_opt = value.(alpha);


And here is the output I get :

julia> optimize!(model);

      COSMO v0.5.0 - A Quadratic Objective Conic Solver
Michael Garstka
University of Oxford, 2017 - 2019


Problem: x ∈ R^{3},
constraints: A ∈ R^{4x3} (6 nnz),
matrix size to factor: 7x7 (19 nnz)
Sets: Nonnegatives of dim: 3
ZeroSet of dim: 1
Settings: ϵ_abs = 1.0e-04, ϵ_rel = 1.0e-04,
ϵ_prim_inf = 1.0e-06, ϵ_dual_inf = 1.0e-04,
ρ = 0.1, σ = 1.0e-6, α = 1.6,
max_iter = 2500,
scaling iter = 10 (on),
check termination every 40 iter,
check infeasibility every 40 iter,
KKT system solver: QDLDL
Setup Time: 1967.41ms

Iter: Objective: Primal Res: Dual Res: Rho:
40 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e-01

Results
Status: Solved
Iterations: 40
Optimal objective: 0.0
Runtime: 2.754s (2754.23ms)

julia> alpha_opt = value.(alpha)
3-element Array{Float64,1}:
0.0
0.0
0.0

I was expecting julia to give me this solution : [2 ,2 ,0].

alpha = zeros(3) doesn’t seem like a solution. The feasible set is given by \alpha = Zu for u \geq 0, where Z = \begin{bmatrix} 1 & 1 \\ 1 & 0 \\ 0 & 1 \end{bmatrix}. The objective by f(\alpha) = \tfrac{1}{2}\alpha^T H \alpha + \alpha^T g, where H = y' .* G .* y and g = -3*ones(3). The reduced subproblem f(Zu) has Hessian
Z^T H Z = \begin{bmatrix} 1 & 1 \\ 1 & 2 \end{bmatrix},
and
Z^T g = \begin{bmatrix} -6 \\ -6 \end{bmatrix}. This implies the solution is u = (6,0)^T, i.e., \alpha = Zu = (6,6,0)^T, and looks like it’s the only one because the problem is convex.

I ran with IPOPT (using NLPModelsIpopt.jl) and found that solution. I’m guessing changing COSMO to IPOPT might also give that same solution.

I’ve also noticed your problem looks like the dual of a SVM. In that case, the -alpha[k] in the sum should be outside the for i, otherwise you’re adding them up for each i. If that change is made, your expected solution is found.

COSMO is a conic solver. I’m surprised it works at all with @NLobjective. I’m guessing that there is a bug in the wrapper that causes the nonlinear objective to be completely ignored (as you see the objective value is reported to be zero).
CC @migarstka

You should use @objective instead of @NLobjective.

2 Likes

Indeed COSMO.jl doesn’t support @NLobjective and should throw a warning when people try to use it in their models. I will fix the bug in the MOI wrapper.
As Miles pointed out, if you change @NLobjective in your example to @objective you should get the correct result.