# Reduced-Order Algorithm for Eigenvalue Assignment of Singularly Perturbed Linear Systems

Authors: Heonjong Yoo; Zoran Gajic; Kyeong-Hwan Lee
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

In this paper, we present an algorithm for eigenvalue assignment of linear singularly perturbed systems in terms of reduced-order slow and fast subproblem matrices. No similar algorithm exists in the literature. First, we present an algorithm for the recursive solution of the singularly perturbed algebraic Sylvester equation used for eigenvalue assignment. Due to the presence of a small singular perturbation parameter that indicates separation of the system variables into slow and fast, the corresponding algebraic Sylvester equation is numerically ill-conditioned. The proposed method for the recursive reduced-order solution of the algebraic Sylvester equations removes ill-conditioning and iteratively obtains the solution in terms of four reduced-order numerically well-conditioned algebraic Sylvester equations corresponding to slow and fast variables. The convergence rate of the proposed algorithm isOε, where ε is a small positive singular perturbation parameter.


---

## Body

## 1. Introduction

The classical method for numerical solution of the Sylvester algebraic equation dates back to reference [1]. Solving the Sylvester algebraic equation numerically is not a simple task [2, 3]. Namely, it was stated in [2, 3] that the algorithm of [1] cannot produce a highly accurate solution. Another method for solving the large-scale Sylvester equations is introduced in [4]. In [4], the authors have shown that researchers have developed some methods for the solution of large-scale Sylvester equations [5–7]. On the other hand, several research studies for solving the Sylvester algebraic equation have gained attention in engineering problem [8–10]. In the image-processing problem, the presentation developed in [8] has shown that the image fusion method applied to a large-scale image facilitates reduction of computational complexity based on the explicit solution of large-scale Sylvester equations. Furthermore, many problems of control theory such as regulator problem [9] and particle swarm theory [10] lead to a Sylvester equation.The aim of our developed algorithm is to solve a large-scale Sylvester equation in order to overcome the numerical ill-conditioning problem of singularly perturbed systems presented in [11]. This leads to reduced-order regular algebraic Sylvester equations [12], combined with the techniques presented in [13, 14] which solves the eigenvalue assignment problem for singularly perturbed linear systems.The general Sylvester equation is defined as(1)TA+MT+N=0.Its unique solutionT exists under the assumption that matrices A and −M have no eigenvalues in common [13].Assumption 1.
λA≠λ−M≠−λM.
Without loss of generality, we will consider the Sylvester equation encountered in the control system design of linear systems:(2)dxtdt=Axt+But,yt=Cxt,where xt∈Rn is the state vector, ut∈Rm is the control input vector, yt∈Rp is the vector of system measurements, and A, B, and C are constant matrices. Forms of the Sylvester algebraic equations that appear in the observer and controller designs are given by(3)ToA+MTo−K¯C=0,K¯=ToK,ATc+TcN−BF¯=0,F¯=FTc,where K stands for the observer feedback gain and F is the system feedback gain. These Sylvester equations were studied in [15]. The system-observer configuration has slow and fast modes since the observer must be much faster than the system, and hence, it represents implicitly a singularly perturbed system. Note that the main difficulty in the numerical solution of algebraic Sylvester equation (3) will come from the fact that the system to be studied has a singularly perturbed structure and not from the system-observer implicit singularly perturbed structure since the controller and the observer eigenvalue assignment problems are done independently using the separation principle. The main difficulty comes from the fact that since the system has slow and fast modes, then the observer must contain very fast modes, which leads to numerical ill-conditioning.


### 1.1. Problem Statement

In this section, we study the Sylvester algebraic equation corresponding to singularly perturbed systems defined by [11] (Chapter 2):(4)x˙1t=A1x1t+A2x2t+B1ut,x1t0=x10,εx˙2t=A3x1t+A4x2t+B2ut,x2t0=x20,yt=C1x1t+C2x2t,where x1t∈Rn1andx2t∈Rn2,n1+n2=n, are, respectively, slow and fast state variables and ε is a small positive singular perturbation parameter. Equation (4) can be obtained from (2) assuming that (2) has eigenvalues clustered into two groups: slow ones closer to the imaginary axis and fast ones farther from the imaginary axis (singularly perturbed structure). In such a case, a similarity transformation converts (2) into (4). The following is the standard assumption used in theory of singular perturbation [11] (Chapter 2).Assumption 2.
The matrixA4 is nonsingular.
We study, without loss of generality, a variant of observer design algebraic Sylvester equation (3) given by(5)TA−AdesT=K¯C,K¯=TK,K=K11εK2.
The matrixAdes is the matrix with the observer-desired closed-loop eigenvalues. The standard observer design assumptions are needed [13].Assumption 3.
.The pair (A, C) is observable, and the pair Ades,K¯ is controllable.
The general existence condition given in Assumption1 and specialized to (5) leads to the following assumption.Assumption 4.
λA∩λAdes=∅.
Having found an invertible solution of (5), then the observer gain is given by K=T−1K¯. Note that Assumption 3 for single-input single-output systems is both sufficient and necessary condition for the existence of an invertible solution of (5). For multi-input multi-output systems, it is only a necessary condition [13], so that a repetitive design algorithm has to be performed until an invertible solution T is obtained (see Section 5). The matrices in (4) and (5) are partitioned as(6)A=A1A21εA31εA4,Ades=As001εAf,C=C1C2,where Ades contains the desired observer closed-loop eigenvalues, that is, λAdes=λA−KC. Note that they are placed far to the left in the complex plane to make the observer asymptotically stable and much faster than the closed-loop system. We have found that the following scaling is appropriate for the solution matrix T:(7)T=T1εT2εT3εT4,which is consistent with the structures of matrices defined in (5) and (6). Namely, the right-hand side of (5) is(8)K¯C=T1K1C1+T2K2C1T1K1C2+T2K2C2εT3K1C1+T4K2C1εT3K1C2+T4K2C2=−Q1Q2Q3Q4=O1O1O1O1.
With the scaling chosen in (7), the left-hand side terms of (5), that is, TA and AdesT, are also both O1. Due to the structure of matrices A and Ades, the singularly perturbed algebraic Sylvester equation defined in (5) is numerically ill-conditioned. To overcome numerical ill-conditioning, we propose a new recursive algorithm for solving (5) in terms of reduced-order well-defined algebraic Sylvester equations. The dual version of (5) used for the system controller design is given by(9)ATc−TcAdesc=BF¯,F¯=FTcB=B11εB2,Adesc=Asc001εAfc.
Multiplying matrices in (9), we can get four algebraic equations for the partitioned matrix Tc as functions of the small parameter ε. It can be found from these equations [16, 17] that the structure of Tc is given by(10)Tc=T1cT2cT3c1εT4c.
Algebraic Sylvester equations (9)-(10) will be solved numerically in terms of reduced-order numerically well-conditioned algebraic Sylvester equations under the standard controller design assumptions [13].Assumption 5.
The pair (A, B) is controllable and the pair Adesc,F¯ is observable.
Moreover, the existence of a unique solution of (9) requires the assumption dual to Assumption 4.Assumption 6.
λA∩λAdesc=∅.


