Translation between SCS status and MOI status

[blegat]

@bodono Could you help me determine the mapping between SCS status and MOI status ?
The termination statuses are described here and result statuses here.
It is especially not easy to determine which of SlowProgress, AlmostSuccess, NumericalError or one of the limit is the reason the solver stopped.
For the result status, it is important to distinguish between FeasiblePoint and NearlyFeasiblePoint. Even if the solution is not optimal, it could be guaranteed to be feasible. Same for infeasibility certificates.
I have a first draft here: https://github.com/JuliaOpt/MathOptInterfaceSCS.jl/blob/master/src/attributes.jl
I have done it quickly and only based it on the names of the C constants so feel free to change it completely :-P

[mlubin]

IIRC, the "Inaccurate" states can also occur when an iteration limit is hit, without any guarantees on the quality of the solution. This doesn't count as AlmostSuccess.

[bodono]

For

    if s in (-3, 1, 2)
        MOI.FeasiblePoint

The SCS status -3 SCS_INDETERMINATE honestly never happens, it's there in the case that SCS converges to all zeros, which is a possible solution to the homogeneous set of equations SCS is solving (though we prove it can theoretically never happen in the paper under normal conditions), but technically speaking if that ever happens it won't be a feasible solution that is returned.

In

    elseif s in (-6, -1)
        MOI.InfeasibilityCertificate
    else
        MOI.InfeasiblePoint
    end

and

    elseif s in (-7, -2)
        MOI.InfeasibilityCertificate
    else
        MOI.InfeasiblePoint
    end

I'm not sure what the difference between these InfeasibilityCertificate and InfeasiblePoint are, but for SCS -6, -1 refers to primal unbounded, dual infeasible, and -2, -7 is primal infeasible, dual unbounded.

Other than that it looks right. One small note - I think it's hard to talk about the difference between feasible and nearly feasible, unless it is up to some defined tolerance. Even just solving a set of linear equations will have some numerical error, e.g. I just ran this python code:

>>> import numpy as np
>>> A = np.random.randn(100, 100)
>>> b = np.random.randn(100)
>>> x = np.linalg.solve(A, b)
>>> np.linalg.norm(A.dot(x) - b)
1.2380685665598168e-12

So the point x which is calculated as A \ b is technically not feasible for the constraint A x == b since there is a residual of 1e-12. Since solving a cone program is likely to include satisfying a linear constraint the same issue of numerical tolerance would apply there. However, if you're willing to take some small residual as still being feasible, then you can plug whatever that tolerance is into SCS and it will find a point with that tolerance given sufficient time (if one exists).

[mlubin]

InfeasiblePoint should never happen in SCS. It's for non-convex solvers that can converge (e.g., accidentally) to infeasible points.

[mlubin]

@bodono, feasibility and near feasibility are all defined with respect to tolerances. A solver should call a point (primal/dual) feasible if it satisfies the pre-defined (primal/dual) feasibility tolerances. Some solvers have a second level of relaxed feasibility tolerances, so if the solution satisfies these relaxed tolerances then it can be called near feasible.

[blegat]

@bodono Note that SCS has some trouble solving linear9 test of MOI. Before SCS 2.0, it used to return SCS_SOLVED but with a low accuracy (it needed a tolerance of 1e-3 to pass the test). Now it returns SCS_SOLVED_INACCURATE with an even lower accuracy (it needs a tolerance of 5e-1). This seems to be caused be large coefficients. See #3

[bodono]

Ok, good to know. I can't find that test anywhere, can you point me to it?

SCS 2.0 is an improvement (in terms of wall clock time) for most problems, but undoubtedly some will do worse. If you want you can switch back to the old behavior for that test by setting acceleration_lookback = 0.

[bodono]

I found the test (I think) https://github.com/JuliaOpt/MathOptInterfaceSCS.jl/blob/master/test/runtests.jl#L13

I'll dig into and see what's happening.

[blegat]

The test is here. You can execute it by replacing the line MOIT.contlineartest(solver, config, ["linear9"]) by MOIT.linear9test(solver, config) in runtests.jl

[bodono]

I put this into cvxpy in python (easier for me access), and SCS can solve it very easily, both SCS 1 and SCS 2, here is the code:

x = cvxpy.Variable(1)
y = cvxpy.Variable(1)

