warmstart SCS

[madeleineudell]

This fix adds an option to warmstart solutions to SCS. Setting the warmstart option to true now doesn't cause a segfault, and actually does warmstart it. Also added a MathProgBase-style setwarmstart! method, in which dual_sol and slack can be set via keyword arguments. Usage examples are in the test file options.jl.

[mlubin]

@madeleineudell, is SCS currently the only solver where you would want to pass this warmstart info?
If so, I think that the most practical approach until JuliaOpt/MathProgBase.jl#64 is worked out would be to leave this out of the MathProgBase interface and just set m.dual_sol and m.slack directly from Convex.jl.

[mlubin]

Also, what does the slack represent? If it is just b- A*x, why is it a separate parameter and not computed directly given primal_sol and the problem data?

[madeleineudell]

It looks like we won't have to mess with the dual/slack variables at all
from Convex; we can just cache the AbstractMathProgModel used the last time
the problem was solved. So we can just use the same setwarmstart!(m, v)
syntax.

In SCS.jl itself, I did make the call signature setwarmstart!(m, v; kwargs...) so you can set the dual_sol and slack if you want to. (No,
I don't know why the slack needs to be set explicitly, but the SCS
documentation says you can't warmstart from one without warmstarting from
all.) Currently my solution is to warmstart any unset solution variable
(primal/dual/slack) as the vector of all zeros (of the appropriate length).
I could set slack = b - A * primal_sol if that seems better. (My guess is
SCS ignores at least one of these three warmstart vectors --- whichever is
updated first in the SCS iterations --- but I don't know which.)

On Mon, Mar 23, 2015 at 9:21 AM, Miles Lubin [email protected]
wrote:

Also, what does the slack represent? If it is just b- A*x, why is it a
separate parameter and not computed directly given primal_sol and the
problem data?

—
Reply to this email directly or view it on GitHub
#22 (comment).

Madeleine Udell
PhD Candidate in Computational and Mathematical Engineering
Stanford University
www.stanford.edu/~udell

[bodono]

The reason you specify all (x,s,y) in SCS for warm-starting rather than just (x,y) is that the same reason SCS returns (x,s,y) instead of just (x,y), which is to give the user more flexibility in how to use the solution. The s variable is not guaranteed to satisfy the equality constraint Ax + s = b exactly, but is guaranteed to be in the cone. So the user can either choose to warm-start using s from a previous run or by passing in b - Ax. If they use s from the previous run and the data is unchanged, then the algorithm will terminate immediately (if the previous run solved the problem and didn't hit max_iters), which is a useful sanity check.

[davidlizeng]

return m just for clarity

[madeleineudell]

Added code to rescale the SDP constraint as in cvxgrp/scs#40. loadconicproblem! rescales the SDP constraints by default, and getconiddual inverts the rescaling. Rescaling can be turned off using the option (:rescale, false).

@bodono does this look good to merge?

[mlubin]

If this passes on travis linux, then we should be good to go. Windows and OS X will definitely fail until we update the binaries.

[madeleineudell]

It passed on the travis linux release, but not the nightly. Is that to be expected?

[mlubin]

No, the error isn't expected:

ERROR: LoadError: LoadError: MethodError: `convert` has no method matching convert(::Type{Ptr{SCS.SCSSettings}}, ::Array{SCS.SCSSettings,1})
This may have arisen from a call to the constructor Ptr{SCS.SCSSettings}(...),
since type constructors fall back to convert methods.
Closest candidates are:
  call{T}(::Type{T}, ::Any)
  convert{T}(::Type{Ptr{T}}, !Matched::UInt64)
  convert{T}(::Type{Ptr{T}}, !Matched::Int64)
  ...
 in call at no file
 in create_scs_data at /home/travis/.julia/v0.4/SCS/src/high_level_wrapper.jl:41
 in create_scs_data at /home/travis/.julia/v0.4/SCS/src/high_level_wrapper.jl:58
 in SCS_solve at /home/travis/.julia/v0.4/SCS/src/high_level_wrapper.jl:130 (repeats 4 times)
 in include at ./boot.jl:253
 in include_from_node1 at ./loading.jl:133
 in anonymous at no file:9
 in include at ./boot.jl:253
 in include_from_node1 at loading.jl:133
 in process_options at ./client.jl:304
 in _start at ./client.jl:404
while loading /home/travis/.julia/v0.4/SCS/test/direct.jl, in expression starting on line 5
while loading /home/travis/.julia/v0.4/SCS/test/runtests.jl, in expression starting on line 7

Have you tried this on 0.4 locally?

[bodono]

It looks like you'll need to merge in #17 first, since I think that PR includes the required API changes for this to work.

[madeleineudell]

I believe I did... take a look at the commit log.

On Tue, Jun 9, 2015 at 5:01 PM, bodono [email protected] wrote:

It looks like you'll need to merge in #17
#17 first, since I think that PR
includes the required API changes for this to work.

—
Reply to this email directly or view it on GitHub
#22 (comment).

Madeleine Udell
PhD Candidate in Computational and Mathematical Engineering
Stanford University
www.stanford.edu/~udell

[madeleineudell]

I fixed the pointer issue @mlubin noted above on .4, but there's now an overflow error (!) I'm not able to reproduce locally:

ERROR: LoadError: LoadError: OverflowError()
 in hcat at abstractarray.jl:641
 in feasible_exponential_conic at /home/travis/.julia/v0.4/SCS/test/direct.jl:73
while loading /home/travis/.julia/v0.4/SCS/test/direct.jl, in expression starting on line 85
while loading /home/travis/.julia/v0.4/SCS/test/runtests.jl, in expression starting on line 7

[bodono]

Line 73 is just a long vector of small ints being read in:

rowval = vec([7 9 18 ... 392 407]');

Could this be a bug in julia about the way it does hcat?

[bodono]

Actually, looks like it might be because of this issue: JuliaLang/julia#11320

We could change the test to make the matrix smaller, if that's confirmed to be the issue.

[bodono]

How close are we to merging this in?

[madeleineudell]

There are two remaining issues I'm aware of:

1. The overflow error above.
2. This branch is still failing a small number of the Convex.jl tests (Test with new symmetric SDP format #14). My current guess is that the problem is in the scaling (8134662) I added to SCS.jl, since Mosek passes all Convex.jl tests. Comparisons to CVX are difficult, since CVX has a sophisticated presolve routine.

[bodono]

If we fix the scaling bug then we should just merge and I will submit a new PR with a smaller test that doesn't cause an overflow. There's no need for that test to be 445 x 148 dimensions. (I don't think I can add commits to this PR for that).

[madeleineudell]

Responses to each of @bodono's three comments:

1. Thanks for finding that bug! All Convex.jl tests pass now with this branch.
2. The MathProgBase conic interface claims to provide a consistent abstraction to conic programming. So the Mosek.jl and SCS.jl implementations of that interfaces should not differ in the scaling they expect. @mlubin this is something that should be clarified on the MathProgBase docs, because I think both interpretations are currently plausible given what's written there.
3. Then let's merge!

[bodono]

Well then should MathProgBase do the scaling? If it's exposing a consistent API for conic programming then it's responsible for the nitty gritty differences between solver APIs.

Who has the power to merge this in?

[madeleineudell]

Just wanted to verify that we have no new errors before merging it in.

[mlubin]

@bodono, SCSMathProgModel shouldn't be considered the SCS wrapper in Julia, it's the MPB abstraction wrapper around the SCS interface. This is where any transformations between what the solver wants and what the interface specifies should take place. @madeleineudell, would you mind putting together a PR to clarify the (lack of) scaling of the input in the MPB docs?

I believe that the SCS_solve function is the "direct" interface, and no transformations are performed there.