Tolerance tuning for inequality constraints

[stephane-caron]

Thanks for pushing SCS to PyPI :grinning: I'm adding support for it to qpsolvers for solving quadratic programs (after QP→SOCP conversion).

The only thing that doesn't work out of the box so far are the default tolerances on inequality constraints. On a sample problem, the solution violates inequalities by about 1e-6, while the other solvers are < 1e-10. I can get around it by setting eps=1e-7 and use_indirect=True (or keep use_indirect=False and set an even lower eps=1e-8), but I don't know how SCS works so this may be inefficient. Would you recommend something else?

For now I've written a comment on it in the doc https://github.com/stephane-caron/qpsolvers/commit/6e12ffe642692814b5f759990156a144fa33df4f, but I'd rather not override solvers' default values if possible.

[bodono]

I'm close to pushing SCS 3.0, which handles both the quadratic x'Px term and the box constraints (lb <= x <= ub) natively, so no conversions required, which should be much easier for qpsolvers to handle and give better results, and handle things like ub = infty. Currently it's sitting in the qp branch. If you want to test it let me know and we can work out any issues together.

[stephane-caron]

Sounds good! I just gave it a whirl it https://github.com/stephane-caron/qpsolvers/pull/43 and it worked out of the box :+1: Performance-wise avoiding the conversion reduces solve times significantly (details in the PR).

Have you documented the API for specifying box constraints? Right now qpsolvers casts all box constraints to linear inequalities, it will be even better if we can avoid this overhead.

[bodono]

Great! Thanks for trying it! The box cone API is as follows:

The cone is determined by two vectors (bl, bu) and defined as K_box = {(t, s): t * bl <= s <= t * bu, t >= 0} (note that the the variable vector is stacked as [t;s], ie, t comes first then s). In most cases the user will also want to add the constraint that t == 1 to the linear equality constraints, which adds a single additional row of all zeros to A and an entry of 1 to b in the appropriate locations. If t=1 then this boils down to the usual constraint that bl <= s <= bu, but without the t variable that set is not a cone and being a cone is required for the math to work.

Recall that the constraints are Ax + s == b with s in K, so depending on your input formulation you might also have to negate A.

The order of cones assumed in the A matrix is now:

1. primal zero / dual free cone
2. linear cone
3. box cone
4. second order cone
5. sd cone
6. primal exponential cone
7. dual exponential cone
8. primal power cone
9. dual power cone

[stephane-caron]

Thank you for these details :slightly_smiling_face: Unfortunately it will be a while before I take the time to look at these box constraints.

Do you have an ETA for the release of SCS 3.0? The results at https://github.com/stephane-caron/qpsolvers/pull/43 are great (and the PR is ready so it wouldn't take me too much time to review and merge it ^.^). Even without the proper box constraints there's already a clear improvement for SCS users.

[bodono]

Hi @stephane-caron , I'm very close to releasing SCS 3.0, just putting on the finishing touches and allowing some time for feedback etc. It should land in the next week or two.

Here is the PR to update core SCS: https://github.com/cvxgrp/scs/pull/169

And here is the PR to update the python interface: https://github.com/bodono/scs-python/pull/50

Also we have a new documentation website that should hopefully make it easier to work with. In fact we should mention qpsolvers on the website somewhere.

Any feedback welcome!

[stephane-caron]

Nice! :smiley: The new Python example is straightforward and helpful.

I'm close to releasing a new version of qpsolvers with SCS 3.0. Feel free to check out the code of the scs_solve_qp function and let me know if there is anything it can do better.

[stephane-caron]

From an API standpoint, it would be nice if SCS would also return the solution to unconstrained problems. The other QP solvers do, and throwing an exception when there is no constraint can become a pain point to some users (for instance those who solve sequences of problems with variable numbers of constraints).

[bodono]

I think it would take some rewriting to support that, but if there are no constraints the solution is just P \ c right? If so then it would be faster for numpy to just compute that in the wrapper.

[stephane-caron]

Sounds good to me. This is done in stephane-caron/qpsolvers@eafa8acbbe2a1ff87f311d72f44b28bc1c276d63.