Robust Positively Invariant Sets for Linear Systems subject to model-uncertainty and disturbances

4th IFAC Nonlinear Model Predictive Control Conference International Federation of Automatic Control Noordwijkerhout, NL. August 23-27, 2012

Robust Positively Invariant Sets for Linear Systems subject to model-uncertainty and disturbances Furqan Tahir ∗ , Imad M. Jaimoukha ∗ ∗

Department of Electrical and Electronic Engineering, Imperial College, London Sw7-2AZ, U.K (e-mail: [email protected], [email protected])

Abstract: In this paper, we propose an algorithm for the computation of Robust Positively Invariant sets for linear discrete-time systems subject to norm-bounded model-uncertainty, additive disturbances and polytopic constraints on the input and state. The invariant set (optimized with respect to a measure of size) and the corresponding controller are computed simultaneously by solving a single LMI optimization problem. The novelty lies in the fact that the proposed scheme explicitly takes account of (norm-bounded) model-uncertainty and does not require any iterative computations or initial estimates of the invariant set or controller. A numerical example, taken from literature, demonstrates the applicability of proposed algorithm. Keywords: Invariant sets; LMI optimization; Model-uncertainty; Farkas Lemma; S-Procedure 1. INTRODUCTION Robust Positively Invariant (RPI) set refers to a bounded state-space region to which the system state can be confined, despite the presence of disturbances/uncertainties, through the application of a state-feedback control law (K). RPI sets are important in the robustness analysis and synthesis of controllers for uncertain systems [Blanchini, 1999]. These sets also play a fundamental role in various Robust Model Predictive Control (RMPC) schemes [Scokaert and Mayne, 1998]. In particular, they can be used to establish stability and recursive feasibility of the RMPC controller [Mayne et al., 2000]. Furthermore, the RPI set is also a suitable target set in robust time-optimal control schemes (see e.g. Bertsekas [1971] and Mayne [1997]). RPI sets have been the subject of extensive research over the past few decades. The two invariant set structures most often considered in the literature are ellipsoidal and polytopic [Kolmanovsky and Gilbert, 1998]. As noted in Blanco et al. [2010], the polytopic RPI sets hold computational advantages over the ellipsoidal sets and are therefore the focus of this paper. Many of the existing methods for computing polytopic RPI sets do not allow for the linear feedback law K to be optimized (see e.g. Gilbert and Tan [1991] and Pluymers et al. [2005]), which in turn has an adverse effect on the size of the computed set. Moreover, in most schemes, computation of polytopic RPI set requires iterative set computations (see e.g. Rakovic et al. [2005], Ong and Gilbert [2006], Rakovic and Kouramas [2006], and Benlaoukli et al. [2009]), which can often be inefficient. An algorithm proposed in Tahir [2010] does enable the computation of RPI sets and controller in one step. However, like most of the above schemes, its applicability is restricted 978-3-902823-07-6/12/$20.00 © 2012 IFAC

213

to systems that only involve bounded disturbances (i.e. no model-uncertainty). An exception to this is the algorithm proposed in Blanco et al. [2010] which does takes account of polytopic model-uncertainties. Here, the RPI set and the corresponding K are computed through an iterative sequence of convex programs. In this paper, we propose an algorithm to compute - in one step - the RPI set and corresponding controller K for systems subject to (norm-bounded) model-uncertainty and additive disturbances with polytopic constraints on the input and state. The RPI set is taken to be a hyperrectangle and is optimized with respect to its perimeter. We show that in the case of unstructured uncertainty, the proposed scheme yields the maximal/minimal hyperrectangle RPI set exactly. Furthermore, we also show that for the case with a fixed K and no uncertainty (only disturbance), the proposed algorithm computes the optimal set by solving a single Linear Program (LP). This paper is organized as follows. Section 2 is concerned with the system description and definition of RPI sets. In section 3, we present the LMI problem for the computation of the optimal RPI set and the corresponding inner controller. In section 4, we present a numerical example. Finally, we conclude in section 5. Notation and background material: The notation we use is fairly standard. R denotes the set of real numbers, Rn denotes the space of n-dimensional (column) vectors whose entries are in R, Rn×m denotes the space of all n × m matrices whose entries are in R and Dn×n denotes the space of diagonal matrices in Rn×n . For A ∈ Rn×m , AT denotes the transpose of A. If A ∈ Rn×n is symmetric, λ(A) denotes the smallest eigenvalue of A and we write A ≥ 0 if λ(A) ≥ 0 and A > 0 if λ(A) > 0. Analogous definitions apply to λ(A), A ≤ 0 and A < 0. For x, y ∈ Rn , x < y (and similarly ≤, > and ≥) is interpreted element-