objective = cvxpy.Maximize(1000 * x + 350 * y)
constraints = [x >= 30, y>=0, x - 1.5 * y >=0, 12 * x + 8 * y <= 1000, 1000*x + 300 * y <= 70000]
prob = cvxpy.Problem(objective, constraints)

prob.solve(solver='SCS', verbose=True)

and the resulting output:

----------------------------------------------------------------------------
	SCS v2.0.0 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012-2017
----------------------------------------------------------------------------
Lin-sys: sparse-indirect, nnz in A = 8, CG tol ~ 1/iter^(2.00)
eps = 1.00e-05, alpha = 1.50, max_iters = 5000, normalize = 1, scale = 1.00
acceleration_lookback = 20, rho_x = 1.00e-03
Variables n = 2, constraints m = 5
Cones:	linear vars: 5
Setup time: 1.09e-02s
----------------------------------------------------------------------------
 Iter | pri res | dua res | rel gap | pri obj | dua obj | kap/tau | time (s)
----------------------------------------------------------------------------
     0| 2.61e+00  5.54e-01  3.91e-01 -2.57e+05 -1.13e+05  3.67e-11  4.46e-01 
   100| 4.61e-03  2.50e-03  2.50e-03 -7.13e+04 -7.17e+04  3.36e-12  4.49e-01 
   200| 2.11e-02  1.81e-02  1.61e-02 -7.04e+04 -7.27e+04  1.06e-11  4.51e-01 
   280| 2.57e-08  5.75e-08  1.15e-08 -7.18e+04 -7.18e+04  2.67e-12  4.54e-01 
----------------------------------------------------------------------------
Status: Solved
Timing: Solve time: 4.54e-01s
	Lin-sys: avg # CG iterations: 1.85, avg solve time: 1.47e-04s
	Cones: avg projection time: 5.85e-08s
	Acceleration: avg step time: 5.80e-06s
----------------------------------------------------------------------------
Error metrics:
dist(s, K) = 1.4083e-15, dist(y, K*) = 0.0000e+00, s'y/|s||y| = 4.0557e-15
|Ax + s - b|_2 / (1 + |b|_2) = 2.5747e-08
|A'y + c|_2 / (1 + |c|_2) = 5.7507e-08
|c'x + b'y| / (1 + |c'x| + |b'y|) = 1.1546e-08
----------------------------------------------------------------------------
c'x = -71818.1842, -b'y = -71818.1859
============================================================================

I was cribbing from the comment in the code you linked, which suggests that SCS is getting the right answer basically

    #   maximize 1000 x + 350 y
    #
    #       s.t.                x >= 30
    #                           y >= 0
    #                 x -   1.5 y >= 0
    #            12   x +   8   y <= 1000
    #            1000 x + 300   y <= 70000
    #
    #   solution: (59.0909, 36.3636)
    #   objv: 71818.1818
    #@test MOI.supportsproblem(solver, MOI.ScalarAffineFunction{Float64},
    #    [
    #        (MOI.ScalarAffineFunction{Float64},MOI.GreaterThan{Float64}),
    #        (MOI.ScalarAffineFunction{Float64},MOI.LessThan{Float64}),
    #        (MOI.SingleVariable,MOI.GreaterThan{Float64})
    #    ]
    #)

I'm not sure what the difference is, but it could be that the input layer to SCS is doing a weird problem transformation that is making it harder to solve.

[bodono]

It's getting x=59.0909071799 and y=36.3636487368, by the way.

[blegat]

I am getting x=59.850517618053175 and y=23.4992567768817 with SCS 2.
Maybe cvxpy is generating a different equivalent formulation. Or maybe it is doing an automatic scaling (e.g. transforming 1000*x + 300 * y <= 70000 into 10x + 3y <= 700).

The output I get is

----------------------------------------------------------------------------
	SCS v2.0.2 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012-2017
----------------------------------------------------------------------------
Lin-sys: sparse-direct, nnz in A = 8
eps = 1.00e-05, alpha = 1.80, max_iters = 5000, normalize = 1, scale = 1.00
acceleration_lookback = 20, rho_x = 1.00e-03
Variables n = 2, constraints m = 5
Cones:	linear vars: 5
Setup time: 6.35e-05s
----------------------------------------------------------------------------
 Iter | pri res | dua res | rel gap | pri obj | dua obj | kap/tau | time (s)
