category: literaturenote citekey: gielerbridgingcontinuoustimelqr2025 title: Bridging Continuous-time LQR and Reinforcement Learning via Gradient Flow of the Bellman Error authors: "Gießler, Armin; Malan, Albertus Johannes; Hohmann, Sören" year: 2025 date: 2025-00-00 2025 doi: 10.1109/CDC57313.2025.11312822 url: "https://arxiv.org/abs/2506.09685v1" zotero_key: 6RJPWE4J zotero_storage: DKCRZZX6 collections: transpordi_planeerimine folder: 001_artiklid firstAuthor: "Gießler, Armin"

status: converted

Bridging Continuous-time LQR and Reinforcement Learning via Gradient Flow of the Bellman Error

Armin Gießler, Albertus Johannes Malan and Soren Hohmann ¨

*Abstract*— In this paper, we present a novel method for computing the optimal feedback gain of the infinite-horizon Linear Quadratic Regulator (LQR) problem via an ordinary differential equation. We introduce a novel continuous-time Bellman error, derived from the Hamilton-Jacobi-Bellman (HJB) equation, which quantifies the suboptimality of stabilizing policies and is parametrized in terms of the feedback gain. We analyze its properties, including its effective domain, smoothness, coerciveness and show the existence of a unique stationary point within the stability region. Furthermore, we derive a closed-form gradient expression of the Bellman error that induces a gradient flow. This converges to the optimal feedback and generates a unique trajectory which exclusively comprises stabilizing feedback policies. Additionally, this work advances interesting connections between LQR theory and Reinforcement Learning (RL) by redefining suboptimality of the Algebraic Riccati Equation (ARE) as a Bellman error, adapting a stateindependent formulation, and leveraging Lyapunov equations to overcome the infinite-horizon challenge. We validate our method in a simulation and compare it to the state of the art.

I. INTRODUCTION

The Linear Quadratic Regulator (LQR) stands as a cornerstone in the realm of control theory. At its core, LQR addresses the fundamental challenge of designing a controller that minimizes a cost function while ensuring system stability. By applying quadratic cost functions to linear systems, LQR is known for its elegant formulation where the optimal solution can be explicitly expressed as a state-feedback control law [1]. Over the last six decades, numerous methods have been developed to solve the LQR problem, the majority of which rely on obtaining the maximum solution1 of the Algebraic Riccati Equation (ARE). These methods include iterative algorithms for discrete-time [2] and continuous-time systems [3], as well as linear-algebra-based approaches that perform eigenvalue decomposition on the Hamiltonian matrix [4]. Additionally, the LQR problem can be transformed into a semidefinite program (SDP) [5], [6], or equivalently expressed as an H² norm minimization problem which can be cast into a convex optimization problem [7], [8].

More recently, Reinforcement Learning (RL) has offered a new perspective on learning optimal controllers for systems with nonlinear, stochastic, or uncertain dynamics [9]. While RL has advanced rapidly with powerful algorithms to search for optimal policies2 , ensuring the stability of these learned policies remains challenging. This has led to a growing intersection between RL and control theory [10]–[12].

Literature review: In contrast to methods relying on the ARE, an alternative approach to solving the LQR problem involves directly minimizing the infinite-horizon LQR cost that is parametrized in the feedback matrix. By using Lyapunov equations to handle the infinite horizon, the LQR problem can be minimized using first-order methods with a closed-form gradient. The gradient expression was initially derived for continuous-time systems in [13], [14] and later adopted for discrete-time systems in [15]. Rigorous analysis of this approach, including properties such as smoothness, coerciveness, and quadratic growth, are presented in [16], [17] for discrete-time systems and [18], [19] for continuoustime systems.

In the context of RL, policy iteration (PI) iteratively determines the optimal policy for continuous-time systems [12]. This method consists of a policy evaluation step, where the continuous-time Bellman equation is used to compute the value function for a stabilizing policy, followed by a policy improvement step to update the control law. PI solves the non-linear Hamilton-Jacobi-Bellman (HJB) equation, which is a sufficient optimality condition3 , by iterations on equations linear in the gradient of the value function. In the LQR case, the HJB equation reduces to the continuous-time ARE (CARE) and PI is exactly Kleinman's Algorithm [3]. Integral RL (IRL) replaces the continuous-time Bellman equation4 with an integral variant, enabling learning without knowledge of the drift dynamics5 or state derivatives [12], [20], [21]. In IRL, the Bellman error, defined on this IRL Bellman equation and parametrized in the value matrix, is minimized during the policy evaluation step.

All the mentioned RL methods learn the optimal value function and derive the policy from it. By directly parametrizing the learning objective in the policy, fewer parameters are generally required, potentially leading to faster convergence and reduced computational complexity in the optimization process. To the best of the authors' knowledge, no existing work parametrizes the HJB equation directly in terms of the feedback gain and formulates a continuous-time Bellman error based on this parametrization.

A. Gießler, A. J. Malan, S. Hohmann are with the Institute of Control Systems, Karlsruhe Insitute of Technology (KIT), 76131, Karlsruhe, Germany. Corresponding author is Armin Gießler, armin.giessler@kit.edu.

¹Maximality is considered in terms of the partial order induced by positive semidefinite matrices.

²Policy is a synonym for control law.

³For input-affine systems and stage costs, which quadratically penalize the input and contain positive definite functions in the states, the HJB equation provides both a necessary and sufficient optimality condition.

⁴ In contrast to the discrete-time Bellman equation, the continuous-time Bellman equation depends on the full system dynamics [12, Sec. 2.6].

⁵The drift dynamics corresponds to the system matrix in the linear case.

Contributions: The main contributions of the paper at hand are:

1) Based on the HJB equation, we formulate a novel Bellman error for a continuous-time LQR, designed to quantify and address the suboptimality of stabilizing policies. The resulting matrix function parametrizes the error directly in terms of the feedback gain matrix.
2) We rigorously analyze the properties of the Bellman error by addressing its effective domain and proving smoothness, coerciveness, and the existence of a unique stationary point.
3) We derive a closed-form expression for the gradient of the Bellman error and establish its well-definedness within the region of stabilizing policies. Furthermore, we show that the corresponding gradient flow generates a unique trajectory of stabilizing feedbacks and converges to the solution of the LQR problem.

