I need to use a Discrete-time Algebraic Riccati Equation (DARE) solver for an embedded controller (with limited processing power) in a research project and sadly, I can't find any implementation of it in C++ that I can understand. I use C++ since it is faster than higher level languages. I found various libraries online that offer Riccati equation or LQR solvers such as the Control Toolbox (https://adrlab.bitbucket.io/ct/v2.3/ct_doc/doc/html/index.html) or Drake (https://drake.mit.edu/) but I can't understand the language and I think they are full of unnecessary functions and things for the simple job of solving DARE. I would like the implementation to use basic open source C++ libraries or header files such as Eigen, Armadillo etc. with a lot of documentation, tutorials, or an active online community using them so that I can learn the syntax on my own. I'm looking for a simple code like Arash's C++ implementation of the Continuous-time Algebraic Riccati Equation (CARE) solver: https://math.stackexchange.com/questions/679989/analytic-solution-to-structured-algebraic-riccati-equation .

To anyone who can help me, I would like to cite you in the bibliography so you can please include a "How to Cite" section. If you don't have a personal C++ code, can you please refer me to a library or something that can help me implement DARE in C++? I attempted to write my own solver but when I read papers of the Riccati solver, I was dumbfounded by the math and terminologies since I am an undergraduate engineering student with only basic knowledge in linear algebra. I'm not sure if this is the right site to ask but I can't find anything else so please help me out.

  • 2
    $\begingroup$ Have you tried implementing the algorithm at en.wikipedia.org/wiki/Algebraic_Riccati_equation#Solution yourself? $\endgroup$ – Kirill Dec 22 '18 at 14:18
  • $\begingroup$ @Kirill, I have not since I thought that the options of implementing it and the theory are complicated. Should I implement the algorithm myself at this point? Can you recommend a reference for a fast DARE solver algorithm that is easy to understand by an amateur in Math like me? Should I implement the one written in Wikipedia? I've read that "State-of-the-art implementations of ARE solvers use a Schur decomposition method. This is what MATLAB is using (SciPy, Octave, and LAPACK also it)." <github.com/RobotLocomotion/drake/issues/1180> Is this true? $\endgroup$ – John Smith Dec 23 '18 at 6:08
  • $\begingroup$ Also, the solution in Wikipedia is for CARE not DARE. $\endgroup$ – John Smith Dec 23 '18 at 6:18
  • $\begingroup$ DARE is at the bottom of the Wikipedia section. It depends on what you need, if you’re stuck with your constraints then a simple algorithm that works for you may well be good enough. $\endgroup$ – Kirill Dec 23 '18 at 12:31
  • 1
    $\begingroup$ @JohnSmith did it work? if yes, a short self-answer would be very useful. $\endgroup$ – Anton Menshov Jun 5 '19 at 22:01

If you want a ten-line solution that is decently fast and stable, you can implement yourself the structured doubling algorithm: set up the coupled iteration

$A_0 = A, G_0 = G = BR^{-1}B^T, H_0 = Q$

While $\frac{\|H_{k+1}-H_k\|}{\|H_{k+1}\|} \geq \varepsilon$:

$\quad \quad A_{k+1} = A_k(I+G_kH_k)^{-1}A_k$

$\quad \quad G_{k+1} = G_k + A_k(I+G_kH_k)^{-1}G_kA_k^T$

$\quad \quad H_{k+1} = H_k + A_k^TH_k(I+G_kH_k)^{-1}A_k$

Return $H_{k+1}$.

The matrices $H_k$ are symmetric, and they converge quadratically fast to the stabilizing solution $X$ of $X = A^TXA - (B^TXA)^T(R+B^TXB)^{-1}B^TXA + Q$.

Reference: http://dx.doi.org/10.1080/00207170410001714988 .


It is an old question but I will answer it. Here is a good library: https://github.com/TakaHoribe/Riccati_Solver

  • $\begingroup$ Could you describe what is in the link? The problem with links is that they tend to disapepar and the answers based on them are less useful. $\endgroup$ – nicoguaro Jun 22 '20 at 16:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.