----------------------------------------------------------------------------
     0| 2.01e+00  5.54e-01  3.11e-01 -2.14e+05 -1.13e+05  0.00e+00  1.98e-05 
   100| 4.47e-02  7.94e-02  6.95e-03 -6.87e+04 -6.97e+04  2.69e-13  1.13e-04 
   200| 2.23e-02  3.94e-02  5.99e-03 -7.03e+04 -6.95e+04  1.61e-12  1.91e-04 
   300| 1.48e-02  8.44e-03  3.43e-03 -7.09e+04 -7.14e+04  2.10e-12  2.67e-04 
   400| 1.31e-02  1.24e-02  4.48e-04 -7.10e+04 -7.10e+04  2.14e-12  3.44e-04 
   500| 5.74e-03  2.20e-02  6.56e-03 -7.12e+04 -7.22e+04  1.81e-12  4.20e-04 
   600| 3.19e-02  1.30e-01  6.58e-02 -7.02e+04 -8.01e+04  2.70e-12  4.95e-04 
   700| 2.18e-03  8.23e-03  4.44e-03 -7.18e+04 -7.24e+04  2.36e-12  5.72e-04 
   800| 1.99e-03  5.27e-03  2.39e-03 -7.18e+04 -7.14e+04  2.26e-12  6.48e-04 
   900| 1.49e-02  1.74e-02  1.75e-03 -7.08e+04 -7.10e+04  1.99e-12  7.24e-04 
  1000| 5.19e-02  7.93e-02  5.53e-02 -6.77e+04 -7.56e+04  2.07e-12  8.01e-04 
  1100| 5.64e-03  1.75e-02  6.81e-03 -7.14e+04 -7.23e+04  1.77e-12  8.77e-04 
  1200| 1.20e-02  1.11e-02  8.30e-03 -7.24e+04 -7.12e+04  1.64e-12  9.53e-04 
  1300| 1.50e-02  1.52e-02  6.85e-04 -7.09e+04 -7.10e+04  1.67e-12  1.03e-03 
  1400| 1.91e-02  4.08e-02  4.72e-03 -7.32e+04 -7.38e+04  9.29e-13  1.11e-03 
  1500| 3.66e-03  9.90e-03  4.66e-03 -7.16e+04 -7.23e+04  1.82e-12  1.18e-03 
  1600| 1.83e-02  3.38e-02  7.10e-03 -7.32e+04 -7.43e+04  1.62e-12  1.26e-03 
  1700| 1.16e-02  6.34e-03  2.44e-03 -7.11e+04 -7.14e+04  1.56e-12  1.33e-03 
  1800| 7.74e-03  5.39e-03  1.38e-03 -7.14e+04 -7.16e+04  1.47e-12  1.43e-03 
  1900| 2.35e-02  1.00e-02  2.84e-04 -7.30e+04 -7.30e+04  1.41e-12  1.51e-03 
  2000| 5.09e-02  3.10e-02  9.66e-03 -7.52e+04 -7.37e+04  1.40e-12  1.58e-03 
  2100| 5.97e-03  3.43e-03  8.55e-04 -7.15e+04 -7.16e+04  7.10e-13  1.66e-03 
  2200| 1.93e-03  4.51e-03  8.04e-05 -7.21e+04 -7.20e+04  6.63e-13  1.74e-03 
  2300| 3.49e-02  4.49e-02  1.18e-02 -7.33e+04 -7.51e+04  6.26e-13  1.81e-03 
  2400| 1.68e-02  3.09e-02  6.10e-03 -7.29e+04 -7.21e+04  3.32e-12  1.89e-03 
  2500| 7.90e-03  1.61e-02  1.60e-03 -7.13e+04 -7.15e+04  2.94e-14  1.96e-03 
  2600| 1.79e-03  2.43e-02  4.40e-03 -7.07e+04 -7.00e+04  1.95e-12  2.04e-03 
  2700| 1.13e-02  5.61e-02  2.03e-02 -7.08e+04 -6.80e+04  3.29e-12  2.12e-03 
  2800| 9.14e-04  7.61e-03  3.04e-03 -7.20e+04 -7.24e+04  3.80e-12  2.19e-03 
  2900| 3.06e-04  5.24e-03  1.87e-03 -7.19e+04 -7.21e+04  3.63e-12  2.27e-03 
  3000| 1.46e-02  1.44e-02  9.69e-04 -7.09e+04 -7.10e+04  2.65e-12  2.35e-03 
  3100| 6.08e-04  9.57e-03  9.35e-04 -7.19e+04 -7.21e+04  2.61e-12  2.42e-03 
  3200| 2.45e-03  4.21e-03  8.90e-05 -7.21e+04 -7.21e+04  1.70e-12  2.50e-03 
  3300| 2.00e-02  1.20e-02  4.43e-03 -7.05e+04 -7.12e+04  1.58e-12  2.57e-03 
  3400| 2.33e-02  1.37e-02  1.56e-02 -7.03e+04 -7.26e+04  1.69e-13  2.65e-03 
  3500| 1.18e-02  3.40e-02  1.16e-02 -7.10e+04 -6.94e+04  6.56e-13  2.73e-03 
  3600| 5.82e-02  5.51e-02  1.74e-02 -6.72e+04 -6.96e+04  7.12e-13  2.80e-03 
  3700| 1.22e-03  2.62e-02  2.18e-03 -7.10e+04 -7.13e+04  5.50e-13  2.88e-03 
  3800| 1.66e-02  1.27e-02  2.56e-03 -7.07e+04 -7.11e+04  5.38e-13  2.95e-03 
  3900| 6.55e-03  4.62e-03  4.81e-03 -7.15e+04 -7.22e+04  1.86e-12  3.03e-03 
  4000| 5.04e-02  4.51e-02  4.10e-02 -7.45e+04 -6.86e+04  1.66e-12  3.11e-03 
  4100| 7.84e-03  1.80e-02  3.12e-03 -7.13e+04 -7.17e+04  1.67e-12  3.19e-03 
  4200| 2.76e-02  3.06e-02  3.42e-03 -7.36e+04 -7.41e+04  1.54e-12  3.26e-03 
  4300| 1.81e-02  1.53e-02  3.60e-03 -7.07e+04 -7.12e+04  3.47e-13  3.34e-03 
  4400| 6.78e-03  1.23e-02  8.77e-03 -7.14e+04 -7.27e+04  3.19e-13  3.44e-03 
  4500| 1.33e-02  2.18e-02  9.75e-04 -7.29e+04 -7.30e+04  3.10e-13  3.55e-03 
  4600| 1.26e-02  7.21e-03  1.00e-02 -7.28e+04 -7.14e+04  2.95e-13  3.65e-03 
  4700| 2.49e-02  5.05e-02  6.73e-03 -7.03e+04 -6.94e+04  3.03e-13  3.76e-03 
  4800| 8.39e-02  3.16e-02  4.23e-02 -7.70e+04 -7.08e+04  2.76e-13  3.86e-03 
  4900| 1.97e-02  9.19e-02  5.14e-02 -7.25e+04 -6.54e+04  2.84e-13  3.97e-03 
  5000| 4.43e-02  5.53e-02  4.98e-04 -6.81e+04 -6.81e+04  2.55e-13  4.06e-03 