Paper Organization: This section ends with some notation and preliminaries. Section II introduces the novel Bellman error formulation and examine its properties. In Section III, the gradient expression and the corresponding gradient flow are presented. Section IV provides simulation results and compares our method to a state-of-the-art approach. Finally, the paper ends with a conclusion in Section V.

Notation: We denote the set of (nonnegative) real numbers as R (R^≥0). The transpose of a matrix A (a vector B) is denoted by A^⊤ (B^⊤). The identity matrix of dimensions n× n is given by In. A symmetric positive definite (semidefinite) matrix A is denoted by A ≻ 0 (A ⪰ 0). The trace of a square matrix A is denoted by tr(A). The spectral (frobenius) norm of a matrix A is denoted by ∥A∥² (∥A∥^F ). The spectrum of a square matrix A, containing the eigenvalues λi(A), is denoted by σ(A). The minimum eigenvalue of A is represented by λmin(A). The vectorization and rank of a matrix A are denoted by vec(A) and rank(A), respectively. The gradient of a scalar function f with respect to a vector x (matrix A), denoted by ∇xf (∇Af), is a column vector (matrix) of partial derivatives. A square matrix A is Hurwitz stable if its eigenvalues have strictly negative real parts.

A. Preliminaries

1) Mathematical Preliminaries: The continuous-time Lyapunov equation is given by

$$AP + PA^{\top} + L = 0, \tag{1}$$

where A, P, L ∈ R ⁿ×ⁿ and P, L are symmetric. If A is Hurwitz stable, the solution of (1) can be written as P = R ^∞ e AtLe^A ^⊤^tdt < ∞.

The vectorization of the matrix product AB, where A ∈ R k×l , B ∈ R ^l×^m, is given by

$$\operatorname{vec}(AB) = (I_m \otimes A) \operatorname{vec}(B) = (B^\top \otimes I_k) \operatorname{vec}(A).$$

The inverse of the vectorization operator vec⁻¹ k×l : R kl → R k×l can be defined as

$$A \mapsto \operatorname{vec}_{k \times l}^{-1}(A) = \left(\operatorname{vec}^{\top}(I_l) \otimes I_k\right) (I_l \otimes A).$$ (2)

The Lyapunov equation (1) can be solved using vectorization6

$$\operatorname{vec}(P) = -\left(I_n \otimes A + A^{\top} \otimes I_n\right)^{-1} \operatorname{vec}(L).$$ (3)

Since the trace operator is a linear operator, tr(A + B) = tr(A) + tr(B) and tr(cA) = c tr(A) hold. Additionally, the trace is invariant under circular shifts and transposition, i.e., tr(ABC) = tr(BCA) = tr(CAB) and tr(A) = tr(A⊤). The trace of a square matrix P A can be calculated as tr(A) = n ⁱ=1 λi(A).

2) Linear Time-Invariant Systems: We consider a continuous-time linear time-invariant (LTI) system

$$\dot{x}(t) = Ax(t) + Bu(t),\tag{4}$$

$$y(t) = Cx(t), (5)$$

where A ∈ R ⁿ×n, B ∈ R ⁿ×m, and C ∈ R ^q×n. The pair (A, B) is stabilizable if and only if, for every λ ∈ σ(A) with Re(λ) ≥ 0, the matrix λIⁿ − A B has full rank. Similarly, the pair (A, C) is detectable if and only if, for every λ ∈ σ(A) with Re(λ) ≥ 0, the matrix λIⁿ − A C has full rank. There exists a time-invariant state feedback controller such that A − BK is Hurwitz stable if and only if the system (A, B) is stabilizable [23, Theorem 14.5]. The set of Hurwitz stabilizing feedback gains is denoted by

$$\mathcal{K} = { K \in \mathbb{R}^{m \times n} : A - BK \text{ is Hurwitz stable} }. \tag{6}$$

The set K is in general open, unbounded and path-connected [24, Section 3].

3) Linear Quadratic Regulator: We briefly recall the standard continuous-time LQR problem which is solved in the remainder of this work [25].

Definition 1 (LQR Problem): The LQR problem is defined as

$$\min_{x,u} \int_0^\infty u^\top(t) R u(t) + x^\top(t) Q x(t) dt,$$ s.t. $\dot{x}(t) = A x(t) + B u(t),$

$$x(0) = x_0,$$ (7)

where Q ⪰ 0, R ≻ 0 are weighting matrices of appropriate dimensions.

If the system (A, B) is stabilizable and the pair (A, ^√ Q) is detectable, the optimal value V ∗ (x0) of the optimization problem is finite and the optimal solution7 u ^∗ = −K^∗x stabilizes the system (A, B). The time-invariant optimal feedback, often referred to as the Kalman gain, is given by K^∗ = R⁻¹B^⊤P ∗ , where the Riccati matrix P ∗ is the positive definite solution of the CARE

$$A^{\top}P + PA - PBR^{-1}B^{\top}P + Q = 0.$$ (8)

The optimal value function of the LQR Problem is given by V ∗ (x0) = x ⊤ ⁰ P ^∗x0. We note that the CARE generally has

⁶Rewriting the equation in this form can be computationally expensive and ill-conditioned. A standard approach to solve the Lyapunov equation is the Bartels-Stewart algorithm [22].

In the following, we omit the explicit time dependency t of variables for simplicity, except where the dependency on t is relevant.

multiple solutions. Under the stabilizability and detectability conditions, there is only one unique positive definite stabilizing solution such that the matrix A−BR−1B⊤P is Hurwitz stable. By substituting K = R−1B⊤P, (8) can be rewritten as

$$A_K^{\top} P + P A_K + Q + K^{\top} R K = 0, \tag{9}$$

where A^K := A − BK. Note that (9) has the form of a Lyapunov equation if K is considered independent of P. We highlight that the optimal feedback cannot be obtained from (9) since there are two unknown matrices P and K. For a stabilizing feedback K, the corresponding value matrix P can be calculated by solving the Lyapunov equation (9).

4) Policy Gradient Flow of LQR: The cost of the LQR problem (7) can be parametrized in K

$$f_K = x_0^{\top} P_K x_0 = \text{tr}(P_K x_0 x_0^{\top}),$$ (10)

where P^K 8 is the solution of (9) for a given (non-optimal) feedback K [18, Section 3.2]. The cost f^K is minimized by the gradient flow [18, Lemma 4.1]

$$\dot{K}(t) = -\nabla_K f_K(t), \tag{11}$$

$$K(0) = K_0 \in \mathcal{K},\tag{12}$$