10.3182/20120823-5-NL-3013.00032

IFAC NMPC'12 Noordwijkerhout, NL. August 23-27, 2012

wise. Given two sets M and V, such that M ⊂ Rn and V ⊂ Rn , the Minkowski (vector) sum is defined by M ⊕ V:= {m + v|m ∈ M, v ∈ V} and the Pontryagin difference is defined by M ∼ V:= {m|m + v ∈ M, ∀v ∈ V}. The identity matrix is denoted by I with the dimension inferred from the context. Let z ∈ Rn and denote the i-th element of z by zi . Then, diag(z) is the diagonal matrix whose (i, i) entry is zi . For square matrices A1 , . . . , Am , diag(A1 , . . . , Am ) denotes the block diagonal matrix whose i-th diagonal block is Ai . In the formulation, we make use the following version of the Farkas Lemma (see e.g. Jonsson [2006]). Lemma 1. Let c ∈ Rn , d ∈ R, A ∈ Rm×n and b ∈ Rm . Suppose there exists yˆ such that Aˆ y < b. Then the following two statements are equivalent: (1) cT y ≤ d ∀y such that Ay ≤ b. (2) ∃µ ∈ Rm such that µ ≥ 0, AT µ = c and bT µ ≤ d. We also use the following lemma, which is based on the results in Sun and Packard [2002]. Lemma 2. Let ∆ be defined as follows: ∆ := { diag(δ1 Iq1 , · · · , δl Iql , ∆l+1 , · · · , ∆l+r ) : δi ∈ R, |δi | ≤ 1 , ∆i ∈ Rqi×qi , k∆i k ≤ 1 } and define the subspaces Σ := { diag(S1 , · · · , Sl , λ1 Iql+1 , · · · , λs Iql+f ) : λj ∈ R, Si = SiT ∈ Rqi×qi } Γ := { diag(G1 , · · · , Gl , 0ql+1 , · · · , 0ql+f ) : Gi = −GTi ∈ Rqi×qi } Let R = RT , F, E, H be matrices of appropriate dimensions. We have det(I − H∆) 6= 0 and the inequality R + F ∆(I − H∆)−1 E + E T (I − ∆T H T )−1 ∆T F T < 0 for every ∆ ∈ ∆ if there exist S ∈ Σ and G ∈ Γ such that S> 0 and  R + E T SE F + E T SH + E T G <0 F T +H T SE + GT E H T SH + H T G + GT H−S 2. PRELIMINARIES We consider the following system

where z > 0. We also consider (non-symmetric) input constraints of the form: u ∈ U := {u ∈ Rnu : u ≤ u ≤ u} (3) where u > 0 and u < 0. An RPI set is defined as follows [Blanchini, 1999]: Definition 3. The set Z ⊂ Rn is a Robust Positively Invariant (RPI) set of system (1) if, with the control law u = Kx ∈ U, it satisfies: (A + Bp ∆Cq )Z ⊕ (Bu + Bp ∆Dqu )KZ ⊕ Bw W ⊆ Z. Thus, if the initial state lies in the set Z, then all subsequent states will be kept within this set by the inner controller u = Kx [Scokaert and Mayne, 1998]. Therefore, the set Z characterizes the evolution of system (1) for all possible disturbances w ∈ W and uncertainties ∆ ∈ ∆ [Kolmanovsky and Gilbert, 1998]. 3. RPI SET FORMULATION As mentioned above, Z is a suitable target set for RMPC schemes since a system state subject to uncertainties/disturbances cannot be regulated to zero. Since the RPI set is not unique, the choice of the size of Z depends on the specific problem being considered. The algorithm we propose in this work can be used the compute both the maximal (hyperrectangle) invariant set - which is a subset of ‘maximal RPI’ set [Perez et al., 2011] - as well as the minimal invariant set along with the corresponding inner controller K for the given input/state constraints. To this end, we consider the following problem to compute the largest possible Z (in terms of set-perimeter) along with the corresponding K: n X ζo =: max zi i=1 KZ ⊆ U (A + Bp ∆Cq )Z ⊕ (Bu + Bp ∆Dqu )KZ ⊕ Bw W ⊆ Z (4)