----------------------------------------------------------------------------
Status: Solved/Inaccurate
Hit max_iters, solution may be inaccurate
Timing: Solve time: 4.06e-03s
	Lin-sys: nnz in L factor: 17, avg solve time: 8.24e-08s
	Cones: avg projection time: 2.68e-08s
	Acceleration: avg step time: 4.62e-07s
----------------------------------------------------------------------------
Error metrics:
dist(s, K) = 1.5636e-13, dist(y, K*) = 0.0000e+00, s'y/|s||y| = 5.5348e-01
primal res: |Ax + s - b|_2 / (1 + |b|_2) = 4.4285e-02
dual res:   |A'y + c|_2 / (1 + |c|_2) = 5.5259e-02
rel gap:    |c'x + b'y| / (1 + |c'x| + |b'y|) = 4.9773e-04
----------------------------------------------------------------------------
c'x = -68075.2575, -b'y = -68143.0578
============================================================================

and the input is (with A both in sparse and dense format)

m = 5
n = 2
A = 5×2 SparseMatrixCSC{Float64,Int64} with 8 stored entries:
  [2, 1]  =  -1.0
  [3, 1]  =  1.5
  [4, 1]  =  8.0
  [5, 1]  =  300.0
  [1, 2]  =  -1.0
  [3, 2]  =  -1.0
  [4, 2]  =  12.0
  [5, 2]  =  1000.0