where

$$\nabla_K f_K = 2 \left( RK - B^\top P_K \right) Y_K, \tag{13}$$

$$Y_K = \int_0^\infty e^{A_K t} x_0 x_0^\top e^{A_K^\top t} \mathrm{d}t. \tag{14}$$

The matrix Y^K can be obtained by solving the matrix Lyapunov equation AKY^K + YKA^⊤ ^K + x0x ⊤ ⁰ = 0. The gradient flow (11) converges to K^∗ if and only if the initial feedback K⁰ stabilizes the system (A, B), i.e., K⁰ ∈ K. The trajectory K(t) generated by (11) always remains in the set of stabilizing feedbacks K.

5) Reinforcement Learning: The LQR problem in Definition 1 is a special case of continuous-time Reinforcement Learning for input-affine systems [12, Chapter 3]. For this work, we focus exclusively on the LQR case and present the relevant formulas. The value function of a stabilizing policy u = µ(x) is defined as

$$V^{\mu}(x_0) = \int_0^{\infty} x^{\top} Q x + \mu^{\top}(x) R \mu(x) dt < \infty,$$ (15)

and can be seen as a cost-to-go from state x⁰ at time instant t under the policy µ(x). By using Leibniz' formula, the infinitesimal version of (15) is

$$0 = \underbrace{x^{\top}Qx + \mu^{\top}(x)R\mu(x) + (\nablax V^{\mu}(x))^{\top}\dot{x}}{=:H(x,\mu(x),\nabla_x V^{\mu})}, \tag{16}$$

where V ^µ(0) = 0. Equation (16) can be considered as a Bellman equation for continuous-time systems9 . The optimal value function V ∗ (x) satisfies the Hamilton-Jacobi-Bellman (HJB) equation

$$0 = \inf_{\mu} H(x, \mu(x), \nabla_x V^*), \tag{17}$$

which reduces to the CARE (8) in the LQR case. The following policy iteration algorithm solves the nonlinear HJB equation (17) by iterating over equations which are linear in the gradient of the value function.

Algorithm 1 (Policy iteration [12, Algorithm 4.1]): Select a stabilizing policy µ0.

1) Policy evaluation step: Solve Bellman equation for V µi

$$0 = x^{\top} Q x + \mu_i^{\top} R \mu_i + (\nabla_x V^{\mu_i}(x))^{\top} (A x + B \mu_i).$$ (18)

2) Policy improvement step: Update the policy using

$$\mu{i+1} = \arg\min{\mu} H(x, \mu(x), \nabla_x V^{\mu_i}),$$ (19)

$$= -\frac{1}{2}R^{-1}B^{\top}\nabla_x V^{\mu_i}.$$ (20)

Algorithm 1 is guaranteed to converge to the optimal policy, but has the drawback that full knowledge of the system dynamics, i.e. (A, B), is needed [26]. For the LQR, Algorithm 1 is equivalent to Kleinman's Algorithm [3], where the two steps are

$$(A - BK_i)^{\top} P_i + P_i (A - BK_i) = -(K_i^{\top} RK_i + Q),$$

$K_{i+1} = R^{-1} B^{\top} P_i.$ (21)

II. PROBLEM FORMULATION AND ANALYTIC PROPERTIES

In this section, we define our Bellman error and study its analytic properties. To ensure that the LQR optimization problem is well-posed and that the optimal solution yields a stabilizing controller, we impose the following assumption.

Assumption 1: The system (A, B) is stabilizable and the pair (A, ^√ Q) is detectable. Additionally, B ̸= 0 and Q ̸= 0.

Throughout the remainder of this paper, we assume that Assumption 1 hold. A policy for a given value function V ^µ(x) = x ^⊤P x, P ≻ 0 can be calculated by (20)10

$$\mu = -R^{-1}B^{\top}Px. \tag{22}$$

By inserting (22) into the HJB equation (17), we obtain

$$0 = x^{\top} \underbrace{(2PA - PBR^{-1}B^{\top}P + Q)}_{-\tilde{M}} = \operatorname{tr}(\tilde{M}xx^{\top}), \qquad (23)$$

$$= x^{\top} \underbrace{\left(\underline{A^{\top}P + PA - PBR^{-1}B^{\top}P + Q}\right)}_{=:M} x = \operatorname{tr}(Mxx^{\top}),$$ (24)

where M = 0 corresponds to the CARE.

Remark 1: Only the symmetric part ¹ (M˜ ^⊤ + M˜ ) = M contributes to the quadratic form and the trace in (23). Therefore, (23) and (24) are equivalent.

It is well known that the optimal feedback gain K of the LQR problem (7) is state-independent. Inspired by the initial

⁸The subscript K of matrices emphasize that the matrix depends on feedback gain K. An exception is the definition A^K = A − BK.

In contrast to the discrete-time Bellman equation, V (xk) = x^⊤ ^k Qx^k + u^⊤ ^k Ru^k + V (xk+1), (16) depends on the system dynamics x˙.

¹⁰In general, there exist infinitely many policies for a value matrix P if P ̸= P ∗.

state-independent formulation of [17], [18], we reformulate (24) as

$$0 = \sum{j=1}^{n} \operatorname{tr}(Mx{j}x_{j}^{\top}), \tag{25}$$

where the states x1, . . . , xⁿ are linearly independent. By using the standard basis vectors as the states and the linearity P of the trace operator, (25) reduces to 0 = tr(M), since n ^j=1 xjx ⊤ ^j = In. The equation (25) holds if P = P ∗ , since M = 0. 11 In order to quantify the suboptimality of a given stabilizing feedback, we make the following definition.

Definition 2 (Bellman error): The Bellman error parametrized in the feedback K is defined as the matrix function e^K : R ^m×ⁿ → R

$$K \mapsto e_K = -\operatorname{tr}(\underbrace{A^{\top}P_K + P_K A - P_K B R^{-1} B^{\top} PK + Q}{=:M_K}),$$ (26)

where the matrix P^K satisfies the Lyapunov equation

$$0 = A_K^{\top} P_K + P_K A_K + Q + K^{\top} R K.$$ (27)

A. Analytic Properties of the Bellman Error

We now investigate for which K ∈ R ^m×ⁿ the Lyapunov equation (27) has a unique solution.

Lemma 1: The Lyapunov equation (27) has a unique solution P^K if and only if A^K and −A^K have no eigenvalue in common, i.e., σ(AK) ∩ σ(−AK) = ∅.