### 1.2. Parallel Algorithm for the Observer Sylvester Equation

The partitioned form of the Sylvester equation given in (5) subject to (6)–(8) is given by(11)T1A1+T2A3−AsT1+Q1=0,T1A2+T2A4−ϵAsT2+Q2=0,εT3A1+T4A3−AfT3+Q3=0,εT3A2+T4A4−AfT4+Q4=0.It can be seen from (11) that this system has two independent sets of linear algebraic equations: the first one for T1 and T2 and the second one for T3 and T4.Settingε=0 in (11), the algebraic equations for zeroth-order approximations of solutions T10, T20, T30, and T40 are obtained as(12)T10A1+T20A3−AsT10+Q1=0,(13)T10A2+T20A4+Q2=0,(14)T40A3−AfT30+Q3=0,(15)T40A4−AfT40+Q4=0.These equations can be solved independently as follows. Unique solutionT40 can be obtained from Sylvester equation (15) under the following assumption.Assumption 7.
λA4≠λAf.
SinceAf, defined in (6), is chosen by the designer as an asymptotically stable matrix, this assumption is easily satisfied. Having obtained T40, from (13) and (14), we can obtain T30andT20 independently as(16)T30=Af−1T40A3+Q3,(17)T20=−Q2+T10A2A4−1.
Substituting (17) into (12) results in(18)T10A0−AsT10+Q0=0,(19)A0=A1−A2A4−1A3,Q0=Q1−Q2A4−1A3.
The unique solutionT10 of algebraic Sylvester equation (18) exists under the following assumption.Assumption 8.
λA0≠λAs.
The solutionsT10,T20,T30, and T40 are Oε close to the exact solutions, that is,(20)T1=T10+εE1,T2=T20+εE2,T3=T30+εE3,T4=T40+εE4.
Now, we show that the values forEi,i=1,2,3,4, can be obtained by running iterations on independent linear reduced-order algebraic equations. Subtracting (12)–(15) from (11) and using (20), we obtain the error equations (after some algebra) in the following form:(21)E1A1−AsE1=−E2A3,(22)E2A4−εAsE2=−E1A2+AsT20,(23)εE3A1−AfE3=−T30A1−E4A3,(24)E4A4−AfE4=−T30A2−εE3A2.
The error equations can be solved iteratively using the fixed-point algorithm, in which the cross-coupling terms multiplied byε are delayed by one iteration. This idea has been used in several algorithms that involve small parameters.Algorithm 1.
(25)E1i+1A0−AsE1i+1=−AsT20A4−1A3−εAsE2iA4−1A3,E2i+1A4−εAsE2i+1=AsT20−E1i+1A2,εE3i+1A1−AfE3i+1=−T30A1−E4iA3,E4i+1A4−AfE4i+1=−T30A2−εE3iA2,E20=0,E30=Af−1T30A1+E40A3,E40A4−AfE40+T30A2=0.
We first solve (22) as(26)E2i+1=AsT20−E1i+1A2+εAsE2iA4−1.
Substituting (26) into (21) gives(27)E1i+1A0−AsE1i+1=−AsT20A4−1A3−εAsE2iA4−1A3.
Equations (26) and (27) have nice forms since the quantity E2 is multiplied by a small parameter ε. Similarly, equations for E3 and E4 can be iteratively solved as(28)εE3i+1A1−AfE3i+1=−T30A1−E4iA3,E4i+1A4−AfE4i+1=−T30A2−εE3iA2.
The following theorem presents the main feature of Algorithm1.Theorem 1.
Under Assumptions7 and 8, Algorithm 1 converges to the exact solution E with the rate of convergence of Oε. The convergence is obtained for sufficiently small values of ε that makes the radius of convergence smaller than 1 in each iteration leading to a contraction mapping. Hence, after i iterations, the solution T is obtained with the accuracy of Oεi, that is,(29)Tji=Tj0+εEji+Oεi,j=1,2,3,4;i=1,2,….Proof of Theorem.
1.
Fori=0, (27) implies(30)E11A0−AsE11=−AsT20A4−1A3−εAsE20A4−1A3=−AsT20A4−1A3.
Note thatE20=0. For i=1, (27) produces(31)E12A0−AsE12=−AsT20A4−1A3−εAsE21A4−1A3.
Subtracting (30) from (31), we have(32)E12−E11A0−AsE12−E11=−εAsE21A4−1A3=Oε.
At this point, we conclude that(33)E12−E11=Oε.
In a similar way, we can write the relationship betweenE13 and E12 as(34)E13−E12A0−AsE13−E12=εAsE21A4−1A3−εAsE22A4−1A3=Oε,which implies that(35)E13−E12=Oε.
Continuing the same procedure, we obtain(36)E1i+1−E1i=Oε.
Now, we work withE2 using (26). For i=0, we have(37)E21=AsT20−E11A2+εAsE20A4−1⇒E21=O1.
Fori=1,(38)E22=AsT20−E12A2+εAsE21A4−1
UsingE20=0 and the result in (33), we get(39)E22−E21=E12−E11A2A4−1+Oε⇒E22−E21=Oε.
Considering (26) for i=2 and using (39), we obtain(40)E23−E22=Oε.
If we keep repeating this process, we conclude that(41)E2i+1−E2i=Oε.
In addition, we have the following relationships. Subtracting (21)-(22) from (25) for i=0 produces(42)E11−E1A0−AsE11−E1=−εAsE20−E2A4−1A3=Oε⇒E11−E1=Oε,E21−E2A4−εAsE21−E2=−E11−E1A2=Oε⇒E21−E2=Oε.
Fori=1, (21)-(22) and (25) produce(43)E12−E1A0−AsE12−E1=−εAsE21−E2A4−1A3=Oε2⇒E12−E1=Oε2,E22−E2A4−εAsE22−E2=−E12−E1A2=Oε2⇒E22−E2=Oε2.
Continuing the same procedure, we have(44)E1i−E1=Oεi,i=1,2,3,…,(45)E2i−E2=Oεi,i=1,2,3,….
Similar procedures applied to (23)-(24) produces(46)E4i+1−E4i=Oε,E4i−E4=Oεi,i=1,2,3,…,(47)E3i+1−E3i=Oε,E3i−E3=Oεi,i=1,2,3,….
Results established in (44)–(47) can be summarized in(48)Eji−Ej=Oεi,j=1,2,3,4;i=1,2,…,which completes the proof of the stated theorem.