where zi denotes the ith element of the polytope z.

xk+1 = Axk + Bu uk + Bw wk + Bp pk ,

(1a)

pk = ∆qk ,

(1b)

qk = Cq xk + Dqu uk ,

(1c)

Remark 4. Note that the constraint (A + Bp ∆Cq )Z ⊕ (Bu +Bp ∆Dqu )KZ ⊕Bw W ⊆ Z ensures that the polytopic set is invariant by construction.

where xk , uk , wk , pk are the state, input, bounded disturbance and uncertainty vectors (respectively) at prediction step k; A is the state transition matrix and Bu , Bw and Bp are the input, disturbance and uncertainty distribution matrices, respectively. We assume that the pair (A,Bu ) is stabilizable and that the state xk is measured. We assume a polytopic (symmetric) disturbance of the form

Theorem 5. Let all definitions and variables be as defined above. Then, there exists an admissible Z and K, i.e. ones satisfying the constraints in (4), if there exist solutions m n i nw ρm , µix ∈ x , ρx , ∈ R , Sij , Si ∈ Σ, Gij , Gi ∈ Γ and µw ∈ R Rn , m ∈ Nu := {1, · · · , nu }, i, j ∈ Nx := {1, 2, . . . , n} to the following LMIs: m ˆT ρm x ≥ 0, ρx + K eum ≥ 0,

wk ∈ W := {w ∈ Rnw : −d ≤ w ≤ d}

T ˆT eTum u − 2eT ρm x − e K eum ≥ 0, ˆ T eum ≥ 0, ρm ≥ 0, ρm − K

(5)

ˆ T eum ≥ 0, −eTum u − 2eT ρm + eT K x

(6)

T µix ≥ 0, µiw ≥ 0, µiw +Bw ei ≥ 0,

(7)

where d > 0. Furthermore, we consider a norm-bounded uncertainty ∆ ∈ ∆, where ∆ is defined above. We consider the polytopic RPI set to be of the form: Z := {z ∈ Rn : −z ≤ x ≤ z} (2) 214

x

x

∀m ∈ Nu ,