Proof: The Lyapunov equation (27) is a special case of the Sylvester equation AX − XB = C. By Sylvester equation Theorem, (27) has a unique solution if and only if σ(AK) ∩ σ(−AK) = ∅ [27, Theorem 2.4.4.1].

The set of feedback matrices that satisfy the condition of Lemma 1 is denoted by

$$\mathcal{K}_{\sigma} := { K \in \mathbb{R}^{n \times m} : \sigma(A_K) \cap \sigma(-A_K) = \emptyset }.$$ (28)

Remark 2: The Lyapunov equation (27) can be solved using vectorization for K ∈ K^σ

$$\operatorname{vec}(P_K) = -\left(I_n \otimes A_K^\top + A_K^\top \otimes I_n\right)^{-1} \operatorname{vec}(Q + K^\top RK). \tag{29}$$

By applying the inverse vectorization operator (2), the Bellman error (26) can be expressed in closed-form for K ∈ Kσ.

Lemma 2: If A^K is Hurwitz stable, then the solution P^K of the Lyapunov equation (27) is unique and positive semidefinite.

Proof: Since A^K is Hurwitz stable, its eigenvalues have strictly negative real parts, ensuring σ(AK) ∩ σ(−AK) = ∅. By Lemma 1, (27) has a unique solution PK. Moreover, the Hurwitz stability of A^K implies that P^K can equivalently be represented as

$$P_K = \int_0^\infty e^{A_K^\top t} \left( Q + K^\top R K \right) e^{A_K t} dt. \tag{30}$$

Since Q ⪰ 0 and R ≻ 0, Q + K⊤RK ⪰ 0. Therefore, the integrand is positive semidefinite for all t ≥ 0, leading to P^K ⪰ 0.

Lemma 3: The matrix M^K of Definition 2 is negative semidefinite for K ∈ Kσ. Furthermore, the image of the matrix function e^K over K^σ satisfies eK(Kσ) ⊆ R≥0.

Proof: In the following, we find an equivalent expression for M^K by using the Lyapunov equation (27). By substituting A^K = A − BK into (27) and solve for A⊤P^K + PKA, we obtain

$$A^{\top} P_K + P_K A = -Q - K^{\top} R K + K^{\top} B^{\top} P_K + P_K B K$$ (31)

Inserting (31) into M^K yields

$$M_K = -K^{\top} R K + K^{\top} B^{\top} P_K + P_K B K - P_K B R^{-1} B^{\top} P_K,$$

= $-\left(K - R^{-1} B^{\top} P_K\right)^{\top} R\left(K - R^{-1} B^{\top} P_K\right).$ (32)

Since R ≻ 0 and P^K is well-defined and finite for K ∈ K^σ (see Lemma 1), the matrix M^K is negative semidefinite [27, Observation 7.1.8]. Moreover, since M^K has only nonpositive eigenvalues, − tr(MK) = − P i λi(MK) ≥ 0 and thus eK(Kσ) ⊆ R^≥⁰ hold.

The following four lemmas are closely related to results presented in [17], [18], where similar conclusions were drawn for the continuous-time and discrete-time LQR cost. We adapt these results to the Bellman error (26). In the next lemma, we study for which feedback K the Bellman error is well-defined and address the effective domain dom(eK) = {K ∈ R ^m×ⁿ : e^K < ∞}.

Lemma 4: The effective domain of the matrix function e^K (26) is given by dom(eK) = Kσ.

Proof: Through Lemma 1, the Lyapunov equation (27) has a unique solution P^K if and only if K ∈ Kσ. When P^K exists, M^K is well-defined as M^K = A^⊤P^K + PKA − PKBR⁻1B^⊤P^K + Q, and consequently e^K = − tr(MK) is finite. If K /∈ Kσ, the Lyapunov eqution does not admit a solution. Specifically, P^K grows unbounded, leading to e^K → ∞. Therefore, e^K is well-defined and finite if and only if K ∈ Kσ.

Remark 3: In contrast to [18], where the set of stabilizing feedback gains K is the interior of dom(fK), the effective domain of e^K contains unstable feedbacks, i.e., K ∈ K^σ \K.

Lemma 5: The Bellman error (26) is real analytic over Kσ, i.e., e^K ∈ C ^ω(Kσ). 12

Proof: For every K ∈ Kσ, the matrix Iⁿ ⊗A^⊤ ^K +A^⊤ ^K ⊗ Iⁿ (see (29)) is invertible, since the eigenvalues of it are {λi(AK) + λ^j (AK) | 1 ≤ i, j ≤ n} and hence nonzero. By Cramer's rule, the map K 7→ P^K (29) is a rational function of polynomials in the entries of K and thus C ^ω. The matrix function e^K is a composition of C ^ω maps

$$K \mapsto P_K \mapsto M_K \mapsto -\operatorname{tr}(M_K)$$ (33)

and thus real analytic over Kσ.

Since the Bellman error e^K is real analytic over Kσ, it is infinitely differentiable and all derivatives are continuous. For the remainder of this section, we analyze e^K over the

¹¹Note that, more generally, any matrix M for which the sum of its eigenvalues is zero satisfies (25), not only M = 0.

¹²The class C^ω denotes smooth functions that can be locally represented by convergent power series.

set of stabilizing feedbacks K ⊂ Kσ. We observe that e^K diverges to infinity either as K approaches the boundary of K13

$$\partial \mathcal{K} \subseteq { K \in \mathbb{R}^{m \times n} : \max_{\lambda \in \sigma(A_K)} \text{Re}(\lambda) = 0 },$$ (34)

or when K ∈ K becomes unbounded.

Lemma 6: The Bellman error e^K is coercive over K in the sense that

$$\lim_{Ki \to K \in \partial \mathcal{K}} e{Ki} = \infty \quad \text{and} \quad \lim{|K|_2 \to \infty, K \in \mathcal{K}} e_K = \infty.$$

Proof: Let {Ki} be a sequence in K converging to a fixed K ∈ ∂K. The corresponding sequence of matrices {PKⁱ } satisfies (29). As Kⁱ → K ∈ ∂K, the matrix Iⁿ ⊗ A^⊤ ^Kⁱ + A^⊤ ^Kⁱ ⊗ Iⁿ has at least one eigenvalue approaching zero and thus ∥ Iⁿ ⊗ A^⊤ ^Kⁱ + A^⊤ ^Kⁱ ⊗ Iⁿ ⁻¹ ∥² → ∞. Since vec(Q + K^⊤ ⁱ RKi) remains bounded for all Kⁱ , it follows that ∥PKⁱ ∥² → ∞ as Kⁱ → K ∈ ∂K. On the other hand, consider ∥K∥² → ∞ with K ∈ K. Taking the spectral norm on the Lyapunov equation (27) yields

