Time-delay systems are ubiquitous in many branches of science and engineering; see the books [1]- [3] for many references and examples. Recently, many theoretical results are developed for time-delay systems, such as input-to-state stability [4], robust H ∞ control [5], and stability analysis of systems with time-varying delay [6], [7]. Examples of timedelay systems are in transportation [8], [9], biological motor control [10], and multi-agent systems [11], [12]. It is thus not surprising that the optimal control problem of time-delay systems has been a fundamentally important, yet challenging, research topic in control theory for several decades. For instance, Eller et al. [13] and Ross et al. [14], [15] proposed solutions to the finite-horizon and infinite-horizon linear quadratic (LQ) optimal control problems of linear time-delay systems, respectively. In these papers, the certain infinitedimensional Riccati equations have to be solved. For this problem, many numerical algorithms have been developed [16]- [18]. However, an accurate model of the system is required for these algorithms, and in reality, it is difficult to derive an exact model due to the complexity of the system and inevitable system uncertainties. Therefore, developing a model-free optimal control approach for time-delay systems is a timely research topic of both theoretical importance and practical relevance. Recent progresses and successes in RL provide an opportunity to advance the state of the art in the area of adaptive optimal control of time-delay systems.
RL is an important branch of machine learning and is aimed at maximizing (or minimizing) the cumulative reward (or cost) through agent-environment interactions. Traditional RL has some fundamental limitations. For example, it often assumes that the environment is depicted by Markov decision processes or discrete-time systems with finite state-action space. Often, the stability aspect of the learned controller by RL is not guaranteed. For many systems described by differential equations, such as autonomous vehicles and quadrupedal robots, the state and action spaces are infinite and the stability of the controller generated by an RL algorithm is innegligible. Therefore, for these safety-critical engineering systems, conventional RL is not directly applicable to learning stable optimal controllers from data, which has motivated the development of ADP [19], [20]. In contrast with conventional RL, the purpose of continuous-time ADP is addressing decision-making problems for dynamical systems described by differential equations, of which both the state and action spaces are continuous. It is theoretically shown that at each iteration of ADP, a stable suboptimal controller with improved performance is obtained. Besides, the sequence of these sub-optimal controllers converges to the optimal one [19].
Recently, ADP techniques are developed for various important classes of linear/nonlinear/periodic dynamical systems and for optimal stabilization, tracking and output regulation problems [19], [21]- [25]. However, a systematic ADP approach to adaptive optimal control of continuous-time time-delay systems is lacking, due to the infinite-dimensional nature of these systems. In [26], although the model-free data-driven control for continuous-time time-delay systems is studied, discretization and/or linearization techniques are applied to transfer the time-delay system to a finite-dimensional delay-free system with augmented states, which leads to an approximate model. In [9], [27]- [32], ADP for discrete-time systems with time delays is studied. Due to the finite dimensionality of discretetime systems with time delays, these proposed ADP methods are not applicable to continuous time-delay systems. In [33], [34], ADP technique is applied for both linear and nonlinear systems with time delays, but the resulting controller cannot achieve optimality [33,Remark 9.1]. Technically, there are several obstacles in the generalization of ADP to time-delay systems. Firstly, for an infinite-dimensional system, optimality properties are hard to analyze, because the corresponding ARE are complex partial differential equations (PDEs). Secondly, stability analysis and controller design for a time-delay system are much more challenging than finite-dimensional systems.
Therefore, the model-free optimal control for a continuoustime time-delay system remains an open problem.
In this paper, in the absence of the precise knowledge of system dynamics, a novel data-driven PI approach for continuous-time linear time-delay systems is proposed based on ADP. The contributions of this paper are as follows. Firstly, inspired by Kleinman's model-based PI algorithm for delayfree linear systems [35], a new model-based PI algorithm is proposed for a class of linear time-delay systems. Given an initial admissible controller, both the stability of the updated sub-optimal controller and the convergence of the algorithm to the (unknown) optimal controller are proved theoretically. It is worth pointing out that due to the infinite dimensionality, both the value function and the control law are functional of the system's state, which in consequence increases the difficulty to design the PI algorithm. Secondly, based on the model-based PI, this paper contributes a data-driven PI approach to adaptive optimal controller design using only the data measured along the trajectories of the system.
The rest of this paper is organized as follows. Section II introduces the class of linear time-delay systems and formulates the adaptive optimal control problem to be addressed in the paper. Section III proposes a model-based PI approach to iteratively solve the LQ optimal control problem for linear time-delay systems. In Section IV, a data-driven PI approach is proposed, and the convergence property of the algorithm is analyzed. Section V illustrates the proposed data-driven PI approach by means of two practical examples. Finally, some concluding remarks are drawn in Section VI.
Notations: In this paper, R (R + ) denotes the set of (nonnegative) real numbers and N + denotes the set of positive integers. | • | denotes the Euclidean norm of a vector or Frobenius norm of a matrix. • ∞ denotes the supremum norm of a function. C 0 (X, Y ) denotes the class of continuous functions from the linear space X to the linear space Y . AC ([-τ, 0], R n ) denotes the class of absolutely continuous functions. df dθ (•) denotes the function which is the derivative of the function f . ⊕ denotes the direct sum. L i ([-τ, 0], R n ) denotes the space of measurable functions for which the ith power of the Euclidean norm is Lebesgue integrable, andL(X, Y ) denote the class of continuous bounded linear operators from X to X and from X to Y , respectively. ⊗ denotes the Kronecker product. vec(A) = a 1 , a
[a] i,j denotes the sub-vector of the vector a comprised of the entries between the ith and jth entries. A † denotes the Moore-Penrose inverse of matrix A.
Consider a linear time-delay system
where τ ∈ R + denotes the constant delay of the system and is assumed to be known,
with the domain of A given by D. Let z 0 = [x (0), x 0 (•)] ∈ D denote the initial state of the system (2). The quadratic performance index is defined as
where
Definition 1 : A control policy u c (x t ) is admissible for system (1) with respect to (3), if system (1) with u c (x t ) is globally asymptotically stable (GAS) at the origin [37,Definition 1.1], and the performance index ( 3) is finite for all z 0 ∈ D.
x(t) is exponentially stabilizable and detectable, which are defined in [36,Definition 5.2.1] and checked by [36,Theorem 5.2.12]. Remark 1. Assumption 1 is a standard prerequisite for LQ optimal control of system (1) to ensure the existence of a unique stabilizing solution [36], [38].
Given the aforementioned assumption, the problem to be studied in this paper can be formulated as follows.
dθ, and without knowing the dynamics of system (1), design a PI-based ADP algorithm to approximate the optimal controller in (4) using only the inputstate data measured along the trajectories of the system.
For a linear system without time delay, i.e. A d = 0 in (1), one can find the optimal solution by solving the ARE as discovered by Kalman [39]. Correspondingly, for system (1), the optimal solution is stated as follows.
Lemma 1 ( [14], [40]). For system (1) with Assumption 1,
is the optimal controller minimizing (3), and the corresponding minimal value functional is
where P * 0 = P * 0 , P * 1 (θ), and P * 2 (θ, ξ) = P * 2 (ξ, θ) for θ, ξ ∈ [-τ, 0] are the unique stabilizing solution to the following PDEs
By [36, Theorem 6.2.7], the time-delay system (1) in closedloop with u * is exponentially stable at the origin.
According to Lemma 1, the optimal controller is obtained by solving (6) . Since ( 6) is nonlinear with respect to P * 0 , P * 1 and P * 2 , it is difficult to solve it directly. In this section, a model-based PI algorithm is proposed to simplify the process of solving (6).
Given an admissible controller
dθ, the model-based PI algorithm for system (1) is proposed as follows.
1) Policy Evaluation: For i ∈ N + , and ξ, θ ∈ [-τ, 0], calculate P 0,i = P 0,i , P 1,i (θ), and P 2,i (θ, ξ) = P 2,i (ξ, θ) by solving the following PDEs,
where
2) Policy Improvement: Update the policy u i+1 by
The policy evaluation step calculates the value functional V i (x 0 ) for the controller u i , which is expressed as
By policy improvement, the value functional is monotonically decreasing (V i+1 (x 0 ) ≤ V i (x 0 )), and converges to the optimal value functional V * (x 0 ). Correspondingly, P 0,i , P 1,i (θ) and P 2,i (ξ, θ) converge to the optimal solutions P * 0 , P * 1 (θ) and P * 2 (ξ, θ), respectively. The convergence of the model-based PI algorithm is rigorously demonstrated in Theorem 1. Before stating Theorem 1, we first introduce Lemma 2 which is instrumental for the proof of Theorem 1. By Lemma 2, for a linear controller u L , if J(x 0 , u L ) is finite, the closed-loop system consisting of (1) and u L is globally exponentially stable.
Lemma 2. Consider system (1) with Assumption 1. If a linear controller u L (x t ) = -Kz(t) satisfies J(x 0 , u L ) < ∞ for any z 0 ∈ D, where K ∈ L(M 2 , R m ), then the closed-loop system with u L is globally exponentially stable at the origin.
Proof. The details of the proof are in [41, Lemma 2].
Theorem 1. Given an admissible control u 1 (x t ), for P 0,i , P 1,i (θ), P 2,i (ξ, θ), and u i+1 (x t ) obtained by solving (7) and (8), and for all i ∈ N + , the following properties hold 1)
Proof. Along the trajectories of (1) driven by u, where u without subscript stands for an arbitrary input, Vi (x t ) is
The detailed derivation of ( 10) is in [41, Equation ( 10)] Property 1) is proved by induction. When i = 1, the admissibility of u 1 (x t ) is given. For i > 1, assume that u i is admissible. When system (1) is driven by u i , by (10), the expression of Vi (
Following the fact that u i is admissible and integrating (11) from 0 to ∞, we have
By (10), along the trajectories of (1) driven by u i+1 ,
Integrating both sides of (13) from 0 to ∞ yields
It follows from ( 14) and Lemma 2 that u i+1 is a globally and exponentially stabilizing controller. By Definition 1, u i+1 is admissible. Via induction, u i is admissible for any i ∈ N + .
Along the trajectories of system (1) driven by u i+1 , by (10),
Since u i+1 is admissible, integrating (15) from 0 to ∞ yields
) is obtained by (14). Furthermore, since V * (x 0 ) = J(x 0 , u * ) is the minimal value functional by Lemma 1, for any i ∈ N + , V * (x 0 ) ≤ V i (x 0 ). Therefore, the proof of 2) is completed. Define P i ∈ L(M 2 ), such that for any z 0 , P i z 0 is
It is easy to check that P i is symmetric [42,Chapter 6], and non-negative [42, Definition 6.3.1], and V i (x 0 ) = z 0 , P i z 0 . Furthermore, according to statement 2), for any i ∈ N + , P * ≤ P i+1 ≤ P i . According to [42, Theorem 6.3.2], there exists P p = P p ≥ 0, such that for all z 0 ∈ M 2 , we have
Therefore, P 0,i , P 1,i (θ) and P 2,i (ξ, θ) pointwisely converge to P 0,p , P 1,p (θ), and P 2,p (ξ, θ), respectively. When P i converges, P 0,p , P 1,p (θ) and P 2,p (ξ, θ) satisfy ( 7) with i replaced by p. K 0,i and K 1,i converge to K 0,p and K 1,p . By the policy improvement step (8), K 0,p and K 1,p satisfy
Substituting ( 18) into (7) with i replaced by p, it is seen that P 0,p , P 1,p (θ) and P 2,p (ξ, θ) solve the PDEs (6). Due to the uniqueness of the solution to (6), P 0,i , P 1,i (θ) and P 2,i (ξ, θ) pointwisely converge to P * 0 , P * 1 (θ) and P * 2 (ξ, θ). Since both P 1,i (θ) and P 2,i (ξ, θ) are continuously differentiable, {P 1,i (θ) : i ∈ N + } and {P 2,i (ξ, θ) : i ∈ N + } are equicontinuous, which leads to the uniform convergence by [43,Chapter 4,Theorem 16]. Hence, 3) is proved.
Notice that although (7) is linear with respect to P 0,i , P 1,i , and P 2,i , since (7) are PDEs, obtaining the analytical solution to ( 7) is still non-trivial. Besides, the accurate knowledge of the system matrices A, A d , and B is required to implement the model-based PI, and in practice due to the complex structure of the system, it is often hard to derive such an accurate model. Therefore, in the next section, a data-driven PI algorithm is proposed to approximate the optimal solution. Remark 2. When A d = 0, (1) is degraded to the normal delay-free systems. According to (7) and (8), we can see that P 1,i (θ) = 0, P 2,i (ξ, θ) = 0, and K 1,i (θ) = 0. As a consequence, (7) and (8) are the same as the model-based PI method in [35]. Therefore, the proposed model-based PI algorithm is a generalization of celebrated Kleinman's PI to linear time-delay systems. Remark 3. In [18], the model-based PI is developed for infinite-dimensional linear systems in the Hilbert space. Although system (1) is one of the infinite-dimensional systems, the concrete expression of PI for linear time-delay systems is not given in [18], and as a consequence, the PI developed in [18] cannot be directly applied to solve the PDEs (6). In this paper, the concrete expression of PI is constructed in (7) and (8), which is one of the major contributions in this paper. Besides, it can be checked that at each iteration, P i defined in (16) satisfies the PI update equations in [18], which is another way to prove the validity of the proposed PI theoretically. Remark 4. As shown in [18], the convergence rate of PI algorithm in the Hilbert space is quadratic, and therefore, the proposed model-based PI for system (1) has the same quadratic convergence rate.
The purpose of this section is to propose a corresponding data-driven PI algorithm that does not require the accurate knowledge of system (1) to solve Problem 1. The input-state trajectories of system (1) is required for the data-driven PI. In other words, the continuous-time trajectories of x(t) and u(t) sampled from system (1) within the interval [t 1 , t L+1 ] is applied to train the control policy. From the RL perspective, u is named behaviour/exploratory policy.
Define v i (t) = u(t)u i (x t ), where u i (x t ) is the value of the control policy u i calculated along the sampled trajectory. By (10), along the trajectories of system (1) driven by the behaviour/exploratory policy u,
Plugging the expressions of u i+1 in ( 8) and V i in ( 9) into (20) yields
As seen in ( 7) and (8), K 1,i (θ) and P 1,i (θ) (P 2,i (ξ, θ)) are continuous functions defined over the set [-τ, 0] ([-τ, 0] 2 ). Next, we use the linear combinations of basis functions to approximate these continuous functions, such that only the weighting matrices of the basis functions should be determined for the function approximation. Let Φ(θ), Λ(ξ, θ), and Ψ(ξ, θ) denote the N -dimensional vectors of linearly independent basis functions. To simplify the notation, we choose the same number of basis functions for Φ, Λ and Ψ. According to the approximation theory [44], the following equations hold
), and e N K,i (θ) ∈ C 0 ([-τ, 0], R mn ) are approximation truncation errors. Therefore, by the uniform approximation theory, as N → ∞, the truncation errors converge uniformly to zero, i.e. for any η > 0, there exists
Therefore, the key idea of the data-driven PI is that W j,i (j = 0, • • • , 3) and U j,i (j = 0, 1) are directly approximated by the data collected from system (1). Define Υ N i as the composite vector of the weighting matrices, i.e.
Let ΥN i be the approximation of Υ N i , and then, the approximations of P j,i (j = 0, 1, 2) can be reconstructed by
where
Furthermore, K0,i+1 and K1,i+1 (θ), the approximations of K 0,i+1 and K 1,i+1 (θ) respectively, can be reconstructed by
where n 6 = n 5 +nm, and n 7 = n 6 +nmN . As a consequence, ûi (x t ), the approximation of u i (x t ), can be expressed as ûi
Based on the approximations in ( 22), we will transfer (21) to a linear equation with respect to ΥN i . Then, the unknown vector Υ N i will be approximated by linear regression, and consequently, P j,i (j = 0, 1, 2) and K j,i+1 (j = 1, 2) can be approximated by ( 25) and (26). In detail, let vi = uûi , ũi = ûiu i . Define the data-constructed matrices
With the help of ( 22) and ( 27), each term in ( 21) is expressed linearly with respect to the weighting matrices
where 1,i (t) 2,i (t), 3,i (t), and ψ i,k are induced by the approximation truncation errors in (23), and ρ 0 i,k and ρ 1 i,k are induced by ũi (their expressions are in [41,Equation (38)]).
With the collected input-state trajectories, define
where 28) and the definitions of M i,k , Y i,k and E i,k in ( 29), ( 21) is finally transferred to a linear equation
Combining (30)
Let Êi be the linear regression error defined as
Assumption 2. Given N ∈ N + , there exist L * ∈ N + and α > 0, such that for all L > L * and i ∈ N + ,
Remark 5. Assumption 2 is reminiscent of the persistent excitation (PE) condition [45], [46]. It is needed to guarantee the uniqueness of the least-square solution to (31), and prove the convergence of the proposed data-driven PI algorithm. As in the literature of ADP-based data-driven control [19], [20], one can fulfill it by means of added exploration noise, such as sinusoidal signals and random noise.
Under Assumption 2, the method of least squares is applied to minimize Ê i Êi , i.e. Ê i Êi is minimized by
With the result of ΥN i in (34), Pj,i (j = 0 • • • 2) and Kj,i (j = 0, 1) can be reconstructed by ( 25) and ( 26) respectively.
The proposed algorithm is shown in Algorithm 1. From (29), M i and Y i are constructed by the input-state trajectory data of system (1). Hence, the system matrices are not involved in the computation of ΥN i . Furthermore, since the behaviour/exploratory policy u is different from the updated policy u i , Algorithm 1 is called off-policy. Remark 6. Due to the property that P 2,i (ξ, θ) = P 2,i (θ, ξ), the diagonal elements of P 2,i satisfy diag(P 2,i (ξ, θ)) = diag(P 2,i (θ, ξ)). Hence, the vector of basis functions Ψ should satisfy Ψ(ξ, θ) = Ψ(θ, ξ) to approximate such functions.
Remark 7. In practice, the integrals in (27) are calculated by Riemann sum, like midpoint, trapezoid, and Simpson's rules. Lemma 3. Under Assumption 2, and given an admissible controller u 1 (x t ) = -K 0,1 x(t) -0 -τ K 1,1 (θ)x t (θ)dθ, for each i ∈ N + and any η > 0, there exists a positive integer
Proof. Define the approximation error as ΥN i = Υ N i -ΥN i . Subtracting (32) from (31) yields
Since Ê i Êi is minimized by the method of least squares,
Furthermore, combining ( 36) and ( 37), we have
Therefore, via Assumption 2, the following inequality holds
Then, the lemma is proved by induction. When i = 1, û1 = u 1 , so ũ1 = 0, and ρ 0 1,k = ρ 1 1,k = ρ 2 1,k = 0. Furthermore, since 1,i (t) 2,i (t), 3,i (t), and ψ i,k are induced by the approximation truncation errors in (22), they converge to zero as N → ∞. Therefore, by the expression of E i,k in (29), for any 1 ≤ k ≤ L, E i,k converges to zero as N → ∞. Consequently, by (39), for i = 1, lim
As a result, the estimation of the weighting matrices in (22) converge to the true values as N → ∞. By the boundedness of the functions Φ(θ), Ψ(ξ, θ), Λ(ξ, θ) on the compact interval θ, ξ ∈ [-τ, 0] and the uniform convergence of the approximation truncation errors in (23), (35) holds for i = 1. Suppose (35) holds for some i -1 > 1. Then, from the second line of (35), it is seen that ũi converges to zero as N → ∞. Since ρ 0 i,k , ρ 1 i,k , and ρ 2 i,k are induced by ũi , they converge to zero as N → ∞. Furthermore, since 1,i (t) 2,i (t), 3,i (t), and ψ i,k are induced by the approximation truncation errors in (22), they converge to zero as N → ∞. Consequently, by the expression of E i,k in (29), for any 1 ≤ k ≤ L, E i,k converges to zero as N → ∞. Therefore, by (39), (40) holds for i. Following the logic of the content below (40), we obtain that (35) holds for i. The proof is completed by induction.
Theorem 2. Given an admissible controller u 1 , for any η > 0, there exist integers i * > 0 and
Proof. The theorem is proven by Theorem 1, Lemma 3, and triangle inequality. See [41,Theorem 2] for details.
By Theorem 2, we see that Kj,i+1 (j = 0, 1) obtained by Algorithm 1 converges to K * j (j = 0, 1) as the iteration step of the algorithm and the number of basis functions tend to infinity. Hence, the proposed data-driven PI solves Problem 1.
Algorithm 1 Data-driven Policy Iteration Calculate ûi (t) = ûi (x t ) along the trajectory of x.
Construct M i and Y i by (29). Collect more data and insert it into M i and Y i . Get ΥN i by solving (34).
11:
Get K0,i+1 and K1,i+1 by (26).
12:
Use ûi (x t ) as the control input.
The proposed data-driven PI algorithm is demonstrated by two practical examples, with regards to regenerative chatter in metal cutting (RCMC) and connected and autonomous vehicles (CAVs) in mixed traffic consisting of autonomous vehicles (AVs) and human-driven vehicles (HDVs).
Consider the example of metal cutting [37, Example 1.1], [47]. The thrust force is proportional to the instantaneous chip thickness ([x(t)] 1 -[x(tτ )] 1 ), leading to the timedelay effect. The model is described by (1) with
and τ = 1.3s. The initial admissible controller is û1 (x t ) = -K 0,1 x(t), with K 0,1 = [1.74, 3.92]. The exploration noise is e(t) = 20 50 i=1 sin ω i t, where ω i is randomly sampled from an independent uniform distribution over [-10, 10]. Q = diag([100, 100]) and R = 1. δ = 10 -3 . For the basis functions, Φ
As shown in Fig. 1, the weights of the basis functions Υ converge after the eighth iteration. In order to inspect the evolution of the performance index, we compare the controllers updated at each iteration for the same initial state x 0 . In Fig. 1, it is seen that the performance index decreases. The responses of the state with the initial controller and the learned ADP controller are compared in Fig. 2. The performance indices are J(x 0 , û1 ) = 5.89 × 10 4 and J(x 0 , û8 ) = 3.03 × 10 4 .
Semi-discretization [48] is applied to discretize (1) into a delay-free system with sampling period Δt = 0.1s. Then, Algorithm 1 is compared with the model-based discrete-time linear quadratic regulator (DLQR) and the discrete-time ADP algorithm in [9] (with the same length of trajectory data). For the same initial state, the performance indices are shown in Table I. The performance index is minimal under Algorithm 1, showing that discretization sacrifices the system performance. Ideally, the performance of the discrete-time ADP is similar as the model-based DLQR. The large deviation between them is induced by the fact that the PE condition for the discretetime ADP is not satisfied. This further illustrates that by semi-discretization, the dramatically increased dimension of the augmented state (26-dimensional) makes the requirements on the sampled data more demanding.
The robustness of Algorithm 1 to measurement noise is evaluated. The measurement of x(t) is disturbed by an independent Gaussian noise ϕ(t) ∼ N (0, 0.2). The result is shown in Fig. 3. Using the noisy data, for the same initial state x 0 , the performance index converges to J = 3.35 × 10 4 . Comparing Fig. 3 with the second figure in Fig. 1, we see that Algorithm 1 can still find a near-optimal solution in the presence of noise.