Full text
Analysis and Control of a Lur’e System with Sector-Bounded, Slope-Restricted Nonlinearities Using Linear and Bilinear Matrix Inequalities Joint Bachelor’s Thesis Degree in Physics Engineering Degree in Mathematics Universitat Politècnica de Catalunya Sandra Wells Cembrano Supervisor: Prof. Richard D. Braatz Co-supervisor: Prof. Carlos A. Ocampo Martínez July 2020
Abstract ENG – The aim of this thesis is to analyze and control a Lur’e system, i.e., a system with a linear time-invariant forward path and a nonlinear static feedback with sector and slope constraints on its nonlinearities. This description encompasses a broad class of systems that commonly arise in a wide range of engineering disciplines. The approach of this work is based on Lyapunov stability theory and uses linear matrix inequalities to propose criteria for the system’s absolute stability and `2and RMS-gains, and bilinear matrix inequalities to propose criteria for state feedback control and state estimation. These types of conditions, in particular linear matrix inequalities, are conveniently treated computationally, as the optimization problems that they give rise to are convex. Numerical results are obtained for several examples, showing a significant improvement with respect to previous conditions in the literature. CAT – L’objectiu d’aquesta tesi és l’anàlisi i el control d’un sistema de Lur’e, és a dir, un sistema amb un camí directe lineal i invariant en el temps i amb una realimentació no lineal estàtica, amb restriccions de sector i de pendent en les seves no-linealitats. Aquesta descripció engloba una àmplia classe de sistemes que apareixen comunament en una àmplia gamma de disciplines d’enginyeria. L’enfocament d’aquest treball es basa en la teoria d’estabilitat de Lyapunov i utilitza desigualtats de matrius lineals per proposar criteris per a l’estabilitat absoluta del sistema i els guanys `2i RMS, i desigualtats de matrius bilineals per proposar criteris per al seu control per realimentació d’estat i estimació d’estat. Aquest tipus de condicions, en particular les desigualtats de matrius lineals, són convenients de tractar computacionalment, ja que els problemes d’optimització a que donen lloc són convexos. S’obtenen resultats numèrics per a diversos exemples, els quals mostren una millora significativa respecte a les condicions anteriors de la literatura. ESP – El objetivo de esta tesis es el análisis y el control de un sistema de Lur’e, es decir, un sistema con un camino directo lineal e invariante en el tiempo y con una realimentación no lineal estática, con restricciones de sector y de pendiente en sus no-linealidades. Esta descripción engloba una amplia clase de sistemas que aparecen comúnmente en una amplia gama de disciplinas de ingeniería. El enfoque de este trabajo se basa en la teoría de estabilidad de Lyapunov y utiliza desigualdades de matrices lineales para proponer criterios para la estabilidad absoluta del sistema y las ganancias `2y RMS, y desigualdades de matrices bilineales para proponer criterios para su control por realimentación de estado y estimación de estado. Este tipo de condiciones, en particular las desigualdades de matrices lineales, son convenientes tratar computacionalmente, ya que los problemas de optimización a que dan lugar son convexos. Se obtienen resultados numéricos para varios ejemplos, los cuales muestran una mejora significativa respecto a las condiciones anteriores de la literatura. Keywords: Nonlinear control, Lyapunov stability, Lur’e systems, LMIs, BMIs. Paraules clau: Control no lineal, estabilitat de Lyapunov, sistemes de Lur’e, LMIs, BMIs. Palabras clave: Control no lineal, estabilidad de Lyapunov, sistemas de Lur’e, LMIs, BMIs. Americal Mathematical Society (AMS) classification: 93C10 1
Acknowledgments I would like to give my biggest thanks to Prof. Richard Braatz for everything he has done as my supervisor. For letting me choose this topic as the focus of my thesis, which not only did I find very exciting, but which also has great potential to improve our current grasp of nonlinear control. For giving me great guidance and advice throughout the development of my thesis, and for always finding opportunities to teach me about applied control, pure math and everything in between. For his generosity in supplementing my funding for my stay in Cambridge. For his flexibility and email availability during the global-pandemic-induced remote working, and for finding me a desk with enhanced access to natural light when my work was still done in person at MIT. It is undoubtedly an honor to have been a part of his group and to now receive his encouragement and support for my professional future. I would like to give huge thanks to soon-to-be-Dr. Anastasia Nikolakopoulou from the Braatz group, for her constant support and encouragement from my very first day. For giving me advice on how to survive at MIT, on how to survive my project and even on how to survive the lockdown. All of this, of course, in addition to answering many questions about my work in a shorter time than I ever took to check my inbox for her reply, regardless of how many conference deadlines she had ahead. I would like to extend my thanks to the whole Braatz group for their warm welcome at my arrival. To Tam, Kaylee, Moo Sun, Rohan, Andy, Patrick, Matt, Andreas, and to all the rest of members who, although I did not get a chance to talk to, contributed to the great atmosphere in the group. Thank you, as well, to Angelique, for all her very kind and fast help with all administrative matters. Thank you, finally, to whoever had the brilliant idea that our group meetings had to have free donuts. I am also very grateful to Prof. Carlos Ocampo Martínez, for all of his comments and suggestions on this written report. In addition, for his fast attention on the bureaucracy matters of this thesis. And last, but definitely not least, for introducing me to Prof. Braatz. This work would not have been possible without him. I would like to give thanks to CFIS for partially funding my stay at MIT. Thanks to their grant, my long lockdown was spent in a very livable room in Cambridge. Thanks to CFIS, as well, it was possible for me to conduct this work at MIT. I would like to show my immense gratitude to my family, for their support, encouragement and motivation not only for this thesis but throughout all of my education. They are behind everything that I achieve. To Paul, for his incomparable company and for everything that we help each other accomplish, and to all my friends from Barcelona, who have stood by me and supported me through all these years. To my friends in Cambridge, too, for providing me with safe and socially-distanced entertainment and city exploration plans during my stay. Finally, I want to thank everyone who has, in every possible way, made the COVID-19 pandemic easier to live through. 2
Contents List of Abbreviations 5 1 Introduction 6 1.1 Mathematical Preliminaries . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 1.1.1 Linear and Bilinear Matrix Inequalities . . . . . . . . . . . . . . . . . . . 6 1.1.2 The S-Procedure for Quadratic Forms . . . . . . . . . . . . . . . . . . . . 8 1.1.3 The Schur Complement Lemma . . . . . . . . . . . . . . . . . . . . . . . 9 1.1.4 Initial Definitions and Notation . . . . . . . . . . . . . . . . . . . . . . . 10 1.2 Literature Review on the Lur’e Problem . . . . . . . . . . . . . . . . . . . . . . 11 1.3 Problem Statement and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . 14 1.3.1 Analytical Approach . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 1.3.2 NumericalApproach ............................. 15 1.4 ThesisOutline..................................... 16 2 Stability Analysis 18 2.1 Analytical LMI Derivation for Stability Analysis . . . . . . . . . . . . . . . . . . 18 2.2 Numerical Results for Stability Analysis . . . . . . . . . . . . . . . . . . . . . . 23 3 Gain Analysis 27 3.1 Analytical LMI Derivation for Gain Analysis . . . . . . . . . . . . . . . . . . . . 27 3.2 Numerical Results for Gain Analysis . . . . . . . . . . . . . . . . . . . . . . . . 33 4 State Feedback Controller Design 35 4.1 Analytical BMI Derivation for State Feedback Controller Design . . . . . . . . . 35 4.2 Illustrative Numerical Example of a State Feedback Controller . . . . . . . . . . 37 5 State Estimator Design 40 5.1 Analytical BMI Derivation for State Estimator Design . . . . . . . . . . . . . . 40 5.2 Illustrative Numerical Example of a State Estimator . . . . . . . . . . . . . . . . 46 6 Conclusions 49 6.1 Contributions ..................................... 49 6.2 Furtherwork ..................................... 50 7 References 51 Appendices 53 3
A MATLAB Scripts for Chapter 2 53 B MATLAB Scripts for Chapter 3 62 C MATLAB Scripts for Chapter 4 71 D MATLAB Scripts for Chapter 5 78 4
List of Abbreviations LTI: Linear time invariant LMI: Linear matrix inequality BMI: Bilinear matrix inequality 5
1 Introduction This thesis was conducted from February 2020 to July 2020 at the Massachusetts Institute of Technology under the supervision of Prof. Richard D. Braatz, as part of the requirements for graduation from CFIS, Universitat Politècnica de Catalunya – BarcelonaTech, with a double Bachelor’s degree in Physics Engineering and in Mathematics. The thesis studies a specific class of nonlinear systems known as Lur’e systems, essentially a combination of a linear time-invariant forward path and a nonlinear static feedback. In particular, the thesis deals with a subclass of these systems, which have sector bounds and slope restrictions on their nonlinearities. Lur’e systems arise very commonly in a wide range of disciplines including chemical engineering, mechanical engineering, and aerospace engineering, among others. The sector and slope constraints on the nonlinearities do not largely reduce the applicability of the problem. A Lur’e system with these constraints can describe broad classes of nonlinear feedback behaviors, included systems described by dynamic neural networks. The neural network model of the nonlinearities will consist of functions such as hyperbolic tangents, which satisfy these constraints. The aim of this thesis is to find conditions to ensure global asymptotic stability, compute input-output gains, and design state feedback control and state estimation for Lur’e systems with sector-bounded, slope-restricted nonlinearities. The approach in this thesis for the study of these systems is based on linear and bilinear matrix inequalities. 1.1 Mathematical Preliminaries 1.1.1 Linear and Bilinear Matrix Inequalities This section summarizes the most important aspects of linear and bilinear matrix inequalities that are used in this thesis, and are obtained from the more detailed tutorial [1]. Another more complete study is available in the book [2]. A linear matrix inequality (LMI) has the form F(x) = F0+ m X i=1 xiFi≺0,(1) where x∈Rmand the matrices Fi∈Rn×n, i = 0,1, . . . , m are symmetric and known. The matrix F(x)is an affine function of the elements of the variable x. The inequality (1) denotes that F(x)is a negative definite matrix, that is, zTF(x)z < 0,∀z6= 0, z ∈Rn. 6
LMIs can be defined analogously for positive definite, negative semidefinite, and positive semidefinite matrices. The first two cases are called strict LMIs, while the last two are called nonstrict LMIs. Linear inequalities, convex quadratic inequalities, matrix norm inequalities, and various constraints from control theory such as Lyapunov and Riccati inequalities can all be written as LMIs. Thus, LMIs are a useful tool for solving a wide variety of optimization and control problems. The LMI (1) is said to be feasible if the set {x|F(x)≺0}is nonempty, and the analogous sets can be defined similarly for the rest of cases. An important property of LMIs is that the set {x|F(x)≺0}is convex, i.e., (1) forms a convex constraint on x. To see this, let x, y be two vectors such that F(x)≺0and F(y)≺0,and let λ∈(0,1).Then, F(λx + (1 −λ)y) = F0+ m X i=1 (λxi+ (1 −λ)yi)Fi =λF0+ (1 −λ)F0+λ m X i=1 xiFi+ (1 −λ) m X i=1 yiFi =λF(x) + (1 −λ)F(y) ≺0. The advantage of formulating control problems in terms of convex optimizations (when possible) is that wide classes of convex optimizations can be solved in polynomial time. Convex optimizations often arise in engineering practice and many can be written as LMIs, which is the strength of using LMI formulations: convex optimizations over LMIs are solvable in polynomial time. A bilinear matrix inequality (BMI) has the form F(x, y) = F0+ m X i=1 xiFi+ l X j=1 yjGj+ m X i=1 l X j=1 xiyjHij ≺0,(2) where Fi, Gj∈Rn×n, i = 0,1, . . . , m, j = 0,1, . . . , l, are symmetric known matrices, Hij ∈ Rn×n, i = 0,1, . . . , m, j = 0,1, . . . , l are known matrices and x∈Rm,y∈Rlare the variables. As for LMIs, BMIs can be negative or positive definite and negative or positive definite, and the same definitions of strict and nonstrict apply. BMI feasibility is defined analogously to that of LMIs. BMIs commonly arise when formulating control design procedures for those uncertain and/or nonlinear systems in which an LMI formulation is not available. In fact, nearly every problem of interest in control can be formulated in terms of optimizations over BMIs. 7
If xis fixed, a BMI is an LMI for yand, therefore, convex in y. Analogously, if yis fixed, a BMI is an LMI for xand convex in x. However, BMIs are not jointly convex in xand y, and control problems expressed as BMIs are not convex optimizations. BMI optimizations are NP-hard and cannot be ensured to be solvable in polynomial time. In more practical terms, this classification implies that algorithms for finding global solutions to optimizations over BMIs are not efficient for large-scale problems. 1.1.2 The S-Procedure for Quadratic Forms The S-procedure is useful for reformulating mathematical structures that commonly arise in Lyapunov control as LMIs. Its statement and proof are extracted from [1] and shown below. Lemma (S-Lemma or S-Procedure): Let fi(x),i= 0, . . . , p be quadratic forms with respect to x∈Rn:fi(x) = xTTix, where Tiare symmetric matrices. If there exist τ1≥ 0, . . . , τp≥0such that f0− p X i=1 τifi≥0,∀x, (3) or, equivalently, T0− p X i=1 τiTi<0, then xTT0x≥0∀xsuch that xTTix≥0, i = 1, . . . , p. Proof: If there exist τ1≥0, . . . , τp≥0such that (3) holds for all xthen (3) also holds for all x such that fi(x)≥0,∀i= 1, . . . , p. Then, for all such x, it must hold that f0(x)≥ p X i=1 τifi(x)≥0 since the summation is over terms that are all nonnegative. The S-procedure also holds, and is proved similarly, for the case where the main inequality is strict. If there exist τ1≥0, . . . , τp≥0such that T0− p X i=1 τiTi0, 8
analysis. The sign definiteness of these integral terms is implied by the sector-bounded property of the nonlinearities. The Lyapunov function is radially unbounded with respect to ¯xk, zero for xk= 0, and strictly positive for all nonzero xk∈Rn. Conservatism is reduced by considering the sector boundedness, which is a global condition on the nonlinearities, and even more reduced by considering the slope restriction, which poses a local condition to their behavior inside the bounded conical sector. These properties are introduced by two means: 1. The local slope restriction is used to find the tightest possible quadratic upper bounds for the integral terms in the Lyapunov function. These quadratic upper bounds lead to LMI conditions ensuring that the Lyapunov function decreases for each k. The tighter that the upper bounds are, the less conservative are the derived conditions. 2. The obtained LMIs need to be satisfied for all values of the state variables that satisfy the sector-boundedness and slope-restriction conditions. This step is introduced through the S-procedure, which is applied to the quadratic forms resulting from the LMI formulation explained above. 1.3.2 Numerical Approach The LMI problems to be solved numerically in this work are LMI feasibility problems, i.e., finding whether a given system of LMIs is feasible. These problems were solved using the LMI Lab [12] which is incorporated into the Matlab Robust Control Toolbox [13, 14]. The commands from this toolbox used in this thesis are: setlmis:Initializes description of LMI system lmivar:Specifies matrix variables in LMI problem newlmi:Attaches identifying tag to LMIs lmiterm:Specifies term content of LMIs getlmis:Calls internal description of LMI system feasp:Computes solution to given system of LMIs In particular, for a given LMI system lmisys, the function call is: [tmin,xfeas] = feasp(lmisys,options,target) The value of tmin indicates whether lmisys is feasible. If negative, lmisys is feasible. If positive, the LMI is infeasible or, if sufficiently small, the system may be feasible but not strictly feasible. For a feasible system, xfeas gives a set of solutions to the variables lmivar that make the problem feasible. 15
For a set of examples, Chapters 2 and 3 solve optimizations over LMIs. Specifically, the maximum or minimum value of a parameter is found such that the given system of LMIs is feasible. In particular, these correspond to the maximum lower bound on the stability margin and the minimum upper bound on the `2or RMS-gain of the examples of Lur’e systems with sector-bounded, slope-restricted nonlinearities. These optimizations are solved by the procedure: 1. Define the system of LMIs. 2. Define the example system. 3. Find a value of the parameter to optimize such that the LMI system is feasible for the example system, and one for which the LMI system is infeasible, using feasp. 4. Use the bisection method to find the maximum or minimum value of the parameter for which tmin is negative, with the previously found feasible and infeasible values as initial search points. The scripts for the described procedure can be found in Appendices A and B. For a given example, Chapters 4 and 5 find whether a system of LMIs is feasible, which is done through the simpler procedure: 1. Define the example system. 2. Define the system of LMIs. 3. Find whether the system of LMIs is feasible for the example system. In these chapters, additionally, the example system is iterated to show its behavior. The scripts for the above procedures and for plotting the system behavior are included in Appendices C and D. 1.4 Thesis Outline Chapter 2 contains the main enabling results by dealing with the stablity analysis of the studied system. The properties of the system’s nonlinearities are used with the aim of finding the tightest possible quadratic upper bounds for the integrals in the Lyapunov function (14), which is the most important contribution of this thesis. The approach was used in this thesis was published in a past Master’s thesis [15]. The previously published bounds on the integrals are disproved by counterexamples, with a simple example being φ(σ) = σand using qk, qk+1 >0. Through the new, corrected quadratic upper bounds presented in this thesis, LMIs are derived as sufficient conditions for the studied system (13) to be absolutely stable. The obtained LMIs are also used to solve the optimization max ξsuch that the system is absolutely stable, which gives 16
a lower bound on the stability/robustness margin of the system. This margin bound is found for a series of numerical examples and compared with the most relevant results in the literature. Chapter 3 considers performance analysis of the Lur’e system. The integral bounds are used again to obtain LMIs for computing upper bounds on the system’s `2-gain and RMS-gain. The tightness of the integral bounds is useful in this chapter as well, to reduce conservatism in the result of this value. The upper bound on the gain is obtained for a few of the numerical example systems analyzed in Chapter 2. Chapter 4 considers the robust stabilization of the system. The bounds and LMIs obtained in Chapter 2 are used to analyze the stability of a system with a proportional controller uk=Kxk. The addition of the new matrix variable Kresults in higher order matrix inequalities, which are reduced to BMIs through the Schur complement lemma. A numerical example is given to illustrate the results. Chapter 5 studies the state estimation of the Lur’e system. The dynamics of the error between the system variables and the estimated variables are analyzed. Because these dynamics are found to have nonlinearities such that φ /∈Φ[0,ξ] sb ∩Φ[0,µ] sr and φ∈Φ[0,µ] sb , new bounds are found for the integrals in the Lyapunov function and, consequently, the Lyapunov function used for the error dynamics system is modified. Further, the corresponding higher order matrix inequalities are found, and again reduced to BMIs through the Schur complement lemma. A numerical example is shown as well to illustrate the results. 17
2 Stability Analysis This chapter derives sufficient conditions for a system of the form (13) to be globally asymptotically stable. The stability or robustness margin is defined accordingly and is then found numerically for a set of examples. The numerical results are compared to the results obtained with the criteria presented in the literature review to quantify the improvement given by the presented criterion. 2.1 Analytical LMI Derivation for Stability Analysis The derivation of this stability criterion requires finding a quadratic upper bound to the variation between sampling instances of the Lyapunov function (14). An essential step for finding this bound is given in Lemma 1 and is later used in the proof of Theorem 1. Lemma 1: Let qk, qk+1 be any two consecutive sampling instances of qand φ∈Φ[0,ξ] sb ∩Φ[0,µ] sr as defined in system (13). Let φk≡φ(qk)and φk+1 ≡φ(qk+1). Then, for each i= 1, . . . , nq, φk,i(qk+1,i −qk,i) + 1 2µi (φk+1,i −φk,i)2 ≤Zqk+1,i qk,i φi(σ)dσ ≤φk+1,i(qk+1,i −qk,i)−1 2µi (φk+1,i −φk,i)2, (15) provided that µi6= 0. For µi= 0, the value of the integral is 0. Proof: The case µi= 0 is trivial. The bounds also clearly hold for qk,i =qk+1,i. For the nontrivial case, first consider the case where qk,i < qk+1,i. Then the local slope restriction property of φgives 0<φi(σ)−φi(ˆσ) σ−ˆσ≤µi,∀σ6= ˆσ∈R=⇒φi(σ)≤min{φk,i +µi(σ−qk,i), φk+1,i}, ∀σ∈[qk,i, qk+1,i]. (16) Let σk∈[qk,i, qk+1,i]be the value of σat which φk,i +µi(σ−qk,i) = φk+1,i, then σk≡qk,i +φk+1,i −φk,i µi and, from (16), φi(σ)≤φk,i +µi(σ−qk,i),∀σ∈[qk,i, σk], φi(σ)≤φk+1,i,∀σ∈[σk, qk+1,i]. 18
The integral in (15) can then be separated into two parts, satisfying: Zqk+1,i qk,i φi(σ)dσ =Zσk qk,i φi(σ)dσ +Zqk+1,i σk φi(σ)dσ ≤Zσk qk,i (φk,i +µi(σ−qk,i)) dσ +Zqk+1,i σk φk+1,idσ =φk+1,i(qk+1,i −qk,i)−1 2µi (φk+1,i −φk,i)2. The lower bound is obtained analogously using that φi(σ)≥max{φk,i, φk+1,i +µi(σ−qk+1,i)},∀σ∈[qk,i, qk+1,i],(17) again implied by the slope restriction, leading to Zqk+1,i qk,i φi(σ)dσ ≥φk,i(qk+1,i −qk,i) + 1 2µi (φk+1,i −φk,i)2. The case for qk,i > qk+1,i consequently holds by using that Zqk+1,i qk,i φi(σ)dσ =−Zqk,i qk+1,i φi(σ)dσ and applying the bounds to integral on the right-hand side. Theorem 1: Given the system (13), a sufficient condition for global asymptotic stability is the existence of a positive semidefinite matrix P=PT∈R(n+np+nq)×(n+np+nq), with a positive definite submatrix P11 =PT 11 ∈Rn×nand diagonal positive semidefinite matrices Q, ˜ Q, T, ˜ T, N ∈ Rnq×nqsuch that G:= G11 G12 G13 GT 12 G22 G23 GT 13 GT 23 G33 ≺0, where G11 =AT(P11 +CTPT 13 +P13C+CTP33C)A−P11 −CTPT 13 −P13C−CTP33C+ +ATCT˜ QXCA −CT˜ QXC, G12 =AT(P11 +CTPT 13 +P13C+CTP33C)B−P12 −CTPT 23 −P13D−CTP33D+ +ATCT˜ QXCB −CT˜ QXD + (CA −C)T˜ Q−CTT+ (CA −C)TN, 19
G13 =ATP12 +ATCTPT 23 +ATP13D+ATCTP33D−(CA −C)TQ+ATCT˜ QXD −ATCT˜ T− −(CA −C)TN, G22 =BT(P11 +CTPT 13 +P13C+CTP33C)B−P22 −DTPT 23 −P23D−DTP33D−QM−1+ +BTCT˜ QXCB −DT˜ QXD +˜ Q(CB −D)+(CB −D)T˜ Q−˜ QM−1−2TX−1−TD− −DTT−2NM−1+N(CB −D) + (CB −D)TN, G23 =BTP12 +BTCTPT 23 +BTP13D+BTCTP33D+M−1Q−(CB −D)TQ+BTCT˜ QXD+ +˜ QD +˜ QM−1−BTCT˜ T−(CB −D)TN+ 2NM−1+ND, G33 =P22 +DTPT 23 +P23D+DTP33D−QD −DTQ−QM−1+DT˜ QXD −˜ QM−1− −2˜ TX−1−˜ TD −DT˜ T−2NM−1−ND −DTN, where M:= diag{µ1, . . . , µnq}and X:= diag{ξ1, . . . , ξnq}. Proof: Given the Lyapunov function (14), a sufficient condition for the global asymptotic stability of the system is for the inequality ∆V(xk)<0,∀k≥0(18) to be satisfied. The variation between the two sampling instances kand k+ 1 is expressed as ∆V(xk) = ζT k(AT aPAa−ET aPEa)ζk+ 2 nq X i=1 Qii Zqk+1,i qk,i φi(σ)dσ + 2 nq X i=1 ˜ Qii Zqk+1,i qk,i [ξiσ−φi(σ)]dσ, (19) where ζk:= xk pk pk+1 , Aa:= A B 0 0 0 I CA CB D , Ea:= I0 0 0I0 C D 0 . In order to find an LMI condition that implies the inequality (18), Lemma 1 is used to find quadratic upper bounds on the two integral terms of ∆V(xk): 2 nq X i=1 Qii Zqk+1,i qk,i φi(σ)dσ ≤2 nq X i=1 Qii φk+1,i(qk+1,i −qk,i)−1 2µi (φk+1,i −φk,i)2=ζT kU1ζk, 20
where M:= diag{µ1, . . . , µnq}<0, and U1=UT 1is defined as U1:= 0 0 −(CA −C)TQ ∗ −QM−1(M−1−(CB −D)T)Q ∗∗−QD −DTQ−QM−1 , and for the second integral term, 2 nq X i=1 ˜ Qii Zqk+1,i qk,i [ξiσ−φi(σ)] dσ = 2 nq X i=1 ˜ Qii (ξi 2(q2 k+1,i −q2 k,i)−Zqk+1,i qk,i φi(σ)dσ) ≤2 nq X i=1 ˜ Qii ξi 2q2 k+1,i −q2 k,i−φk,i(qk+1,i −qk,i)−1 2µi (φk+1,i −φk,i)2=ζT kU2ζk, where X:= diag{ξ1, . . . , ξnq}<0, and U2=UT 2is defined as U2:= ATCT˜ QXCA −CT˜ QXC ATCT˜ QXCB −CT˜ QXD+ +(CA −C)T˜ QATCT˜ QXD ∗ BTCT˜ QXCB −DT˜ QXD+ +˜ Q(CB −D)+ +(CB −D)T˜ Q−˜ QM−1 BTCT˜ QXD +˜ QD +QM−1 ∗ ∗ DT˜ QXD −˜ QM−1 . Thus, ∆V(xk)≤ζT k(AT aPAa−ET aPEa+U1+U2)ζk,∀k≥0, which implies that ζT k(AT aPAa−ET aPEa+U1+U2)ζk<0 =⇒∆V(xk)<0.(20) Therefore, a sufficient condition for the global asymptotic stability of the system is for the left-hand side of (20) to be satisfied for all ζkthat satisfy the sector-boundedness and sloperestriction conditions on φ. These conditions on ζkare introduced via the S-procedure. First, the LMI form of the conditions is found. For the sector boundedness: φ∈Φ[0,ξ] sb ⇐⇒ φk,i[ξ−1 iφk,i −qk,i]≤0, i = 1, . . . , nq,∀k≥0.(21) A useful notation for using the S-procedure with condition (40) is nq X i=1 2τiφk,i[ξ−1 iφk,i −qk,i] = ζT kS1ζk≤0, where T:= diag{τ1, . . . , τnq}<0and S1=ST 1is defined as S1:= 0CTT0 ∗2TX−1+TD +DTT0 ∗ ∗ 0 . 21
Similarly, for the next sampling instance, nq X i=1 2˜τiφk+1,i[ξ−1 iφk+1,i −qk+1,i] = ζT kS2ζk≤0, where ˜ T:= diag{˜τ1,...,˜τnq}<0, and S2=ST 2is defined as S2:= 0 0 ATCT˜ T ∗0BTCT˜ T ∗ ∗ 2˜ TX−1+˜ TD +DT˜ T . For the slope restriction: φ∈Φ[0,µ] sr ⇐⇒ 0≤φk+1,i −φk,i qk+1,i −qk,i ≤µi, i = 1, . . . , nq ⇐⇒ (φk+1,i −φk,i)µ−1 i(φk+1,i −φk,i)−(qk+1,i −qk,i)≤0. (22) A useful notation for using the S-procedure with condition (22) is nq X i=1 2νi(φk+1,i −φk,i)µ−1 i(φk+1,i −φk,i)−(qk+1,i −qk,i)=ζT kS3ζk≤0, where N:= diag{ν1, . . . , νnq}<0, and S3=ST 3is defined as S3:= 0−(CA −C)TN(CA −C)TN ∗2NM−1−N(CB −D)−(CB −D)TN(CB −D)TN−N(2M−1+D) ∗ ∗ 2NM−1+ND +DTN . Finally, applying the S-procedure gives that, if the LMI G:= AT aPAa−ET aPEa+U1+U2− S1−S2−S3≺0is feasible, then ∆V(xk)<0is satisfied ∀k≥0, and the system is globally asymptotically stable. The sufficient condition for the stability of system (13) given by Theorem 1 can be used to find a lower bound for the stability/robustness margin for each nonlinear input pi, which is the maximum value of ξisuch that the sufficient condition for stability is satisfied. Since stability depends on the nqcomponents of the nonlinear input vector p, a simplified case for which to find the robustness margin is that in which all components ξiof the vector ξ are equal, thus finding the maximum sector boundedness restriction imposed on all nonlinear input components such that the sufficient condition for stability is satisfied. 22
Thus, the lower bound on the robustness margin is the value of ξthat solves the optimization: max ξ subject to G≺0 Q, ˜ Q, T, ˜ T, N <0 P=PT<0, P11 0 with all matrices as defined in Theorem 1. Because the condition in Theorem 1 is not a necessary condition, this margin bound can be conservative with respect to the real robustness margin. 2.2 Numerical Results for Stability Analysis The sufficient condition for global asymptotic stability imposed by Theorem 1 is used to find lower bounds on the robustness margin of several examples of systems of the form (13). All examples consider sector bounds and slope restrictions with ξiand µiare taken to be the same for all nonlinearities φi. In addition, µis taken to be linearly dependent on ξ. Theorem 1 deals with systems of the form (13) that can either have D= 0 or D6= 0. All examples in this section have D= 0 to allow comparison with results from other criteria in the literature. Example 1: G(z) = −0.5z+ 0.1 (z2−z+ 0.89) (z+ 0.1), µ = 2ξ. Example 2: A= 0.2948 0 0 0 0 0 0.4568 0 0 0 0 0 0.0226 0 0 0 0 0 0.3801 0 0 0 0 0 −0.3270 , B = −1.1878 0.2341 −2.2023 0.0215 0.9863 −1.0039 −0.5186 −0.9471 0.3274 −0.3744 , C=−1.1859 1.4725 −1.2173 −1.1283 −0.2611 −1.0559 0.0557 −0.0412 −1.3493 0.9535 , D = 02×2 , µ =ξ. 23
Example 3: A= 0.0469 −0.3992 −0.0835 0.3902 −0.5363 −0.2744 0.4378 −1.3576 0.4651 , B = −0.5673 −0.2785 0.1155 −0.0649 −2.1849 −0.5976 , C=0.3587 −1.0802 −0.6802 −1.3833 −1.0677 1.1497 , D = 02×2, µ =ξ. Example 4: A= 0.4030 0 0 0−0.1502 0 0 0 −0.1502 , B = −0.2494 0.2542 −0.2036 , C=0.9894 0.6649 0.4339 , D = 0, µ = 2ξ. Example 5: A= 0.4783 0 0 0 0 0.7871 0 0 0 0 0.7871 1 0 0 0 0.7871 , B = −1.5174 1.2181 0.2496 −0.5181 , C=0.8457 −2.0885 1.2190 0.1683 , D = 0, µ = 2ξ. Example 6: A= 0.5359 0 0 0 0 0 0 0 0 0 0.9417 0 0 0 0 0 0 0 0 0 0.9802 0 0 0 0 0 0 0 0 0 0.5777 0 0 0 0 0 0 0 0 0 −0.1227 0 0 0 0 0 0 0 0 0 −0.0034 0 0 0 0 0 0 0 0 0 −0.5721 0 0 0 0 0 0 0 0 0 0.2870 0 0 0 0 0 0 0 0 0 −0.3599 , 24
L2,24 =BT pCT q˜ QXDqp +˜ QDqp +M−1˜ Q, L2,25 =BT pCT q˜ QXDqw +˜ QDqw, L2,33 =BT wCT q˜ QXCqBw−DT qw ˜ QXDqw, L2,34 =BT wCT q˜ QXDqp, L2,35 =BT wCT q˜ QXDqw, L2,44 =DT qp ˜ QXDqp −˜ QM−1, L2,45 =DT qp ˜ QXDqw, L2,55 =DT qw ˜ QXDqw. The terms zT kzkand wT kwkfrom the inequality (24) can also be expressed in matrix form as zT kzk=ζT kZζk, wT kwk=ζT kWζk, where Z=ZTand W=WTare defined as Z:= CT zCzCT zDzp CT zDzw 0 0 ∗DT zpDzp DT zpDzw 0 0 ∗ ∗ DT zwDzw 0 0 ∗ ∗ ∗ 0 0 ∗ ∗ ∗ ∗ 0 , W := 0 0 0 0 0 ∗0 0 0 0 ∗ ∗ I0 0 ∗ ∗ ∗ 0 0 ∗ ∗ ∗ ∗ 0 . Since ∆V(xk) + zT kzk−γ2wT kwk≤ζT k(AT aPAa−ET aPEa+L1+L2+Z−γ2W)ζk,∀k≥0, then ζT k(AT aPAa−ET aPEa+L1+L2+Z−γ2W)ζk≤0,∀k≥0 =⇒ =⇒∆V(xk) + zT kzk−γ2wT kwk≤0,∀k≥0.(25) Thus, the upper bound on the `2-gain of the system is found as the minimum value of γsuch that the left-hand side of (25) holds, given that ζksafisfies the sector-boundedness and sloperestriction conditions on φ. These conditions on ζkwill be introduced again via the S-procedure. The matrix-form of the conditions is ζT kS1ζk≤0and ζT kS2ζk≤0 for the sector boundedness at the sampling instances kand k+ 1, respectively, and ζT kS3ζk≤0 31
for the slope restriction. The matrices S1=ST 1,S2=ST 2and S3=ST 3are obtained analogously to the proof of Theorem 1 and, for the input-output formulation system (23), are defined as S1:= 0CT qT0 0 0 ∗2TX−1+ +TDqp +DT qpTTDqw 0 0 ∗ ∗ 0 0 0 ∗ ∗ ∗ 0 0 ∗ ∗ ∗ ∗ 0 , S2:= 000 ATCT q˜ T0 ∗0 0 BT pCT q˜ T0 ∗ ∗ 0BT wCT q˜ T0 ∗ ∗ ∗ 2˜ TX−1+ +˜ TDqp +DT qp ˜ T˜ TDqw ∗ ∗ ∗ ∗ 0 , and S3:= 0−(CqA−Cq)TN0 (CqA−Cq)TN0 ∗ 2NM−1− −N(CqBp−Dqp)− −(CqBp−Dqp)TN −N(CqBw−Dqw) −2M−1N+ +(CqBp−Dqp)TN− −NDqp −NDqw ∗ ∗ 0 (CqBw−Dqw)TN0 ∗ ∗ ∗ 2NM−1+ +NDqp +DT qpNNDqw ∗ ∗ ∗ ∗ 0 , where, again, T, ˜ T, N ∈Rnq×nqare diagonal positive semidefinite matrices. Finally, applying the S-procedure gives that an upper bound on the `2-gain of the system is the minimum value of γsuch that the LMI H:= AT aPAa−ET aPEa+L1+L2+Z−γ2W− S1−S2−S340holds. The RMS-gain of the system (23) is defined by sup kwkRMS 6=0 kzkRMS kwkRMS , where the RMS-norm of the input wis defined by kwkRMS ≡v u u tlim sup K→∞ K X k=0 wT kwk and analogously for the output z. 32
The same optimization in Theorem 2 can be used to find the RMS-gain of the system, which follows since inequality (24) also implies that K X k=0 zT kzk≤γ2 K X k=0 wT kwk∀K≥0(since V≥0)⇐⇒ ⇐⇒ 1 K K X k=0 zT kzk≤γ21 K K X k=0 wT kwk∀K≥0 =⇒ =⇒lim sup K→∞ 1 K K X k=0 zT kzk≤γ2lim sup K→∞ 1 K K X k=0 wT kwk⇐⇒ ⇐⇒ kzk2 RMS ≤γ2kwk2 RMS =⇒kzkRMS kwkRMS ≤γ∀kwkRMS 6= 0. Thus, the upper bound on the `2-gain and RMS-gain obtained through Theorem 2 are the same. As with the robustness margin bound derived from Theorem 1, since condition (24) is sufficient but not necessary for the `2-gain and the RMS-gain to be smaller or equal to γ, the results might also present conservatism with respect to the true value of the gain. 3.2 Numerical Results for Gain Analysis The optimization in Theorem 2 is solved for an extended version of Example 2 from Section 2.2, i.e., A= 0.2948 0 0 0 0 0 0.4568 0 0 0 0 0 0.0226 0 0 0 0 0 0.3801 0 0 0 0 0 −0.3270 , Bp= −1.1878 0.2341 −2.2023 0.0215 0.9863 −1.0039 −0.5186 −0.9471 0.3274 −0.3744 , Bw= 1 1 1 1 1 , Cq=−1.1859 1.4725 −1.2173 −1.1283 −0.2611 −1.0559 0.0557 −0.0412 −1.3493 0.9535 , Dqp = 02×2, Dqw =0 0, Cz=10000, Dzp = 0, Dzw = 1, µ =ξ, by using the LMI solver tools from the MATLAB Robust Control Toolbox and by implementing the bisection method to find the lowest value of γ2for which the LMI constraints are feasible, as explained in Section 1.3.2. This procedure gives that the `2-gain, and equivalently an RMS-gain, for the system is less than or equal to 2.42751 for ξ= 0.01 and 2.50171 for ξ= 0.1. 33
The same method was applied to an extended version of Example 3 from Section 2.2, i.e., A= 0.0469 −0.3992 −0.0835 0.3902 −0.5363 −0.2744 0.4378 −1.3576 0.4651 , Bp= −0.5673 −0.2785 0.1155 −0.0649 −2.1849 −0.5976 , Bw= 0.0962 −0.5829 −0.0482 0.4739 −1.1274 1.1238 , Cq=0.3587 −1.0802 −0.6802 −1.3833 −1.0677 1.1497 , Dqp = 02×2, Dqw =0.1562 0.4342 0.5472 0.0356 , Cz=0.9792 0.1112 −0.8091 0.6970 1.3471 −0.0023 , Dzp =0.0010 −0.7238 1.2356 0.2360 , Dzw =0.5474 0.0242 0.2762 0.0486 , µ=ξ. The `2-gain, and RMS-gain, is obtained to be less than or equal to 4.76056 for ξ= 0.01 and 5.41419 for ξ= 0.1. 34
4 State Feedback Controller Design This chapter derives sufficient BMI conditions for a proportional state feedback controller to globally asymptotically stabilize a nominally unstable system of the form (13). The conditions are based on the bounds from the stability analysis and, thus, conservatism is correspondingly reduced. 4.1 Analytical BMI Derivation for State Feedback Controller Design The system controlled by state feedback has the form xk+1 =Axk+Bpk+Buuk= (A+BuK)xk+Bpk qk=Cxk+Dpk uk=Kxk pk=−φ(qk) (26) where u∈Rnuis the input control variable, Bu∈Rn×nu,K∈Rnu×nis the state feedback controller matrix, and the rest of the variables and matrices are as defined for the system (13). Theorem 3: A sufficient condition for the controller matrix Kto globally asymptotically stabilize the system (26) is the existence of a positive definite matrix P=PT∈ R(n+np+nq)×(n+np+nq), with a positive definite submatrix P11 =PT 11 ∈Rn×n, a diagonal positive definite matrix ˜ Q∈Rnq×nq, and diagonal positive semidefinite matrices Q, T, ˜ T, N ∈Rnq×nq that satisfy the bilinear matrix inequality JT=J:= J1J2J3 JT 2−P0 JT 30J4 ≺0,(27) 35
where J1:= −P11 −CTPT 13− −P13C−CTP33C− −CT˜ QXC −P12 −CTPT 23 −P13D− −CTP33D+ +ATCT˜ QXCB+ +KTBT uCT˜ QXCB− −CT˜ QXD +ATCT˜ Q+ +KTBT uCT˜ Q−CT˜ Q− −CTT+ATCTN+ +KTBT uCTN−CTN −ATCTQ−KTBT uCTQ+ +CTQ+ATCT˜ QXD+ +KTBT uCT˜ QXD− −ATCT˜ T−KTBT uCT˜ T− −ATCTN−KTBT uCTN+ +CTN ∗ −P22 −DTPT 23 −P23D− −DTP33D−QM−1+ +BTCT˜ QXCB− −DT˜ QXD +˜ Q(CB −D)+ +(CB −D)T˜ Q−˜ QM−1− −2TX−1−TD −DTT− −2NM−1+N(CB −D)+ +(CB −D)TN M−1Q−(CB −D)TQ+ +BTCT˜ QXD+ +˜ QD +˜ QM−1−BTCT˜ T− −(CB −D)TN+ +2NM−1+ND ∗ ∗ −QD −DTQ−QM−1+ +DT˜ QXD −˜ QM−1− −2˜ TX−1−˜ TD −DT˜ T− −2NM−1−ND −DTN , J2:= ATP11 +KTBT uP11+ +ATCTPT 13+ +KTBT uCTPT 13 ATP12 +KTBT uP12+ +ATCTPT 23+ +KTBT uCTPT 23 ATP13 +KTBT uP13+ +ATCTP33+ +KTBT uCTP33 BTP11 +BTCTPT 13 BTP12 +BTCTPT 23 BTP13 +BTCTP33 PT 12 +DTPT 13 P22 +DTPT 23 P23 +DTP33 , J3:= ATCT˜ QX +KTBT uCT˜ QX 0 0 , J4:= −˜ QX. Proof: Let ˜ A:= A+BuK. Then the stability of system (26) is equivalent to the stability of the system (13) by substituting the original matrix Aby ˜ A. Thus, the controller matrix Kstabilizes the system if ˜ G:= ˜ AT aP˜ Aa−ET aPEa+˜ U1+˜ U2−S1−˜ S2−˜ S3≺0, where all matrices are defined in Theorem 1, and tildes have been used to indicate matrices where matrix Ais present and therefore replaced by ˜ A. A first glance at this inequality shows that the term 36
˜ AT aP˜ Aawill give rise to the submatrix variables in Pbeing pre- & postmultiplied by matrix variable K. To address these higher order terms, the Schur complement lemma is applied as ˜ G≺0 P0⇐⇒ −ET aPEa+˜ U1+˜ U2−S1−˜ S2−˜ S3˜ AT aP P˜ Aa−P≺0. In order for the Schur complement lemma to hold, matrix Pmust be positive definite rather than merely positive semidefinite. This matrix inequality is not yet bilinear, given that the trilinear term ˜ ATCT˜ QXC ˜ Ais in the submatrix ˜ U2,11. Let V:= ˜ ATCT˜ QXC ˜ A0 0 0 0 0 0 0 0 and let ˜ ˜ U2:= ˜ U2−V. The Schur complement lemma is then applied once more to give "−ET aPEa+˜ U1+˜ ˜ U2−S1−˜ S2−˜ S3˜ AT aP P˜ Aa−P#−−V0 0 0 ≺0 ˜ QX 0 ⇐⇒ J≺0. Note that here, again, using the Schur complement lemma imposes that ˜ Qand Xbe positive definite. 4.2 Illustrative Numerical Example of a State Feedback Controller The BMI obtained in Theorem 3 has the advantage that it is an LMI for fixed K. This property is used for obtaining numerical results showing the stabilization of a system through a matrix K. In particular, the numerical results illustate how Theorem 3 holds for an example for which Kis known, which is implemented in LMI software as explained in Section 1.3.2. The system with the form (26) defined by the matrices A= 0.8−0.25 0 1 1 0 0 0 0 0 0.2 0.3 0 0 1 0 , B = 0 0 1 0 , Bu= 1 0 0 0 , Cq=0.8−0.5 0 1 ,(28) with K=0 0 0 0 and with ξ=µ= 2 is not globally asymptotically stable. Unstable behavior of the system can be seen in Figure 1. 37
Figure 1: First 500 iterations of the system (28) without control. The initial values of the four state variables are set to random values between 0 and 104. Some points are labeled to show that the absolute values of the state variables are many (and increasing) orders of magnitude larger than their initial values. Theorem 3 holds for this example since the BMI (27) is infeasible for system (28) with K=0 0 0 0 . The BMI is feasible with the state feedback controller K=0 0 −1−1, implying that the controlled system is globally asymptotically stable. This stable behavior is shown in Figure 2, again as predicted by Theorem 3. 38
(a) (b) Figure 2: (a) First 500 iterations of the system (28) with the state feedback controller K=0 0 −1−1. The initial values of the four state variables are set to random values between 0 and 104. (b) Close-up of Figure (a), which confirms that the values of the state variables converge to zero. 39
5 State Estimator Design In many real systems for which the dynamics are known, their variables cannot be measured. It can be useful in practice to be able to generate a state estimate ˆxk, and the output of a state estimator can be fed to a state feedback controller to collectively produce an output feedback controller. Obtaining an estimator system that converges to the true values of the variables involves analyzing the dynamics of the estimation error variables, to ensure that they absolutely converge to zero. The nonlinearities in the error dynamics system do not have a local slope restriction and, thus, the tight integral bounds obtained in the stability analysis are not applicable to this chapter. New bounds for the Lyapunov function variation integrals are developed based on the sector boundedness of the nonlinearities. 5.1 Analytical BMI Derivation for State Estimator Design The design problem is to find an estimator matrix L∈Rn×nyto estimate the state of the system xk+1 =Axk+Bppk+Buuk qk=Cqxk+Dqppk+Dquuk yk=Cyxk+Dyuuk pk=−φ(qk) (29) through the state estimator system ˆxk+1 =Aˆxk+Bpˆpk+Buuk+L(ˆyk−yk) ˆqk=Cqˆxk+Dqp ˆpk+Dquuk ˆyk=Cyˆxk+Dyuuk ˆpk=−φ(ˆqk) (30) with estimation error dynamics ˆxk+1 −xk+1 = (A+LCy)(ˆxk−xk) + Bp(ˆpk−pk) ˆqk−qk=Cq(ˆxk−xk) + Dqp(ˆpk−pk) ˆyk−yk=Cy(ˆxk−xk) ˆpk−pk=−(φ(ˆqk)−φ(qk)) ≡ −f(ˆqk−qk;qk) (31) where all variables and matrices in the system (29) are as defined in the system (23) with the control variable u∈Rnubeing analogous to the input vector w, and all variables with hats in the system (30) are the estimates of the corresponding variables in the system (29). Considering the dynamics of the error (31) in place of the dynamics of the original system (13) leads to another Lur’e problem where now the nonlinearity f /∈Φ[0,ξ] sb ∩Φ[0,µ] sr . In fact, 40
Figure 3: First 500 iterations of the error dynamics system (41) with zero estimator matrix L. The initial values of the four state variables are set to random values between 0 and 104. Some points are labeled to show that the absolute values of the state variables are many (and increasing) orders of magnitude higher than their initial values. Theorem 4 holds for this example since the BMI (37) is infeasible for system (41) with L=0 0 0 0 T. The BMI is found to be feasible for the estimator matrix L=−1 0 0 0 T, implying that the error dynamics system is globally asymptotically stable. This stable behavior is shown in Figure 4, again as predicted by Theorem 4. 47
(a) (b) Figure 4: (a) First 500 iterations of the system (41) with the estimator matrix L=−1 0 0 0 T. The initial values of the four error state variables are set to random values between 0 and 104. (b) Close-up of (a), which confirms that the values of the error state variables converge to zero. 48
6 Conclusions 6.1 Contributions The main contributions of this thesis are: •Rigorous mathematical formulations are developed for new criteria for the absolute stability, `2and RMS-gain, state feedback controller design, and state estimator design of a discrete-time Lur’e system with sector-bounded, slope-restricted nonlinearities. •Tight bounds are found for the variation in the integrals in the used modified Lur’ePostnikov Lyapunov function, based on the local slope restriction on the system’s nonlinearities. These bounds lead to a significant decrease in the conservatism of the conditions for Lyapunov stability of the system. The reduced conservatism is illustrated by obtaining the results for the lower bound on the robustness margin of six numerical examples and comparing them to the most relevant criteria in the literature. •The same bounds are applied in the `2and RMS-gain analysis and for state feedback controller design, and therefore conservatism is significantly decreased for those problems as well. The upper bound on the gain is obtained for illustrative numerical examples and the behavior of a nominally unstable system stabilized by state feedback is shown. •The state estimator error dynamics are analyzed and the resulting nonlinearities are characterized, leading to a modification of the previously used Lyapunov function. Bounds on the variation of the integrals in the Lyapunov function are derived based on the sector boundedness of the nonlinearities. A sufficient condition for the stability of the error dynamics system is provided and illustrated through a numerical example. •The criteria provided for the Lyapunov stability and gain of the system are based on LMIs and therefore lead to convex optimizations that are highly computationally efficient. The criteria provided for the state feedback controller design and state estimator design are based on BMIs, which can be solved with software that is readily available. Although BMIs are not as efficient for problems of large dimensions, the resulting problems are LMIs for fixed Kor for fixed L. •A software implementation in MATLAB is provided for the set of numerical examples of the four LMI criteria developed in the thesis. The provided scripts can be easily adapted to any other example of a Lur’e system with sector-bounded, slope-restricted nonlinearities. The bisection method scripts used to solve the feasibility optimization problems in the absolute stability and gain analysis can also be easily adapted to different desired accuracies in the results. 49
6.2 Further work Natural continuations of this work might include: •Studying the state feedback controller design for the case with Dqu 6= 0, and attempting to reduce the resulting matrix inequalities to BMIs. Although the addition of this term makes the reduction more complicated, the term appears in some systems including for some types of dynamic neural networks. •Designing an estimator-based output feedback controller for the studied system. This approach would enable the control of systems with unmeasured variables. •Implementing the optimizations for the design of a state feedback controller and a state estimator using a BMI solver and applying to several examples. •Exploring how close iteratively solving the LMI optimizations of each of the variables in the BMI formulations converges to the global optimum, for a range of numerical examples. 50
7 References [1] J. G. VanAntwerp and R. D. Braatz. A tutorial on linear and bilinear matrix inequalities. Journal of Process Control, 10(4):363–385, 2000. [2] S. Boyd, L. El Ghaoui, E. Feron, and V. Balakrishnan. Linear Matrix Inequalities in System and Control Theory, volume 15 of Studies in Applied Mathematics. SIAM, Philadelphia, PA, June 1994. [3] A. I. Lur’e and V. N. Postnikov. On the theory of stability of control systems. Applied Mathematics and Mechanics, 8(3):246–248, 1944. [4] I. W. Sandberg. A frequency-domain condition for the stability of feedback systems containing a single time-varying nonlinear element. The Bell System Technical Journal, 43(4):1601–1608, 1964. [5] R. Brockett. The status of stability theory for deterministic systems. IEEE Transactions on Automatic Control, 11(3):596–606, 1966. [6] W. M. Haddad, J. P. How, S. R. Hall, and D. S. Bernstein. Extensions of mixed-µbounds to monotonic and odd monotonic nonlinearities using absolute stability theory: part I. In Proceedings of the 31st IEEE Conference on Decision and Control, pages 2813–2819 vol.3, 1992. [7] Y. Z. Tsypkin. Fundamentals of the theory of non-linear pluse control systems. IFAC Proceedings Volumes, 1(2):172 – 180, 1963. 2nd International IFAC Congress on Automatic and Remote Control: Theory, Basel, Switzerland, 1963. [8] E. Jury and B. Lee. On the stability of a certain class of nonlinear sampled-data systems. IEEE Transactions on Automatic Control, 9(1):51–61, 1964. [9] V. Kapila and W. M. Haddad. A multivariable extension of the Tsypkin criterion using a Lyapunov-function approach. IEEE Transactions on Automatic Control, 41(1):149–152, 1996. [10] P. Park and S. W. Kim. A revisited Tsypkin criterion for discrete-time nonlinear Lur’e systems with monotonic sector-restrictions. Automatica, 34(11):1417–1420, 1998. [11] S.M. Lee and J. H. Park. Robust stabilization of discrete-time nonlinear Lur’e systems with sector and slope restricted nonlinearities. Applied Mathematics and Computation, 200(1):429–436, 2008. [12] P. Gahinet and A. Nemirovskii. LMI Lab: A package for manipulating and solving LMIs. 1993. [13] G. Balas, R. Cihang, A. Packard, and M. Safonov. Robust Control ToolboxTM 3 user’s guide. January 2012. 51
[14] P. Gahinet, A. Nemirovski, A.J. Laub, and M. Chilali. LMI Control Toolbox user’s guide. The MathWorks, Inc., May 1995. [15] K.-K. K. Kim. Robust control for systems with sector-bounded, slope-restricted, and odd monotonic nonlinearities using linear matrix inequalities. Master’s thesis, University of Illinois at Urbana-Champaign, Illinois, USA, 2009. 52
Appendices A MATLAB Scripts for Chapter 2 To use these scripts, for each example i, run the code findmargin_examplei.m, which calls the scripts lmis_stability.m and bisection.m. findmargin_example1.m 1% G(z) = (−0.5z + 0.1) / (z^3 −z^2 + 0.89z + 0.1z^2 −0.1z + 0.089) = 2% = (−0.5z + 0.1) / (z^3 −0.9z^2 + 0.79z + 0.089) 3% 4% ==> b0 = 0, b1 = 0, b2 = −0.5, b3 = 0.1 5% (a0 = 1), a1 = −0.9, a2 = 0.79, a3 = 0.089 6% 7% Using the Observable Canonical Form 8 9clear variables 10 11 % Previous results for bisection: 12 % xi = 2 feasible, xi = 3 infeasible 13 xi_a = 2; xi_b = 3; 14 15 % Definition of example 16 global A B C D n nq 17 18 A = [ 0.9 1 0; 19 -0.79 0 1; 20 -0.089 0 0 ]; 21 B = [ 0 ; 22 -0.5; 23 0.1 ]; 24 C = [1 0 0]; 25 D = 0; 26 27 % Dimensions of the system 28 n = 3; 29 nq = 1; 30 31 mu_ct = 2; % Linear dependence of mu on xi 32 nq_vec = [1 0]; % Parameter for matrix variables 33 34 %% 35 53
36 margin_example1 = bisection ( @lmis_stability , xi_a , xi_b , mu_ct , nq_vec ) findmargin_example2.m 1clear variables 2 3% Previous results for bisection: 4% xi = 0.1 feasible, xi = 0.5 infeasible 5xi_a = 0.1; xi_b = 0.5; 6 7% Definition of example 8global A B C D n nq 9 10 A = diag ([0.2948 0.4568 0.0226 0.3801 -0.3270]); 11 B = [ -1.1878 0.2341; 12 -2.2023 0.0215; 13 0.9863 -1.0039; 14 -0.5186 -0.9471; 15 0.3274 -0.3744 ]; 16 C = [ -1.1859 1.4725 -1.2173 -1.1283 -0.2611; 17 -1.0559 0.0557 -0.0412 -1.3493 0.9535 ]; 18 D = zeros(2); 19 20 % Dimensions of the system 21 n = 5; 22 nq = 2; 23 24 mu_ct = 1; % Linear dependence of mu on xi 25 nq_vec = [1 0; 1 0]; % Parameter for matrix variables 26 27 %% 28 29 margin_example2 = bisection ( @lmis_stability , xi_a , xi_b , mu_ct , nq_vec ) findmargin_example3.m 1%clear variables 2 3% Previous results for bisection: 4% xi = 0.1 feasible, xi = 0.5 infeasible 5xi_a = 0.1; xi_b = 1; 6 7% Definition of example 54
8global A B C D n nq 9 10 A = [ 0.0469 -0.3992 -0.0835; 11 0.3902 -0.5363 -0.2744; 12 0.4378 -1.3576 0.4651 ]; 13 B = [-0.5673 -0.2785; 14 0.1155 -0.0649; 15 -2.1849 -0.5976 ]; 16 C = [ 0.3587 -1.0802 -0.6802; 17 -1.3833 -1.0677 1.1497 ]; 18 D = zeros(2); 19 20 % Dimensions of the system 21 n = 3; 22 nq = 2; 23 24 mu_ct = 1; % Linear dependence of mu on xi 25 nq_vec = [1 0; 1 0]; % Parameter for matrix variables 26 27 %% 28 29 margin_example3 = bisection ( @lmis_stability , xi_a , xi_b , mu_ct , nq_vec ) findmargin_example4.m 1clear variables 2 3% Previous results for bisection: 4% xi = 1 feasible, xi = 3 infeasible 5xi_a = 1; xi_b = 3; 6 7% Definition of example 8global A B C D n nq 9 10 A = diag ([0.4030 -0.1502 -0.1502]); 11 B = [ -0.2494; 12 0.2542; 13 -0.2036 ]; 14 C = [0.9894 0.6649 0.4339]; 15 D = 0; 16 17 % Dimensions of the system 18 n = 3; 19 nq = 1; 20 55
21 mu_ct = 2; % Linear dependence of mu on xi 22 nq_vec = [1 0]; % Parameter for matrix variables 23 24 %% 25 26 margin_example4 = bisection ( @lmis_stability , xi_a , xi_b , mu_ct , nq_vec ) findmargin_example5.m 1clear variables 2 3% Previous results for bisection: 4% xi = 0.01 feasible, xi = 0.1 infeasible 5xi_a = 0.01; xi_b = 0.1; 6 7% Definition of example 8global A B C D n nq 9 10 A = diag ([0.4783 0.7871 0.7871 0.7871]); A(3 ,4) = 1; 11 B = [ -1.5174; 12 1.2181; 13 0.2496; 14 -0.5181 ]; 15 C = [0.8457 -2.0885 1.2190 0.1683]; 16 D = 0; 17 18 % Dimensions of the system 19 n = 4; 20 nq = 1; 21 22 mu_ct = 2; % Linear dependence of mu on xi 23 nq_vec = [1 0]; % Parameter for matrix variables 24 25 %% 26 27 margin_example5 = bisection ( @lmis_stability , xi_a , xi_b , mu_ct , nq_vec ) findmargin_example6.m 1clear variables 2 3% Previous results for bisection: 4% xi = 0.1 feasible, xi = 1 infeasible 56
39 gain_example2 = bisection_semidef ( @lmis_gain , gamma_a , gamma_b , ... 40 mu_ct , nq_vec ) findgain_example3.m 1clear variables 2 3% Previous results for bisection: 4% gamma = 0.1 infeasible, gamma = 100 feasible 5gamma_a = 0.1; gamma_b = 100; 6 7% Definition of example 8global A Bp Bw Cq Cz Dqp Dqw Dzp Dzw n nq nw 9 10 A = [ 0.0469 -0.3992 -0.0835; 11 0.3902 -0.5363 -0.2744; 12 0.4378 -1.3576 0.4651 ]; 13 Bp = [ -0.5673 -0.2785; 14 0.1155 -0.0649; 15 -2.1849 -0.5976 ]; 16 Bw = [ 0.0962 -0.5829; 17 -0.0482 0.4739; 18 -1.1274 1.1238 ]; 19 Cq = [ 0.3587 -1.0802 -0.6802; 20 -1.3833 -1.0677 1.1497 ]; 21 Dqw = [ 0.1562 0.4342; 22 0.5472 0.0356 ]; 23 Dqp = zeros(2); 24 Cz = [ 0.9792 0.1112 -0.8091; 25 0.6970 1.3471 -0.0023 ]; 26 Dzp = [ 0.0010 -0.7238; 27 1.2356 0.2360 ]; 28 Dzw = [ 0.5474 0.0242; 29 0.2762 0.0486 ]; 30 31 32 % Dimensions of the system 33 n = 3; 34 nq = 2; 35 nw = 2; 36 37 mu_ct = 1; % Linear dependence of mu on xi 38 nq_vec = [1 0; 1 0]; % Parameter for matrix variables 39 40 %% 63
41 42 gain_example3 = bisection_semidef ( @lmis_gain , gamma_a , gamma_b , ... 43 mu_ct , nq_vec ) lmis_gain.m 1function tmin = lmis_gain(gamma, mu_ct , nq_vec ) 2 3poszero = 1e -12; % Numerical zero for nonstrict LMIs 4negzero = -poszero ; % For positive semidefinite LMIs 5 6global A Bp Bw Cq Cz Dqp Dqw Dzp Dzw n nq nw 7 8% Definition of M and X: change xi according to example 9xi = 0.01; 10 % xi = 0.1; 11 mu = mu_ct * xi; 12 13 X = xi * eye (nq ); 14 M = mu * eye (nq ); 15 Xi = inv (X); 16 Mi = inv (M); 17 18 %% 19 setlmis ([]); 20 P11 = lmivar (1, [n 1]); 21 P12 = lmivar (2 , [n nq ]); 22 P13 = lmivar (2 , [n nq ]); 23 P14 = lmivar (2 , [n nw ]); 24 P22 = lmivar (1, [nq 1]); 25 P23 = lmivar (2 , [nq nq ]); 26 P24 = lmivar (2 , [nq nw ]); 27 P33 = lmivar (1, [nq 1]); 28 P34 = lmivar (2 , [nq nw ]); 29 P44 = lmivar (1, [nw 1]); 30 31 Q = lmivar (1, nq_vec ); 32 tQ = lmivar (1, nq_vec ); 33 T = lmivar (1, nq_vec ); 34 tT = lmivar (1, nq_vec ); 35 N = lmivar (1, nq_vec ); 36 37 38 H = newlmi ; 39 64
40 lmiterm ([H 1 1 P11], A’, A); 41 lmiterm ([H 1 1 P13 ], A’, Cq*A, ’s’); 42 lmiterm ([H 1 1 P33 ], A ’*Cq ’, Cq*A); 43 lmiterm ([H 1 1 P11 ], -1, 1); 44 lmiterm ([H 1 1 P13], -1, Cq , ’s’); 45 lmiterm ([H 1 1 P33 ], -Cq ’, Cq ); 46 lmiterm ([H 1 1 tQ], A ’*Cq ’, X*Cq*A); 47 lmiterm ([H 1 1 tQ], -Cq ’, X*Cq ); 48 lmiterm ([H 1 1 0], Cz ’* Cz ); 49 50 lmiterm ([H 1 2 P11 ], A’, Bp ); 51 lmiterm ([H 1 2 P13 ], A’, Cq*Bp ); 52 lmiterm ([H 1 2 -P13], A’*Cq ’, Bp ); 53 lmiterm ([H 1 2 P33 ], A ’*Cq ’, Cq*Bp ); 54 lmiterm ([H 1 2 P12 ], -1, 1); 55 lmiterm ([H 1 2 -P23], -Cq ’, 1); 56 lmiterm ([H 1 2 P13 ], -1, Dqp ); 57 lmiterm ([H 1 2 P33 ], -Cq ’, Dqp ); 58 lmiterm ([H 1 2 tQ], A ’*Cq ’, X*Cq*Bp ); 59 lmiterm ([H 1 2 tQ], Cq ’, X*Dqp ); 60 lmiterm ([H 1 2 tQ], (Cq*A - Cq)’, 1); 61 lmiterm ([H 1 2 0], Cz ’* Dzp ); 62 lmiterm ([H 1 2 T], -Cq ’, 1); 63 lmiterm ([H 1 2 N], (Cq*A - Cq)’, 1); 64 65 lmiterm ([H 1 3 P11 ], A’, Bw ); 66 lmiterm ([H 1 3 P13 ], A’, Cq*Bw ); 67 lmiterm ([H 1 3 -P13], A’*Cq ’, Bw ); 68 lmiterm ([H 1 3 P33 ], A ’*Cq ’, Cq*Bw ); 69 lmiterm ([H 1 3 P13 ], -1, Dqw ); 70 lmiterm ([H 1 3 P33 ], -Cq ’, Dqw ); 71 lmiterm ([H 1 3 P14 ], -1, 1); 72 lmiterm ([H 1 3 P34 ], -Cq ’, 1); 73 lmiterm ([H 1 3 tQ], A ’*Cq ’, X*Cq*Bw ); 74 lmiterm ([H 1 3 tQ], -Cq ’, X*Dqw ); 75 lmiterm ([H 1 3 0], Cz ’* Dzw ); 76 77 lmiterm ([H 1 4 P12 ], A’, 1); 78 lmiterm ([H 1 4 -P23], A ’*Cq ’, 1); 79 lmiterm ([H 1 4 P13 ], A’, Dqp ); 80 lmiterm ([H 1 4 P33 ], A ’*Cq ’, Dqp ); 81 lmiterm ([H 1 4 Q], -(Cq*A - Cq)’, 1); 82 lmiterm ([H 1 4 tQ], A ’*Cq ’, X* Dqw ); 83 lmiterm ([H 1 4 tT], -A ’*Cq ’, 1); 84 lmiterm ([H 1 4 N], -(Cq*A - Cq)’, 1); 85 65
86 lmiterm ([H 1 5 P13 ], A’, Dqw ); 87 lmiterm ([H 1 5 P33 ], A ’*Cq ’, Dqw ); 88 lmiterm ([H 1 5 P14 ], A’, 1); 89 lmiterm ([H 1 5 P34 ], A ’*Cq ’, 1); 90 lmiterm ([H 1 5 tQ], A ’*Cq ’, X* Dqw ); 91 92 lmiterm ([H 2 2 P11 ], Bp ’, Bp ); 93 lmiterm ([H 2 2 P13 ], Bp ’, Cq*Bp , ’s’); 94 lmiterm ([H 2 2 P33 ], Bp ’*Cq ’, Cq*Bp ); 95 lmiterm ([H 2 2 P22 ], -1, 1); 96 lmiterm ([H 2 2 P23 ], -1, Dqp , ’s’); 97 lmiterm ([H 2 2 P33 ], -Dqp ’, Dqp ); 98 lmiterm ([H 2 2 Q], -1, Mi ); 99 lmiterm ([H 2 2 tQ], Bp ’*Cq ’, X*Cq*Bp); 100 lmiterm ([H 2 2 tQ], -Dqp ’, X*Dqp ); 101 lmiterm ([H 2 2 tQ], 1, Cq*Bp - Dqp , ’s’); 102 lmiterm ([H 2 2 tQ], -1, Mi ); 103 lmiterm ([H 2 2 0] , Dzp ’* Dzp ); 104 lmiterm ([H 2 2 T], -2, Xi ); 105 lmiterm ([H 2 2 T], -1, Dqp , ’s’); 106 lmiterm ([H 2 2 N], -2, Mi ); 107 lmiterm ([H 2 2 N], 1, Cq*Bp - Dqp , ’s’); 108 109 110 lmiterm ([H 2 3 P11 ], Bp ’, Bw ); 111 lmiterm ([H 2 3 -P13], Bp ’*Cq ’, Bw ); 112 lmiterm ([H 2 3 P13 ], Bp ’, Cq*Bw ); 113 lmiterm ([H 2 3 P33 ], Bp ’*Cq ’, Cq*Bw ); 114 lmiterm ([H 2 3 P23 ], -1, Dqw ); 115 lmiterm ([H 2 3 P33 ], -Dqp ’, Dqw ); 116 lmiterm ([H 2 3 P24 ], -1, 1); 117 lmiterm ([H 2 3 P34 ], -Dqp ’, 1); 118 lmiterm ([H 2 3 tQ], Bp ’*Cq ’, X*Cq*Bw); 119 lmiterm ([H 2 3 tQ], -Dqp ’, X*Dqw ); 120 lmiterm ([H 2 3 tQ], 1, Cq*Bw - Dqw ); 121 lmiterm ([H 2 3 0] , Dzp ’* Dzw ); 122 lmiterm ([H 2 3 T], -1, Dqw ); 123 lmiterm ([H 2 3 N], 1, Cq*Bw - Dqw ); 124 125 126 lmiterm ([H 2 4 P12 ], Bp ’, 1); 127 lmiterm ([H 2 4 -P23], Bp ’*Cq ’, 1); 128 lmiterm ([H 2 4 P13 ], Bp ’, Dqp ); 129 lmiterm ([H 2 4 P33 ], Bp ’*Cq ’, Dqp ); 130 lmiterm ([H 2 4 Q], Mi , 1); 131 lmiterm ([H 2 4 Q], -(Cq*Bp - Dqp)’, 1); 66
132 lmiterm ([H 2 4 tQ], Bp ’*Cq ’, X*Dqp ); 133 lmiterm ([H 2 4 tQ], 1, Dqp ); 134 lmiterm ([H 2 4 tQ], Mi , 1); 135 lmiterm ([H 2 4 tT], -Bp ’*Cq ’, 1); 136 lmiterm ([H 2 4 N], 2*Mi , 1); 137 lmiterm ([H 2 4 N], -(Cq*Bp - Dqp)’, 1); 138 lmiterm ([H 2 4 N], 1, Dqp ); 139 140 lmiterm ([H 2 5 P13 ], Bp ’, Dqw ); 141 lmiterm ([H 2 5 P33 ], Bp ’*Cq ’, Dqw ); 142 lmiterm ([H 2 5 P14 ], Bp ’, 1); 143 lmiterm ([H 2 5 P34 ], Bp ’*Cq ’, 1); 144 lmiterm ([H 2 5 tQ], Bp ’*Cq ’, X*Dqw ); 145 lmiterm ([H 2 5 tQ], 1, Dqw ); 146 lmiterm ([H 2 5 N], 1, Dqw ); 147 148 lmiterm ([H 3 3 P11 ], Bw ’, Bw ); 149 lmiterm ([H 3 3 P13 ], Bw ’, Cq*Bw , ’s’); 150 lmiterm ([H 3 3 P33 ], Bw ’*Cq ’, Cq*Bw ); 151 lmiterm ([H 3 3 P33 ], -Dqw ’, Dqw ); 152 lmiterm ([H 3 3 P34], -Dqw ’, 1, ’s’); 153 lmiterm ([H 3 3 P44 ], -1, 1); 154 lmiterm ([H 3 3 tQ], Bw ’*Cq ’, X*Cq*Bw); 155 lmiterm ([H 3 3 tQ], -Dqw ’, X*Dqw ); 156 lmiterm ([H 3 3 0] , Dzw ’* Dzw ); 157 lmiterm ([H 3 3 0], -gamma^2); 158 159 lmiterm ([H 3 4 P12 ], Bw ’, 1); 160 lmiterm ([H 3 4 -P23], Bw ’*Cq ’, 1); 161 lmiterm ([H 3 4 P13 ], Bw ’, Dqp ); 162 lmiterm ([H 3 4 P33 ], Bw ’*Cq ’, Dqp ); 163 lmiterm ([H 3 4 Q], -(Cq*Bw - Dqw)’, 1); 164 lmiterm ([H 3 4 tQ], Bw ’*Cq ’, X*Dqp ); 165 lmiterm ([H 3 4 tT], -Bw ’*Cq ’, 1); 166 lmiterm ([H 3 4 N], -(Cq*Bw - Dqw)’, 1); 167 168 lmiterm ([H 3 5 P13 ], Bw ’, Dqw ); 169 lmiterm ([H 3 5 P33 ], Bw ’*Cq ’, Dqw ); 170 lmiterm ([H 3 5 P14 ], Bw ’, 1); 171 lmiterm ([H 3 5 P34 ], Bw ’*Cq ’, 1); 172 lmiterm ([H 3 5 tQ], Bw ’*Cq ’, X*Dqw ); 173 174 lmiterm ([H 4 4 P22 ], 1, 1); 175 lmiterm ([H 4 4 P23 ], 1, Dqp , ’s’); 176 lmiterm ([H 4 4 P33 ], Dqp ’, Dqp ); 177 lmiterm ([H 4 4 Q], -1, Dqp , ’s’); 67
178 lmiterm ([H 4 4 Q], -1, Mi ); 179 lmiterm ([H 4 4 tQ], Dqp ’, X*Dqp ); 180 lmiterm ([H 4 4 tQ], -1, Mi ); 181 lmiterm ([H 4 4 tT], -1, Xi ); 182 lmiterm ([H 4 4 tT], -1, Dqp , ’s’); 183 lmiterm ([H 4 4 N], -2, Mi ); 184 lmiterm ([H 4 4 N], -1, Dqp , ’s’); 185 186 lmiterm ([H 4 5 P23 ], 1, Dqw ); 187 lmiterm ([H 4 5 P33 ], Dqp ’, Dqw ); 188 lmiterm ([H 4 5 P24 ], 1, 1); 189 lmiterm ([H 4 5 P34 ], Dqp ’, 1); 190 lmiterm ([H 4 5 Q], -1, Dqw ); 191 lmiterm ([H 4 5 tQ], Dqp ’, X*Dqw ); 192 lmiterm ([H 4 5 tT], -1, Dqw ); 193 lmiterm ([H 4 5 N], -1, Dqw ); 194 195 lmiterm ([H 5 5 P33 ], Dqw ’, Dqw ); 196 lmiterm ([H 5 5 P34], Dqw ’, 1, ’s’); 197 lmiterm ([H 5 5 P44 ], 1, 1); 198 lmiterm ([H 5 5 tQ], Dqw ’, X*Dqw ); 199 200 % semidefiniteness 201 lmiterm ([-H 1 1 0] , poszero ); 202 lmiterm ([-H 2 2 0] , poszero ); 203 lmiterm ([-H 3 3 0] , poszero ); 204 lmiterm ([-H 4 4 0] , poszero ); 205 lmiterm ([-H 5 5 0] , poszero ); 206 207 208 P11l = newlmi ; 209 lmiterm ([- P11l 1 1 P11], 1, 1); 210 211 Pl = newlmi ; 212 lmiterm ([- Pl 1 1 P11], 1, 1); 213 lmiterm ([- Pl 1 2 P12], 1, 1); 214 lmiterm ([- Pl 1 3 P13], 1, 1); 215 lmiterm ([- Pl 1 4 P14], 1, 1); 216 lmiterm ([- Pl 2 2 P22], 1, 1); 217 lmiterm ([- Pl 2 3 P23], 1, 1); 218 lmiterm ([- Pl 2 4 P24], 1, 1); 219 lmiterm ([- Pl 3 3 P33], 1, 1); 220 lmiterm ([- Pl 3 4 P34], 1, 1); 221 lmiterm ([- Pl 4 4 P44], 1, 1); 222 % semidefiniteness 223 lmiterm ([ Pl 1 1 0], negzero ); 68
224 lmiterm ([ Pl 2 2 0], negzero ); 225 lmiterm ([ Pl 3 3 0], negzero ); 226 lmiterm ([ Pl 4 4 0], negzero ); 227 228 Ql = newlmi ; 229 lmiterm ([- Ql 1 1 Q], 1, 1); 230 % semidefiniteness 231 lmiterm ([ Ql 1 1 0], negzero ); 232 tQl = newlmi ; 233 lmiterm ([- tQl 1 1 tQ], 1, 1); 234 % semidefiniteness 235 lmiterm ([ tQl 1 1 0], negzero ); 236 Tl = newlmi ; 237 lmiterm ([- Tl 1 1 T], 1, 1); 238 % semidefiniteness 239 lmiterm ([ Tl 1 1 0], negzero ); 240 tTl = newlmi ; 241 lmiterm ([- tTl 1 1 tT], 1, 1); 242 % semidefiniteness 243 lmiterm ([ tTl 1 1 0], negzero ); 244 Nl = newlmi ; 245 lmiterm ([- Nl 1 1 N], 1, 1); 246 % semidefiniteness 247 lmiterm ([ Nl 1 1 0], negzero ); 248 249 lmisys = getlmis ; 250 251 %% 252 target = []; options = zeros(1 ,5); 253 options (2) = 300; 254 options (3) = 10; 255 [tmin ,xfeas ] = feasp (lmisys , options , target ); 256 end bisection_semidef.m 1function p = bisection_semidef (f, a, b, mu_ct , nq_vec ) 2 3% The goal is not to find a zero for tmin, 4% but a value at which it is smaller than TOL. 5TOL = 1e -10; 6% The result of feasp then yields "may be feasible 7% but not strictly feasible". 8 9dif_tol = 1e -6; % tolerance for iteration difference 69
10 max_it = 100; % max number of iterations 11 12 fa = f(a, mu_ct , nq_vec ); 13 fb = f(b, mu_ct , nq_vec ); 14 15 if (fa >= TOL && fb >= TOL) || (fa < TOL && fb < TOL) 16 disp (’Wrong choice ’) 17 else 18 it = 1; 19 dif_p = 1; 20 p = (a + b )/2; 21 fp = f(p, mu_ct , nq_vec ); 22 while (it < max_it && dif_p >= dif_tol ) 23 it = it + 1; 24 fa = f(a, mu_ct , nq_vec ); 25 if (fa < TOL && fp >= TOL) || (fa >= TOL && fp < TOL ) 26 b = p; 27 else 28 a = p; 29 end 30 p_old = p; 31 p = (a + b )/2; 32 p_new = p; 33 dif_p = abs( p_new - p_old ); 34 fp = f(p, mu_ct , nq_vec ); 35 if it == max_it ; ... 36 disp (’Stopped because bisection it = max_it ’); end 37 if dif_p < dif_tol ; ... 38 disp (’Stopped because bisection dif_p < dif_tol ’); end 39 end 40 if fp > TOL; p = b; end 41 % Note: this gives a feasible bound 42 end 43 end 70
C MATLAB Scripts for Chapter 4 The code controller_feasibility.m solves the feasibility problem for the given example. The code controller_simulation.m plots the behavior of the example. The controller matrix Khas to be changed within the definition of the example in both scripts. controller_feasibility.m 1clear variables 2 3poszero = 1e -12; % Numerical zero for nonstrict LMIs 4negzero = -poszero ; % For positive semidefinite LMIs 5 6% Dimensions of the system 7n = 4; 8nq = 1; 9nu = 1; 10 11 % Definition of example 12 A = [0.8 -0.25 0 1; 13 1 0 0 0; 14 0 0 0.2 0.3; 15 0 0 1 0 ]; 16 B = [0; 17 0; 18 1; 19 0 ]; 20 Bu = [1; 21 0; 22 0; 23 0]; 24 C = [0.8 -0.5 0 1]; 25 D = 0; 26 K = zeros(1, n); 27 % K = [0 0 −1−1]; % Select K for example 28 29 mu_ct = 1; % Linear dependence of mu on xi 30 nq_vec = [1 0]; % Parameter for matrix variables 31 32 % Definition of M and X: change xi according to example 33 xi = 2; 34 mu = mu_ct * xi; 35 36 X = xi * eye (nq ); 37 M = mu * eye (nq ); 71
38 Xi = inv (X); 39 Mi = inv (M); 40 41 %% 42 setlmis ([]); 43 P11 = lmivar (1, [n 1]); 44 P12 = lmivar (2 , [n nq ]); 45 P13 = lmivar (2 , [n nq ]); 46 P22 = lmivar (1, [nq 1]); 47 P23 = lmivar (2 , [nq nq ]); 48 P33 = lmivar (1, [nq 1]); 49 50 Q = lmivar (1, nq_vec ); 51 tQ = lmivar (1, nq_vec ); 52 T = lmivar (1, nq_vec ); 53 tT = lmivar (1, nq_vec ); 54 N = lmivar (1, nq_vec ); 55 56 57 J = newlmi ; 58 59 lmiterm ([J 1 1 P11 ], -1, 1); 60 lmiterm ([J 1 1 P13 ], -1, C, ’s’); 61 lmiterm ([J 1 1 P33 ], -C’, C); 62 lmiterm ([J 1 1 tQ], -C’, X*C); 63 64 lmiterm ([J 1 2 P12 ], -1, 1); 65 lmiterm ([J 1 2 -P23], -C’, 1); 66 lmiterm ([J 1 2 P13], -1, D); 67 lmiterm ([J 1 2 P33 ], -C’, D); 68 lmiterm ([J 1 2 tQ], A ’*C’, X*C*B); 69 lmiterm ([J 1 2 tQ], K ’*Bu ’*C’, X*C*B); 70 lmiterm ([J 1 2 tQ], -C’, X*D); 71 lmiterm ([J 1 2 tQ], A ’*C’, 1); 72 lmiterm ([J 1 2 tQ], K ’*Bu ’*C’, 1); 73 lmiterm ([J 1 2 tQ], -C’, 1); 74 lmiterm ([J 1 2 T], -C’, 1); 75 lmiterm ([J 1 2 N], A’*C’, 1); 76 lmiterm ([J 1 2 N], K’*Bu ’*C’, 1); 77 lmiterm ([J 1 2 N], -C’, 1); 78 79 lmiterm ([J 1 3 Q], -A ’*C’, 1); 80 lmiterm ([J 1 3 Q], -K ’*Bu ’*C’, 1); 81 lmiterm ([J 1 3 Q], C’, 1); 82 lmiterm ([J 1 3 tQ], A ’*C’, X*D); 83 lmiterm ([J 1 3 tQ], K ’*Bu ’*C’, X*D); 72
38 %% 39 setlmis ([]); 40 P11 = lmivar (1, [n 1]); 41 P12 = lmivar (2 , [n nq ]); 42 P13 = lmivar (2 , [n nq ]); 43 P22 = lmivar (1, [nq 1]); 44 P23 = lmivar (2 , [nq nq ]); 45 P33 = lmivar (1, [nq 1]); 46 47 Q = lmivar (1, nq_vec ); 48 T = lmivar (1, nq_vec ); 49 tT = lmivar (1, nq_vec ); 50 51 52 R = newlmi ; 53 54 lmiterm ([R 1 1 P11 ], -1, 1); 55 lmiterm ([R 1 1 P13], -1, Cq , ’s’); 56 lmiterm ([R 1 1 P33 ], -Cq ’, Cq ); 57 58 lmiterm ([R 1 2 P12 ], -1, 1); 59 lmiterm ([R 1 2 -P23], -Cq ’, 1); 60 lmiterm ([R 1 2 P13 ], -1, Dqp ); 61 lmiterm ([R 1 2 P33 ], -Cq ’, Dqp ); 62 lmiterm ([R 1 2 Q], A’*Cq ’, M*Cq*Bp); 63 lmiterm ([R 1 2 Q], Cy ’*L’*Cq ’, M*Cq*Bp ); 64 lmiterm ([R 1 2 T], -Cq ’, 1); 65 66 lmiterm ([R 1 3 Q], A’*Cq ’, M*Dqp ); 67 lmiterm ([R 1 3 Q], Cy ’*L’*Cq ’, M* Dqp ); 68 lmiterm ([R 1 3 tT], -A ’*Cq ’, 1); 69 lmiterm ([R 1 3 tT], -Cy ’*L’*Cq ’, 1); 70 71 lmiterm ([R 1 4 P11 ], A’, 1); 72 lmiterm ([R 1 4 P11 ], Cy ’*L’, 1); 73 lmiterm ([R 1 4 -P13], A ’*Cq ’, 1); 74 lmiterm ([R 1 4 -P13], Cy ’*L’*Cq ’, 1); 75 76 lmiterm ([R 1 5 P12 ], A’, 1); 77 lmiterm ([R 1 5 P12 ], Cy ’*L’, 1); 78 lmiterm ([R 1 5 -P23], A ’*Cq ’, 1); 79 lmiterm ([R 1 5 -P23], Cy ’*L’*Cq ’, 1); 80 81 lmiterm ([R 1 6 P13 ], A’, 1); 82 lmiterm ([R 1 6 P13 ], Cy ’*L’, 1); 83 lmiterm ([R 1 6 P33 ], A ’*Cq ’, 1); 79
84 lmiterm ([R 1 6 P33 ], Cy ’*L’*Cq ’, 1); 85 86 lmiterm ([R 1 7 Q], A’*Cq ’, M); 87 lmiterm ([R 1 7 Q], Cy ’*L’*Cq ’, M); 88 89 lmiterm ([R 2 2 P22 ], -1, 1); 90 lmiterm ([R 2 2 P23 ], -1, Dqp , ’s’); 91 lmiterm ([R 2 2 P33 ], -Cq ’, Dqp ); 92 lmiterm ([R 2 2 Q], Bp ’*Cq ’, M*Cq*Bp ); 93 lmiterm ([R 2 2 T], -2, Mi ); 94 lmiterm ([R 2 2 T], -1, Dqp , ’s’); 95 96 lmiterm ([R 2 3 Q], Bp ’*Cq ’, M*Dqp ); 97 lmiterm ([R 2 3 tT], -Bp ’*Cq ’, 1); 98 99 lmiterm ([R 2 4 P11 ], Bp ’, 1); 100 lmiterm ([R 2 4 -P13], Bp ’*Cq ’, 1); 101 102 lmiterm ([R 2 5 P12 ], Bp ’, 1); 103 lmiterm ([R 2 5 -P23], Bp ’*Cq ’, 1); 104 105 lmiterm ([R 2 6 P13 ], Bp ’, 1); 106 lmiterm ([R 2 6 P33 ], Bp ’*Cq ’, 1); 107 108 lmiterm ([R 3 3 Q], Dqp ’, M* Dqp ); 109 lmiterm ([R 3 3 tT], -2, Mi ); 110 lmiterm ([R 3 3 tT], -1, Dqp , ’s’); 111 112 lmiterm ([R 3 4 -P12], 1, 1); 113 lmiterm ([R 3 4 -P13], Dqp ’, 1); 114 115 lmiterm ([R 3 5 P22 ], 1, 1); 116 lmiterm ([R 3 5 -P23], Dqp ’, 1); 117 118 lmiterm ([R 3 6 P23 ], 1, 1); 119 lmiterm ([R 3 6 P33 ], Dqp ’, 1); 120 121 lmiterm ([R 4 4 P11 ], -1, 1); 122 123 lmiterm ([R 4 5 P12 ], -1, 1); 124 125 lmiterm ([R 4 6 P13 ], -1, 1); 126 127 lmiterm ([R 5 5 P22 ], -1, 1); 128 129 lmiterm ([R 5 6 P23 ], -1, 1); 80
130 131 lmiterm ([R 6 6 P33 ], -1, 1); 132 133 lmiterm ([R 7 7 Q], -1, M); 134 135 136 P11l = newlmi ; 137 lmiterm ([- P11l 1 1 P11], 1, 1); 138 139 Pl = newlmi ; 140 lmiterm ([- Pl 1 1 P11], 1, 1); 141 lmiterm ([- Pl 1 2 P12], 1, 1); 142 lmiterm ([- Pl 1 3 P13], 1, 1); 143 lmiterm ([- Pl 2 2 P22], 1, 1); 144 lmiterm ([- Pl 2 3 P23], 1, 1); 145 lmiterm ([- Pl 3 3 P33], 1, 1); 146 147 Ql = newlmi ; 148 lmiterm ([- Ql 1 1 Q], 1, 1); 149 Tl = newlmi ; 150 lmiterm ([- Tl 1 1 T], 1, 1); 151 % semidefiniteness 152 lmiterm ([ Tl 1 1 0], negzero ); 153 tTl = newlmi ; 154 lmiterm ([- tTl 1 1 tT], 1, 1); 155 % semidefiniteness 156 lmiterm ([ tTl 1 1 0], negzero ); 157 158 lmisys = getlmis ; 159 160 %% 161 target = []; options = zeros(1 ,5); 162 options (2) = 300; 163 options (3) = 10; 164 feasp ( lmisys , options , target ) estimator_simulation.m 1clear variables 2 3% Dimensions of the system 4n = 4; 5nq = 1; 6nu = 1; 7ny = 1; 81
8 9% Assign value of xi 10 constant = 2; 11 12 % Initial values of state variables 13 x = 1e4* rand (n ,1); 14 15 % Number of simulation iterations 16 n_it = 499; 17 18 % Definition of example 19 A = [0.8 -0.25 0 1; 20 1 0 0 0; 21 0 0 0.2 0.3; 22 0 0 1 0 ]; 23 Bp = [0; 24 0; 25 1; 26 0 ]; 27 Cy = [0 0 1 1]; 28 Cq = [0.8 -0.5 0 1]; 29 Dqp = zeros(1); 30 L = zeros(1, n)’; 31 % L = [−1 0 0 0]’; % Select L for example 32 33 q = Cq*x; 34 p = - constant * q; 35 36 for i = 1: n_it 37 x_next = (A + L*Cy )*x(:, end) + Bp*p(:, end ); 38 q_next = Cq* x_next ; 39 p_next = -constant * q_next ; 40 41 x = [x x_next ]; 42 q = [q q_next ]; 43 p = [p p_next ]; 44 end 45 %% 46 figure (1); 47 plot (x (1 ,:)); hold on; 48 plot (x (2 ,:)); plot(x(3 ,:)); plot (x(4 ,:)); hold off; 49 set (gca,’fontsize’, 14); 50 legend (’ x _ {k ,1} - x_{k ,1} ’,’ x _ {k ,2} - x_{k ,2} ’,... 51 ’ x _ {k ,3} - x_{k ,3} ’,’ x _ {k ,4} - x_{k ,4} ’); 52 title(’Trajectory of error state variables ’); 53 xlabel (’k’); ylabel (’ x _ k - x_k ’); 82
54 55 figure (2); 56 plot (x (1 ,:)); hold on; 57 plot (x (2 ,:)); plot(x(3 ,:)); plot (x(4 ,:)); hold off; 58 set (gca,’fontsize’, 14); 59 legend (’ x _ {k ,1} - x_{k ,1} ’,’ x _ {k ,2} - x_{k ,2} ’,... 60 ’ x _ {k ,3} - x_{k ,3} ’,’ x _ {k ,4} - x_{k ,4} ’); 61 ylim ([ -1e-4 1e -4]); 62 title(’Trajectory of error state variables ’); 63 xlabel (’k’); ylabel (’ x _ k - x_k ’); 83