osqp / OSQP.jl

Julia interface for OSQP: The Operator Splitting QP Solver
https://osqp.org/
Other
68 stars 25 forks source link

OSQP violates constraints for a quadratic programming problem #47

Open JackDevine opened 5 years ago

JackDevine commented 5 years ago

Hi all,

Thanks for the hard work put into OSQP.

When using OSQP.jl from JuMP.jl to solve a quadratic programming problem, I find that OSQP does not satisfy the constraints that I require.

The first bit of code is just some code to make some random labeled data, where each data point is a 2D vector in $[-1,1]\times[-1,1]$. The points are separated by a straight line. You can ignore this code, you just need it to be able to run the example that I give

function create_line(line_start,line_end)
    m = (line_end[2]-line_start[2])/(line_end[1]-line_start[1])
    c = line_start[2]-m*line_start[1]

    [m,c], x -> m*x+c
end

function classify(line::Function,point::AbstractArray)
    ifelse(line(point[1]) > point[2],-1,1)
end

function classify(weights::AbstractArray,point::AbstractArray)
    sign(weights'*point)
end

function create_data!(data,labels)
    N = size(data)[1]
    data_valid = false
    while !data_valid
        data .= 2rand(N,2).-1
        line_start = 2rand(1,2).-1
        line_end = 2rand(1,2).-1

        global line_coefficients,line = create_line(line_start,line_end)
        for i in 1:N
            labels[i] = classify(line,data[i,:])
        end

        !(all(labels .== 1) || all(labels .== -1)) && (data_valid = true)
    end

    data,labels,line
end

The next bit of code is the important bit:

using LinearAlgebra
using JuMP
using OSQP

N = 100
data = Array{Float64}(undef,N,2)
labels = Array{Int64}(undef,N)
data,labels,line = create_data!(data,labels)

m = Model(with_optimizer(OSQP.Optimizer))
# m = Model(with_optimizer(Ipopt.Optimizer,print_level=0))
@variable(m,α[1:N])
Q = [labels[i]*labels[j]*dot(data[i,:],data[j,:]) for i in 1:N, j in 1:N]
@objective(m,Min,0.5α'*Q*α-ones(N)'*α)
@constraint(m,[i=1:N],0.0 <= α[i])
@constraint(m,dot(labels,α) == 0)

@time optimize!(m)

αsol = value.(α)
minimum(αsol)

Upon running this, I usually find that the minimum value of α is much less than zero, which is a violation of the constraint that I gave. If I use Ipopt and uncomment the commented line, then the constraint will be satisfied.

julia> versioninfo()
Julia Version 1.1.0
Commit 80516ca202 (2019-01-21 21:24 UTC)
Platform Info:
  OS: Linux (x86_64-pc-linux-gnu)
  CPU: Intel(R) Core(TM) i7-9700K CPU @ 3.60GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, skylake)

(v1.1) pkg> st OSQP
    Status `~/.julia/environments/v1.1/Project.toml`
  [ab2f91bb] OSQP v0.5.1

(v1.1) pkg> st JuMP
    Status `~/.julia/environments/v1.1/Project.toml`
  [f6369f11] ForwardDiff v0.10.3
  [4076af6c] JuMP v0.19.0
bstellato commented 5 years ago

This is a problem with adaptive_rho. If you disable it, OSQP converges without problems. I have recreated your script in python and opened an issue in the main repository about this: https://github.com/oxfordcontrol/osqp/issues/151

odow commented 2 years ago

This seems like it can be closed in favor of the upstream issue. It isn't a bug in OSQP.jl.