## 1.1. Problem Statement

In this section, we study the Sylvester algebraic equation corresponding to singularly perturbed systems defined by [11] (Chapter 2):(4)x˙1t=A1x1t+A2x2t+B1ut,x1t0=x10,εx˙2t=A3x1t+A4x2t+B2ut,x2t0=x20,yt=C1x1t+C2x2t,where x1t∈Rn1andx2t∈Rn2,n1+n2=n, are, respectively, slow and fast state variables and ε is a small positive singular perturbation parameter. Equation (4) can be obtained from (2) assuming that (2) has eigenvalues clustered into two groups: slow ones closer to the imaginary axis and fast ones farther from the imaginary axis (singularly perturbed structure). In such a case, a similarity transformation converts (2) into (4). The following is the standard assumption used in theory of singular perturbation [11] (Chapter 2).Assumption 2.
The matrixA4 is nonsingular.
We study, without loss of generality, a variant of observer design algebraic Sylvester equation (3) given by(5)TA−AdesT=K¯C,K¯=TK,K=K11εK2.
The matrixAdes is the matrix with the observer-desired closed-loop eigenvalues. The standard observer design assumptions are needed [13].Assumption 3.
.The pair (A, C) is observable, and the pair Ades,K¯ is controllable.
The general existence condition given in Assumption1 and specialized to (5) leads to the following assumption.Assumption 4.
λA∩λAdes=∅.
Having found an invertible solution of (5), then the observer gain is given by K=T−1K¯. Note that Assumption 3 for single-input single-output systems is both sufficient and necessary condition for the existence of an invertible solution of (5). For multi-input multi-output systems, it is only a necessary condition [13], so that a repetitive design algorithm has to be performed until an invertible solution T is obtained (see Section 5). The matrices in (4) and (5) are partitioned as(6)A=A1A21εA31εA4,Ades=As001εAf,C=C1C2,where Ades contains the desired observer closed-loop eigenvalues, that is, λAdes=λA−KC. Note that they are placed far to the left in the complex plane to make the observer asymptotically stable and much faster than the closed-loop system. We have found that the following scaling is appropriate for the solution matrix T:(7)T=T1εT2εT3εT4,which is consistent with the structures of matrices defined in (5) and (6). Namely, the right-hand side of (5) is(8)K¯C=T1K1C1+T2K2C1T1K1C2+T2K2C2εT3K1C1+T4K2C1εT3K1C2+T4K2C2=−Q1Q2Q3Q4=O1O1O1O1.
With the scaling chosen in (7), the left-hand side terms of (5), that is, TA and AdesT, are also both O1. Due to the structure of matrices A and Ades, the singularly perturbed algebraic Sylvester equation defined in (5) is numerically ill-conditioned. To overcome numerical ill-conditioning, we propose a new recursive algorithm for solving (5) in terms of reduced-order well-defined algebraic Sylvester equations. The dual version of (5) used for the system controller design is given by(9)ATc−TcAdesc=BF¯,F¯=FTcB=B11εB2,Adesc=Asc001εAfc.
Multiplying matrices in (9), we can get four algebraic equations for the partitioned matrix Tc as functions of the small parameter ε. It can be found from these equations [16, 17] that the structure of Tc is given by(10)Tc=T1cT2cT3c1εT4c.
Algebraic Sylvester equations (9)-(10) will be solved numerically in terms of reduced-order numerically well-conditioned algebraic Sylvester equations under the standard controller design assumptions [13].Assumption 5.
The pair (A, B) is controllable and the pair Adesc,F¯ is observable.
Moreover, the existence of a unique solution of (9) requires the assumption dual to Assumption 4.Assumption 6.
λA∩λAdesc=∅.


