OSQP giving an inaccurate result for a simple MPC problem

[robin...@gmail.com]

I've been trying to use OSQP (through the Julia interface) for some simple Model-Predictive Control applications, but I'm getting a lot of suboptimal answers from the solver. I was concerned that there might be a bug in the Julia MathOptInterface layer, but I've now found a simple case that gives a very wrong answer even when using OSQP from Julia and Python without any high-level interface. The complete problem setup (in Julia) is as follows:

```
P = spzeros(4, 4)
P[1, 1] = 20

q = zeros(4)

A = spzeros(5, 4)
A[2, 1] = -100
A[1, 2] = 10
A[2, 2] = 1
A[3, 2] = -100
A[3, 3] = 1
A[1, 4] = -0.01
A[4, 4] = 1
A[5, 4] = 1

l = [0., -100, 0, -100, -Inf]
u = [0., -100, 0, Inf, 100]

m = OSQP.Model()
OSQP.setup!(m; P=P, A=A, l=l, u=u)
results = OSQP.solve!(m)
results.x
```

The expected solution is [0.999, -0.1, -10.0, -100.0], which is exactly what I get using Gurobi. OSQP, on the other hand, returns [0.999999, -7.93976e-5, -0.00793976, -0.0793976], which is far from the optimal solution.

The complete solver output is:

```
-----------------------------------------------------------------
OSQP v0.4.0 - Operator Splitting QP Solver
(c) Bartolomeo Stellato, Goran Banjac
University of Oxford - Stanford University 2018
-----------------------------------------------------------------
problem: variables n = 4, constraints m = 5
nnz(P) + nnz(A) = 9
settings: linear system solver = qdldl,
eps_abs = 1.0e-03, eps_rel = 1.0e-03,
eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
rho = 1.00e-01 (adaptive),
sigma = 1.00e-06, alpha = 1.60, max_iter = 4000
check_termination: on (interval 25),
scaling: on, scaled_termination: off
warm start: on, polish: off

iter objective pri res dua res rho time
1 0.0000e+00 1.00e+02 1.00e+04 1.00e-01 4.46e-05s
50 1.0000e+01 1.58e-07 2.00e-04 1.00e-01 7.30e-05s

status: solved
number of iterations: 50
optimal objective: 10.0000
run time: 7.80e-05s
optimal rho estimate: 3.98e-04
```

Gurobi, on the other hand, reports:

```
Academic license - for non-commercial use only
Optimize a model with 10 rows, 4 columns and 16 nonzeros
Model has 1 quadratic objective term
Coefficient statistics:
Matrix range [1e-02, 1e+02]
Objective range [0e+00, 0e+00]
QObjective range [2e+01, 2e+01]
Bounds range [0e+00, 0e+00]
RHS range [1e+02, 1e+02]
Presolve removed 10 rows and 4 columns
Presolve time: 0.00s
Presolve: All rows and columns removed

Barrier solved model in 0 iterations and 0.00 seconds
Optimal objective 9.98001000e+00
```

note that OSQP is reporting an optimal objective of 10.0000, while Gurobi was able to find an objective of 9.98001 (the true optimal value).

I have tried the exact same problem using the OSQP Python bindings and gotten the exact same result as I got in Julia, so I don't think this is a problem with the Julia interface.

These results were generated with OSQP v0.4.0, but I've had similar problems with 0.3.x versions.

Do you have any suggestions? Is this a bug in OSQP?

[robin...@gmail.com]

By the way, I tried removing the infinities and adjusting the scaling of the problem a little bit:

```
P = spzeros(4, 4)
P[1, 1] = 20

q = zeros(4)

A = spzeros(5, 4)

A[1, 2] = 10

A[1, 4] = -0.1
A[2, 1] = -10

A[2, 2] = 1
A[3, 2] = -10

A[3, 3] = 1

A[4, 4] = 1
A[5, 4] = 1

l = [0., -10, 0, -100, -100]
u = [0., -10, 0, 100, 100]
```

This still results in a wrong answer from OSQP: [0.995844, -0.0392951, -0.392951, -3.93114] with an objective value of 9.9171, vs. Gurobi's [0.9, -1.0, -10.0, -100.0] with an objective value of 8.1000.

[Bartolomeo Stellato]

It might be for two reasons:

• The cost function is not strictly convex. I haven't inspected the problem carefully, but there might be multiple solutions with the same cost function value.
  
  
• This behavior can be linked to the tolerances. The default ones for OSQP are less accurate than interior point solvers. For that small problem, the reported solution satisfies the tolerances eps_abs=1e-03 and eps_rel=1e-03. If you set the tolerances eps_abs=1e-06, eps_rel=1e-06 you should get exactly the solution that Gurobi reports. 
  However, this might require a large number of iterations for some badly scaled problems.

Do these differences affect the closed-loop behavior? If yes how much? I believe the feedback should take care of most of it. Also, in some applications the linearized models are already quite far from the true ones and it is not necessary to solve the QPs to high precision.

As a side, in the future we are planning to implement some acceleration schemes to converge to high accuracy solutions. It will take some time though.

Bartolomeo