$$||(A - BK)^{\top} P_K + P_K (A - BK)||_2 = || - Q - K^{\top} RK||_2.$$ (35)

The left-hand side can be upper bounded with the triangle inequality and submultiplicativity as

$$|(A - BK)^{\top} P_K + P_K (A - BK)|_2$$ (36)

$$\leq 2|A|_2|P_K|_2 + 2|B|_2|K|_2|P_K|_2 \tag{37}$$

and a lower bound for the right-hand side can be obtained by using the reverse triangle inequality14 and the inequality |∥A∥ − ∥B∥| ≥ ∥B∥ − ∥A∥

$$|| - Q - K^{\top} R K ||_2 \ge ||| - Q ||_2 - || K^{\top} R K ||_2 |$$ (38)

$$\geq |K^{\top} R K|_2 - | - Q|_2 \tag{39}$$

$$\geq \lambda_{\min}(R) ||K||_2^2 - ||Q||_2. \tag{40}$$

Combining these bounds result in

$$\begin{split} 2\left(|A|{2} + |B|{2}|K|{2}\right)|P{K}|{2} &\geq \lambda{\min}(R)|K|{2}^{2} - |Q|{2}, \ \Leftrightarrow |P{K}|{2} &\geq \frac{\lambda{\min}(R)|K|{2}^{2} - |Q|{2}}{2|A|{2} + 2|B|{2}|K|{2}}. \end{split}$$

Therefore ∥PK∥² → ∞ as ∥K∥² → ∞ with K ∈ K. Next, we observe that

$$\begin{split} |M_K|_2 &\geq ||A^\top P_K + P_K A + Q|_2 - |PBR^{-1}B^\top P|_2| \ &\geq |P_K B R^{-1}B^\top P_K|_2 - |A^\top P_K + P_K A + Q|_2 \ &\geq |P_K B R^{-1}B^\top P_K|_2 - |2P_K A|_2 - |Q|2 \ &\geq \lambda{\min}^+ (BR^{-1}B^\top) |P_K|_2^2 - 2|A|_2 |P_K|_2 - |Q|_2, \end{split}$$

¹³For LTI systems with m = 1, the equality in (34) holds. However, for m > 1, establishing whether equality holds or only a subset relation applies becomes significantly more intricate (see [24, Prop. 3.5]).

where λ + min(·) denotes the smalles non-zero eigenvalue of a matrix. For sufficiently large ∥PK∥2, the quadratic term dominates the other two terms, leading to ∥MK∥² → ∞ as ∥PK∥² → ∞. Since M^K ⪯ 0 (see Lemma 3), its eigenvalues are non-positive, ensuring e^K = − tr(MK) ≥ 0. Given that B ̸= 0 (see Assumption 1), it follows that e^K → ∞ as ∥PK∥² → ∞. Combining these results, we conclude that e^K → ∞ as Kⁱ → K ∈ ∂K or ∥K∥² → ∞.

With coercivity established, we now analyze the compactness of the sublevel sets of e^K over K.

Lemma 7: For every α > 0, the sublevel set of the Bellman error e^K over K,

$$L_{\alpha} = { K \in \mathcal{K} : e_K \le \alpha },\tag{41}$$

is compact.

Proof: To prove that the sublevel set L^α is compact, we establish that L^α is both closed and bounded. By Lemma 5, e^K is real analytic over K^σ and hence continuous over K ⊂ Kσ. The sublevel set L^α is the preimage of the closed interval (−∞, α] under eK. Continuity of e^K ensures that L^α is closed in K [28, Theorem 18.1]. Lemma 6 establish that e^K → ∞ as ∥K∥² → ∞ or K → ∂K ∈ K. Suppose now that L^α is unbounded. Then there exists a sequence {Kn} ⊂ L^α with ∥Kn∥² → ∞ or Kⁿ → K ∈ ∂K. By coercivity, eKⁿ → ∞, which contradicts eKⁿ ≤ α. Thus, L^α is bounded. Since L^α is both closed and bounded, it is compact.

III. GRADIENT FLOW OF THE BELLMAN ERROR

In this section, we study the behavior of the gradient flow of the Bellman error. Since the Bellman error e^K is real analytic over K (see Lemma 5), the gradient is well-defined at all points within this domain. The closed from expression of the gradient is provided in the following proposition.

Proposition 1: For K ∈ K, the gradient of the Bellman error e^K is given by16

$$\nabla_K e_K = -4 \left( RK - B^\top P_K \right) X_K, \tag{42}$$

where

$$0 = A_K X_K + X_K A_K^{\top} + \frac{1}{2} \left( \tilde{A}_K + \tilde{A}_K^{\top} \right),$$ (43)

$$\tilde{A}_K := A - BR^{-1}B^{\top}P_K. \tag{44}$$

Proof: The Bellman error is a composition of C ^ω maps, K 7→ P^K 7→ M^K 7→ − tr(MK). The differential of (27) with respect to K is given by

$$0 = dA_K^{\mathsf{T}} P_K + A_K^{\mathsf{T}} dP_K$$

$$+ dP_K A_K + P_K dA_K + dK^{\mathsf{T}} RK + K^{\mathsf{T}} RdK \quad (45)$$

$$= (-BdK)^{\mathsf{T}} P_K + A_K^{\mathsf{T}} dP_K$$

$$+ dP_K A_K + P_K (-BdK) + dK^{\mathsf{T}} RK + K^{\mathsf{T}} RdK$$

$$= A_K^{\mathsf{T}} dP_K + dP_K A_K^{\mathsf{T}} + \underbrace{dK^{\mathsf{T}} (RK - B^{\mathsf{T}} PK)}{=:U}$$

$$+ (K^{\mathsf{T}} R - P_K B) dK \quad (46)$$

¹⁶We note that the gradient is given in denominator layout notation such that ∇Ke^K ∈ Rm×ⁿ has the same dimensions as K.

¹⁴The reverse triangle inequality state that ∥A − B∥ ≥ |∥A∥ − ∥B∥|.