## 1.2. Parallel Algorithm for the Observer Sylvester Equation

The partitioned form of the Sylvester equation given in (5) subject to (6)–(8) is given by(11)T1A1+T2A3−AsT1+Q1=0,T1A2+T2A4−ϵAsT2+Q2=0,εT3A1+T4A3−AfT3+Q3=0,εT3A2+T4A4−AfT4+Q4=0.It can be seen from (11) that this system has two independent sets of linear algebraic equations: the first one for T1 and T2 and the second one for T3 and T4.Settingε=0 in (11), the algebraic equations for zeroth-order approximations of solutions T10, T20, T30, and T40 are obtained as(12)T10A1+T20A3−AsT10+Q1=0,(13)T10A2+T20A4+Q2=0,(14)T40A3−AfT30+Q3=0,(15)T40A4−AfT40+Q4=0.These equations can be solved independently as follows. Unique solutionT40 can be obtained from Sylvester equation (15) under the following assumption.Assumption 7.
λA4≠λAf.
SinceAf, defined in (6), is chosen by the designer as an asymptotically stable matrix, this assumption is easily satisfied. Having obtained T40, from (13) and (14), we can obtain T30andT20 independently as(16)T30=Af−1T40A3+Q3,(17)T20=−Q2+T10A2A4−1.
Substituting (17) into (12) results in(18)T10A0−AsT10+Q0=0,(19)A0=A1−A2A4−1A3,Q0=Q1−Q2A4−1A3.
The unique solutionT10 of algebraic Sylvester equation (18) exists under the following assumption.Assumption 8.
λA0≠λAs.
The solutionsT10,T20,T30, and T40 are Oε close to the exact solutions, that is,(20)T1=T10+εE1,T2=T20+εE2,T3=T30+εE3,T4=T40+εE4.
Now, we show that the values forEi,i=1,2,3,4, can be obtained by running iterations on independent linear reduced-order algebraic equations. Subtracting (12)–(15) from (11) and using (20), we obtain the error equations (after some algebra) in the following form:(21)E1A1−AsE1=−E2A3,(22)E2A4−εAsE2=−E1A2+AsT20,(23)εE3A1−AfE3=−T30A1−E4A3,(24)E4A4−AfE4=−T30A2−εE3A2.
The error equations can be solved iteratively using the fixed-point algorithm, in which the cross-coupling terms multiplied byε are delayed by one iteration. This idea has been used in several algorithms that involve small parameters.Algorithm 1.
(25)E1i+1A0−AsE1i+1=−AsT20A4−1A3−εAsE2iA4−1A3,E2i+1A4−εAsE2i+1=AsT20−E1i+1A2,εE3i+1A1−AfE3i+1=−T30A1−E4iA3,E4i+1A4−AfE4i+1=−T30A2−εE3iA2,E20=0,E30=Af−1T30A1+E40A3,E40A4−AfE40+T30A2=0.
We first solve (22) as(26)E2i+1=AsT20−E1i+1A2+εAsE2iA4−1.
Substituting (26) into (21) gives(27)E1i+1A0−AsE1i+1=−AsT20A4−1A3−εAsE2iA4−1A3.
Equations (26) and (27) have nice forms since the quantity E2 is multiplied by a small parameter ε. Similarly, equations for E3 and E4 can be iteratively solved as(28)εE3i+1A1−AfE3i+1=−T30A1−E4iA3,E4i+1A4−AfE4i+1=−T30A2−εE3iA2.
The following theorem presents the main feature of Algorithm1.Theorem 1.
Under Assumptions7 and 8, Algorithm 1 converges to the exact solution E with the rate of convergence of Oε. The convergence is obtained for sufficiently small values of ε that makes the radius of convergence smaller than 1 in each iteration leading to a contraction mapping. Hence, after i iterations, the solution T is obtained with the accuracy of Oεi, that is,(29)Tji=Tj0+εEji+Oεi,j=1,2,3,4;i=1,2,….Proof of Theorem.
1.
Fori=0, (27) implies(30)E11A0−AsE11=−AsT20A4−1A3−εAsE20A4−1A3=−AsT20A4−1A3.
Note thatE20=0. For i=1, (27) produces(31)E12A0−AsE12=−AsT20A4−1A3−εAsE21A4−1A3.
Subtracting (30) from (31), we have(32)E12−E11A0−AsE12−E11=−εAsE21A4−1A3=Oε.
At this point, we conclude that(33)E12−E11=Oε.
In a similar way, we can write the relationship betweenE13 and E12 as(34)E13−E12A0−AsE13−E12=εAsE21A4−1A3−εAsE22A4−1A3=Oε,which implies that(35)E13−E12=Oε.
Continuing the same procedure, we obtain(36)E1i+1−E1i=Oε.
Now, we work withE2 using (26). For i=0, we have(37)E21=AsT20−E11A2+εAsE20A4−1⇒E21=O1.
Fori=1,(38)E22=AsT20−E12A2+εAsE21A4−1
UsingE20=0 and the result in (33), we get(39)E22−E21=E12−E11A2A4−1+Oε⇒E22−E21=Oε.
Considering (26) for i=2 and using (39), we obtain(40)E23−E22=Oε.
If we keep repeating this process, we conclude that(41)E2i+1−E2i=Oε.
In addition, we have the following relationships. Subtracting (21)-(22) from (25) for i=0 produces(42)E11−E1A0−AsE11−E1=−εAsE20−E2A4−1A3=Oε⇒E11−E1=Oε,E21−E2A4−εAsE21−E2=−E11−E1A2=Oε⇒E21−E2=Oε.
Fori=1, (21)-(22) and (25) produce(43)E12−E1A0−AsE12−E1=−εAsE21−E2A4−1A3=Oε2⇒E12−E1=Oε2,E22−E2A4−εAsE22−E2=−E12−E1A2=Oε2⇒E22−E2=Oε2.
Continuing the same procedure, we have(44)E1i−E1=Oεi,i=1,2,3,…,(45)E2i−E2=Oεi,i=1,2,3,….
Similar procedures applied to (23)-(24) produces(46)E4i+1−E4i=Oε,E4i−E4=Oεi,i=1,2,3,…,(47)E3i+1−E3i=Oε,E3i−E3=Oεi,i=1,2,3,….
Results established in (44)–(47) can be summarized in(48)Eji−Ej=Oεi,j=1,2,3,4;i=1,2,…,which completes the proof of the stated theorem.


