Question on 4 different functions for solving the KKT system of OCP QP

[syangav]

Hello,

Thanks for open-sourcing such great code!

I noticed that there are 4 different functions who all do some kinds of Riccati recursion in HPIPM. They are:

1. OCP_QP_FACT_SOLVE_KKT_UNCONSTR. This is the easiest to understand, both by name and by content. It is doing Riccati recursion for unconstrained LQR problem. Depending on square_root_alg, it is either using the classical version (no Cholesky on the cost-to-go matrix P) or the square root version (Cholesky on cost-to-go matrix P).
2. OCP_QP_FACT_SOLVE_KKT_STEP. This is used when lq_fact = 0 and the initial phase of lq_fact = 1.
3. OCP_QP_FACT_LQ_SOLVE_KKT_STEP. This is used when lq_fact = 2 and will be switched to if lq_fact = 1 and some inaccuracy happens. Thanks for your explanantion here.
4. OCP_QP_SOLVE_KKT_STEP. This seems to be the most widely used one among the 4 in the file x_ocp_qp_ipm.c.

To be honest, option 2, 3, 4 are too long to understand, so it would be great if you can elaborate the differences between them and their respective usage scenarios. I can tell that all of them are solving the perturbed KKT system with the inequality constraint matrix C, current guess of inequality Lagrangian multipliers \Lambda and current guess of slack variables T.

The perturbed KKT system still has the OCP QP KKT structure hence can be solved efficiently by Riccati recursion. If I were the code developer, I would add the extra term to the Hessian matrix H (the functions COND_SLACKS_FACT_SOLVE, COND_SLACKS_SOLVE and EXPAND_SLACKS seem to do the job) and then simply call OCP_QP_FACT_SOLVE_KKT_UNCONSTR . May I ask why you did not choose to do so?

Thanks for your time in advance!

Best regards,
Yang

[giaf]

Hi,

about the 4 different functions,

1. Correct.
2. KKT system factorization and solution by means of Riccati recursion performed on the KKT system after the Lagrange multipliers and the relative slack variables have been eliminated (i.e. the Hessian updated) as in the figure above.
3. Same as 2. but it uses an array algorithm based on LQ factorization for improved numerical stability, at the expense of more costly computations. This only supports the square root variant of Riccati (thus needing the Hessian to be positive definite), and call back 2. for the Classical Riccati algorithm.
4. This is only the solution without the factorization, so to be performed after 2 when the corrector step in the Mehrotra's IPM is computed.

Updating the Hessian within the factorization improves performance for small matrices (less routines calls, better data locality, etc ...), compared to doing as you suggest and updating the Hessian before calling the unconstrained Riccati.