¹⁵The inequality holds universally due to the property of absolute values |x| ≥ ±x for all x ∈ R. Specifically, letting x = ∥A∥ − ∥B∥, we have |x| ≥ −x, which simplifies to |∥A∥−∥B∥| ≥ ∥B∥−∥A∥. Equality occurs when ∥B∥ ≥ ∥A∥.

$$= A_K^{\mathsf{T}} dP_K + dP_K A_K^{\mathsf{T}} + U + U^{\mathsf{T}}. \tag{47}$$

Since A^K is Hurwitz stable for every K ∈ K, the solution of the Lyapunov equation (47) is unique and can be expressed as

$$dP_K = \int_0^\infty e^{A_K^\top t} (U + U^\top) e^{A_K t} dt.$$ (48)

The differential of M^K 26 is given by

$$dM_K = A^{\top} dP_K + dP_K A - dP_K B R^{-1} B^{\top} P_K - P_K B R^{-1} B^{\top} dP_K.$$ (49)

Putting (49) into the differential of Bellman error yields

$$de_K = -\operatorname{tr}(dM_K) \tag{50}$$

$$= -2\operatorname{tr}(dP_K A) + 2\operatorname{tr}(dP_K B R^{-1} B^{\top} P_K)$$ (51)

$$= -2\operatorname{tr}\left(dP_K(A - BR^{-1}B^{\top}P_K)\right) \tag{52}$$

$$= -2\operatorname{tr}\left(\int{0}^{\infty} e^{A{K}^{\top}t}(U+U^{\top})e^{A{K}t}dt\tilde{A}{K}\right)$$ (53)

$$= -4\operatorname{tr}\left(\int_0^\infty e^{A_K^\top t} U e^{A_K t} \mathrm{d}t \tilde{A}_K\right). \tag{54}$$

By using the cyclic property and the linearity of the trace operator17, we obtain

$$de_K = -4\operatorname{tr}\left(U\int_0^\infty e^{A_K t} \tilde{A}_K e^{A_K^{\mathsf{T}} t} dt\right). \tag{55}$$

Since the trace of a skew-symmetric matrix is always zero, only the symmetric part of the integral in (55) contributes to the trace. The symmetric part of the matrix R ^∞ e ^AKtA˜Ke A ⊤ ^Ktdt is equal to the solution X^K of Lyapunov equation

$$A_K X_K + X_K A_K^{\top} + \frac{1}{2} \left( \tilde{A}_K + \tilde{A}_K^{\top} \right) = 0.$$ (56)

Combining everything yields

$$de_K = -4\operatorname{tr}\left(dK^{\top}\left(RK - B^{\top}P_K\right)X_K\right). \tag{57}$$

By using the relation df = tr(dZW) ⇒ ∇Zf = W [29], we obtain (42).

Remark 4: Note that the gradient (42) only exists for K ∈ K. Otherwise, for K ∈ K^σ \ K, the integral in (55) diverges because A^K is not Hurwitz stable.

Theorem 1: For K ∈ K, the minimization problem of the Bellman error (26), i.e.,

$$\min_{K \in \mathcal{K}} e_K, \tag{58}$$

has a unique global minimizer K^∗ = R⁻¹B^⊤P ∗ , where P ∗ is the unique positive definite solution of the CARE (8). Furthermore, K^∗ is the only stationary point of e^K on K.

Proof: The Bellman error e^K is continuous (see Lemma 5) and has nonempty compact sublevel sets on K (see Lemma 7), given that (A, B) is stabilizable (see Assumption 1). By the extreme value theorem [30, Theorem 4.15.], a minimum exists on these sublevel sets. Since e^K is not constant, the minimum lies within Lα. To identify the stationary point, we prove that X^K is nonsingular and solve ∇Ke^K = 0 with respect to K. First, we show that A˜^⊤ ^K ⁺A˜^K is nonsingular. By substituting A^K = A˜^K +BR−1B⊤P^K − BK into (27), we obtain

$$\tilde{A}{K}^{\top}P{K} + P{K}\tilde{A}{K} + 2P{K}BR^{-1}B^{\top}P{K} - K^{\top}B^{\top}P{K} - P{K}BK$$

$$= -Q - K^{\top}RK,$$

$$\Leftrightarrow \tilde{A}{K}^{\top}P{K} + P{K}\tilde{A}{K} = -Q - L^{\top}RL - (K+L)^{\top}R(K+L),$$ (59)

where L = R−1B^T PK. The right-hand side of (59) is negative semidefinite, and nonzero due to Q ̸= 0 (see Assumption 1). Therefore,

$$v^{T} \left( \tilde{A}{K}^{\top} P{K} + P{K} \tilde{A}{K} \right) v \le 0 \tag{60}$$

has to hold for all vectors v ∈ R ⁿ and A˜^⊤ ^KP^K ⁺PKA˜^K ̸= 0. Assume for a contradiction that A˜^⊤ ^K ⁺ ^A˜^K is singular. Then there exists a nonzero vector v such that (A˜^⊤ ^K ⁺ ^A˜K)^v = 0. Substituting A˜Kv = −A˜^⊤ ^Kv into (60) leads to

$$v^{T}(\underbrace{\tilde{A}{K}^{\top}P{K} - P{K}\tilde{A}{K}^{\top}})v \leq 0.$$

$$=:N$$ (61)

Now observe that

$$\operatorname{tr}(N) = \operatorname{tr}(\tilde{A}_K^{\top} P_K) - \operatorname{tr}(P_K \tilde{A}_K^{\top}) \tag{62}$$

$$= \operatorname{tr}(\tilde{A}_K^{\top} P_K) - \operatorname{tr}(\tilde{A}_K^{\top} P_K) = 0.$$ (63)

Since the sum of the eigenvalues of N is zero, N cannot be negative semidefinite and nonzero, leading to a contradiction. Hence, A˜^⊤ ^K ⁺ ^A˜^K is nonsingular. Next we show that ^X^K is nonsingular. Assume X^K is singular, then there exists a nonzero vector v such that XKv = 0. Pre- and postmultiplying (43) by v ^⊤ and v, respectively, yields

$$0 = v^T \left( A_K X_K + X_K A_K^\top + \frac{1}{2} \left( \tilde{A}_K + \tilde{A}_K^\top \right) \right) v \quad (64)$$

$$= v^T \frac{1}{2} \left( \tilde{A}_K + \tilde{A}_K^\top \right) v. \tag{65}$$