## 2. Parallel Algorithm for Controller Sylvester Equation

The controller design algebraic Sylvester equation defined in (9)-(10) can be partitioned as(49)A1T1c+A2T3c−T1cAsc+Q1c=0,εA1T2c+A2T4c−T2cAfc+Q2c=0,A3T1c+A4T3c−εT3cAsc+Q3c=0,εA3T2c+A4T4c−T4cAfc+Q4c=0,(50)Qc=BF¯=BFTc=−Q1c1εQ2c1εQ3c1ε2Q4c=O1O1εO1εO1ε2.Settingε=0 in (49), the zeroth-order approximations T1c0,T2c0,T3c0, and T4c0 can be obtained as follows:(51)A1T1c0+A2T3c0−T1c0Asc+Q1c=0,(52)A2T4c0−T2c0Afc+Q2c=0,(53)A3T1c0+A4T3c0+Q3c=0,(54)A4T4c0−T4c0Afc+Q4c=0.The unique solutionT4c0 can be found from algebraic Sylvester equation (54) under the following assumption, which is easily satisfied since Afc is an arbitrary stable closed-loop matrix.Assumption 9.
λA4≠λAfc.
From (52)-(53), we can obtain T2c0 and T3c0 as(55)T2c0=A2T4c0+Q2cAfc−1,(56)T3c0=−A4−1A3T1c0+Q3c.
Substituting (56) into (51) results in(57)A0T1c0−T1c0Asc+Q0c=0,(58)A0=A1−A2A4−1A3,Q0c=Q1c−A2A4−1Q3c.
The solutionT10 can be obtained by solving algebraic Sylvester equation (57) under the following assumption.Assumption 10.
λA0≠λAsc.
SinceAsc is chosen by the designer, this assumption is easily satisfied. We define the approximation errors as(59)T1c=T1c0+εE1c,T2c=T2c0+εE2c,T3c=T3c0+εE3c,T4c=T4c0+εE4c.
Subtracting (51)–(54) from (49) and using (59), we obtain the error equations in the following form:(60)A1E1c−E1cAsc=−A2E3c,E2cAfc−εA1E2c=−A1T2c0+A2E4c,A4E3c−εE3cAsc=−A3E1c−T3c0Asc,E4cAfc−A4E4c=A3T2c0−εA3E2c.
These error equations can be solved using the fixed-point algorithm, dual to Algorithm1, as follows.Algorithm 2.
(61)A0E1ci+1−E1ci+1Asc=−εA2A4−1E3ciAsc−A2A4−1T3c0Asc,A4E3ci+1−εE3ci+1Asc=−A3E1ci+T3c0Asc,E2ci+1Afc−εA1E2ci+1=A1T2c0+A2E4ci,E4ci+1Afc−A4E4ci+1=A3T2c0+εA3E2ci,E3c0=0,E2c0=0,A0E1c0−E1c0Asc+A2A4T3c0Asc=0.
The convergence proof of Algorithm2 can be done via the dual arguments used in Algorithm 1. Similarly, we can state the corresponding theorem dual to Theorem 1.Theorem 2.
Under Assumptions5, 9, and 10, Algorithm 2 converges for sufficiently small values of ε with the rate of Oε to the sought solution Tj,j=1,2,3,4; that is, after i iterations, we have(62)Tjci=Tjc0+εEjci+Oεi,j=1,2,3,4;i=1,2,3,….
The proof of Theorem2 parallels the one of Theorem 1.


## 3. Observer and Controller Designs

The general design of an observer and a controller using the Sylvester approach is presented in [13]. We will exploit two-time scale property so that the design is done in terms of reduced-order problems. The observer design procedure for the system defined in (2) has the following steps [13].(63)z˙t=Adeszt+TBut+K¯yt,x^t=T−1zt.The state feedback controller for the system defined in (2) can be obtained using the steps.Comment 1.
according to our experience, we need only one repetition to obtain invertibleTc. According to [13], the pair A,B controllable and the pair Ades,F¯ observable are both the necessary and sufficient condition for invertibility of Tc in the case of single-input single-output systems. For multiple-input multiple-output systems, this condition is only sufficient.
The feedback system is given by [13](64)x˙t=A−BFxt,λA−BF=λdesiredsystem.Comment 2.
well-conditioning with respect to matrix inversion of matricesT and Tc can be established by using the formula for matrix inversion of partitioned matrices [13]. For matrix T defined in (7), we obtain(65)T−1=T1εT2εT3εT4−1=T1−1+OεOε2Oε21εI+Oε−1.
Hence, theT matrix is well conditioned with respect to matrix inversion if matrix T1 is well conditioned with respect to matrix inversion. Similarly, for matrix inversion of Tc, the matrix T4c must be well conditioned with respect to the matrix inversion.Comment 3.
systems (63) and (64) have the structure of linear singularly perturbed systems.


