Approximate Dynamic Programming For Linear Systems with State and Input Constraints
Ankush Chakrabarty Rien Quirynen Claus Danielson Weinan Gao

	
	Selective Inference via Marginal Screening for High Dimensional Classification
    Yuta Umezu         Ichiro Takeuchi 

    Received: date / Accepted: date
==============================================================================


	Enforcing state and input constraints during reinforcement learning (RL) in continuous state spaces is an open but crucial problem which remains a roadblock to using RL in safety-critical applications. This paper leverages invariant sets to update control policies within an approximate dynamic programming (ADP) framework that guarantees constraint satisfaction for all time and converges to the optimal policy (in a linear quadratic regulator sense) asymptotically. An algorithm for implementing the proposed constrained ADP approach in a data-driven manner is provided. The potential of this formalism is demonstrated via numerical examples.
	
	
	

§ INTRODUCTION

	Combining optimal control theory and reinforcement learning (RL) has yielded many excellent algorithms for generating control policies that imbue the closed-loop system with a desired level of performance in spite of unmodeled dynamics or modeling uncertainties <cit.>. Specifically, adp (sometimes also referred to as adaptive dynamic programming), a modern embodiment of RL <cit.> applied to continuous state and action spaces has gained traction for its ability to provide tractable solutions (in spite of the curse of dimensionality) to optimal control problems via function approximation and iterative updates of control policies and value functions <cit.>.
	
	There are two main classes of adp algorithms: policy iteration and value iteration <cit.>. A policy iteration algorithm for discrete-time linear systems was formulated in <cit.> that leverages Q-functions proposed in <cit.>, enabling control policy design without a complete system description. This methodology has been extended to continuous-time systems <cit.>, ℋ_2 and ℋ_∞ formulations <cit.>, tracking <cit.>, output regulation <cit.>, and  game-theoretic settings <cit.>. To reiterate, a particularly beneficial feature of this class of iterative methods is that control policies generated by policy  iteration converge to the optimal control policy with data obtained by  exciting the system dynamics, in spite of incomplete model knowledge <cit.>. 
	While optimality is important for certifying performance in a control system, often times the more critical concern is safety. A key aspect of safe control design is the ability of the system to respect both state and input constraints. To the best of our knowledge, this critical problem remains an open challenge in the context of ADP (and RL at large) in continuous state and action spaces.
	
	In this paper, we modify the classical policy iteration algorithm to incorporate safety through constraint satisfaction. The key idea is to compute control policies and associated constraint admissible invariant sets that ensure the system states and control inputs never violate design constraints. In the spirit of ADP, these policies and invariant sets are computed iteratively, and the sequence of policies are guaranteed to converge asymptotically to the optimal constraint-satisfying policy, provided that the system is sufficiently excited. The use of invariant sets to incorporate safety in learning/adaptive control algorithms via constraint handling has been done in model-based control design, such as model predictive control (MPC) <cit.>, but its application to data-driven or model-free RL methods is relatively unexplored. A recent paper <cit.> is a noteworthy exception, although  our method is distinct from this work in that we do not compute a model of the system using the data obtained during operation; that is, our method is a direct data-driven approach, as defined in <cit.>.
		
	The main contributions of this paper are: (i) we extend classical policy iteration in continuous state-action spaces to enforce state and input constraints; (ii) we provide a data-driven variant of this constrained policy iteration algorithm with unknown state matrix; and, (iii) we provide new sufficient conditions for safety (via constraint satisfaction), stability, and convergence of the policies generated by our proposed algorithm to the optimal constrained control policy.
	
	The rest of the paper is organized as follows. In Section <ref> and <ref>, we present our notation and describe the problem statement in formal terms. We also discuss the standard (unconstrained) policy iteration algorithm. Our proposed alterations to the standard policy iteration algorithm for enabling state and input constraint satisfaction is discussed in Section <ref>. Theoretical performance certificates such as safety, stability, and algorithm convergence are provided in Section <ref>. Numerical examples including a 2D illustrative example and a 5D example are provided in Section <ref> to illustrate the potential of this method. Conclusions are drawn in Section <ref>.
	