Since A˜^K + A˜^⊤ ^K is nonsingular, (65) implies v = 0, which contradicts the assumption that v is nonzero. Hence, X^K is nonsingular. With X^K nonsingular, the equation ∇Ke^K = 0 implies (RK − B^⊤PK) = 0. By setting (RK − B^⊤PK) = 0, we obtain K^∗ = R⁻¹B^⊤P ∗ , where P ∗ is the unique positive definite solution of (8). Therefore, K^∗ is the unique minimizer of e^K within K.

In the next theorem, we propose an optimization algorithm to solve the optimization problem of Theorem 1.

Theorem 2: The gradient flow

$$\dot{K} = -\beta \nabla_K e_K, \quad K(0) = K_0 \in \mathcal{K}, \tag{66}$$

where β > 0, solves the optimization problem of Theorem 1 and converges to the optimal feedback K^∗ . The trajectory K(t) remains within K for all t ≥ 0.

Proof: Let e^K be the Lyapunov functional. Since M^K is negative semidefinite, e^K = − tr(MK) = 0 if and only if M^K = 0, which implies K = K^∗ . Additionally, e^K > 0 for

¹⁷By combining both properties, the equality tr(R R A(x)B(x)dt) = tr(A(x)B(x))dt = R tr(B(x)A(x))dt = tr(R B(x)A(x)dt) hold.

K ̸= K^∗ over K because the codomain of e^K is R≥⁰ (see Lemma 3). The time derivative of e^K is given by

$$\dot{e}_K = \frac{d}{dt}e_K = \operatorname{tr}\left(\left(\nabla_K e_K\right)^\top \dot{K}\right) = -\alpha |\nabla_K e_K|_F^2 \le 0.$$ (67)

Since e˙^K < 0 for all K ∈ K \ {K^∗}, e^K is a proper Lyapunov function over K and the equilibrium point K^∗ is asymptotically stable. Additionally, the trajectory K(t) for t > 0 remains in the sublevel set LeK⁰ for any K⁰ ∈ K, because e^K is monotonically decreasing over time (see (67)). Hence, the region of attraction of K^∗ is K.

Remark 5: The gradient flow of Theorem 2 can be interpreted as a continuous-time analogue of a PI algorithm (see Algorithm 1), specifically Kleinman's Algorithm (21). For a given policy K, solving the Lyapunov equation (27) with respect to P^K can be seen as the policy evaluation step, where the value matrix is determined for the given policy. Simultaneously, the gradient flow updates the policy by moving in the direction of the negative gradient, thereby reducing the Bellman error and improving the policy. This corresponds to the policy improvement step of the PI algorithm. By using the gradient flow, both steps occur concurrently in a continuous manner.

Remark 6: Since the Bellman error e^K is directly parametrized in the feedback matrix (see Remark 2), the gradient flow of Theorem 2 can be interpreted as a policy gradient method or actor-only approach. Our gradient flow is a model-based method, requiring exact knowledge of the system (A, B).

IV. SIMULATION RESULTS

In this section, we illustrate the matrix function e^K and its gradient flow, and compare our method with the gradient flow of the LQR cost (11) of [18]. To visualize the results, we consider a system with dimensions n = 2 and m = 1 given by

$$A = \begin{bmatrix} -2 & 1\ 0 & -1 \end{bmatrix}, \quad B = \begin{bmatrix} 1\ 1 \end{bmatrix}, \quad K = \begin{bmatrix} K_1 & K_2 \end{bmatrix}, \quad (68)$$

and Q = I2, R = 2. The system (68) is stabilizable and the pair (A, ^√ Q) is detectable. By using the Routh-Hurwitz stabilty criterion, the set of stabilizing feedbacks is given by

$$\mathcal{K} = { K \in \mathbb{R}^{1 \times 2} : K_2 > -K_1 - 1 }. \tag{69}$$

The Bellman error e^K and the LQR cost f^K are explicitly given by the rational functions18

$$\begin{split} e_K = & \left( K_1^6 + 4K_1^5K_2 + 12K_1^5 + 7K_1^4K_2^2 + 34K_1^4K_2 + 49K_1^4 \right. \ & + 8K_1^3K_2^3 + 40K_1^3K_2^2 + 84K_1^3K_2 + 72K_1^3 \ & + 7K_1^2K_2^4 + 36K_1^2K_2^3 + 58K_1^2K_2^2 + 32K_1^2K_2 + 29K_1^2 \right. \ & + 4K_1K_2^5 + 28K_1K_2^4 + 60K_1K_2^3 + 16K_1K_2^2 \ & - 52K_1K_2 - 8K_1 + K_2^6 + 10K_2^5 + 37K_2^4 \ & + 56K_2^3 + 17K_2^2 - 22K_2 + 5 \right) / \ & \left. \left( 2\left(K_1^2 + 2K_1K_2 + 4K_1 + K_2^2 + 4K_2 + 3\right)^2 \right), \end{split}$$

¹⁸The explicit expressions are obtained by using the vectorization of the Lyapunov equations (27) and (43).

Fig. 1. Contour plot of the Bellman error e^K

Fig. 2. Contour plot of the LQR cost f^K − fK^∗ 19

$$\begin{split} f_K = & 2 \big( K_1^3 + 2 K_1^2 K_2 + 5 K_1^2 + 2 K_1 K_2^2 \ & + 4 K_1 K_2 + 4 K_1 + 2 K_2^3 + 7 K_2^2 + 2 K_2 + 5 \big) / \ & \big( 2 \left( K_1^2 + 2 K_1 K_2 + 4 K_1 + K_2^2 + 4 K_2 + 3 \right) \big). \end{split}$$

The Bellman error e^K has a degree of 6 and 4 in the numerator and denominator, respectively, whereas the LQR cost f^K has degrees of 3 and 2. The denominators of both functions are the same, except that one is squared, while the numerator polynomials differ. Notably, the stability boundary K² = −K¹ − 1 (see (69)) corresponds to a root of the denominators of both e^K and fK. The contour plots of e^K and f^K − fK^∗ are shown in Figures 2 and 1, where the same 28 level sets are depicted. While the sublevel sets in the region K are similar, they exhibit different shapes. Both functions approach infinity by moving towards the stability boundary or as ∥K∥² → ∞, K ∈ K. Both matrix functions have a zero value at K^∗ and are non-injective on K, i.e., there are infinitely many policies for a given function value.

