Continuous-Time Algebraic Riccati Equations

• Continuous-Time Algebraic Riccati Equations are nonlinear matrix equations central to infinite-horizon optimal control, yielding the unique stabilizing feedback matrix.
• They blend operator-theoretic foundations with diverse computational methods, including Schur-based factorizations, Krylov subspace iterations, and ADI-type iterations.
• Extensions to singular, stochastic, indefinite, and multilinear tensor cases demonstrate their broad applicability in robust controller design, model reduction, and optimization.

Continuous-time algebraic Riccati equations (CAREs) comprise a central class of nonlinear matrix equations framing the optimal control and filtering of linear dynamical systems. The canonical CARE arises in infinite-horizon linear-quadratic regulator (LQR) problems, as well as in observer and robust controller synthesis, model reduction, and inverse problems. The solution theory of CAREs combines operator-theoretic foundations (Hamiltonian structures, stabilizability/detectability) with rich computational methodology ranging from direct factorization (Schur, sign-function) to structural iterative algorithms (ADI, Krylov, doubling, Newton–Kleinman) and polynomial optimization relaxations. Generalizations encompass singular, stochastic, indefinite, and multilinear tensor variants.

1. Mathematical Formulation, Core Properties, and Key Existence Criteria

The standard continuous-time CARE seeks a symmetric matrix $P=P^T \in \mathbb{R}^{n\times n}$ solving

$A^T P + P A - P B R^{-1} B^T P + Q = 0,$

where $A \in \mathbb{R}^{n \times n}$, $B \in \mathbb{R}^{n \times m}$, $Q=Q^T \succeq 0$, $R=R^T \succ 0$ (2002.04246, Kao et al., 2020, Zhang et al., 2024). This equation governs the optimal state-feedback gain $K = R^{-1} B^T P$ for systems $\dot{x} = A x + B u$, minimizing the infinite-horizon quadratic cost $J = \int_0^\infty [ x^T Q x + u^T R u ] dt$.

A matrix $P \succeq 0$ is called stabilizing if the closed-loop $A - B R^{-1} B^T P$ is Hurwitz. Classical results (Kucera, Kano, Lancaster–Rodman) establish equivalence between the existence of a unique maximal stabilizing solution $P \succeq 0$ and the dual conditions:

• $(A,B R^{1/2})$ is stabilizable.
• $(Q^{1/2},A)$ is detectable.
• The Hamiltonian matrix

$$H = \begin{pmatrix} A & -B R^{-1} B^T \ -Q & -A^T \end{pmatrix}$$

has no purely imaginary eigenvalues and exactly $n$ eigenvalues in the left half-plane (Zhang et al., 2024, Kao et al., 2020, 2002.04246).

The impulse-free singular LQR generalization accommodates $R \geq 0$ (possibly singular) and a cross-weight $S \in \mathbb{R}^{n \times m}$, producing the GCARE: $$A^T X + X A - (X B + S) R^\dagger (B^T X + S^T) + Q = 0,$$ with the constraint $\ker R \subseteq \ker(X B + S)$ (Ferrante et al., 2013).

2. Unified Framework: Differential–Algebraic Connection and Discrete–Continuous Limits

Finite-horizon problems yield the differential Riccati equation (DRE): $$\dot{P}(t) + A^T P(t) + P(t) A - P(t) B R^{-1} B^T P(t) + Q = 0, \quad P(T) = H,$$ whose steady-state ($T \to \infty$) is exactly the algebraic CARE (2002.04246).

A comprehensive Riccati hierarchy connects continuous-time and sampled-data (discrete-time) regulators: as the time-horizon $T\to\infty$ and sampling step $\Delta \to 0$, the sampled-data algebraic Riccati equation (SD-ARE) converges to the permanent CARE (P-ARE). All limits (sample–continuous, finite–infinite horizon) commute: $$\text{SD-DRE} \xrightarrow{\Delta\to 0} \text{P-DRE} \xrightarrow{T\to\infty} \begin{array}{ll} \text{SD-ARE} \xrightarrow{\Delta\to 0} \text{P-ARE} \end{array}$$ (2002.04246).

3. Computational Algorithms: Direct, Iterative, and Low-Rank Factorizations

Dense problems ($n\lesssim10^3$) are frequently resolved by Schur-based methods (eigendecomposition of Hamiltonian pencils), yielding $O(n^3)$ complexity (Kao et al., 2020, Haqiri et al., 2015). For large-scale, sparse or low-rank instances, structure-exploiting iterative methods dominate:

Projection Methods: Extended/rational/tangential Krylov subspace approaches construct a basis $Q_j$ (block-Krylov or rational Arnoldi) and solve small projected CAREs for $Y_j$, lifting to $X_j=Q_j Y_j Q_j^T$ (Benner et al., 2018, Bertram et al., 2023, Bertram et al., 2023). The block rational Krylov/Arnoldi framework unifies RADI, qADI, Cayley, and invariant subspace approaches—each yields the same sequence $(X_j)$ for shared shifts (Bertram et al., 2023, Bertram et al., 2023).

ADI-Type Iterations (RADI): Low-rank ADI variants apply successive shifted solves, updating factors $Z_k$ and $Y_k$ via rank-$p$ corrections, with explicit residual monitoring via $R_k$ (Benner et al., 2015, Zhang et al., 2024). Sparse linear algebra and Sherman–Morrison–Woodbury formula are recurrent.