§ NOTATION

	We denote by ℝ the set of real numbers, ℝ_+ as the set of positive reals, and ℕ as the set of natural numbers.
	For every v∈ℝ^n, we denote v=√(v^⊤ v), where v^⊤ is the transpose of v. The sup-norm is defined as v_∞≜sup_t∈ℝv(t). We denote by σ(P) and σ(P) as the smallest and largest singular value of a square, symmetric matrix P, respectively. The symbol ≻(≺) indicates positive (negative) definiteness and A≻ B implies A-B≻ 0 for A,B of appropriate dimensions. Similarly, ≽ (≼) implies positive (negative) semi-definiteness. The operator norm is denoted P and is defined as the maximum singular value of P, vec(P) denotes the column-wise vectorization of P, and ⊗ denotes the Kronecker product. We parameterize an ellipsoid = {x: x^⊤ P_0 x ≤ρ} using a scalar ρ>0 and a matrix P_0≻ 0.
	
	

§ MOTIVATION

	In this section, we begin by describing a general approximate dynamic programming formulation for solving the unconstrained discrete-time LQR problem.
	

 §.§ Problem Statement

	We consider discrete-time linear systems of the form
	
    x_t+1 = Ax_t + Bu_t,

	where t∈ℝ is the time index, x∈𝕏⊂ℝ^n is the state of the system, u∈𝕌⊂ℝ^m is the control input, and x_t_0 is a known initial state of the system. We assume the admissible state and input constraints sets 𝕏 and 𝕌 are polytopic, and therefore, can be represented as
	
    𝒳' = {[ x; u ]∈ℝ^n+m: c_i^⊤ x + d_i^⊤ u≤ 1},

	for i=1,…, r, where r is the total number of state and input constraints and c_i∈ℝ^n and d_i∈ℝ^m.
	The sets 𝕏⊂ℝ^n and 𝕌⊂ℝ^m are known, compact, convex, and contain the origin in their interiors. Note that with any fixed control policy K, the constraint set described in (<ref>) is equivalent to the set
	
    𝕏' = {x∈ℝ^n: (c_i^⊤ + d_i^⊤ K) x ≤ 1},

	for i=1,…, r.
	
	The inequalities (<ref>) define a polytopic admissible state and input constraint set. Note that c_i=0 implies that the ith constraint is an input constraint, and d_i=0 implies that it is a state constraint.
	
	The system matrix A and input matrix B have appropriate dimensions. 
	We make the following assumption on our knowledge of the system; these are standard assumptions in policy iteration methods.
	
		The matrix A is unknown, the matrix B is known, and the pair (A, B) is stabilizable. Furthermore, there exists a known control gain K_0 such that u=K_0x is a stabilizing control policy for the system (<ref>).
	
	While the knowledge of the input matrix B is not needed in approaches like Q-learning <cit.>, it is fairly standard for policy improvement in policy iteration methods, even with function approximators <cit.>. From a practical perspective, it is not uncommon for a designer to have knowledge of input channels and channel gains that represent the elements of the B matrix.
	
	Our objective is to design an optimal control policy  such that the state-feedback controller u= x stabilizes the partially known system (<ref>) while minimizing a cost functional
	
    V := ∑_t=0^∞ x_t^⊤ Q x_t + u_t^⊤ R u_t

	where Q≽ 0 and R≻ 0 are user-defined symmetric matrices, with the pair (A,Q^1/2) being observable. The main contribution of this paper is to derive controller gains that stabilize the system (<ref>) while strictly enforcing state and input constraints.
	
	

 §.§ Overview of optimal control for  discrete-time LQR

	Let the value function be defined as
	
    V_t(x_t, u_t) := ∑_k=t^∞ x_k^⊤ Q x_k + u_k^⊤ R u_k.

	Clearly, V_t satisfies the recurrence relation
	
    V_t(x_t, u_t) = x_t^⊤ Q x_t + u_t^⊤ R u_t + V_t+1(x_t+1, u_t+1).

	We know from optimal control theory that the optimization problem
	
    V_∞(x_t) := min_u V_t(x_t, u_t)

	is solved in order to obtain the optimal control action 
	
    u_∞ := min_u V_t(x_t, u_t)

	for each time instant t≥ t_0.
	For discrete-time linear time-invariant systems of the form (<ref>), we know that the value function V_t is quadratic in the state <cit.>. Therefore, solving (<ref>) is equivalent to finding a symmetric matrix ≻ 0 that satisfies the equation
	
    A^⊤ A -  + Q - A^⊤ B (R+ B^⊤ B)^-1 B^⊤ A= 0.

	Upon solving for , the optimal unconstrained discrete-time LQR policy generated by solving (<ref>) is given by
	
    = -(R+ B^⊤ B)^-1B^⊤ A.


	Since by assumption, the model A is unknown, one cannot directly compute   from (<ref>) or  from (<ref>). Instead, we resort to adp, an iterative method for `learning' the optimal control policy (<ref>) by using on-line data without knowing a full model of the system (<ref>). A popular embodiment of adp is policy iteration, wherein an initial stabilizing control policy K_0 is iteratively improved using operational data, that is, without full model information. The sequence of control policies converges asymptotically to the optimal control policy  defined in (<ref>). The key steps of policy iteration without constraints are described next.
	
	

 §.§ Unconstrained policy iteration

	Let K_t be the t-th policy iterate, where t ∈ℕ. Policy iteration has two key steps: policy evaluation and policy improvement. We begin by describing the steps in model-based policy iteration and subsequently demonstrate how to perform the same steps in a data-driven manner.
	

  §.§.§ Model-based policy evaluation

	In the policy evaluation step, the value function parameter P_t+1≻ 0 is estimated with the control gain K_t using the relation
	
    (A+BK_t)^⊤ P_t+1(A+BK_t) - P_t+1 + Q + K_t^⊤ R K_t = 0.

	Note that (<ref>) can be derived from (<ref>) when V_t = x_t^⊤ P_t x_t and replacing u_t=K_t x_t and x_t+1=(A+BK_t)x_t.
	

  §.§.§ Model-based policy improvement

	Upon updating the value function via (<ref>), one needs to update the corresponding control policy. This is done by computing the new controller gain via
	
    K_t+1 = -(R+B^⊤ P_t+1 B)^-1 B^⊤ P_t+1A.

	This equation is reminiscent of the optimal control policy equation (<ref>); in fact, the unique stationary point of the system of equations (<ref>)–(<ref>) is at P_t = P_∞ and K_t = K_∞ as demonstrated in <cit.>.
	
	This model-based implementation can be performed in a data-driven manner, described next.
	

  §.§.§ Data-driven policy evaluation

	We assume that policy iteration is performed a discrete-time instances t_i where 
	
    𝒯={t_i}_i=1^∞

	denotes the set of all policy iteration times. The minimum number of data-points obtain between policy iterations [t_i,t_i+1] is given by 
	
    N = inf_i∈ℕ {t_i+1 - t_i| t_i, t_i+1∈𝒯},

	that is, N denotes the minimum number of data points contained within any learning cycle. In a model based implementation, 𝒯=ℕ.
	
	At each learning time instant t_i∈𝒯, one can rewrite (<ref>) as
	
    x_t^⊤ P^+ x_t = x_t^⊤ Q x_t + u^⊤_t R u_t +x_t+1^⊤ P^+ x_t+1,

	for every t∈{t_i+1, t_i+2,…, t_i+1}, with P^+ representing the updated value function matrix.
	Assuming that the state and input data is available to us, and that Q and R are known, we can rewrite (<ref>) as 
	
    Δ_xxvec(P^+) = [ x_t_i+1^⊤ Q x_t_i +1 + u_t_i+1^⊤ R u_t_i+1; x_t_i+2^⊤ Q x_t_i +2 + u_t_i+2^⊤ R u_t_i+2;                                          ⋮;                                           ;  x_t_i+1^⊤ Q x_t_i+1 + u_t_i+1^⊤ R u_t_i+1 ],

	where
	
    Δ_xx = [         x_t_i⊗ x_t_i - x_t_i+1⊗ x_t_i+1;                                       ⋮; x_t_i+1⊗ x_t_i+1 - x_t_i+1+1⊗ x_t_i+1+1 ].


	Under well-known persistence of excitation conditions <cit.>, one can solve (<ref>) as a (regularized) least squares problem subject to the constraint that P^+≻ 0 to obtain P^+ without knowing A or B. For the time instants t∈𝒯 when the learning occurs, the new value function matrix P_t+1 is set to P^+ obtained by solving (<ref>). For other time instants between learning time instants, that is t∉𝒯, the value function matrix obtained in the previous learning cycle is utilized, that is, P_t+1 := P_t.
	

  §.§.§ Data-driven policy improvement

	Since the control policy is restricted to be linear in this paper, finding an optimal policy is tantamount to finding the minimizer K_t+1 of the optimization problem
	
    min_K∑_t=t_i+1^t_i+1(x_t^⊤ K^⊤ R K x_t + x_t^⊤ Q x_t + x_t^⊤ (A+BK)^⊤ P_t+1 (A+BK)x_t),

	where t_i, t_i+1∈𝒯.
	This is a quadratic optimization problem in K because {x_t}, Q, R, and P_t+1 are all known quantities in the window {t_i+1, t_i+2,…, t_i+1}. 
	Note that K_t+1 can be updated recursively within each learning window t_i≤ t ≤ t_i+1 using P_t+1 for these time instants. Since (<ref>) is a quadratic problem, using Newton-type iterative solvers are expected to yield quick convergence; in this case, in one step.
	

§ CONSTRAINED ADP

	In this section, we elucidate upon how to use invariant sets to generate new control policies that are both stabilizing and constraint satisfying. We also propose an algorithm for implementing a constrained ADP in a data-driven manner.
	
	We begin with the following definition.
		
		
		A non-empty set ℰ within the admissible state space 𝕏 is a cais for the closed-loop system (<ref>) under a control law u = Kx if, for every initial condition x_t_0∈ℰ, all subsequent states x_t∈ℰ and inputs Kx_t∈𝕌 for all t≥ t_0.
	

According to Assumption 1, the ADP iteration is initialized with a stabilizing linear controller K_0.  This stabilizing controller renders a subset of the state-space invariant while satisfying state and input constraints. In particular, there exists an ellipsoidal region 

    = {x: x^⊤ P_0 x ≤ρ},
 
such that ⊂𝕏 and K_0 ⊂𝕌. We assume that the value function matrix P_0 defining the initial CAIS ellipsoid  is known. This is encapsulated formally herein.

[Constrained ADP]

There exists a symmetric positive definite matrix P_0 such that ⊂𝕏 is a CAIS for the closed-loop system (1) under the initial control policy u=K_0 x, and K_0x ∈𝕌 for all x ∈.


	

 §.§ Model-based constrained policy iteration

	

  §.§.§ Model-based constrained policy evaluation

	Let 
	
    := (A+BK_t)^⊤ P(A+BK_t) - P + Q + K_t^⊤ R K_t.

	In order to implement constrained model-based policy evaluation (that is, obtain P_t+1 from K_t and P_t), we need to solve the following semi-definite programming problem:
	
	
	
    P_t+1,ρ_t+1 = _P, ρ
    subject to:   
        (A+BK_t)^⊤ P(A+BK_t) - λ P ≼ 0
        x_t^⊤ P x_t ≤ρ
        (c_k^⊤+d_k^⊤ K_t)^⊤ρ (c_k^⊤+d_k^⊤ K_t)≼ P
       α_1𝐈≼ P ≼α_2 𝐈
       ρ> 0

	
	for some α_1, α_2 > 0 and 
	
    λ<(α_1/α_2)^2/N.

	Here, k∈{1,…, r}. Note that ensuring this problem is convex involves fixing the scalars α_1 and α_2, and pre-computing λ using (<ref>).
		
	The rationale behind (<ref>) can be explained as follows. Since (<ref>) ensures that P≻ 0, this constraint, along with the objective (<ref>), is equivalent to (<ref>), which is identical to the unconstrained policy evaluation step. Therefore, constraint satisfaction is made possible by equipping the constraints (<ref>)–(<ref>) and (<ref>).
	
	The inequality (<ref>) ensures that the value function is contractive, and therefore, non-increasing for every t≥ t_0. To see this, we multiply (<ref>) by x^⊤ and x from the left and right, respectively, which yields
	
    x_t+1^⊤ P x_t+1 - x_t^⊤ P x_t  ≤ -(1-λ)   x_t^⊤ P x_t  < 0,

	for any t, since 0<λ<1. This is a key ingredient to ensure that the updated control policies will provide stability certificates for the closed-loop system. The two inequalities (<ref>) and (<ref>) enforce that the state and input constraints with the current policy are satisfied in spite of the value function update, given the current state x_t. The condition (<ref>) ensures that the value function matrix P is positive definite, and the positive scalar ρ allows the selection of sub- and super-level sets of the Lyapunov function. More details about how these conditions relate to theoretical properties of the proposed constrained ADP algorithm are provided in Section <ref>.
	

  §.§.§ Model-based constrained policy improvement

	Unlike unconstrained policy iteration, adding state and input constraints could result in nonlinear optimal control policies. In this paper, we restrict ourselves to design linear control policies of the form u=Kx, and hence, our optimal policy improvement step is analogous to the unconstrained case (<ref>), that is,
	
    K^⋆_t+1 = -(R+B^⊤ P_t+1 B)^-1 B^⊤ P_t+1A.

	
	In spite of parameterizing via linear control policies, our controller is actually nonlinear since K_t depends on P_t which depends on the states through (<ref>).
	
	
	We adopt a backtracking strategy in order to update the current constrained policy K_t to a new constrained policy K_t+1 that is as close as possible to the unconstrained policy K^⋆_t+1 in (<ref>) that enforces state and input constraints. A simplified version of this backtracking strategy is outlined in Algorithm <ref>. 
	
	
	
	Note that the conditional statement in step 2 can be implemented efficiently based on a Cholesky factorization to check whether this particular symmetric matrix is positive definite.
	
	
	A particular benefit of our proposed method is that it enables both expansion, contraction, and rotation of the constraint admissible invariant sets. This is important in reference tracking for instance where a more aggressive controller is required when the state is near the boundary of the state constraints. This could also be useful for applying this approach to nonlinear systems where (A,B) is a local linear approximation of the globally nonlinear dynamics. Our approach allows the ellipsoidal invariant sets to adapt its size and shape based on the local vector field. For example, suppose ℰ_ denote the CAIS that is associated with the constrained optimal control policy  and optimal value function defined by . Also suppose that we have an initial admissible policy K_0 whose associated CAIS   is contained within ℰ_. Then our proposed method will generate a sequence of ℰ_P_t such that these invariant sets will expand, contract, and rotate as necessary until the sequence of invariant sets {ℰ_P_t} converges to the optimal ℰ_.

	

 §.§ Data-driven constrained policy iteration

	
	In order to obtain a data-driven implementation of the constrained ADP method, one needs to gather a sequence of state-input data points {,,} and control policies {} which will be used to update the value function matrix and control policies at the learning time instants defined by 𝒯 in (<ref>). Given the discrete-time system dynamics in (<ref>), the relation between these data points is given by
	
    = A  + B  = A  + B ( +),

	where  represents a known exploration noise signal that ensures the system (<ref>) is persistently excited; see <cit.>. To arrive at a more compact notation, let us define 
    :=  - B   and   :=
 such that 
	
    = A  + B  = (A+B).
 
	
	

  §.§.§ Data-driven constrained policy evaluation

	
	Consider the i-th learning cycle, occuring at the time instant t_i∈𝒯. Let
	
    := ^⊤ P  - ^⊤ P  +  Q  + ^⊤ R .

	The data-driven analogue of the constrained policy evaluation step discussed in the previous section is given by the following semi-definite program (SDP) with α_1 and α_2 fixed:
	 
		
    , ρ_t+1 :=
    		min_ρ, P  1/2 ∑_t=0^-1()^2 - λ_ρρ
    subject to:   
       ^⊤ P  - λ^⊤ P ≤ 0 
        x_t_i+1^⊤ P x_t_i+1≤ρ
        
    		(c_k^⊤+d_k^⊤)^⊤ρ (c_k^⊤+d_k^⊤)≼ P
    		
       α_1𝐈≼ P ≼α_2 𝐈
       ρ > 0,

	
	for t ∈{t_i+1,t_i+2,…,t_i+1} and k ∈{1,…,r}.
	Note that the final four inequalities in (<ref>) are exactly the set of inequalities presented in (<ref>) with the model information replaced by state and input data. Also, replacing  in (<ref>) with (A+B) using  equality (<ref>) shows that it is equivalent to the inequality (<ref>).
	
	
	
	

  §.§.§ Data-driven constrained policy improvement

	Once a value function is found whose sub-level set is constraint admissible, the corresponding policy K_t+1 is to be computed. If A and B are known, this step would be easy: indeed, one could utilize Eq. (<ref>) to this end.
	
	However, since only B is known (by assumption), we resort to a data-driven iterative update methodology for generating the new policy. 
	
	Given the current policy , we gather another batch of measurements {,,,}_t=t_i+1,…,t_i+1 where a new policy K̅_t is the optimizer of the least squares problem
 	
    Kmin     1/2∑_t=t_i+1^t_i+1^⊤(K^⊤ R K + (A+BK)^⊤ (A+BK)).

	
	The problem (<ref>) can be solved in a data-driven manner efficiently using a real-time recursive least squares (RLS) implementation <cit.> 
	 
		
		H_t+1    = H_t + ^⊤⊗(R+B^⊤B),  

		g_t+1    = ⊗(R + B^⊤),  

		vec(K̅_t+1)    = vec(K̅_t) - β_t  H_t+1^-1   g_t+1, 
		
	
	for t=,…,-1. Note that (<ref>) is solved without knowledge of A using the updates. Also,  the starting value for the Hessian matrix is chosen as the identity matrix ρ and ρ > 0 to ensure non-singularity. The step size β_t is typically equal to one, even though a smaller step β_t≤ 1 can be chosen, e.g., based on the backtracking procedure in Algorithm <ref> in order to impose the affine state and input constraints in (<ref>) for each updated control policy K̅_t+1. The Hessian matrix in (<ref>) can be reset to H = q ≻ 0 whenever a new value function is obtained from solving the SDP in (<ref>). Note that (<ref>) corresponds to a rank-m matrix update, where m denotes the number of control inputs. Therefore, its matrix inverse H_t+1^-1 can be updated efficiently using the Sherman-Morrison formula, for example, in the form of m rank-one updates.
	











	
	

  §.§.§ Algorithm Implementation: Pseudocode

	
	Algorithm <ref> provides a detailed description of our proposed approach for data-driven constrained adaptive dynamic programming for linear systems. The general procedure corresponds to the sequence of high-level steps:
	
		
		
  * We require an initial stabilizing policy K̅_0 and a corresponding constraint admissible invariant set (CAIS) ; see Assumptions 1 and <ref>.
		
  * Obtain a sequence of at least  data points {,,,} while the system is persistently excited and compute a new ellipsoidal set defined by the matrix  by solving the least squares SDP in (<ref>).
		
  * At each time step, perform the policy improvement step to compute K̅_t+1 based on the real-time recursive least squares method as described in (<ref>), in combination with the backtracking procedure of Algorithm <ref> to enforce state and input constraints.
		
  * If the policy improvement has converged based on the condition g_t≤ϵ, return to step (ii).
		
	
	
	
	
	
	

 §.§ Remarks on computational complexity

	
	Semidefinite programs (SDP) of the form (<ref>) are convex optimization problems that can be solved in polynomial time, for example, using interior point methods (IPMs). However, in general, standard implementations of IPMs for solving SDPs have a computational complexity Ø(n^6) when solving for n × n matrix variables and a memory complexity of Ø(n^4) <cit.>. Instead, the per iteration complexity and memory requirements for first order optimization algorithms such as, e.g., the alternating direction method of multipliers (ADMM) can be much smaller, even though they typically require more iterations in practice <cit.>. Note that, instead, the policy improvement steps are computationally cheap because both the low-rank update techniques for the Hessian matrix (<ref>) and the matrix-vector multiplication in (<ref>) can be performed with a complexity Ø(n^2m^2) that scales quadratically with the dimensions of the policy matrix K.
		
	The policy evaluation step, based on the SDP solution in (<ref>), could be computed also using a recursive least squares type implementation or in a receding horizon or sliding window manner. However, given the computational complexity of treating linear matrix inequalities in the SDP formulation, a batch-type approach as in Algorithm <ref> would typically be preferred for real-time feasible control applications under strict timing requirements. Additionally, it is important to note that the SDP solution in Algorithm <ref> is not necessarily required to be real-time feasible, unlike the recursive least squares based policy improvement step in (<ref>), which is computationally cheap. 
	
	
	

§ CONSTRAINT SATISFACTION,  STABILITY, AND ALGORITHM CONVERGENCE

	We present theoretical guarantees for our proposed constrained policy iteration. For the data-driven case, we adhere to the standard assumption that the system is persistently excited. The following theorem demonstrates constraint enforcement and stability guarantees of the closed-loop system.
	
	
		Suppose Assumptions <ref> and <ref> hold. Then the system (<ref>) in closed-loop with the time-varying controller u_t = K_t x_t has the following properties:
		
			
  * The constraints x_t ∈𝕏 and u_t ∈𝕌 are satisfied for all t ∈ℕ.
			
  * The closed-loop system is asymptotically stable.
		
	
	
	
	
	(i) Consider the ellipsoid ℰ_P_t+1^ρ.
	The inequality (<ref>) yields
	
    1/ρ P_t+1≽ (c_i^⊤ + d_i^⊤ K)^⊤ (c_i^⊤ + d_i^⊤ K).

	For x ∈ℰ_P_t+1^ρ we have x^⊤ P x ≤ρ. Thus, 
    x^⊤ (c_i^⊤ + d_i^⊤ K)^⊤ (c_i^⊤ + d_i^⊤ K) x ≤ 1
 which implies x ∈𝕏' since (c_i^⊤ + d_i^⊤ K) x ≤ 1. Thus, ℰ_P_t+1^ρ⊆𝕏', which implies that state and input constraints are satisfied for all states inside the ellipsoid ℰ_P_t+1^ρ. 
	
	Note that the state x_t is contained in this ellipsoid while the t-th controller u_t = K_t x_t is active. This can be inferred from (<ref>) and (<ref>), which (respectively) imply that the ellipsoid ℰ_P_t+1^ρ is positive-invariant and that the initial state is contained in the ellipsoid when the controller is first engaged.
			
	(ii) Since the closed-loop system is a switched system, we will use the concept of multiple Lyapunov functions to prove stability. Consider the set of Lyapunov functions 
    V_P(x) = x^⊤ P x
 for all P that satisfies (<ref>). We will show that the t-th controller u_t = K_t x_t decreases all of these Lyapunov functions[The Lyapunov functions do not necessarily decrease monotonically for all t as long as the function values are decreasing at each learning instant t_i.] over the time period t ∈{t_i+1, t_i+2,…, t_i+1} for the i-th learning cycle in which it was engaged (recall t_i, t_i+1∈𝒯 defined in (<ref>)). Note
	
    V_P(x_t_i+1)    ≤α_2/α_1 V_P_t+1(x_t_i+1) 
       ≤α_2/α_1λ^N V_P_t+1(x_t_i) 
       < α_1/α_2 V_P_t+1(x_t_i) 
       ≤ V_P(x_t_i),

	where the first and last inequalities are a consequence of (<ref>). The second inequality is a consequence of (<ref>) and (<ref>), along with Algorithm 1 which uses a convex combination of K_t and K_t+1^⋆, both of which are guaranteed to contract the value function by λ. The third inequality is a consequence of the condition (<ref>). Since V_P(x_t_i+1)<V_P(x_t_i), each Lyapunov function V_P(x) converges to zero. As (<ref>) ensures that these Lyapunov functions are positive-definite, we get x_t→ 0 as t→∞, which concludes the proof.
	
	
	Previous stability results for approximate dynamic programming rely on the tacit assumption that the learning converges after a finite number of batch iterations (typically one). In other words, the adaptive controller only works because it stops adapting. In contrast, for constraint satisfaction, the controller may need to continually adapt since the set of active constraints will change as the state evolves. This necessitates the development of a more involved set of conditions to ensure that feedback control loop and the learning loop do not destabilize each other. 
	
	
	
	Suppose Assumptions <ref> and <ref> hold. Let α_1 ≤σ(), α_2 ≥σ(), and 
    λ≥σ( I - ^-1/2 (Q+K_∞^⊤ R K_∞) ^-1/2).
 Under the iteration (<ref>) and (<ref>), the value P_t and policy K_t converge to the LQR cost-to-go  and controller gain . That is,
	
    lim_t→∞ P_t =   and  lim_t→∞ K_t = K_∞.

		
	
	
	If feasible, the LQR cost-to-go  with controller gain  will be the optimal solution of (<ref>). By the assumptions on α_1,α_2, the LQR value  satisfies (<ref>). Note that
	
    σ( I - ^-1/2 (Q+^⊤ R ) ^-1/2) 
           = sup_z  z^⊤ (I - ^-1/2 (Q+^⊤ R ) ^-1/2) z /z^⊤ z
           = sup_x  x^⊤ ( - (Q+^⊤ R )) x /x^⊤ x
	
	where z = ^-1/2 x. Thus, the LQR satisfies
	
    (A+B)^⊤ P_∞ (A+B) 
            =  - (Q + ^⊤ R ) 
           ≤σ( I - ^-1/2 (Q+^⊤ R ) ^-1/2) 
           ≤λ.

	Thus, the LQR value  and policy  satisfy (<ref>).
		
	Finally, we show that (<ref>) and (<ref>) are satisfied at some finite time T≥ t_0. Since, the state 𝕏 and input 𝕌 constraints contain the origin in their interiors, there exists ρ > 0 such that the ellipsoidal region ℰ_P_∞^ρ of the LQR cost-to-go x^⊤ x satisfies (<ref>). 
		
	Furthermore, since the closed-loop system is asymptotically stable and ℰ_P_t^ρ is a CAIS for every t≥ t_0, there exists a finite time T ≥ t_0 such that x_t ∈ℰ_P_t^ρ for all t≥ T. Thus, (<ref>) will be satisfied by the LQR controller after time T. As for all t≥ T, the state and input constraints are automatically satisfied (and are therefore, inactive), one can use the same arguments as in classical model-based policy iteration <cit.> to conclude the proof.
	
	
	In order to solve SDP (<ref>) using least squares methods, one need to collect a sequence of states such that the matrix Δ_x̅x̅ (obtained by replacing x in (<ref>) with x̅) has full column rank. This full rank condition is like the condition of persistence excitation (PE) in adaptive control theory. In order to satisfy this full rank condition, we add an exploration noise ν_t into the input to excite the system as in <cit.>. As the exploration noise ν_t goes to zero, the solution to (<ref>) will converge to the solution to model-based constrained policy evaluation problem (<ref>).  
	
	

§ NUMERICAL EXAMPLE

	

 §.§ Linear system with two states, one control input

	We randomly generate controllable systems of the form (<ref>) to test the proposed algorithm. A particular realization of these randomly generated systems,
	A = [[  1.1387  0.0491; -0.8680  0.9679 ]],    B = [[ -0.5507;  0.0758 ]]
	is investigated to illustrate constraint satisfaction and stability of the algorithm. Of course, A is unknown (and unstable), B is known, and it is verified that (A,B) is a controllable pair. The admissible state space is given by 𝕏={x∈ℝ^2: x_∞≤ 1}, and the operational cost is parameterized by Q=I_2 and R=0.5.
	
	For learning, the window length is fixed at N=8 samples (𝒯 = {8, 16, 24,…}), and the regularization parameter for policy updating is given by ρ_K=10^-4. Persistence excitation is ensured by generating uniformly distributed noise bounded within [-0.02, 0.02]. An initial policy is generated that satisfies state constraints using the randomly chosen cost matrices that are distinct from Q and R, and an initial condition is generated randomly on the boundary of the initial domain of attraction. Therefore, the initial state is ensured to be within 𝕏 but sufficiently far from the origin to require non-trivial control for stabilization.

	
	The results of the constrained policy iteration algorithm are illustrated in Fig. <ref>. In Fig. <ref>[A], a sequence of ellipsoids generated by our proposed algorithm is presented. Note that the ellipsoids generated in subsequent learning cycles after the first (the orange elongated ellipsoid) are not mere sub- or super-level sets of the initial ellipsoid; instead, the policy iterator allows for contractions and expansions on both x_1 and x_2 axes until the true policy is learned. As evident from subplots [B] and [C], state constraints are not violated throughout the learning procedure. The subplots [D, E] demonstrate the convergence of a sub-optimal initial control policy at t=0 to the true and optimal LQR policy at around t=24, after three learning cycles.
	
	

 §.§ Higher-dimensional linear system

	
		
	In addition, let us illustrate the performance for a dynamic system with 5 states and 2 control inputs. For this purpose, we have randomly generated 50 unique dynamic systems (<ref>) that are constructed to be (slightly) unstable but controllable, such that an initial stabilizing policy can be obtained. Similar to before, the system matrix A is unknown but the input matrix B is known. The admissible state and input space, respectively, is given by 𝕏={x∈ℝ^5: x_∞≤ 1} and 𝕌={u∈ℝ^2: u_∞≤ 1}, and the cost matrices read as Q=I_5 and R=0.5 I_2. Figure <ref> presents the closed-loop state and input trajectories when applying the proposed constrained ADP implementation (see Algorithm <ref>) to each of these generated test problems. Note that the learning window has been chosen to be equal to N=30 samples; indicating that 𝒯 = {30, 60,…}.
	
	From Figure <ref>, it can be observed that the 5-dimensional dynamic system is stabilized (the small perturbations in subplot [A] and [B] are due to the exploratory noise) in all of the generated test cases and both the state and input constraints are respected at all time. In addition, in most of the cases, the optimal policy is obtained relatively quickly in an amount of time that corresponds to 2 learning window lengths, i.e., 60 time steps in Figure <ref>. The policy error, computed using the matrix 2-norm of the difference between the current and optimal policy, generally decreases over time for all cases, under the necessary condition for persistence of excitation.	
	

§ CONCLUDING REMARKS

	In this paper, we provide a methodology for implementing constraint satisfying policy iteration for continuous-time, continuous-state systems via invariant sets. Benefits of our approach include computational tractability, and safety guarantees through constraint satisfaction. In future work, we will extend this framework to nonlinear systems.
	
	
	IEEEtran
	