## 4. Simulation Results

Consider a4th−order system with the matrices A,B, and C taken from [11] (Chapter 3, 1999) as follows:(66)A=00.400000000.345000−5.2400−4.65002.6200000−10.0000,B=00010,C=10000010.The pairA,C is observable, and we can proceed with the observer design algorithm. The designer decides to place observer eigenvalues at the desired location by choosing matrices As and Af. In the following, we will design a controller with the desired closed-loop eigenvalues placed at −0.2,−0.3,−7,−8. Note that A,B is controllable.


### 4.1. Observer Design Algorithm1

We choose the observer eigenvalues such that it is roughly ten times faster than the closed-loop system. Consequently, we chooseAdes as(67)Ades=As001εAf=−500000−600000−5000000−600.We chooseK¯ as(68)K¯=13243546,so that Ades,K¯ is controllable, as required in Step 2 of Algorithm 1.Algorithm 1: (Observer design).
Step 1: choose Ades in (6) such that Assumption 4 is satisfied.Step 2: guess K¯ in (5) and (7) such that Assumption 3 is satisfied.Step 3: solve TA−AdesT=K¯C using Algorithm 1.Step 4: if T−1 does not exist, go back to Step 2 and guess another K¯ and repeat the process until T−1 exists. The observer structure for zt is given as in [13]:The matrixQ defined in (8) is given as(69)Q=−K¯C=−10−30−20−40−30−50−40−60.The zeroth-order approximationsT10,T20,T30, and T40 are obtained as(70)T10=0.0200000−0.06830390.0333333−0.0758390,T20=−6.5022900−1.7035999−8.6584182−2.2685055,T30=0.06000000.0010570.06666660.0008801,T40=0.1009387−0.00053970.1007810−0.0004475.Performing iterations, we obtain the sixteen decimal digits accuracy afteri=50:(71)E150=00.7507135000.81924891,E250=71.6329673716.6030405993.8071279822.30657343,E350=10−3×0−0.48007370−0.4444867,E450=10−5×−0.70332390.0037606−0.48428400.0021505.The corresponding iterative solutionT^ and the exact solution T (obtained by using the MATLAB function lyap to solve the full-order Sylvester equation) are given by(72)T^=0.02000000.00676730.0661006−0.00432950.03333330.00608580.0722294−0.00378480.00600000.00010090.0100938−0.00005390.00666660.00008350.0100780−0.0000447=T.The differenceT−T^ is(73)T−T^=1.882474211376878×10−16.The solutionT is invertible in the first run of Algorithm 1 (see Section 5). The corresponding observer gain K=T^−1K¯ is(74)K=105×0.01122900.0101385−0.7819285−0.9259744−0.00234090.0007244−1.2083683−1.3963744.Checking the corresponding observer closed loop eigenvalues, we have(75)λA−KC=−49.9999999999540−60.0000000000346−500.000000001148−599.9999999999814,which with the accuracy of O10−12 is close to the chosen desired eigenvalues of the matrix Ades. The general Sylvester equation for the observed case is given as in (63).


### 4.2. Result of Comparison Algorithm1 with Algorithm 432