Figure 3 shows the residuals ^∥K(t)−K∗∥^F ∥K(0)−K∗∥^F plotted over time for gradient flows of the LQR cost and Bellman error, with K(0) = -0 0 . For the third gradient flow, which is a gradient flow of the Bellman error, the learning rate β is chosen such that the initial adaptation is similar to gradient flow of the LQR cost. All three gradient flows exhibit linear convergence to the optimal feedback K^∗ . For

¹⁹Although the integral of the LQR cost function diverges for unstable feedbacks, positive and finite costs can be obtained in the instability region K² < −K¹ − 1 by using the vectorization of the Lyapunov equations.

Fig. 3. Plot of the residuals ^∥K(t)−K∗∥^F ∥K(0)−K∗∥^F for the gradient flow of the LQR cost f^K and the Bellman error e^K

the given system (68), the gradient flows of the Bellman error e^K converge faster than that of the LQR cost fK. Due to numerical limitations, the residuals of the gradient flows stabilize at different values.

V. CONCLUSION

In this paper, we introduced a novel continuous-time Bellman error for the LQR problem, parametrized by the feedback gain. We analyzed its properties, including its effective domain, smoothness, and coerciveness, and showed the existence of a unique stationary point. Furthermore, we derived a closed-form gradient expression of the Bellman error and established its well-definedness within the stabilizing policy region. The corresponding gradient flow solves the LQR problem and converges to the optimal feedback. Additionally, this work provided novel insights on the LQR problem by bridging RL perspectives on the suboptimality of the CARE with first-order approaches in LQR theory. We validated our method through simulation and demonstrated that, in a specific example, our approach converges faster than the state-of-the-art approach.

REFERENCES

[1] R. E. Kalman, "Contributions to the theory of optimal control," Bol. Soc. Mat. Mex., no. 1, 1960.
[2] G. Hewer, "An iterative technique for the computation of the steady state gains for the discrete optimal regulator," IEEE Transactions on Automatic Control, no. 4, 1971.
[3] D. Kleinman, "On an iterative technique for Riccati equation computations," IEEE Trans. Autom. Control, no. 1, 1968.
[4] P. Lancaster and L. Rodman, Algebraic Riccati Equations. New York, USA: Oxford University Press Inc., 1995.
[5] V. Balakrishnan and L. Vandenberghe, "Semidefinite programming duality and linear time-invariant systems," IEEE Transactions on Automatic Control, no. 1, 2003.
[6] D. Yao, S. Zhang, and X. Y. Zhou, "LQ control via semidefinite programming," in Proceedings of the 38th IEEE Conference on Decision and Control, 1999.
[7] T. Chen and B. A. Francis, Optimal Sampled-Data Control Systems. Springer London, 1995.
[8] E. Feron et al., "Numerical Methods for H² Related Problems," in IEEE American Control Conference, 1992.
[9] R. S. Sutton and A. G. Barto, Reinforcement Learning: An Introduction, 2nd ed. The MIT Press, 2018.
[10] B. Kiumarsi et al., "Optimal and Autonomous Control Using Reinforcement Learning: A Survey," IEEE Transactions on Neural Networks and Learning Systems, no. 6, 2018.
[11] B. Recht, "A Tour of Reinforcement Learning: The View from Continuous Control," Annual Review of Control, Robotics, and Autonomous Systems, no. 1, 2019.
[12] D. Vrabie, K. Vamvoudakis, and F. L. Lewis, Optimal Adaptive Control and Differential Games by Reinforcement Learning Principles. London, UK: The Institution of Engineering and Technology, 2013.
[13] W. Levine and M. Athans, "On the determination of the optimal constant output feedback gains for linear multivariable systems," IEEE Trans. Autom. Control, no. 1, 1970.
[14] J. Geromel and J. Bernussou, "Optimal decentralized control of dynamic systems," Automatica, no. 5, 1982.
[15] K. Martensson, "Gradient Methods for Large-Scale and ˚ Distributed Linear Quadratic Control," Ph.D. dissertation, Lund University, Lund, 2012.
[16] M. Fazel et al., "Global Convergence of Policy Gradient Methods for the Linear Quadratic Regulator," in Proceedings of the 35th ICML, 2018.
[17] J. Bu et al., LQR through the Lens of First Order Methods: Discrete-time Case, 2019. arXiv: 1907 . 08921 [eess.SY].
[18] J. Bu, A. Mesbahi, and M. Mesbahi, Policy Gradient-based Algorithms for Continuous-time Linear Quadratic Control, 2020. arXiv: 2006.09178 [eess.SY].
[19] J. Bu, A. Mesbahi, and M. Mesbahi, "LQR via First Order Flows," in American Control Conference, 2020.
[20] D. Vrabie et al., "Adaptive optimal control for continuoustime linear systems based on policy iteration," Automatica, no. 2, 2009.
[21] D. Vrabie and F. Lewis, "Neural network approach to continuous-time direct adaptive optimal control for partially unknown nonlinear systems," Neural Networks, no. 3, 2009.
[22] R. H. Bartels and G. W. Stewart, "Solution of the Matrix Equation AX + XB = C," Communications of the ACM, no. 9, 1972.
[23] J. P. Hespanha, Linear Systems Theory, 2nd ed. Princeton, USA: Princeton University Press, 2018.
[24] J. Bu, A. Mesbahi, and M. Mesbahi, On Topological and Metrical Properties of Stabilizing Feedback Gains: the MIMO Case, 2019. arXiv: 1904.02737 [cs.SY].
[25] B. D. O. Anderson and J. B. Moore, Optimal Control: Linear Quadratic methods. Englewood Cliffs, USA: Prentice-Hall, Inc., 1990.
[26] R. J. Leake and R.-W. Liu, "Construction of Suboptimal Control Sequences," SIAM Journal on Control, no. 1, 1967.
[27] R. A. Horn and C. R. Johnson, Matrix Analysis, 2nd ed. Cambridge University Press, 2012.
[28] J. Munkres, Topology, 2nd ed. Harlow, UK: Pearson Education Limited, 2014.
[29] T. Minka, "Old and New Matrix Algebra Useful for Statistics," 2000.
[30] W. Rudin, Principles of Mathematical Analysis, 3rd ed. New York, USA: McGraw-Hill, 1964.

gielerbridgingcontinuoustimelqr2025.md 51 KB История Исходник