Lyapunov & Riccati equation solvers

[jamestjsp]

Summary

Implement core matrix equation solvers that are prerequisites for controller design, observer design, system norms, gramians, and balanced realization.

Functions

• Lyap(A, Q) — solve continuous Lyapunov equation A*X + X*A' + Q = 0
• DLyap(A, Q) — solve discrete Lyapunov equation A*X*A' - X + Q = 0
• Care(A, B, Q, R) — solve continuous algebraic Riccati equation A'*X + X*A - X*B*R⁻¹*B'*X + Q = 0
• Dare(A, B, Q, R) — solve discrete algebraic Riccati equation A'*X*A - X - A'*X*B*(R+B'*X*B)⁻¹*B'*X*A + Q = 0

Implementation Notes

• Use Schur decomposition (real Schur form) for Lyapunov solvers (Bartels-Stewart algorithm)
• Use generalized eigenvalue approach or iterative Newton for Riccati solvers
• SLICOT routines SB03MD (Lyapunov) and SB02OD (Riccati) are the gold standard reference
• Consider optional cross-term S for generalized Riccati: Care(A, B, Q, R, S)

Blocked by

None — this is a foundation issue.

Blocks

• Gramians & state-space analysis
• System norms (H2, H∞)
• Controller design (LQR/LQE)
• Observer & estimator design (Kalman)
• Balanced realization & advanced model reduction