A = 5×2 Array{Float64,2}:
   0.0    -1.0
  -1.0     0.0
   1.5    -1.0
   8.0    12.0
 300.0  1000.0
b = [-30.0, 0.0, 0.0, 1000.0, 70000.0]
c = [-350.0, -1000.0]
f = 0
l = 5

If I scale the problem to

   maximize 10 x + 3.5 y

       s.t.                x >= 30
                           y >= 0
                 x -   1.5 y >= 0
            3   x +    2  y <= 250
            10 x +   3   y <= 700

I still get the same status but I get x = 58.89854797801036 and y = 38.460053472974366 which is closer.

[bodono]

Here is what SCS is getting from cvxpy:

A:
[[   -1.      0. ]
 [    0.     -1. ]
 [   -1.      1.5]
 [   12.      8. ]
 [ 1000.    300. ]]

b: [ -3.00000000e+01  -0.00000000e+00  -0.00000000e+00   1.00000000e+03 7.00000000e+04]

c: [-1000.  -350.]

cone: {'s': [], 'q': [], 'f': 0, 'ep': 0, 'l': 5}

almost exactly the same, just with the columns of A and c flipped, which shouldn't make a difference. I'll keep investigating.

[bodono]

Something fishy is going on, because if I manually put in the input values you have, I get the following:

----------------------------------------------------------------------------
	SCS v2.0.2 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012-2017
----------------------------------------------------------------------------
Lin-sys: sparse-indirect, nnz in A = 8, CG tol ~ 1/iter^(2.00)
eps = 1.00e-05, alpha = 1.50, max_iters = 5000, normalize = 1, scale = 1.00
acceleration_lookback = 20, rho_x = 1.00e-03
Variables n = 2, constraints m = 5
Cones:	linear vars: 5
Setup time: 1.42e-04s
----------------------------------------------------------------------------
 Iter | pri res | dua res | rel gap | pri obj | dua obj | kap/tau | time (s)
----------------------------------------------------------------------------
     0| 2.61e+00  5.54e-01  3.91e-01 -2.57e+05 -1.13e+05  3.67e-11  1.52e-02
   100| 4.61e-03  2.50e-03  2.50e-03 -7.13e+04 -7.17e+04  1.55e-12  1.54e-02
   200| 1.07e-02  1.45e-02  3.04e-04 -7.21e+04 -7.22e+04  3.77e-12  1.56e-02
   300| 1.31e-03  1.25e-02  5.91e-03 -7.15e+04 -7.24e+04  2.05e-13  1.57e-02
   380| 7.57e-07  1.51e-06  9.04e-07 -7.18e+04 -7.18e+04  7.26e-13  1.58e-02
----------------------------------------------------------------------------
Status: Solved
Timing: Solve time: 1.58e-02s
	Lin-sys: avg # CG iterations: 1.98, avg solve time: 3.28e-07s
	Cones: avg projection time: 4.00e-08s
	Acceleration: avg step time: 7.64e-07s
----------------------------------------------------------------------------
Error metrics:
dist(s, K) = 3.9345e-14, dist(y, K*) = 0.0000e+00, s'y/|s||y| = 1.9830e-15
primal res: |Ax + s - b|_2 / (1 + |b|_2) = 7.5666e-07
dual res:   |A'y + c|_2 / (1 + |c|_2) = 1.5141e-06
rel gap:    |c'x + b'y| / (1 + |c'x| + |b'y|) = 9.0419e-07
----------------------------------------------------------------------------
c'x = -71818.2074, -b'y = -71818.0775
============================================================================

[bodono]

Ah, I found it, it's the alpha = 1.80, by default it's alpha = 1.50 which is more conservative. If I switch to alpha = 1.80 in python with your input parameters I get the same bad performance.

[bodono]

By the way, the alpha term is just a type of overstep parameter, so the update takes a larger step by that factor and in theory it should converge for all choices of \alpha \in (0,2), but in practice being too aggressive can hurt performance, which we see here. I think the default of 1.5 should be propagated to SCS.jl.