--
You received this message because you are subscribed to the Google Groups "OSQP" group.
To unsubscribe from this group and stop receiving emails from it, send an email to osqp+uns...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/osqp/c3efb18e-b732-4c7e-8327-525c01e5e230%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Robin Deits]

> The cost function is not strictly convex. I haven't inspected the problem carefully, but there might be multiple solutions with the same cost function value.

That's definitely true, but it's a real part of the problem I'm trying to solve. i can get OSQP and Gurobi to converge to similar primal solutions, but only by adding 0.1 to the entire diagonal of P, which is so much regularization that the answer I get is nowhere near optimal for the un-regularized problem.

> If you set the tolerances eps_abs=1e-06, eps_rel=1e-06 you should get exactly the solution that Gurobi reports. 

Yes, that works! I have experimented with this in the past, and it does frequently result in very long solve times, but it's nice to be able to get the right answer.

> Do these differences affect the closed-loop behavior? If yes how much?

Yeah, they really do. In the original problem I posted, the fourth variable is the control action. I get a control action of -0.079 from OSQP vs. a control action of -100.0 from Gurobi. So the control action is wrong by more than 3 orders of magnitude. Since this *is* the feedback controller, there's no hope of stabilizing the closed loop.

> Also, in some applications the linearized models are already quite far from the true ones and it is not necessary to solve the QPs to high precision.

Definitely true. I'd be fine with significant inaccuracy in the value of the primal solution (1-10% would probably be fine). But in this case I'm getting a primal solution which is off by more than a factor of 1000, which doesn't work.

I'll try tightening the optimality tolerances, but my experience has been that this tends to make the solver too slow to use in feedback. But I'll give it a try.

To unsubscribe from this group and stop receiving emails from it, send an email to osqp+unsubscribe@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/osqp/c3efb18e-b732-4c7e-8327-525c01e5e230%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--

Robin L. H. Deits

Robot Locomotion Group @CSAIL

Massachusetts Institute of Technology

[Paul Goulart]

This is unrelated to OSQP, but I think it is generally not a good idea to leave the control input in your MPC problem completely unpenalised like this.   In the unconstrained LQR case what you are doing amounts to setting R = 0, which will typically lead to difficulties in even ensuring that an optimal controller exists.   

Related to this is that, even if you get the exact optimal solution at each time step, your controller is likely to exhibit chatter and related undesirable behaviours.

As a practical matter within the solver, you are saying that you have absolutely no preference about what control input is chosen so long as the constraints are satisfied.   You should therefore not be surprised if a slightly suboptimal solution differs greatly from the optimal one.   I would suggest adding at least a small quadratic (or even 1-norm) penalty on the control input.

To unsubscribe from this group and stop receiving emails from it, send an email to osqp+uns...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/osqp/CACV469GaM6_Ja_8KO_PqW5bEOY%2BeoTTXAx%3D2%3DO%3DYZuTe_mUD%2BA%40mail.gmail.com.

[robin...@gmail.com]

> osqp+uns...@googlegroups.com.
>
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/osqp/c3efb18e-b732-4c7e-8327-525c01e5e230%40googlegroups.com.
>
> For more options, visit
> https://groups.google.com/d/optout.
>
>
>
>
>
>
>
>
>
>
>
> --
>
>
> Robin L. H. Deits
>
> Robot Locomotion Group @CSAIL
>
>
>
>
> Massachusetts Institute of Technology
>
>
>
>
>
>
>
>
> --
>
> You received this message because you are subscribed to the Google Groups "OSQP" group.
>
> To unsubscribe from this group and stop receiving emails from it, send an email to
> osqp+uns...@googlegroups.com.
>
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/osqp/CACV469GaM6_Ja_8KO_PqW5bEOY%2BeoTTXAx%3D2%3DO%3DYZuTe_mUD%2BA%40mail.gmail.com.
>
> For more options, visit https://groups.google.com/d/optout.

Yeah, understood. This isn't the actual MPC problem I'm trying to solve, merely a minimal example that showed the solver behavior. In practice, for the problem specified above, I have to put a very large penalty (1e-1) on control action to get OSQP and Gurobi to agree with the default tolerances, at which point I've penalized control so much that the controller itself is ineffective.

[Paul Goulart]

Yeah, understood. This isn't the actual MPC problem I'm trying to solve, merely a minimal example that showed the solver behavior. In practice, for the problem specified above, I have to put a very large penalty (1e-1) on control action to get OSQP and Gurobi to agree with the default tolerances, at which point I've penalized control so much that the controller itself is ineffective. 

Could you supply us with a larger size example that exhibits the same behaviour?  We are trying to collect a data set of problems for which first order methods seem to have difficulty.

I was able to rescale your problem manually such that it converges in about 50 iterations to the right solution even with the tolerances set to 1e-3, but it is hard to say whether the method I used is robust or if I just got lucky.

-- 
You received this message because you are subscribed to the Google Groups "OSQP" group.
To unsubscribe from this group and stop receiving emails from it, send an email to osqp+uns...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/osqp/5e15464a-4db8-4d69-8845-ce485415e7df%40googlegroups.com.