Simple-Robotics / proxsuite

The Advanced Proximal Optimization Toolbox
BSD 2-Clause "Simplified" License
399 stars 48 forks source link

Detect infeasible QP problem before solving it #291

Open singya opened 8 months ago

singya commented 8 months ago

Hello,

I am building an MPC controller that uses ProxQP to get the optimal control actions.

It is about a mobile robot that avoid obstacles. I generate at each iteration (every 0.1 seconds) a warm-starting obtained by receding the previous result.

Obstacles have some degree of uncertainty on its identification so that the resulting constraints moves resulting sometimes on infeasible problems. When an infeasible problem arrives, the solver takes a long time (several minutes) to end up with an anyway infeasible solution. In those cases, I actually prefer to take the warm-start as the optimal solution, hoping for future iterations to be feasible, or eventually bringing the robot to a stop.

Currently, I have set a maximum number of iterations for the problem to finish early when the solution gets too complex:

results = proxsuite.proxqp.sparse.solve(
            H=H,
            g=g,
            C=C,
            l=l,
            u=u,
            x=warm_start,
            max_iter = 10
            rho = 1e-8
        )

But this also takes much longer (~0.3 seconds) for a solution that I won't really use. If the results.info.status is not proxsuite.proxqp.PROXQP_SOLVED, I use instead the warm-start. It is also interesting to mention that I have never obtained the status of PROXQP_PRIMAL_INFEASIBLE or PROXQP_DUAL_INFEASIBLE.

So, I was wondering if it would be possible to stop the solver as soon as it recognizes that the problem is infeasible, or if there exist some simple feasibility check that I could do before calling the method solve.

In the QP problem that I am considering, all equality constraints are implicit. So, I do not have matrices A and b. The matrices H, g, C, l and u for an infeasible example can be found in this link:

https://drive.google.com/file/d/1_p4KUv4KTVvmG2fM0cdPNp7mky_XKJYL/view?usp=sharing

And it can be tested using the following script:

import pickle
import proxsuite

with open("matrices.pickle", "rb") as f:
    matrices = pickle.load(f)

results = proxsuite.proxqp.sparse.solve(
    H=matrices["H"],
    g=matrices["g"],
    C=matrices["C"],
    l=matrices["l"],
    u=matrices["u"],
    x=matrices["x"],
    max_iter = 1000,
    rho = 1e-8
)

print(results.info.status)
whtqh commented 8 months ago

I got the same problem, bro.

fabinsch commented 8 months ago

Hey @singya,

currently, our default settings for infeasibility detection are very low, see here

T eps_primal_inf = 1.E-12,
T eps_dual_inf = 1.E-12,

and do not correspond to what's written in our doc. We will fix this in #292 soon.

For now, if you want to change these values the best is to do the following

import pickle
import proxsuite

with open("matrices.pickle", "rb") as f:
    matrices = pickle.load(f)

n = matrices["H"].shape[0]
n_eq = 0
n_in = matrices["C"].shape[0]

qp = proxsuite.proxqp.sparse.QP(n, n_eq, n_in)
qp.settings.verbose = False
qp.settings.eps_primal_inf = 1e-4
qp.settings.eps_dual_inf = 1e-4
qp.settings.max_iter = 1000

qp.init(
    H=matrices["H"],
    g=matrices["g"],
    A=None,
    b=None,
    C=matrices["C"],
    l=matrices["l"],
    u=matrices["u"],
    rho=1e-8,
)

# warm start x, y and z
qp.solve(matrices["x"], None, None)
print(f"{qp.results.info.status=}")
print(f"{qp.results.info.iter=}")

This setup still results in QPSolverOutput.PROXQP_MAX_ITER_REACHED. Only I set eps_primal_inf and eps_dual_inf to 1e-1 I get QPSolverOutput.PROXQP_DUAL_INFEASIBLE. I hope this is already helpful to you.