IFAC NMPC'12 Noordwijkerhout, NL. August 23-27, 2012

 T i  ˆ T ei + eTi Bp Sij BpT ei ? −ej µx − eTj (AZd + Bu K)  ≤ 0 1 ˆ j + GTij BpT ei − (Cq Zd + Dqu K)e −Sij 2 (8)   ˆ µix , µiw ) Li (Zd , K, ?  ≤ 0 (9)  1 ˆ + GT B T ei − Si − (Cq Zd + Dqu K)e i p 2 ∀j ∈ Nx , ∀i ∈ Nx , T ˆ µix , µiw ) := −eT Zd e + dT (2µiw + Bw where Li (Zd , K, ei ) + i T T i T T ˆ e (2µx + (AZd + Bu K) ei ) + ei Bp Si Bp ei .

Furthermore, ei denotes the ith column of the n × n identity matrix, eum denotes the mth column of the nu ×nu identity matrix, e is the n-dimensional vector of ones, ˆ := KZd . Zd := diag(z) > 0 so that z = Zd e, and K

eTum Kx ≤ eTum u , ∀x ∈ Z

(10)

eTum Kx

(11)

≥

, ∀x ∈ Z

for all m ∈ Nu .

z T (2µix + (A + Bu K)T ei + (Cq + Dqu K)T∆T BpT ei )+ · · · T + dT (2µiw + Bw ei ) − eTi z ≤ 0. (15c) Pre-multiplying the inequalities in (15b) by Zd and using ˆ := KZd yields the set the re-definitions µix := Zd µix and K of inequalities: T µix ≥ 0 , µiw ≥ 0 , µiw + Bw ei ≥ 0 , i T ˆ ˆ T ∆T BpT ei ≥ 0 , µx + (AZd + Bu K) ei + (Cq Zd + Dqu K)

e

T

(2µix +

(16) T T T T ˆ ˆ (AZd +Bu K) ei +(Cq Zd +Dqu K) ∆ Bp ei )+· · · T + dT (2µiw +Bw ei )−eTi Zd e ≤ 0 (17)

Note that inequality conditions (16) and (17) are both necessary and sufficient for the invariance constraint in (4). Now inequality (16) can be written in the form:

By applying Lemma 1, it can be shown that (10) is satisfied m n if and only if there exist ρm x , ρx ∈ R such that m m m T ρm x ≥ 0 , ρx ≥ 0 , ρx = ρx + K eum , T m eTum u − z T ρm x − z ρx ≥ 0.

By eliminating ρm x from the above, we obtain the following equivalent conditions m T ρm x ≥ 0 , ρx + K eum ≥ 0, T T eTum u − 2z T ρm x − z K eum ≥ 0.

(12)

Note that (12) is nonlinear in z and K. By defining Zd := diag(z) > 0 so that z = Zd e, pre-multiplying the first and second inequality by Zd and introducing the m ˆ re-definitions ρm x := Zd ρx and K := KZd results in the inequalities (5), which are linear in all the variables. Analogous to the above procedure, using Lemma 1 on (11) followed by the linearization results in (6). Note that inequality conditions (5) and (6) are both necessary and sufficient for the input constraints in (4). Now, since the sets Z, ∆ and W are symmetric, the invariance constraint in (4) can be written as: eTi [(A + Bu K)x + Bp ∆(Cq + Dqu K)x + Bw w]

(15a)

µix ≥ 0 , µix +(A + Bu K)Tei +(Cq + Dqu K)T∆T BpT ei ≥ 0, (15b)

∀i ∈ Nx .

Proof. The input constraints in (4) can be written as: eTum u

T µiw ≥ 0 , µiw + Bw ei ≥ 0 ,

≤

eTi z

Rij + Fj ∆(I − H∆)−1 Ei + EiT (I −∆TH T )−1∆TFjT ≤ 0 (18) with     ˆ T ei BpT ei −eTj µixk − eTj (AZd + Bu K) Rij Ei  :=  1 ˆ T Fj H − eTj (Cq Zd + Dqu K) 0 2 for all j ∈ Nx , and i ∈ Nx , where ej denotes the jth column of the n × n identity matrix. Applying Lemma 2 on (18) yields the LMI (8). Similarly, using Lemma 2 on inequality (17) yields the LMI (9). This completes the proof. 2 Now let us define a variable ζ such that: n X ζ− eTi Zd ei ≤ 0

(19)

i=1

Then, it follows from Theorem 5 that an upper bound on ζo , call it ζ o , and the optimal Z and K can be computed by solving the following LMI problem: ζ o= max{ζ : (5, 6, 7, 8, 9, 19) are satisfied for some µiw∈ Rnw, m n µix ∈ Rn , ρm x , ρx ∈ R , Sij , Si ∈ Σ, Gij , Gi ∈ Γ,

m ∈ Nu , and i, j ∈ Nx }.

(20)

(13)

for all x ∈ Z, ∆ ∈ ∆, w ∈ W, and i ∈ Nx . It follows from Lemma 1 that (13) is satisfied if and only if there exist µix , µix ∈ Rn , and µiw , µiw ∈ Rnw such that µix ≥ 0 , µix ≥ 0 , µix = µix + (A + Bu K)T ei + (Cq + Dqu K)T ∆T BpT ei , T µiw ≥ 0 , µiw ≥ 0 , µiw = µiw + Bw ei ,

z T µix + z T µix + dT (µiw + µiw ) ≤ eTi z. Eliminating µix and µiw from the above inequalities yields: 215

Remark 6. The above problem can be modified to compute the smallest (hyperrectangle) invariant set. It is easy to verify that for this case, (20) becomes a minimization Pn problem and the inequality (19) changes to: ζ − i=1 eTi Zd ei ≥ 0. Remark 7. If the model-uncertainty is unstructured (that is, if ∆ := {∆ ∈ Rq×q , k∆k ≤ 1 }) then the matrix variables from Lemma 2 simply reduce to G = 0 and S = λI for some scalar λ > 0. Note that in this case, the LMI conditions (8) and (9) become both necessary and

IFAC NMPC'12 Noordwijkerhout, NL. August 23-27, 2012

sufficient for the invariance constraint in (4). It follows that for the unstructured uncertainty, problem (20) computes the maximal/minimal (hyperrectangle) invariant set and the optimal feedback gain K exactly (i.e. ζ o = ζ0 ). Remark 8. Note that we can easily incorporate state constraints in the proposed RPI set formulation by including, in problem (20), an LMI condition of the form: Zd e ≤ zc . Remark 9. The proposed algorithm can be modified to optimize the volume of the RPI set (instead of perimeter). In this case, we can simply replace (19) by: 

n X

eTi Zd ei

  i=1  n X  T e1 Zd e1 − eTi Zd ei   i=2   ζ

eT1 Zd e1 − n X

n X



eTi Zd ei

ζ

i=2

eTi Zd ei

0

i=1

n X

0

eTi Zd ei

    ≥ 0    

i=1

Remark 10. Suppose there is no uncertainty (i.e. ∆ = 0). nT Moreover, let y = [z T µ1T µ1T · · · µnT w · · · µw x x 1T nu T 1T nu T T n(nw +n+2nu +1) ρx · · · ρx ρx · · · ρx ] ∈ R and define a column vector c whose first n elements are 1 while the rest are 0 i.e. c := [1, . . . , 1, 0, . . . , 0]T ∈ Rn(nw +n+2nu +1) . Then, it is easy to verify that for the case with a fixed K, problem (20) can be transformed into a Linear Program (LP) which optimizes the objective function cT y subject to constraints of the form Ey ≤ G. 4. NUMERICAL EXAMPLE We consider an uncertainty-added version of the example in Rakovic et al. [2005]. In particular, we have dynamics (1) with:       11 1 10 A= , Bu = , Bw = , 01 1 01       0.2 0 11 1 Bp , Cq = , Dqu = , 0 0.2 01 1 Here the disturbance is given as w ∈ W where  W :=

w ∈ R2 :

.



−0.5 −0.5



 ≤w≤

0.5 0.5



Furthermore, we consider model-uncertainty ∆ ∈ ∆ where ∆ := {diag(δ1 , δ2 ) : δi ∈ R, |δi | ≤ 1}. Finally, we consider input constraints u ∈ U := {u ∈ R : −3 ≤ u ≤ 3}. Solving the LMI problem (20) yields the following (largest) RPI set Z and controller K:  Z=

x ∈ R2 :



−3.269 −2.038



 ≤x≤

K = −[0.294

3.269 2.038

 ,

1]

Now solving (20) as a minimization problem yields the following Z and controller K: 216

Fig. 1. Simulation result with wk = d cos(k), ∆ = diag(1,1)      −0.5 0.5 2 Z= x∈R : ≤x≤ , (21) −1.3 1.3 K = −[1

1].

Figure 1 shows the simulation results for RPI set (21). We note that even with the initial state on the hyperrectangle vertex (circled) and the worst-case persistent uncertainty along with a sinusoidal disturbance, the computed controller K manages to keep the state within the RPI set (dashed). In the absence of uncertainty (i.e. Bp = 0) and with T disturbance such that d = [1 1] , solving (20) as a minimization problem yields the RPI set with polytope T z = [1 2] and K = [−1 −1]. Now, for this case, the minimal RPI (mRPI) set is given by: F∞ = ⊕∞ i=0 (A + Bu K)i Bw W [Kolmanovsky and Gilbert, 1998]. Since F∞ invloves Minkowski’s sum of infinite many terms, therefore it is generally impossible to compute unless the system dynamics are nilpotent [Mayne, 1997]. It is interesting to note that for the considered example, our algorithm yields a K (above) for which, the closed-loop dynamics are nilpotent i.e. (A + Bu K)i = 0 ∀ i ≥ 2. In fact, the set F∞ (computed through the evaluation of the Minkowski sum) comes out to be exactly the same as the Z obtained through our algorithm above. In other words, for this example, our algorithm yields the mRPI set exactly. 5. CONCLUSION In this paper, we have proposed an algorithm for the one-step computation of polytopic RPI sets along with their corresponding linear state-feedback gain for systems involving norm-bounded model-uncertainties and additive disturbances as well as constraints on the state and the input. In the literature, the problem of computing a suitable RPI set generally requires iterative computations and is mostly restricted to systems involving only disturbances (i.e. no uncertainty). Furthermore, many of the existing methods require the gain K to be fixed which has an adverse effect

IFAC NMPC'12 Noordwijkerhout, NL. August 23-27, 2012

on the size of the computed set. The proposed algorithm, however, takes account of model-uncertainty and computes both an optimal (smallest or largest perimeter-wise) RPI set and the corresponding control law K through a single LMI optimization problem. We have shown that in the case of unstructured uncertainty, the algorithm computes the optimal Z and K exactly. Also, for the case with a fixed K (and no uncertainty), the optimal (hyperrectangle) RPI set can be computed through a single Linear Program. REFERENCES Blanchini, F. (1999), “Set Invariance in Control,” Automatica, 35, 1747–1767. Scokaert, P.O.M., and Mayne, D.Q. (1998), “Min-max Feedback Model Predictive Control for Constrained Linear Systems”, IEEE Transactions on Automatic Control, 43, 1136–1142. Mayne, D.Q., Rawlings, J.B., Rao, C.V., and Scokaert, P.O.M. (2000), “Constrained Model Predictive Control: Stability and Optimality,” Automatica, 36, 789–814. Bertsekas, D.P., and Rhodes, I.B. (1971), “On the minimax reachability of target sets and target Tubes,” Automatica, 7, 233-247. Mayne, D.Q. and Schroeder, W.R. (1997), “Robust timeoptimal control of constrained linear systems,” Automatica, 33, 2103-2118. Kolmanovsky, I., and Gilbert, E.G. (1998), “Theory and Computation of Disturbance Invariant Sets for Discretetime Linear Systems,” Mathematical Problems in Engineering, 4, 317–367. Blanco, T.B., Cannon, M., and De Moor, B. (2010), “On efficient computation of low-complexity controlled invariant sets for uncertain linear systems,” International Journal of Control, 83, 1339-1346. Gilbert, E and Tan, K. (1991), “Linear systems with state and control constraints: The theory and application of maximal output admissible sets,” IEEE Trans. on Automatic Control, 36, 1008-1020. Pluymers, B., Rossiter, J. Suykens, J., and De Moor, B. (2005) “Interpolation based MPC for LPV systems using Polyhedral invariant sets,” in Proc. of the American Control Conf., 810-815. Rakovic, S.V., Kerrigan, E.C., Kouramas, K.I., and Mayne, D.Q. (2005), “Invariant approximations of the minimal robust positively invariant set,” IEEE Transactions on Automatic Control, 50, 406–410. Ong, C.-J., and Gilbert, E. G., (2006) “The minimal disturbance invariant set: Outer approximations via its partial sums,” Automatica, 42, 1563-1568. Rakovic, S. V., Kouramas, K. I. (2006), “The minimal robust positively invariant set for linear discrete time systems: Approximation methods and control applications,” in Proc. of the 45th IEEE Conf. on Decision and Control, San Diego, USA. Benlaoukli, H., Olaru, S., Marinkovic, S. and Boucher, P. (2009), “Invariant set constructions for feasible reference tracking,” in Proc. of the IEEE International Conf. on Control and Automation, Christchurch, New Zealand. Tahir, F. (2010), “Efficient computation of Robust Positively Invariant Sets with Linear State-feedback Gain as a Variable of Optimization,” in Proc. of the 7th Interna217

tional Conference on Electrical Engineering, Computing Science and Automatic Control, Mexico. Jonsson, U.T. (2006), “A Lecture on the S-procedure,” Lecture notes, Division of Optimization and Systems Theory, Royal Institute of Technology, Stockholm, Sweden. Sun, K., and Packard, A. (2002), “Robust H2 and H∞ Filters for Uncertain LFT Systems,” in Proceedings of the 41st IEEE Conference on Decision and Control, Las Vegas, Nevada, USA. P´erez, E., Ari˜ no, C., Blasco, F.X., and Mart´ınez, M.A. (2011), “Maximal Closed Loop Admissible Set for Linear Systems with Non-convex Polyhedral Constraints”, Journal of Process Control, 21, 529–537.