To compare our proposed Algorithm1 against the existing design algorithm, we implement the Algorithm 432 presented in [1]. Applying the Schur transformation to equation (5), we have(76)Ades′T′+T′A′=Q′,where(77)Ades′=UT−AdesU,T′=UTTV,A′=VTAV,Q′=UTQV,where U and V are the Schur transformation in order to construct a lower triangular form for Ades and an upper triangular form for A given in (5) and used in [1]. The formula for the solution T in (5), based on the relationship in (77), is given in(78)T=U−TT′V−1.Using the MATLAB Schur function and the similarity transformation given in (77), we have(79)U=0001001001001000,V=1.000000000000.62739100.778704300−0.77870430.627391000001.0000000,Ades′=−6000000−5000000−600000−50,A′=00.25095640.311481700−0.42820665.5850000−2.040205400−4.22179331.6437644000−10.0000000,Q′=−4.00000004.6722262−3.76434610−3.00000003.8935218−3.13695510−2.00000003.1148175−2.50956400−1.00000002.3361131−1.88217300.Partitioning (76) results in(80)Ades11′0Ades21′Ades22′T11′T12′T21′T22′+T11′T12′T21′T22′A11′A12′0A22′=Q11′Q12′Q21′Q22′,where dimensions of Ades11′, Ades22′, A11′, A22′, T11′, T22′, and Q11′, Q22′ are, respectively, Ades11′∈Rn1, Ades22′∈Rn2, A11′∈Rn1, A22′∈Rn2, T11′∈Rn1, T22′∈Rn2, and Q11′∈Rn1, Q22′∈Rn2 with the corresponding dimensions of Ades12′, Ades21′, A12′, A21′, T12′, T21′, and Q12′,Q21′.From equation (80) and the condition Ades21′=0, we have(81)Ades11′T11′+T11′A11′=Q11′,(82)Ades22′T21′+T21′A11′=Q21′,(83)Ades22′T22′+T22′A22′=Q22′−T21′A12′,(84)Ades11′T12′+T12′A22′=Q12′−T11′A12′.Solving reduced-order Sylvester equations (81)–(84), we have(85)Ades11′=−60000−500,Ades22′=−6000−50,A11′=00.25095640−0.4282066,A12′=0.311481705.5850000−2.0402054,A22′=−4.22179331.64376440−10.0000000,Q11′=−4.00000004.6722262−3.00000003.8935218,Q12′=−3.76434610−3.13695510,Q21′=−2.00000003.1148175−1.00000002.3361131,Q22′=−2.50956400−1.88217300.The corresponding solution ofT′ is given by(86)T′=T11′T12′T21′T22′=0.0066666−0.00777870.00616160.00004260.0060000−0.00777730.00613890.00005080.0333333−0.05140730.03476750.00231470.0200000−0.04622590.03006590.0023955.Original solutionT using the inverse transformation of (78) is given by(87)T=0.0200000−0.00558920.05485950.00239550.0333333−0.00517880.06184390.00231470.0060000−0.00009900.00990770.00005080.0066666−0.00008220.00992300.0000426.The difference betweenT and T^ of our proposed Algorithm 1 is(88)T^−T=0.0244206.Since they are different, we conclude that our proposed Algorithm1 has a better accuracy regarding the solution of the general Sylvester equation.Algorithm 2(Controller design).
Similarly, we design a controller for the same system using the algorithm for solving the controller algebraic Sylvester equation fromSection 4. We choose(89)Adesc=Asc001εAfc=−0.20000−0.30000−70000−8.
The zeroth-order approximationsT1c0,T2c0,T3c0, and T4c0 are obtained as(90)T1c0=−10.2973342−20.43701555.148667115.3277616,T2c0=001.46528871.5177425,T3c0=−2.9847345−13.32848835.00000007.0000000,T4c0=−2.9730496−3.51940292.66666664.5000000.
Using the proposed algorithm, we obtain afteri=50:(91)E1c50=2.477925231.3924384−1.2389626−23.5443288,E2c50=−0.7544016−0.7109162−1.4508588−0.9591003,E3c50=0.718239220.47332941.02040812.1649484,E4c50=2.94377152.224000600.
The iterative solutionT^c and the exact solution Tc are given by(92)T^c=abcd=Tc,where(93)a=−10.0495417−17.29777165.024770812.9733287,b=−0.0754401−0.07109161.32020281.4218325,c=−2.9129106−11.28115545.10204087.2164948,d=−26.7867249−32.970029126.666666645.0000000.
Their difference is(94)Tc−T^c=1.761132598848588×10−14.
The solutionTc is invertible in the first run of Algorithm 2 (see Section 5). The controller gain F=F¯T^c−1 is given as(95)F=0.929306331.027256330.431286250.08500000.
CheckingλA−BF, we have(96)λA−BF=−0.199999999999998−0.300000000000001−7.000000000000027−7.999999999999988,which has produced the desired eigenvalues with the accuracy of O10−13.Algorithm 2: (Controller design).
Step 1: choose Adesc in (9) such that Assumption 6 is satisfied.Step 2: guess F¯ in (9) and (10) such that Assumption 5 is satisfied.Step 3: solve ATc−TAdesc=BF¯ using Algorithm 2.Step 4: if Tc−1 exists, then F=F¯Tc−1, otherwise go back to Step 2, guess another F¯, and repeat the process.


## 4.1. Observer Design Algorithm1

We choose the observer eigenvalues such that it is roughly ten times faster than the closed-loop system. Consequently, we chooseAdes as(67)Ades=As001εAf=−500000−600000−5000000−600.We chooseK¯ as(68)K¯=13243546,so that Ades,K¯ is controllable, as required in Step 2 of Algorithm 1.Algorithm 1: (Observer design).
Step 1: choose Ades in (6) such that Assumption 4 is satisfied.Step 2: guess K¯ in (5) and (7) such that Assumption 3 is satisfied.Step 3: solve TA−AdesT=K¯C using Algorithm 1.Step 4: if T−1 does not exist, go back to Step 2 and guess another K¯ and repeat the process until T−1 exists. The observer structure for zt is given as in [13]:The matrixQ defined in (8) is given as(69)Q=−K¯C=−10−30−20−40−30−50−40−60.The zeroth-order approximationsT10,T20,T30, and T40 are obtained as(70)T10=0.0200000−0.06830390.0333333−0.0758390,T20=−6.5022900−1.7035999−8.6584182−2.2685055,T30=0.06000000.0010570.06666660.0008801,T40=0.1009387−0.00053970.1007810−0.0004475.Performing iterations, we obtain the sixteen decimal digits accuracy afteri=50:(71)E150=00.7507135000.81924891,E250=71.6329673716.6030405993.8071279822.30657343,E350=10−3×0−0.48007370−0.4444867,E450=10−5×−0.70332390.0037606−0.48428400.0021505.The corresponding iterative solutionT^ and the exact solution T (obtained by using the MATLAB function lyap to solve the full-order Sylvester equation) are given by(72)T^=0.02000000.00676730.0661006−0.00432950.03333330.00608580.0722294−0.00378480.00600000.00010090.0100938−0.00005390.00666660.00008350.0100780−0.0000447=T.The differenceT−T^ is(73)T−T^=1.882474211376878×10−16.The solutionT is invertible in the first run of Algorithm 1 (see Section 5). The corresponding observer gain K=T^−1K¯ is(74)K=105×0.01122900.0101385−0.7819285−0.9259744−0.00234090.0007244−1.2083683−1.3963744.Checking the corresponding observer closed loop eigenvalues, we have(75)λA−KC=−49.9999999999540−60.0000000000346−500.000000001148−599.9999999999814,which with the accuracy of O10−12 is close to the chosen desired eigenvalues of the matrix Ades. The general Sylvester equation for the observed case is given as in (63).


## 4.2. Result of Comparison Algorithm1 with Algorithm 432