Doubling Algorithms: Structure-preserving doubling (SDA), alternating direction doubling (ADDA), and decoupled/low-rank adaptions exhibit quadratic or doubly-exponential convergence via symplectic pencil transformations. Low-rank truncation (dSDA$_\mathrm{t}$, R-ADDA) with rank-revealing SVD restricts factor growth (Guo et al., 2020, Zhang et al., 2024).

Newton–Kleinman and Lyapunov Embedding: Newton iteration transforms CARE into a Lyapunov equation for the next iterate. Sparse/low-rank ADI solvers compute Lyapunov steps. Indefinite LDL$^T$ factorizations (LDL–Newton) facilitate general (possibly indefinite) CAREs (Saak et al., 2024, Li et al., 2022).

4. Polynomial Optimization, Verification, and Tensor Extensions

Polynomial optimization via Lasserre's SDP hierarchy offers a direct method for existence/nonexistence certificates of PSD CARE solutions: two convex polynomial programs encode the Cholesky-factor or direct $P$ formulation, relaxed to a sequence of SDPs whose feasibility indicates solution existence (Zhang et al., 2024). The flat-extension criterion allows exact recovery of the global minimizer for finite $k$.

Verification methods (Krawczyk intervals, ADI-inspired fixed-point maps) compute guaranteed enclosures for the stabilizing solution $X_s$, employing interval arithmetic and wrapping-reduction techniques and leveraging preconditioned Hamiltonian decompositions for computation and uniqueness detection (Haqiri et al., 2015).

CAREs generalize to tensor equations (ARTEs) via the Einstein product in even-order tensors. Hamiltonian/symplectic tensor structures enable Schur–Hamiltonian decomposition and symplectic tensor SVD; existence and uniqueness are characterized by multilinear generalizations of stabilizability/detectability. Newton iteration and bounded-real/small-gain theorems extend robust control theory to the multilinear domain (Wang et al., 2024).

5. Extensions: Stochastic, Indefinite, and Singular Problems

Stochastic continuous-time algebraic Riccati equations (SCARE) emerge in optimal control of systems with multiplicative or additive noise: $$A^T X + X A + Q + \sum_{i} A_0^{iT} X A_0^i - (...)\Bigl(R + ...\Bigr)^{-1}(...) = 0$$ (Guo et al., 2024, Huang et al., 2023, Huang et al., 2024). RADI-type and FP-doubling algorithms, combined with SDA, are effective for large-scale SCAREs, exhibiting monotonic and robust convergence properties.

Indefinite quadratic forms (arising e.g. in $H_\infty$ design) are handled by splitting and forward-shifting techniques, reducing indefiniteness to monotonic sequences of definite CAREs, and maintaining low-rank factorization throughout (Benner et al., 2021, Saak et al., 2024). Theoretical convergence is preserved under standard operator-theoretic assumptions, and numerical evidence confirms applicability in both dense and high-dimensional sparse cases.

Singular LQ problems (with $R$ positive semi-definite or singular) necessitate generalized pseudo-inverse variants (GCARE). The impulse-free optimal control law and complementary projector $G=I-R^\dagger R$ distinguish regular (non-impulsive) solutions, with the closed-loop dynamics remaining stable or non-expansive in the singular case (Ferrante et al., 2013).

6. Numerical Performance, Scalability, and Solver Selection

Comprehensive benchmarks indicate that RADI-type solvers, block rational Krylov methods (RKSM, PGKSM), and low-rank doubling algorithms (R-ADDA, dSDA$_\mathrm{t}$) achieve competitive scalability, with memory and compute cost scaling linearly in problem size $n$ for numerically low-rank $X$. Newton–Kleinman (LDL$^T$) achieves high accuracy and quadratic convergence; projection methods are best when subspace dimension remains modest and direct sparse factorizations are feasible (Benner et al., 2018, Zhang et al., 2024, Bertram et al., 2023, Guo et al., 2020).

RADI implementation is particularly favorable for applications needing stabilizing feedback only (memory $O(n(p+m))$), whereas projection methods deliver full low-rank $X$ with cubic cost in the subspace rank. Monotonic algorithms admit robust initializations, and hybrid schemes (FP-care SDA $\to$ Newton/Lyap SDA) balance global convergence and local speed, especially crucial in stochastic or ill-posed control instances (Huang et al., 2023, Huang et al., 2024).

SDP-based polynomial optimization guarantees detection of nonexistence. Automatic differentiation frameworks now support direct differentiation of CARE solvers in both forward and reverse mode—a key enabler for inverse control and learning applications (Kao et al., 2020).

7. Research Directions and Open Problems

Recent advances include generalization to stochastic and multilinear tensor systems, rank-adaptive and memory-efficient implementations (R-ADDA, dSDA$_\mathrm{t}$), and verification techniques via interval/Krawczyk and polynomial optimization hierarchies. Open problems concern the a priori sharpening of existence/nonexistence criteria (especially for singular/indefinite cases), extension of verified solvers to discrete-time, higher-index descriptor, and non-symmetric Riccati equations, and further integration with automatic differentiation and operator learning pipelines (Zhang et al., 2024, Huang et al., 2023, Haqiri et al., 2015).

The continuous-time algebraic Riccati equation remains a rich area at the confluence of control theory, numerical analysis, operator theory, and optimization, with ongoing work addressing more challenging classes (indefinite, singular, stochastic, multilinear), faster/robust solvers, and deeper links to computational convexity and learning frameworks.