To compare our proposed Algorithm1 against the existing design algorithm, we implement the Algorithm 432 presented in [1]. Applying the Schur transformation to equation (5), we have(76)Ades′T′+T′A′=Q′,where(77)Ades′=UT−AdesU,T′=UTTV,A′=VTAV,Q′=UTQV,where U and V are the Schur transformation in order to construct a lower triangular form for Ades and an upper triangular form for A given in (5) and used in [1]. The formula for the solution T in (5), based on the relationship in (77), is given in(78)T=U−TT′V−1.Using the MATLAB Schur function and the similarity transformation given in (77), we have(79)U=0001001001001000,V=1.000000000000.62739100.778704300−0.77870430.627391000001.0000000,Ades′=−6000000−5000000−600000−50,A′=00.25095640.311481700−0.42820665.5850000−2.040205400−4.22179331.6437644000−10.0000000,Q′=−4.00000004.6722262−3.76434610−3.00000003.8935218−3.13695510−2.00000003.1148175−2.50956400−1.00000002.3361131−1.88217300.Partitioning (76) results in(80)Ades11′0Ades21′Ades22′T11′T12′T21′T22′+T11′T12′T21′T22′A11′A12′0A22′=Q11′Q12′Q21′Q22′,where dimensions of Ades11′, Ades22′, A11′, A22′, T11′, T22′, and Q11′, Q22′ are, respectively, Ades11′∈Rn1, Ades22′∈Rn2, A11′∈Rn1, A22′∈Rn2, T11′∈Rn1, T22′∈Rn2, and Q11′∈Rn1, Q22′∈Rn2 with the corresponding dimensions of Ades12′, Ades21′, A12′, A21′, T12′, T21′, and Q12′,Q21′.From equation (80) and the condition Ades21′=0, we have(81)Ades11′T11′+T11′A11′=Q11′,(82)Ades22′T21′+T21′A11′=Q21′,(83)Ades22′T22′+T22′A22′=Q22′−T21′A12′,(84)Ades11′T12′+T12′A22′=Q12′−T11′A12′.Solving reduced-order Sylvester equations (81)–(84), we have(85)Ades11′=−60000−500,Ades22′=−6000−50,A11′=00.25095640−0.4282066,A12′=0.311481705.5850000−2.0402054,A22′=−4.22179331.64376440−10.0000000,Q11′=−4.00000004.6722262−3.00000003.8935218,Q12′=−3.76434610−3.13695510,Q21′=−2.00000003.1148175−1.00000002.3361131,Q22′=−2.50956400−1.88217300.The corresponding solution ofT′ is given by(86)T′=T11′T12′T21′T22′=0.0066666−0.00777870.00616160.00004260.0060000−0.00777730.00613890.00005080.0333333−0.05140730.03476750.00231470.0200000−0.04622590.03006590.0023955.Original solutionT using the inverse transformation of (78) is given by(87)T=0.0200000−0.00558920.05485950.00239550.0333333−0.00517880.06184390.00231470.0060000−0.00009900.00990770.00005080.0066666−0.00008220.00992300.0000426.The difference betweenT and T^ of our proposed Algorithm 1 is(88)T^−T=0.0244206.Since they are different, we conclude that our proposed Algorithm1 has a better accuracy regarding the solution of the general Sylvester equation.Algorithm 2(Controller design).
Similarly, we design a controller for the same system using the algorithm for solving the controller algebraic Sylvester equation fromSection 4. We choose(89)Adesc=Asc001εAfc=−0.20000−0.30000−70000−8.
The zeroth-order approximationsT1c0,T2c0,T3c0, and T4c0 are obtained as(90)T1c0=−10.2973342−20.43701555.148667115.3277616,T2c0=001.46528871.5177425,T3c0=−2.9847345−13.32848835.00000007.0000000,T4c0=−2.9730496−3.51940292.66666664.5000000.
Using the proposed algorithm, we obtain afteri=50:(91)E1c50=2.477925231.3924384−1.2389626−23.5443288,E2c50=−0.7544016−0.7109162−1.4508588−0.9591003,E3c50=0.718239220.47332941.02040812.1649484,E4c50=2.94377152.224000600.
The iterative solutionT^c and the exact solution Tc are given by(92)T^c=abcd=Tc,where(93)a=−10.0495417−17.29777165.024770812.9733287,b=−0.0754401−0.07109161.32020281.4218325,c=−2.9129106−11.28115545.10204087.2164948,d=−26.7867249−32.970029126.666666645.0000000.
Their difference is(94)Tc−T^c=1.761132598848588×10−14.
The solutionTc is invertible in the first run of Algorithm 2 (see Section 5). The controller gain F=F¯T^c−1 is given as(95)F=0.929306331.027256330.431286250.08500000.
CheckingλA−BF, we have(96)λA−BF=−0.199999999999998−0.300000000000001−7.000000000000027−7.999999999999988,which has produced the desired eigenvalues with the accuracy of O10−13.Algorithm 2: (Controller design).
Step 1: choose Adesc in (9) such that Assumption 6 is satisfied.Step 2: guess F¯ in (9) and (10) such that Assumption 5 is satisfied.Step 3: solve ATc−TAdesc=BF¯ using Algorithm 2.Step 4: if Tc−1 exists, then F=F¯Tc−1, otherwise go back to Step 2, guess another F¯, and repeat the process.


## 5. Conclusions

It was shown that the numerically ill-conditioned Sylvester algebraic equation for singularly perturbed systems can be decomposed into four lower-order well-conditioned Sylvester equations. The recursive fixed-point-type methods were utilized in order to obtain numerical solutions for such lower-order algebraic Sylvester equations. The corresponding observer and controller design algorithms for assignment of observer and controller closed-loop eigenvalues in terms of reduced-order slow and fast subproblems were presented.



---

*Source: 3948564-2020-05-30.xml*