A quasi-LMI approach to computing stabilizing parameter ranges of multi-loop PID controllers

Journal of Process Control 17 (2007) 59–72 www.elsevier.com/locate/jprocont

A quasi-LMI approach to computing stabilizing parameter ranges of multi-loop PID controllers Qing-Guo Wang *, Chong Lin, Zhen Ye, Guilin Wen, Yong He, Chang Chieh Hang Department of Electrical and Computer Engineering, National University of Singapore, Singapore 119260, Singapore Received 1 February 2006; received in revised form 28 July 2006; accepted 9 August 2006

Abstract This paper addresses the problem of determining the parameter ranges of multi-loop PID (proportional-integral-derivative) controllers which stabilize a given process. An effective computational scheme is established by converting the considered problem to a quasi-LMI problem connected with robust stability test. The descriptor model approach is employed together with linearly parameter-dependent Lyapunov function method. Examples are given for illustration. The results are believed to facilitate real time tuning of multi-loop PID controllers for practical applications.  2006 Elsevier Ltd. All rights reserved. Keywords: PID control; Uncertain systems; Robust stability

1. Introduction Proportional-integral-derivative (PID) controllers have dominated industrial applications for more than fifty years because of their simplicity in controller structure, robustness to modeling errors and disturbances, and the availability of numerous tuning methods [1,6]. Stability analysis of single-input single-output (SISO) PID systems is straightforward. Usually, the Nyquist stability theorem is utilized with help of the Nyquist curve of the open loop transfer function. For multi-input multi-output (MIMO) systems, the generalized Nyquist stability theorem was addressed by Rosenbrock [40], MacFarlane and Belletrutti [31], Nwokah [35] and Morari [34]; and effectively unified by Nwokah et al. [36]. The relevant tools such as characteristic loci, Nyquist arrays and Gershgorin bands are developed to help MIMO system analysis and design in frequency domain, which is similar to the SISO case in nature but

*

Corresponding author. Tel.: +65 65162282; fax: +65 67791103. E-mail address: [email protected] (Q.-G. Wang).

0959-1524/$ - see front matter  2006 Elsevier Ltd. All rights reserved. doi:10.1016/j.jprocont.2006.08.006

not as convenient as their counterparts in SISO case, owing to their complexity. As for PID controller tuning, most of the existing works have been on SISO systems, such as Ziegler–Nichols like method [18], pole-placement design [15], root-locus based methods [42,52], frequency response methods [20,13,53], and optimization techniques [43,2]. Comparatively, much fewer methods are available for MIMO PID tuning. Among the most known is the biggest log modulus tuning (BLT) method [30], which is viewed as a direct extension of classical Ziegler–Nichols method to the MIMO case. Palmor et al. [37] proposed a decentralized PID controller design for two-input two-output (TITO) systems, in which the desired critical point is used to tune the PID controller by Zieler–Nichols rule or its modifications. Loh et al. [29] presented a PID tuning method for multi-loop systems based on relay feedback experiments. Zhang et al. [55] presented a technique for multi-loop PI controller design to achieve dominant pole-placement for TITO processes. Wang et al. [51] proposed a tuning method for fully cross-coupled multi-variable PID controller from decentralized relay feedback test, which first estimates the frequency response matrix at zero and oscillation

60

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

frequencies from the test and then designs the diagonal elements of multi-variable PID controller independent of offdiagonal ones. The internal model control (IMC) was also generalized to MIMO systems [12]. Although great progress on PID control has been made recently, some fundamental issues remain to be addressed for better understanding and applications of PID controllers, especially for MIMO case. The very first task at the outset of PID control is to get a stabilizing PID controller for a given process; If possible, it would be desirable to find their parameter regions for stabilizing a given process. This problem is of great importance, both theoretically and practically, and also related to the problem of stability margins or robustness. Unfortunately, most of existing methods mentioned above can only determine some values (but not ranges) of stabilizing PID parameters. For SISO systems, the gain and phase margins are well defined and can easily be determined graphically or numerically. The characterization of the set of all stabilizing PID controllers is developed for a SISO delay-free linear time invariant (LTI) plant in Datta et al. [10], and for a SISO LTI plant with time-delay in Silva et al. [44], based on the HermiteBiehler Theorem, its extensions, and some optimization techniques. However, it is pointed out [48,50] that their methods are unlikely to be extended to the MIMO case. In the context of MIMO PID systems, not much work has been done. Safonov and Athans [41] proposed a singular value approach to multi-loop stability analysis, where the sufficient condition of stability and some characterization of frequency-dependent gain and phase margins for multi-loop systems is developed. But their criterion is conservative. Morari [34] introduced the concept of integral controllability, that is, for any k such that 0 < k 6 k* (k* > 0), the feedback system with the open loop as G(s)k/s is stable. He also gave the necessary and sufficient conditions of integral controllability for MIMO systems but did not tell how to determine k*. On the other hand, Yaniv [54] developed a control method to meet some stability margins which are defined loop by loop like a single variable system. Li and Bruce Lee [24] showed that the H1 norm of a sensitivity function matrix for a stable multivariable closed-loop system is related to some common gain and phase margins for all the loops. Ho et al. [21] defined gain and phase margins and use them for multi-variable control system design assuming that the process is diagonally dominant or made so. Such definitions of gain and phase margins based on Gershgorin bands or other frequency domain techniques are more or less conservative, which brings some limitation of their applications. Doyle [11] developed the l-analysis, which is utilized as an effective tool for robust stabilizing analysis in multi-variable feedback control [45]. As a method in frequency domain, the l-analysis treats system uncertainties as complex valued. But the parameters of PID controllers are all real. Thus, when l-analysis is used to determine the stabilizing ranges of PID controllers, conservativeness is inevitable, and a detailed analysis on this is shown in the example of

Section 6. In summary, it can been concluded that there seems neither satisfactory definition for MIMO gain and phase margins, nor effective technique for determining them so far. To our best knowledge, no results are available to find the stabilizing PID ranges for MIMO processes. It should be noted that recent developments in the timedomain approaches to MIMO PID control is appealing [56,27,17]. The basic idea in such approaches is to transform MIMO PID control system to an equivalent static output feedback (SOF) system for which plenty of results can be used. Though the static output feedback stabilizability is still hard to solve, Lyapunov-like conditions [47] and the solution of some linear quadratic (LQ) control problem [46,22] have been developed to enable stability analysis and stabilization. Bernussou et al. [4] showed how to convert an LQ problem in a new parameter space such that the resulting equivalent problem is convex. Boyd et al. [5] showed how to convert control design problems to a class of convex programming problems with linear objective function and constraints expressed in terms of linear matrix inequalities (LMI). Cao et al. [7] proposed an iterative LMI approach for static output feedback stabilization, and sufficient LMI conditions for such a control problem were given by Crusius and Trofino [8]. It seems that time domain approach with help of the LMI-like tools opens an new direction to analysis and design of MIMO PID control systems and makes it possible to give better results than classical frequency domain methods mentioned above. In this paper, we investigate a linear MIMO plant under a diagonal or block-diagonal PID control structure using time-domain approach to determine the PID stabilizing ranges as well as the gain margins. We will modify our recent descriptor model approach to transform the problem into a robust stability problem for a linear polytopic system. In this way, a detailed scheme in descriptor version is provided for the robust stability test and an effective procedure is given to find the parameter ranges of PID controllers. The scheme incorporates a relaxed LMI technique which not only effectively solves the considered PID problem, but also leads to better results than the existing methods for special cases of standard polytopic systems [16,38,39,19]. The present procedure is a kind of quasiLMI based convex computation which can be fulfilled through LMI-Toolbox [14]. Notation. Rn denotes the n-dimensional real Euclidean space; the superscript ‘T’ stands for the matrix transpose; W > 0 (W P 0) means that W is real, symmetric and positive-definite (positive-semidefinite).

2. Problem formulation To illustrate mutual dependence of loop gains which stabilizes a coupled multi-variable system under decentralized control, let us consider a 2 · 2 system with the transfer function matrix:

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

2

1 6sþ 1 GðsÞ ¼ 6 4 3 sþ1

3 2 s þ 17 7: 4 5 sþ1

10 8 6 4

A decentralized proportional controller K(s) = diag{k1, k2} is applied to it in the unity negative feedback configuration, as shown in Fig. 1. It follows [40,49] that the characteristic equation of the closed-loop system is

2 0 –2

P c ðsÞ ¼ P G ðsÞP K ðsÞ det½I þ GðsÞKðsÞ

–4

¼ s2 þ ðk 1 þ 4k 2 þ 2Þs þ ðk 1 þ 4k 2 þ 1  2k 1 k 2 Þ ¼ 0;

–6

ð1Þ

–8

where PG(s) and PK(s) are the pole polynomials of G(s) and K(s), respectively. The closed-loop system is stable if and only if all the roots of Pc(s) have negative real parts, or  k 1 þ 4k 2 þ 2 > 0; ð2Þ k 1 þ 4k 2 þ 1  2k 1 k 2 > 0: The solution of (2) is 8 k 2 >  14 k 1  12 ; > > < k 2 < 2ðkk11þ1 ; if k 1 > 2; 2Þ > > : k > k1 þ1 ; if k < 2:

g11

y1

+ +

g12

g21 +

k2

u2

g22

+

+

–

Fig. 1. Block diagram of TITO system.

–0.4

–0.6

–0.8

–1 –0.2

0

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

Re Fig. 3. Nyquist curve of g~22 for k1 = 1.

Then, the Nyquist stability theorem for SISO systems can be applied to determine the stabilizing range of k2, if k1 is specified. For example, when k1 = 1, the Nyquist curve of g~22 ðjxÞ is depicted as Fig. 3. Since the open-loop

u1

10

–0.2

y 2 ðsÞ k 1 g12 g21 4s þ 4  2k 1 ¼ g22  : ¼ u2 ðsÞ 1 þ k 1 g11 ðs þ 1Þðs þ 1 þ k 1 Þ

k1

5

ð3Þ

Thus, the equivalent open-loop transfer function which k2 stabilizes when the first loop is closed with gain k1 is

–

0

0

This stabilizing range is drawn as the shaded region in Fig. 2. For example, when k1 = 1, we get from the figure or (3) that k2 > 3/4. Alternately, we may look at the system loop by loop. From Fig. 1, it is straightforward [32] to see that   k 1 g12 g21 y 2 ¼ g22 u2 þ g21 u1 ¼ g22  ð4Þ u2 : 1 þ k 1 g11

+

–5

Fig. 2. Stabilization region of (k1, k2).

1

2ðk 1 2Þ

g~22 ðsÞ ¼

–10 –10

Im

2

61

y2

system g~22 has no poles in the RHP, the closed-loop system is stable if and only if the Nyquist curve of g~22 does not encircle the point (1/k2, 0), or k2 > 3/4, which is the same as before. In general, the characteristic equation for the second equivalent loop is 1 þ k 2 g~22 ðsÞ ¼ 0, which gives exactly (1). Similarly, we may work with the first equivalent open-loop transfer function and lead to the same results. It is seen from this example that the stabilizing range of one loop gain depends on the value of other loop’s gain. This range can be computed with the SISO method for the equivalent SISO plant derived from the given MIMO system with all other loops closed with the fixed loop gains kj, j 5 i. If k1 is fixed at some value, the stabilizing range for k2 is uniquely determined. For instance, k1 = 1 yields k2 > 3/4, and graphically such a stabilizing range for k2 is between two intersection points of line k1 = 1 with the lower and upper boundaries of the shaded (stabilizing) region of Fig. 2. Note that loop 1 may have some uncertainties on its parameters and/or k1 needs to be tuned or

62

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

de-tuned separately. When k1 or loop 1 has some changes, the previous stabilizing range for k2 may not be stabilizing any more. Such results are not very useful in the context of MIMO gain margins and their applications as they are sensitive to other loops’ gains. Therefore, it is more practical and useful to prescribe a range for k1 when determining the stabilizing range for k2. In general, if k1 varies in some range, the stabilizing range for k2 can be uniquely determined. For instance, k1 2 [1,2] yields k2 2 [3/4, + 1). Graphically such a stabilizing region for both k1 and k2 is a rectangle with length k1 from 1 to 2 and width k2 from 3/4 to +1. When the range of k1 changes, so does the stabilizing range of k2. For instance, {(k1, k2)|k1 2 [3,4], k2 2 [5/4,5/4]} gives another stabilizing rectangle for k1 and k2. In view of the above observations, we are motivated to find such stabilizing rectangles and formulate the problem as follows. Consider an m · m square plant G(s) with n-dimensional state-space realization: x_ ðtÞ ¼ AxðtÞ þ BuðtÞ; yðtÞ ¼ CxðtÞ;

ð5Þ

where x 2 Rn is the state, y 2 Rm is the output, B and C are real constant matrices with appropriate dimensions. This paper focuses on the following form of PID controllers: U(s) = K(s)E(s), where e(t) = r(t)  y(t), r(t) is the set point and K2 þ K 3s KðsÞ ¼ K 1 þ s ¼ diagfk 11 I 11 ; . . . ; k 1r1 I 1r1 g 1 þ diagfk 21 I 21 ; . . . ; k 2r2 I 2r2 g s þ s diagfk 31 I 31 ; . . . ; k 3r3 I 3r3 g;

ð6Þ

where k1i, k2j and k3l are scalars to be determined, I1i, I2j and I3l are identity matrices dimensionsP m1i, m2j and Pr1 with P r2 r3 m3l, respectively, and i¼1 m1i ¼ j¼1 m2j ¼ l¼1 m3l ¼ m. Since our concern in this paper is stabilization, r(t) has no effect and can be ignored. The controller in (6) can be re-written in time domain as Z t uðtÞ ¼ K 1 yðtÞ  K 2 yðhÞ dh  K 3 y_ ðtÞ 0

:¼ 

r1 X i¼1

k 1i I 1i yðtÞ

r2 X

Z k 2i I 2i

i¼1

t

yðhÞ dh 

0

r3 X

Note that in Problem 1, the controller of the form  KðsÞ ¼ k 1 þ ks2 þ k 3 s I m corresponds to the specially chosen r1 = r2 = r3 = 1; the controller of the form KðsÞ ¼ diagfk 1i g þ 1s diagfk 2i g þ s diagfk 3i g corresponds to the specially chosen r1 = r2 = r3 = m (or m1i = m2j = m3l = 1). It is worth mentioning that the gain margins for MIMO systems can readily be defined and obtained as by-products of solutions to Problem 1. Consider the example again in the special case where K(s) = kI2, or k1 and k2 are equal to each other. Then, it follows from pffiffiffiffiffi pffiffiffiffiffi (3) that the stabilizing range is k 2 ½ð5  33Þ=4; ð5 þ 33Þ=4. Graphically, such a stabilizing range is obtained as the straight line, BD, in Fig. 2, where B and D are two intersection points of line k1 = k2 with the boundary of shaded region. BD is uniquely determined. In general, for an m · m square plant in (5) under the decentralized proportional controller form in (7) with the common gain for all loops, K(s) = kIm, suppose that the solution to Problem 1 is k 2 ½k; k:

ð9Þ

This stabilizing range is uniquely determined and called as the common gain margin of the system. Graphically, such a stabilizing range (9) is the largest line segment of k1 = k2 =    = km available in the stabilizing region for ki, i = 1, 2, . . . , m. Consider now the decentralized proportional controller, K(s) = diag{k1, k2, . . . , km}, with probably different gains for different loops. Suppose that the solution to Problem 1 is k i 2 ½k i ; k i ; i ¼ 1; 2; . . . ; m: ð10Þ Then, the closed-loop remains stable even when the gain for the ith loop, ki, varies between ki and k i , provided that other loop gains, kj, j = 1, 2, . . . , m, j 5 i, are (arbitrary but) also within their respective ranges. ½k i ; k i  is called the gain margin for the ith loop, subject to other loops’ gain margins within ½k j ; k j , j = 1, 2, . . . , m, j 5 i. Note that the margins so defined allow variations/uncertainties of other loops’ gains within ½k j ; k j , which facilitates their use in loop tunings. Actually, such ranges are usually quite large for normal stable plants. The formulation here truly reflects the distinct feature of a MIMO system from the SISO case that the stabilizing range for a loop depends on other loops’ gains in general.

_ k 3i I 3i yðtÞ;

i¼1

3. The proposed approach ð7Þ

where I mi ¼ diagf0; . . . ; 0; I mi ; 0; . . . ; 0g 2 Rmm ; m ¼ 1; 2; 3;

i ¼ 1; 2; . . . ; rm :

ð8Þ

The problem considered in this paper is as follows. Problem 1. For a plant (5) under the controller (7), find the ranges of scalars k1i, k2j and k3l, i = 1, . . . , r1, j = 1, . . . , r2, l = 1, . . . , r3, such that the closed-loop system is stable for all allowable k1i, k2j and k3l in these ranges.

We will transform the closed-loop system into a descriptor form analogous to that of Lin et al. [27,28]. It should be pointed out that the descriptor model in Lin et al. [27] is obtained an augmented state xðtÞ ¼ R t by introducing T ½xT ðtÞ; 0 xT ðhÞ dh; x_ T ðtÞ . This brings conservatism since the resulting design is only applicable to a narrow class of systems with matrix C being of full column rank. To overcome such a drawback, in thisR paper, we replace the t augmented state by xðtÞ ¼ ½xT ðtÞ; 0 y T ðhÞ dh; y_ T ðtÞT . The new output remains the same as in Lin et al. [27], i.e.,

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

Rt y ðtÞ ¼ ½y T ðtÞ; 0 y T ðhÞ dh; y_ T ðtÞT . Noticing the fact that y_ ðtÞ ¼ CAxðtÞ þ CBuðtÞ, system (5) with (7) is then transformed into the following descriptor control system: Ex_ ðtÞ ¼ AxðtÞ þ BuðtÞ; ð11Þ

y ðtÞ ¼ CxðtÞ; uðtÞ ¼ Ky ðtÞ; where 2

In

6 E¼40 0 2 C 6 C¼40 0

0

0

2

3

A

0

0

3

6 7 7 0 5; A ¼ 4 C 0 0 5; 0 CA 0 I m 3 0 7 0 5; K ¼ ½ K 1 K 2 K 3 : Im

Im 0 0 Im 0

2

B

3

6 7 B ¼ 4 0 5; CB

In view of the diagonal structure of PID controller (7), the closed-loop system of (11) is rewritten as Ex_ ðtÞ ¼ ðA  BKCÞxðtÞ ¼

A

r1 X

k 1i A1i 

i¼1

r2 X

k 2j A2j 

j¼1

r3 X

6 A1i ¼ 4 2

BI 1i C

0

0 CBI 1i C

0 0

0 0 6 A3l ¼ 4 0 0 0 0

0

3

7 0 5; 0 3

BI 3l 7 0 5; CBI 3l

0

k 3l A3l xðtÞ

BI 2j

6 A2j ¼ 4 0 0 0 CBI 2j

0

3

7 0 5; 0

r1 X

k 1i A1i 

i¼1

r2 X

k 2i A2i 

i¼1

r3 X

k 3i A3i :

ð14Þ

i¼1

The task now is to compute the perturbation ranges for scalars k mi such that ðE; Acl Þ remains admissible. To this end, we specify the lower and upper bounds for k mi as upp blow mi and bmi , respectively, i.e., k mi 2 ½blow ; bupp ; mi mi

m ¼ 1; 2; 3; i ¼ 1; 2; . . . ; rm : blow i

ð15Þ

bupp i

For brevity, relabel them as and with i = 1, 2, . . . , upp r1 + r2 + r3. Let r0 = r1 + r2 + r3 and b ¼ ½blow 1 ; b1 ; . . . ; low upp br0 ; br0 . Then, Acl is equivalently recast as a matrix polytope with r ¼ 2r0 vertices denoted by Aj ðbÞ 2 Rðnþ2mÞðnþ2mÞ , ( r X Acl 2 AðaÞ : AðaÞ ¼ aj Aj ðbÞ; j¼1

)

ð16Þ

aj ¼ 1; aj P 0; j ¼ 1; 2; . . . ; r :

j¼1

l¼1

2

Acl ¼ A0cl 

!

ð12Þ 2

techniques available [7,56,27]. A specific procedure is provided in Section 5 to fulfill this step. Next, set k mi ¼ k mi  k 0 , m = 1, 2, 3, i = 1, 2, . . . , rm. Then, mi

r X

:¼ AclxðtÞ; where

63

By the work of Masubuchi et al. [33], it is known that a nominal pair (E, A) is admissible if and only if there exists a matrix P such that PTA + ATP < 0 with PTE = ETP P 0. Therefore, the pair ðE; Acl Þ is robustly admissible if and only if there exists a parameter-dependent Lyapunov matrix P(a) such that T

P ðaÞ E ¼ ET P ðaÞ P 0; T

ð17Þ

T

P ðaÞ AðaÞ þ AðaÞ P ðaÞ < 0: ð13Þ

for i = 1, 2, . . . , r1, j = 1, 2, . . . , r2 and l = 1, 2, . . . , r3. A descriptor system of the form (12) is called admissible if the system, or say, the pair ðE; Acl Þ, is regular, impulse-free and stable. Please refer to Dai [9] and Masubuchi et al. [33] for detailed definitions. So far, Problem 1 has been converted to the following problem. Problem 2. Find the ranges of scalars k1i, k2j and k3l, i = 1, . . . , r1, j = 1, . . . , r2, l = 1, . . . , r3, such that the closedloop system (12) is admissible for all allowable k1i, k2j and k3l in these ranges. To solve Problem 2, one could adopt the structured singular value method (i.e., l-analysis) as presented in Lin et al. [26,25]. Noticing the fact that the l-analysis may produce conservative results due to the requirement of common perturbation bounds, we next suggest an alternative method in the polytopic context. To address Problem 2, the first step k 0mi are such Pr1 is0to findP r2 0 0 0 that Pr3 0ðE; Acl Þ with Acl ¼ A  i¼1 k 1i A1i  i¼1 k 2i A2i  i¼1 k 3i A3i is admissible. This step can be done by standard

ð18Þ

An alternative result which is equivalent to the above criterion is easy to be established as follows. Lemma 3.1. The pair ðE; AðaÞÞ is robustly admissible if and only if there exist parameter-dependent matrices P(a), F(a) and H(a) such that T

P ðaÞ E ¼ ET P ðaÞ P 0; " T T F ðaÞAðaÞ þ AðaÞ F ðaÞ T T P ðaÞ  F ðaÞ þ H ðaÞ AðaÞ

H H ðaÞ  H ðaÞT

#

ð19Þ < 0: ð20Þ

Here and in the sequel, an ellipsis w denotes a block induced by symmetry. Proof. The proof is parallel to that for standard systems in Geromel et al. [16] and Peaucelle et al. [38]. h Using Lemma 3.1, we have the following LMI-based result. Proposition 3.1. The pair ðE; Acl Þ is robustly admissible if there exist matrices Pj, Fj, Hj and Xjl with X jj ¼ X Tjj , l 6 j, j; l ¼ 1; 2; . . . ; r, such that P Tj E ¼ ET P j P 0;

ð21Þ

64

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

Hjl þ Hlj < X jl þ X Tjl ; j ¼ 1; 2; . . . ; r; l 6 j; 2 3 X 11 H    H 6X 7 6 21 X 22    H 7 6 . .. 7 .. .. 6 . 7 6 0; 4 . . 5 . . X r1 X r2    X rr

ð22Þ

Xjl ¼ ð23Þ

Hjl ¼

T

#

F j Al ðbÞ þ Al ðbÞ F Tj

H

P j  F Tj þ H Tj Al ðbÞ

H j  H Tj

:

Proof. Let the parameter-dependent matrices P(a), F(a) and H(a) be P ðaÞ ¼

r X

aj P j ;

F ðaÞ ¼

r X

j¼1

H ðaÞ ¼

r X

aj F j ;

j¼1

ð24Þ

aj H j :

j¼1

If conditions (22) and (23) are true, substituting (24) into the matrix of (20), yields " # F ðaÞAðaÞ þ AðaÞT F ðaÞT H T

T

P ðaÞ  F ðaÞ þ H ðaÞ AðaÞ H ðaÞ  H ðaÞ ¼

r X r X j¼1

<

r X

aj al Hjl ¼

l¼1

a2j X jj þ

j¼1

r X

a2j Hjj þ

r X

r X

T

aj al ðHjl þ Hlj Þ

l
j¼1

T

F j Al ðbÞ þ Al ðbÞ F Tj

H

Z j E þ LY j  F Tj þ H Tj Al ðbÞ H j  H Tj

# :

aj al ðX jl þ X Tjl Þ

Step 1. Find a setPof scalars P k 0mi such thatPðE; A0cl Þ with r r2 r3 0 1 A0cl ¼ A  i¼1 k 1i A1i  i¼1 k 02i A2i  i¼1 k 03i A3i is admissible. Step 2. Find the maximum b0 P 0 such that LMIs (23) and (25) are feasible for b = [b0, b0, . . . , b0, b0]. Step 3. Find blow 6 b0 such that LMIs (23) and (25) are 1 feasible for b ¼ ½blow 1 ; b0 ; . . . ; b0 ; b0 . Step 4. Find bupp P b such that LMIs (23) and (25) are 0 1 feasible for b ¼ ½blow ; bupp 1 ; b0 ; b0 ; . . . ; b0 ; b0 . 1 Step 5. Repeat Steps 3 and 4 such that LMIs (23) and (25) upp low upp are feasible for b ¼ ½blow 1 ; b1 ; . . . ; br0 ; br0 . Step 6. Calculate the range of kmi from (15) by upp k mi ¼ k 0mi þ ½blow mi ; bmi . The above procedure is a theoretic summary on determining stabilizing ranges of PID parameters. Yet, it leaves for many practical implementation problems, such as how to find the initial stabilizing parameter k 0mi and how to make the stabilizing ranges as large as possible. Also, it should be pointed out that different solutions may be obtained if the parameters of PID controller are reordered. To obtain a reasonable stabilizing ranges as large as possible, we add some modifications to Procedure 3.1 from the practical point of view, which will be described in detail in Section 5 and summarized as Algorithm 5.1.

l
2

X 11

6 6 . ¼ ½a1 I;    ; ar I 6 .. 4 X r1

32

3



H

..

76 7 .. 76 .. 7 6 . 7 6 0: . 7 54 5

.

   X rr

a1 I

4. Special cases

ar I

This completes the proof from Lemma 3.1. h We remark that, when E ¼ I (i.e., for standard systems), Proposition 3.1 reduces to the robust stability test for Acl to be Hurwitz. In this case, the present method is less conservative than those given in Geromel et al. [16], Peaucelle et al. [38], Ramos and Peres [39] because the conditions in (21)–(23) reduce to those in these papers when setting Xjl = 0. Proposition 3.1 provides a quasi-LMI condition to search for b. Based on Proposition 3.1, we present an LMI-based algorithm to compute ranges of PID controller gains. Note that (21) and (22) can be combined to a single T LMI. Let L ¼ ½0; 0; I n  2 Rð2nþmÞn . Then, similar to Lin et al. [28], (22) with (21) is equivalent to the following LMI for additional matrices Zj > 0 and Y j 2 Rnð2nþmÞ : Xjl þ Xlj < X jl þ

"

Procedure 3.1

where "

where

X Tjl ;

j ¼ 1; 2; . . . ; r; l 6 j;

ð25Þ

For three special cases of PID control, namely, P, PD and PI control, their transformed state-space representations are different, which leads to different LMI conditions. Hence, in this section, we would like to give such representations and conditions for these three special cases for easy reference and applications. 4.1. Proportional control In this special case, K2 = 0 and K3 = 0 in (6). Then, x_ ðtÞ ¼ AxðtÞ þ BuðtÞ; yðtÞ ¼ CxðtÞ;

ð26Þ

uðtÞ ¼ K 1 yðtÞ: or, rewritten as x_ ðtÞ ¼

A

r1 X

! k 1i A1i xðtÞ :¼ Acl xðtÞ;

ð27Þ

i¼1

where A1i ¼ BI 1i C;

i ¼ 1; 2; . . . ; r1 :

ð28Þ

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

As processed before, equivalently recast Acl as a matrix polytope with r ¼ 2r1 vertices denoted by Aj ðbÞ 2 Rnn ,

b b AðaÞ : AðaÞ ¼

^r X

b j ðbÞ; aj A

j¼1

( Acl 2

( b cl 2 A

65

AðaÞ : AðaÞ ¼

r X

^r X

aj Aj ðbÞ;

j¼1 r X

)

aj ¼ 1; aj P 0; j ¼ 1; 2; . . . ; ^r :

ð35Þ

j¼1

) ð29Þ

Therefore, in the PD control case, we have the following result.

Therefore, in the proportional control case, we have the following result.

b cl Þ is robustly admissible if b A Proposition 4.2. The pair ð E; there exist matrices Pj, Fj, Hj and Xjl with X jj ¼ X Tjj , l 6 j, j; l ¼ 1; 2; . . . ; ^r, such that

aj ¼ 1; aj P 0; j ¼ 1; 2; . . . ; r :

j¼1

Proposition 4.1. The polytope Acl is robustly stable if there exist matrices Pj > 0, Fj, Hj and Xjl with X jj ¼ X Tjj , l 6 j, j, l = 1, 2, . . . , r, such that Hjl þ Hlj < X jl þ X Tjl ; j ¼ 1; 2; . . . ; r; l 6 j; 2 3 X 11    H 6 7 6 .. .. 7 .. 6 0; 6 . . 7 . 4 5 X r1



ð30Þ

ð31Þ

P Tj E ¼ ET P j P 0; b jl þ H b lj < X jl þ H 2 3 X 11    H 6 . .. 7 .. 6 . 7 . 5 6 0; . 4 . X ^r1    X ^r^r " b jl ¼ H

" Hjl ¼

T

F j Al ðbÞ þ Al ðbÞ F Tj

H

P j  F Tj þ H Tj Al ðbÞ

H j  H Tj

# :

In this special case, K2 = 0 in (6). Let ^xðtÞ ¼ ½xT ðtÞ; y_ T T T ðtÞ and ^y ðtÞ ¼ ½y T ðtÞ; y_ T ðtÞ . Then we have the following control system: " # " # " # In 0 B A 0 ^xðtÞ þ ^x_ ðtÞ ¼ uðtÞ; 0 0 CB CA I m " # C 0 ð32Þ ^xðtÞ; ^y ðtÞ ¼ 0 Im

b A

i¼1

b 1i  k 1i A

r3 X

! b 3i ^xðtÞ :¼ A b cl^xðtÞ; k 3i A

#

b l ðbÞ þ A b l ðbÞT F T F jA j

H

b l ðbÞ P j  F Tj þ H Tj A

H j  H Tj

or, rewritten as e A

r1 X

e 1i  k 1i A

i¼1

or, rewritten as b ^x_ ðtÞ ¼ E

ð38Þ

:

this special case, K3 = R0 in (6). Let ~xðtÞ ¼ ½xT ðtÞ; R t In t T T T y ðhÞ dh and e y ðtÞ ¼ ½y T ðtÞ; 0 y T ðhÞ dh . Then we have 0 the following control system:   A 0 B _~xðtÞ ¼ ~xðtÞ þ uðtÞ; C 0 0  C 0 ð39Þ ~y ðtÞ ¼ ~xðtÞ; 0 Im uðtÞ ¼ ½ K 1 K 2 ~y ðtÞ;

~x_ ðtÞ ¼

K 3 ^y ðtÞ;

r1 X

ð37Þ

4.3. PI control

4.2. PD control

uðtÞ ¼ ½ K 1

j ¼ 1; 2; . . . ; ^r; l 6 j;

where

X rr

where

ð36Þ X Tjl ;

ð33Þ

where  A e A¼ C

0 0

;

r2 X

! e 2i ~xðtÞ :¼ A e cl~xðtÞ; k 2i A

ð40Þ

i¼1

" e 1i ¼ A

BI 1i C

0

0

0

# ;

" e 2j ¼ A

0

BI 2j

0

0

# ;

i¼1

ð41Þ where   I 0 A 0 b¼ b¼ n E ; ; A 0 0 CA I m " # " # C 0 BI 0 BI 1i 3l b 1i ¼ b 3l ¼ A ; ; A CBI 1i C 0 0 CBI 3l

ð34Þ

for i = 1, 2, . . . , r1 and l = 1, 2, . . . , r3. As proceeded before, b cl can be equivalently recast as a matrix polytope with A b j ðbÞ 2RðnþmÞðnþmÞ , ^r ¼ 2r1 þr3 vertices denoted by A

. Assume that for i = 1, 2, . . . , r1 and j = 1, 2, . . . , r2P P a set of e0 ¼ A e  r1 k 0 A e 1i  r2 k 0 A e scalars k 0mi are such that A i¼1 1i i¼1 2i 2i cl is Hurwitz stable. Set ~k mi ¼ k mi  k 0mi , m = 1, 2, i = 1, 2, . . . , rm. Then, e cl ¼ A e0  A cl

r1 X i¼1

~k 1i A e 1i 

r2 X

~k 2i A e 2i :

ð42Þ

i¼1

The following task is to compute the perturbation ranges e cl remains stable. Now, we specify for scalars ~k mi such that A

66

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

upp the lower and upper bounds for e k mi as blow mi and bmi , respeclow upp tively, and relabel them as bi and bi with i = 1, 2, . . . , e cl is equivalently recast as a matrix polyr1 + r2. Then, A e j ðbÞ 2RðnþmÞðnþmÞ , tope with ~r ¼ 2r1 þr2 vertices denoted by A ( ~r X e cl 2 AðaÞ e e e j ðbÞ; A : AðaÞ ¼ aj A j¼1 ~r X

)

aj ¼ 1; aj P 0; j ¼ 1; 2; . . . ; ~r :

ð43Þ

j¼1

Therefore, in the PI control case, we have the following result. e cl is robustly stable if there Proposition 4.3. The polytope A exist matrices Pj > 0, Fj, Hj and Xjl with X jj ¼ X Tjj ; l 6 j; j; l ¼ 1; 2; . . . ; ~r, such that e jl þ H e lj < X jl þ X T ; H jl 2 3 X 11    H 6 . .. 7 .. 6 . 7 . 5 6 0; . 4 . X ~r1    X ~r~r

j ¼ 1; 2; . . . ; ~r; l 6 j;

ð44Þ ð45Þ

where " e jl ¼ H

e l ðbÞ F T e l ðbÞ þ A F jA j

H

e l ðbÞ P j  F Tj þ H Tj A

H j  H Tj

T

# :

5. Computational algorithm In Step 1 of Procedure 3.1, scalars k 0mi ; m ¼ 1; 2; 3; i ¼ 1; 2; . . . ; rm , are determined such that the closed-loop system is stable. The selection of k 0mi will determine the location of the stabilizing range of kmi obtained in later steps. From a practical point of view, a control engineer would like the origin to be contained in the stabilizing range of kmi to facilitate control tuning if possible (this is the case if the plant is stable). The reason is that the open-loop corresponds to a zero gain controller, or the origin in the parameter space of kmi. To have a closed-loop control, a practising engineer will typically gradually increase loop gains from zero, and will have great choice of such gains and beneficial loop performance if he or she is given a sufficiently large stabilizing range containing the origin. In this section, we modify Procedure 3.1 such that the initial settings and the subsequent search for the desired stabilizing ranges are carried out in a systematic way, and the largest stabilizing range containing the origin and other ranges of interests are obtained if they are not empty. To illustrate our ideas and resulting modifications, consider again the example in Section 2, where a decentralized proportional controller, K(s) = diag{k1, k2} is applied to stabilize the plant with the transfer function matrix: " # GðsÞ ¼

1 sþ1 3 sþ1

2 sþ1 4 sþ1

in the unity negative feedback configuration. We proceed as follows. (i) Start from the simplest common gain controller, K(s) = kI2 (or k1 = k2 = k) to stabilize G(s), where k is a scalar. Let k0 be a stabilizing point. Since this plant, G(s), is stable, the origin (k0 = 0) is already a stabilizing point. (ii) Let k ¼ k  k 0 . By Barmish [3] (see Proposition 5.1 below), the stabilizing range of k is calculated as k 2 ð0:1861; 2:6861Þ. Thus, the stabilizing range in terms of k is k ¼ k 0 þ k 2 ð0:1861; 2:6861Þ. Graphically, such a stabilizing range is the straight line, BD, shown in Fig. 2. (iii) The mid-point of the above stabilizing range of k is (2.6861  0.1861)/2 = 1.25. One may reset k 01 ¼ k 02 ¼ 1:25 and calculate b0 = 1.4361 by Step 2 of Procedure 3.1 such that LMIs (30) and (31) are feasible for b = [b0, b0, . . . , b0, b0]. Thus, the initial stabilizing range of the independent gain controller, K(s) = diag{k1, k2}, is k i 2 ½k 0i  b0 ; k 0i þ b0  ¼ ½0:1861; 2:6861, i = 1,2, which contains the origin, k1 = k2 = 0. Note that so calculated k1 and k2 can now vary mutually independently within this range while the closedloop remains stable. Graphically, this range is the square, ABCD, shown in Fig. 2. One may wish to shift the range to reflect different scaling and/or importance of different gains. (iv) Before Step 3 of Procedure 3.1 is applied, we arrange k1 and k2 in decreasing order of their importance, i.e., if k1 needs to be as large as possible (most important), then resize its stabilizing range firstly, and so on. As we pointed out in Section 2 that the stabilizing range of k1 usually has effects on the stabilizing range of k2, our ordering of ki helps to obtain a stabilizing range as large as possible at the desired location. Suppose that k1 is more important than k2, one wishes to shift the initial stabilizing range, ki 2 [0.1861, 2.6861], i = 1, 2, to the new location as desired. Note that B (or D) lies on the boundary of the stabilizing region (see Fig. 2), the upper bound (or the lower bound) of ki cannot be extended. Hence, to obtain a stabilizing range as large as possible, we shift the initial stabilizing ranges of ki from the point A or C, that is, choose initial values ðk 01 ; k 02 Þ ¼ ½k 01  b0 ; k 02 þ b0  ¼ ð0:1861; 2:6861Þ or ðk 01 ; k 02 Þ ¼ ½k 01 þ b0 ; k 02  b0  ¼ ð2:6861; 0:1861Þ. Thus, LMIs (30) and (31) are feasible for b ¼ ½k 01  b0 ; k 01 þ b0 ; k 02  b0 ; k 02 þ b0   ½k 01 ; k 01 ; k 02 ; k 02  ¼ ½0; 2:8722; 2:8722; 0 or [2.8722,0,0,2.8722]. (v) Since our range shift starts from the point A or C, the corner of the square ABCD, some ‘‘0’’ items will appear in b, which implies that the upper bound or the lower bound of ki may be possibly extended so that they should be resized a priori. Moreover, b should also be relaxed as b* = ab, where a 2 (0, 1) with a = 0.5 by default. Then, the stabilizing range of k1 and k2 at such a shifted position can be calcu-

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

lated by Steps 3–6 of Procedure 3.1. For example, suppose that ðk 01 ; k 02 Þ ¼ ð0:1861; 2:6861Þ, k1 is more important than k2. Let a = 0.5 and b* = ab = [0,1.4361,1.4361,0]. We firstly find the lower bound of k1, secondly the upper bound of k2, thirdly the upper bound of k1 and finally the lower bound of k2. The resulting stabilizing range of k1 and k2 is k 1 2 ½7; 1:2561;

k 2 2 ½0:8117; 1:25:

Such a range shifting to the different location will lead to a larger stabilizing range if it exists there. One sees from the above example that our modifications to Procedure 3.1 consist of four steps. Firstly, find a common gain controller which can stabilize the plant. In the case of stable plants, the default gain is zero; Secondly, determine, by some formula, the stabilizing range of the common gain controller containing the above stabilizing point; Thirdly, compute the initial stabilizing range for independent gain controller by Step 2 of Procedure 3.1; Finally, shift it to the desired location, resize it, and compute the new range with Steps 3 to 6 of Procedure 3.1. The same technique is now applied to general cases as follows. P control. Start with the common gain controller, K(s) = kIm, or u(t) = ky(t), where k is a scalar. By (27), the closed-loop system becomes x_ ðtÞ ¼ ðA  kBCÞxðtÞ ¼ Acl xðtÞ; which is a regular system. Let k ¼ k  k 0 . Then Acl ¼ A0cl  kBC;

ð46Þ

where A0cl ¼ A  k 0 BC, and k0 is to stabilize the plant. If the plant is stable, take k0 = 0; otherwise, find such a stabilizing non-zero k0 [46]. If G(s) cannot be stabilized by any common gain controller, we have to find a controller in the general form of (6), i.e., k 0mi , m = 1, 2, 3, i = 1, 2, . . . , rm. This is the stabilization problem by static output feedback, which can be solved by a few standard techniques [7,56,27]. þ Denote by k min and kmax , respectively, the minimum negative eigenvalue and the maximum positive eigenvalue of a square matrix (set as zero if none). For regular stable systems, the stabilizing range of k is calculated from the following formula. Proposition 5.1 [3]. The matrix Acl given in (46) with a stable A0cl and an uncertain k remains robustly stable if k 2 ðk min ; k max Þ;

where k min ¼

k max ¼

k min



kþ max

1 ðA0cl



 In þ In 

A0cl Þ1 ððBCÞ

 In þ In 

A0cl Þ1 ððBCÞ

1 ðA0cl

;  I n þ I n  ðBCÞÞ :  I n þ I n  ðBCÞÞ

Here, ‘’ denotes the Kronecker product.

k 2 2 ½1:25; 202:6861:

Similarly, if the alternative ðk 01 ; k 02 Þ ¼ ð2:6861, 0:1861Þ is used, the stabilizing ranges become k 1 2 ½1:2468; 4;

67

ð47Þ

Consider now the independent gain controller. Reset k 0mi ! k 0mi þ ðk min þ k max Þ=2 and calculate the b0 with Step 2 of Procedure 3.1, i.e., find the maximum b0 P 0 such that LMIs (30) and (31) are feasible for b = [b0, b0, . . . , b0, b0]. This yields the mutually independent stabilizing range of kmi as k mi 2 ½k 0mi  b0 ; k 0mi þ b0 ; m ¼ 1; ;2; 3; i ¼ 1; 2; . . . ; m. Graphically, these initial stabilizing ranges form an m-dimension super cube. Next, to get the stabilizing range at any desired location, start the range shift from the corners of the above cube, i.e., choose k 0mi ¼ k 0mi  b0 or k 0mi ¼ k 0mi þ b0 . As the sequence of resizing the stabilizing range for each loop is important, arrange kmi in decreasing order of their importance, that is, if k11 needs to be as large as possible (most important), take it at the first place in the sequence of kmi, and so on. Suppose that kmi is arranged in decreasing order of their 0 ; . . . ; b0 ; importance as [k11, . . . , k3m], then, b ¼ ½b011 ; b 11 3m 0 0 0 0 0 0    b3m  with bmi ¼ k mi  b0  k mi and bmi ¼ k mi þ b0  k 0mi . Since we start from the corner of the super cube, some ‘‘0’’ items may appear in b and they should be tuned a priori. Moreover, b should also be relaxed as b* = ab, where a 2 (0, 1) with a = 0.5 by default. All these measures guarantee that our search will succeed in getting the stabilizing range of kmi at the desired location. Finally, following Steps 3 to 6 of Procedure 3.1, the stabilizing range of kmi is actually determined. PI control. Start with the common gain controller, K(s) = k(1 + 1/s)Im, or uðtÞ ¼ k½I m I m ~y ðtÞ, where k is a scalar. By (40), the closed-loop system becomes e  kH e cl~xðtÞ; e Þ~xðtÞ ¼ A ~x_ ðtÞ ¼ ð A ð48Þ which is a regular system, where  r1 r2 X X BC B e 1i þ e 2i ¼ e ¼ H A A : 0 0 i¼1 i¼1 Let e k ¼ k  k 0 . Then e cl ¼ A e 0  k H e; A cl e0 ¼ A e  k0 H e , and k0 is to stabilize the plant. Note where A cl that unlike the P control case, even if the plant is stable, we cannot take k0 as zero because an integrator is present here. This is the problem of so called integral controllability, and the following proposition gives the criterion for the existence of such a k0. Proposition 5.2. Consider the plant G(s) under the PI control K(s) = k0(1 + 1/s)Im. Suppose that G(s) is strictly proper, then

68

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

e0 ¼ A e  k0 H e is stable for some k0 if: (i) A cl 0 (a) k > 0 and all the eigenvalues of G(0) lie in the open right half complex plane; or, (b) k0 < 0 and all the eigenvalues of G(0) lie in the open left half complex plane. e0 ¼ A e  k0 H e is unstable for any k0 if: (ii) A cl 0 (a) k > 0 and the number of eigenvalues of G(0) in open left half complex plane is odd; or, (b) k0 < 0 and the number of eigenvalues of G(0) in open right half complex plane is odd. Proof. Let G*(s) = (s + 1)G(s). The results follow directly from Morari [34]. h If G(s) cannot be stabilized by any common gain controller, we have to find k 0mi , m = 1, 2, 3, i = 1, 2, . . . , rm, to stabilize the plant instead. The same techniques mentioned in the case of P control can be still applied to find such k 0mi . Once k0 (or k 0mi ) is determined, Barmish’s formula becomes applicable to calculate the stabilizing range of e k by replace . After that, follow the same procedure as in ing BC with H the case of P control. PD control. Start with the common gain controller, KðsÞ ¼ kð1 þ sÞI m ; or uðtÞ ¼ k½I m I m ^y ðtÞ, where k is a scalar. By (33), the closed-loop system becomes

where k1 ¼

1 ; ðCBÞ k min

k2 ¼

1 ; kþ ðCBÞ max

1 1 ; k4 ¼ þ ; 1 1 ðT T Þ ðT k k 2 min max 1 1 T 2Þ 2 3 A1  I n þ I n  A1 I n  A2 A2  I n 6 7 T1 ¼ 6 I n  A3 I n  A4 0 7 4 5; k3 ¼

2 6 T2 ¼ 6 4

A3  I n

A4  I n

0

ðBCÞ  I n þ I n  ðBCÞ

I n  ðBÞ

ðBÞ  I n

I n  ðCBCÞ

I n  ðCBÞ

0

ðCBCÞ  I n

0

ðCBÞ  I n

3 7 7: 5

Here, ‘’ denotes the Kronecker product. Once ð^k min ; ^k max Þ is found as above, let k mi ¼ k 0mi þ ^k 2 ðk 0 þ ^k min ; k 0 þ ^k max Þ. Then, follow the same procedure mi mi as in the case of P control. PID control. Start with the common gain controller, K(s) = k(1 + 1/s + s)Im, or uðtÞ ¼ k ½ I m I m I m y ðtÞ, where k is a scalar. By (12), the closed-loop system becomes  ÞxðtÞ ¼ AclxðtÞ; Ex_ ðtÞ ¼ ðA  k H

ð51Þ

which is a singular system, where

b  kH b cl^xðtÞ; b ^x_ ðtÞ ¼ ð A b Þ^xðtÞ ¼ A E which is a singular system, where b ¼ H

r1 X

b 1i þ A

i¼1

r3 X

b 3i ¼ A

i¼1



BC

B

CBC

CB

: Let k ¼ k  k 0 . Then

Let ^k ¼ k  k 0 . Then

; Acl ¼ A0cl  ^k H

b cl ¼ A b 0  ^k H b; A cl

ð49Þ

where

where " b  k0 H b0 ¼ A b ¼ A cl " :¼

A1

A2

A3

A4

#

A  k 0 BC

k 0 B

CA  k 0 CBC

I m  k 0 CB

#

b 0 Þ is admissible. If the b A and k0 is such that the pair ð E; cl 0 plant is stable, take k = 0; otherwise, use the same method as in the case of P control to determine a non-zero k0 or b0 Þ b A k 0mi ; m ¼ 1; 2; 3; i ¼ 1; 2; . . . ; m, such that the pair ð E; cl is admissible. For singular stable systems, the stabilizing range of ^k can be calculated from the following formula. Proposition 5.3 [23]. The largest interval of ^k such that the b 0  ^k H b A b g remains admissible is given by pair f E; cl ^k 2 ð^k min ; ^k max Þ ¼ ðk 1 ; k 2 Þ \ ðk 3 ; k 4 Þ;

ð50Þ

and k0 is such that the pair ðE; A0cl Þ is admissible. Note that unlike the PD control case, even if the plant is stable, we cannot take k0 as zero because an integrator is present here. If G(s) cannot be stabilized by any common gain controller, we have to find k 0mi ; m ¼ 1; 2; 3; i ¼ 1; 2; . . . ; rm , to stabilize the plant instead. Once k0 (or k 0mi ) is determined, Lee et al.’s formula becomes applicable to calculate the stabilizing range of k by replacing BC with H1, B with H2, CBC with H3, and CB with H4, respectively. After that, follow the same procedure as in the case of P control. The above development is summarized as follows.

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

Algorithm 5.1 Step 1. Find a common gain controller, K(s), to stabilize the plant, G(s). If K(s)G(s) is stable, take k0 = 0; otherwise, use any standard technique [7,56,27] to find the scalar k0. Let k 0mi ¼ k 0 . If G(s) cannot be stabilized by any common gain controller, find a controller in the general form of (6), i.e., k 0mi ; m ¼ 1; 2; 3; i ¼ 1; 2; . . . ; rm . Step 2. Let k ¼ k  k 0 or k mi ¼ k mi  k 0mi . Calculate the stabilizing ranges of k or k mi as ðk min ; k max Þ by formula of Barmish [3] for P/PI control or Lee et al. [23] for PD/PID control. Step 3. Reset k 0mi ! k 0mi þ ðk min þ k max Þ=2 and find the maximum b0 P 0 such that LMIs (23) and (25) are feasible for b = [b0, b0, . . . , b0, b0]. Obtain the mutually independent stabilizing range of kmi as k mi 2 ½k 0mi  b0 ; k 0mi þ b0 . Step 4. Arrange kmi in decreasing order of their importance and choose initial values k 0mi ¼ k 0mi  b0 or k 0 ¼ k 0 þ b0 . Thus, LMIs (23) and (25) are still mi mi 0 ; . . . ; b0 ; b 0 , where b0 ¼ feasible for b ¼ ½b0m1 ; b m1 mm mm mi 0 ¼ k 0 þ b  k 0 . k 0mi  b0  k 0mi and b 0 mi mi mi Step 5. Relax b as b* = ab, a 2 (0, 1) with a = 0.5 by  ¼ 0), find blow 6 0 (or default. If bmi ¼ 0 (or b mi mi upp bmi P 0) such that LMIs (23) and (25) are feasible for i = 1, 2, . . . , m.  6¼ 0), find blow 6 ab0 (or Step 6. If bmi 6¼ 0 (or b mi mi mi upp 0  bmi P abmi ) such that LMIs (23) and (25) are still feasible for i = 1, 2, . . . , m. Step 7. Calculate the range of kmi from (15) by upp k mi ¼ k 0mi þ ½blow mi ; bmi .

Consider a process with transfer function [6] 2 3 s1 4 6 ðs þ 1Þðs þ 3Þ s þ 3 7 7: GðsÞ ¼ 6 4 1 3 5 sþ2 sþ2 Its state-space minimum realization is of the form (5) with 2 3 1:255 0:2006 0:1523 6 7 A¼6 0:465 0:8761 7 4 1:452 5; 5:496

4:108

0:1458 6 B¼6 4 0:113

0:03811

0:4033

0:7228

" C¼

0:3773

4:28 3

7 0:06613 7 5;

Case 1. P control. Consider a common gain controller, K(s) = kI2, to stabilize the plant, G(s). Since G(s) is stable, take k0 = 0. Let k ¼ k  k 0 . By Barmish’s formula in Proposition 5.1, we compute the stabilizing range of k as k 2 ½0:5522; 1:5513. Thus, the stabilizing range of k is obtained as k ¼ k 0 þ k 2 ½0:5522; 1:5513: Reset k 01 ¼ k 02 ¼ k 0 þ ð0:5522 þ 1:5513Þ=2 ¼ 0:4995 and calculate b0 = 1.0518. Then, the stabilizing range with mutually independent gains of ki is ki 2 [0.5522, 1.5513], i = 1, 2. Suppose that k1 is more important than k2 and choose k 0 ¼ 0:5522 and k 0 ¼ 1:5513 as initial values. Then, LMIs 1 2 (30) and (31) are still feasible for b = [0.5522, 1.5513, 0.5522, 1.5513]  [0.5522, 0.5522, 1.5513, 1.5513] = [0, 2.1035, 2.1035, 0]. Let a = 0.5 and relax b as b* = ab. The sequence of range shifting is as follows: firstly find the lower bound of k1, secondly the upper bound of k2, thirdly the upper bound of k1 and finally the lower bound of k2. Fix the stabilizing range of k1 as k1 2 [1.6556,1.2] and compute the stabilizing range of k2 as upp ½blow 2 ; b2  ¼ ½2:00355; 100;

which yields the stabilizing proportional controller gain ranges as k 1 2 ½1:6556; 1:2;

k 2 2 ½0:45225; 100:

ð52Þ

7:907 8:879

upp ½blow 1 ; b1  ¼ ½3:2436; 1:94905;

which yields the stabilizing proportional controller gain ranges as k 1 2 ½3:7958; 1:39685;

k 2 2 ½0:3529; 4:0967:

ð53Þ

Comparison with Ho’s method. Ho et al. [21] gave a definition for loop’s gain margin of MIMO systems based on Gershgorin bands under the assumption that the plant is diagonally dominant. For this example, the gain margins of each loop are computed by Ho’s method as k 1 2 ½1:6556; 1:2;

k 2 2 ½0:3529; 4:0967:

One sees that when the stabilizing range for one loop is the same, the range for other loop is much more conservative by Ho’s method than ours. Furthermore, the common gain margin by Ho’s method is k 2 ½1:6556; 1:2 \ ½0:3529; 4:0967 ¼ ½0:3529; 1:2;

4:831

# :

5:273

We now use the proposed algorithm to compute the stabilizing ranges of different types of controllers. Suppose the largest available range of parameters is ±100.

If the stabilizing range of k2 is fixed to k2 2 [0.3529, 4.0967] and the stabilizing range of k1 is calculated as

6. An example

2

69

3:06

which is also more conservative than ours k 2 [0.5522, 1.5513]. Note that our method goes beyond P-control and finds the stabilizing parameter ranges for PI, PD and

70

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

PID controlllers, which is not possible by Ho’s method or Gershgorin’s theorem. Comparison with the l-analysis. Let K = diag{k1, k2} and e i Þ, where k 0 are the nominal stabilizing k i ¼ k 0i ð1 þ wi D i e gains and D i are the parameter uncertainties scaled by weights wi, i = 1, 2. To get the maximum gain ranges from the l-analysis, one may proceed as follows with knowledge of our result in (52). Set k 0i as the mid-point of the stabilizing range of ki given by (52) as k 01 ¼ ð1:6556 þ 1:2Þ=2 ¼ 0:2278 and k 02 ¼ ð0:45225 þ 100Þ=2 ¼ 49:7739. If w1 = w2, the l-analysis yields l = 1.0179, which results in e i j < 1=l ¼ 0:9825 and the the allowable perturbation j D stabilizing ranges of ki as k 1 2 ½0:4516; 0:0040;

k 2 2 ½0:8710; 98:6768:

ð54Þ

Adjusting wi will lead to different stabilizing ranges of ki from which the least conservative one is obtained as k 1 2 ½1:2494; 0:7938;

k 2 2 ½13:8371; 85:7107;

ð55Þ

which is more conservative than those in (52). Similarly, if we set k 0i with knowledge of our ranges in (53), the least conservative stabilizing ranges of ki are obtained as k 1 2 ½1:9143; 0:4847;

k 2 2 ½1:2594; 2:4844;

ð56Þ

which is still conservative than ours. It should be pointed out that all the above calculations with l-analysis have made use of the known stabilizing ranges of ki obtained by our method, and that the stabilizing ranges of ki depend not only on the weights, but also on the nominal gains. If such prior knowledge about the stabilizing ranges of ki is unknown, one has to start with some stabilizing gains determined by users, which are unlikely to be the mid-point of the actual (but unknown yet) stabilizing ranges, and the results from l-analysis will certainly become more conservative. For example, If k 0i deviate from the mid-point of the ranges in (52), say, k 01 ¼ ð1:6556Þ 3=4 þ 1:2=4 ¼ 0:9417 and k 02 ¼ ð0:45225Þ  3=4 þ 100=4 ¼ 24:6608, then the least conservative stabilizing ranges of ki in this case are k 1 2 ½1:2494; 0:6340;

k 2 2 ½13:8347; 35:4869;

ð57Þ

which is even more conservative than those in (55) indeed. It is concluded that the l-analysis gives conservative results for computing stabilizing controller gains. One reason is that the parameters of PID controllers are all real, while the l-analysis treats all systems uncertainties as complex valued. The other reason is that the actual stabilizing ranges are generally not symmetric with respect to the nominal value while the allowable perturbations in the lanalysis is always symmetric with respect to the nominal stabilizing value. The proposed method do not suffer such disadvantages and thus produce much stronger stabilizing results. Besides, the computational complexity of computing l has a combinatoric growth with the number of

parameters involved. Although practical algorithms are possible in such a case, they are very time consuming. Case 2. PI control. Consider the common gain controller, K(s) = k(1 + 1/s)I2, to stabilize the plant, G(s). Since " Gð0Þ ¼

 13 1 2

4 3 3 2

#

with eigenvalues k1 = 0.6442 < 0 and k2 = 1.8109 > 0, by Proposition 5.2, no k exists to stabilize the plant. Hence, let KðsÞ ¼ diagfk 1 ; k 2 g þ 1s diagfk 3 ; k 4 g. One sees that k 01 ¼ 1, k 02 ¼ 5, k 03 ¼ 1 and k 04 ¼ 5 can stabilize the plant. By Algorithm 5.1, we compute the stabilizing PI controller ranges as k 1 2 ½15:4516; 1:0542;

k 2 2 ½4:9458; 100;

k 3 2 ½1:0545; 0:0001;

k 4 2 ½4:9451; 5:0550:

Additionally, if we choose another stabilizer as k 01 ¼ 1, k 02 ¼ 1, k 03 ¼ 0:1 and k 04 ¼ 0:1, the stabilizing PI ranges become k 1 2 ½1:9482; 0:7263; k 3 2 ½0:1; 0:05;

k 2 2 ½0:1479; 100;

k 4 2 ½0:1; 1:1035:

Case 3. PD control. Consider the common gain controller, K(s) = k(1 + s)I2, to stabilize the plant, G(s). Since G(s) is stable, take k0 = 0. Let ^k ¼ k  k 0 . By Lee et al.’s formula in Proposition 5.3, we compute the stabilizing range of ^k as ^k 2 ð0:2361; 4:2389Þ \ ð0:2361; 1:5515Þ ¼ ð0:2361; 1:5515Þ. Thus, the stabilizing range of k is obtained as k ¼ k 0 þ ^k 2 ð0:2361; 1:5515Þ. Suppose that K(s) = (k1 + k2s)I2. Algorithm 5.1 then yields the stabilizing PD ranges as k 1 2 ½0:5522; 0:6578;

k 2 2 ½0:6577; 2:0822:

To find other possible stabilizing ranges of k1 and k2, take k0 = 5 and k0 = 5, respectively. It is easy to check that k0 in both cases can stabilize G(s). Then, we compute the stabilizing PD controller ranges as k 1 2 ½1:5515; 100;

k 2 2 ½4:2389; 100;

k 1 2 ½100; 0:5683;

k 0 ¼ 5;

k 2 2 ½100; 0:2361;

k 0 ¼ 5:

Case 4. PID control. Consider the common gain controller, K(s) = k(1 + 1/s + s)I2, to stabilize the plant, G(s). By the standard techniques [7,56,27], we obtain k0 = 5. Let k ¼ k  k 0 . By Lee et al.’s formula in Proposition 5.3, we compute the stabilizing range of k as k 2 ð0:7611; þ1Þ \ ð2:2222; 1:7203Þ ¼ ð0:7611; 1:7203Þ. Thus, the stabilizing range of k is obtained as k ¼ k 0 þ k 2 ð4:2389; 6:7203Þ.

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

Suppose that K(s) = (k1 + k2/s + k3s)I2. After descriptor transformation, we have the following closed-loop system of the form (12) Ex_ ðtÞ ¼ ðA  k 1 A1  k 2 A2  k 3 A3 ÞxðtÞ :¼ AclxðtÞ; where

2 3 A 0 0 I5 0 6 7 E¼ ; A ¼ 4 C 0 0 5; 0 0 CA 0 I 2 2 3 2 0 B 0 0 0 B 6 7 6 A2 ¼ 4 0 0 0 5 ; A3 ¼ 4 0 0 0 

0 CB 0

ð58Þ

2

3 BC 0 0 6 7 A1 ¼ 4 0 0 0 5 ; CBC 0 0 3 7 5:

0 0 CB

k 01

For ¼ k 02 ¼ k 03 ¼ k 0 , the pair ðE; A0cl Þ is A0cl ¼ A  k 01 A1  k 02 A2  k 03 A3 . To find k i ¼ k i  k 0 ; i ¼ 1; 2; 3; such that i

admissible where the ranges of

Ex_ ðtÞ ¼ AclxðtÞ ¼ ðA0cl  k 1 A1  k 2 A2  k 3 A3 ÞxðtÞ;

ð59Þ

upp is robustly admissible, we let k i 2 ½blow i ; bi . Then, Acl is equivalently recast as a matrix polytope with 8 vertices low low A1 ðbÞ ¼ A0cl  blow 1 A 1  b2 A 2  b 3 A 3 ; low low A2 ðbÞ ¼ A0cl  bupp 1 A 1  b 2 A 2  b3 A 3 ; upp low A3 ðbÞ ¼ A0cl  blow 1 A 1  b2 A 2  b 3 A 3 ;

upp upp A5 ðbÞ ¼ A0cl  blow 1 A 1  b2 A 2  b 3 A 3 ; upp low A6 ðbÞ ¼ A0cl  bupp 1 A 1  b 2 A 2  b3 A 3 ; upp low A7 ðbÞ ¼ A0cl  bupp 1 A 1  b 2 A 2  b3 A 3 ; upp upp A8 ðbÞ ¼ A0cl  bupp 1 A 1  b 2 A 2  b3 A 3 :

Reset k 01 ¼ k 02 ¼ k 03 ¼ k 0 þ ð0:7611 þ 1:7302Þ=2 ¼ 5:4796 and calculate b0 = 1.1342, the stabilizing range with the mutually independent gains of ki is ki 2 [4.3454, 6.6138], i = 1,2,3. Suppose that the importance of ki is in order of k1, k2 and k3. Let a = 0.5 and choose initial values k 0 ¼ k 0 ¼ k 0 ¼ 6:6138. By Algorithm 5.1, we compute 1 2 3 upp ½blow 1 ; b1  ¼ ½5:0794; 100; upp ½blow 2 ; b2  ¼ ½1:1535; 100; upp ½blow 3 ; b3  ¼ ½1:3316; 100;

which yields the stabilizing PID ranges as k 1 2 ½1:5344; 100;

k 2 2 ½5:4603; 100;

troller parameter ranges. Numerical examples have been given to illustrate the use of the present procedure. It has been seen that the stabilizing ranges obtainable from our procedure is large and sufficient for practical tuning purpose. It should be pointed out that our algorithm provides stabilizing ranges in ‘‘sufficiency’’ sense only. In other words, the ranges are not necessary since there might be some PID settings which are not in the computed ranges but yet stabilize the plant. Besides, if the plant has time delay, one may use Pade approximation for time delay so as to apply our procedure. For instance, consider the example in Section 6, and suppose that a time delay L = 0.4 is present at the second loop, that is 2 3 s1 4 0:4s e 6 ðs þ 1Þðs þ 3Þ s þ 3 7 7: GðsÞ ¼ 6 4 1 3 0:4s 5 e sþ2 sþ2 Using the Pade approximation eLs  (1  Ls/2)/(1 + Ls/2), G(s) is approximated as 2 3 s1 4ð1  0:2sÞ 6 ðs þ 1Þðs þ 3Þ ðs þ 3Þð1 þ 0:2sÞ 7 b 7: GðsÞ ¼6 4 1 3ð1  0:2sÞ 5 sþ2

upp low A4 ðbÞ ¼ A0cl  blow 1 A 1  b2 A 2  b 3 A 3 ;

k 3 2 ½5:2821; 100:

7. Conclusions The problem of determining the parameter ranges of stabilizing multi-loop PID controllers has been investigated in this paper. A detailed scheme has been proposed using the descriptor model approach. Linearly parameter-dependent technique and convex optimization method have been employed to establish basic criteria for computing the con-

71

ðs þ 2Þð1 þ 0:2sÞ

For a proportional control K(s) = diag{k1, k2}, the stabilizing region for k1 and k2 are calculated as: k1 2 [1.7174, 0.1079] and k2 2 [0.5580, 0.8910]. The time-delay case for our problem will lead to a different system description. The feedback of delay output gives rise to a more complicated state equation, for which the stabilizing ranges of PID parameters may not be transformed into a polytopic problem, whereas the technique used in this paper is suitable for a polytopic problem. One needs to find a totally different technique to solve the delay problem, which could be a future study. References ˚ stro¨m, T. Hagglund, PID Controllers: Theory, Design, and [1] K.J. A Tuning, second ed., Instrument Society of America, Research Triangle Park, NC, 1995. [2] K.J. Astrom, H. Panagopoulos, T. Hagglund, Design of PI controllers based on non-convex optimization, Automatica 34 (5) (1998) 585–601. [3] B.R. Barmish, New Tools for Robustness of Linear Systems, Macmillan, New York, 1994. [4] J. Bernussou, P.L.D. Peres, J.C. Geromel, A linear programming oriented procedure for quadratic stabilization of uncertain systems, Systems and Control Letters 13 (1) (1989) 65–72. [5] S.P. Boyd, L. El Ghaoui, E. Feron, V. Balakrishnan, Linear Matrix inequalities in system and control theory, Siam Studies in Applied Mathematics, Philadelphia, PA, 1994. [6] G.F. Bryant, L.F. Yeung, Multivariable Control System Design Techniques: Dominance and Direct Methods, John Wiley & Sons, Chichester, England, 1996. [7] Yong-Yan Cao, James Lam, You-Xiam Sun, Static output feedback stabilization: an ILMI approach, Automatica 34 (12) (1998) 1641– 1645.

72

Q.-G. Wang et al. / Journal of Process Control 17 (2007) 59–72

[8] Ce´sar A.R. Crusius, Alexandre Trofino, Sufficient LMI conditions for output feedback control problems, IEEE Transactions on Automatic Control 44 (5) (1999) 1053–1657. [9] L. Dai, Singular Control Systems, Springer-Verlag, Berlin, 1989. [10] Aniruddha Datta, Ming-Tzu Ho, Shankar P. Bhattacharyya, Structure and Synthesis of PID Controllers, Springer-Verlag, London, 2000. [11] J. Doyle, Analysis of feedback systems with structured uncertainty, IEE Proceedings, Part D 129 (6) (1982) 242–250. [12] C.G. Economou, M. Morari, Internal model control 6: multiloop design, Industrial and Engineering Chemistry Process Design and Development 25 (1986) 411–419. [13] H.W. Fung, Qing-Guo Wang, T.H. Lee, PI tuning in terms of gain and phase margins, Automatica 34 (9) (1998) 1145–1149. [14] P. Gahinet, A. Nemirovski, A.J. Laub, M. Chilali, LMI Control Toolbox for Use with Matlab, Math Works, Natick, MA, 1995. [15] P.J. Gawthrop, Self-tuning PID controllers: algorithm and implementation, IEEE Transactions on Automatic Control 31 (1986) 201–209. [16] J.C. Geromel, M.C. de Oliveira, L. Hsu, LMI characterization of structural and robust stability, Linear Algebra and its Applications 285 (1998) 69–80. [17] L. Guo, H. Wang, PID controller design for output PDFs of stochastic systems using linear matrix inequalities, IEEE Transactions on Systems, Man and Cybernetics – Part B: Cybernetics 35 (1) (2005) 65–71. [18] C.C. Hang, K.K. Tan, S.L. Ong, A comparative study of controller tuning formulae, IEE Proceedings, Part D 138 (2) (1979) 111–118. [19] Y. He, M. Wu, J.-H. She, Improved bounded-rela-lemma representation and H1 control of systems with polytopic uncertainties, IEEE Transactions on Circuits and Systems-II: Express Briefs 52 (7) (2005) 380–383. [20] W.K. Ho, C.C. Hang, L.S. Cao, Tuning of PID controller based on gain and phase margin specifications, Automatica 31 (3) (1995) 497– 502. [21] Weng Khuen Ho, Tong Heng Lee, Oon P. Gan, Tuning of multiloop PID controllers based on gain and phase margin specifications, Industrial and Engineering Chemistry Research 36 (6) (1997) 2231– 2238. [22] V. Kucera, C.E. de Souza, A necessary and sufficient condition for output feedback stabilizability, Automatica 31 (9) (1995) 1357–1359. [23] Li Lee, Chun-Hsiung Fang, Jer-Guang Hsieh, Exact unidirectional perturbation bounds for robustness of uncertain generalized statespace systems: Continuous-time case, Automatica 33 (10) (1997) 1923–1927. [24] Yanlin Li, E. Bruce Lee, Stability robustness characterization and related issues for control systems design, Automatica 29 (2) (1993) 479–484. [25] C. Lin, J.-L. Wang, C.B. Soh, Maximum bounds for robust stability of linear uncertain descriptor systems with structured perturbations, International Journal of Systems Science 34 (7) (2003) 463–467. [26] C. Lin, J. Lam, J.-L. Wang, G.H. Yang, Analysis on robust stability for interval descriptor systems, Systems and Control Letters 42 (2001) 267–278. [27] C. Lin, Q.-G. Wang, T.H. Lee, An improvement on multivariable PID controller design via iterative LMI approach, Automatica 40 (2004) 519–525. [28] C. Lin, Q.-G. Wang, T.H. Lee, Robust normalization and stabilization of uncertain descriptor systems with norm-bounded perturbations, IEEE Transactions on Automatic Control 50 (4) (2005) 515–520. [29] A.P. Loh, C.C. Hang, C.X. Quek, V.U. Vasnani, Autotuning of multiloop proportional-integral controllers using relay feedback, Industrial and Engineering Chemistry Research 32 (1993) 1102–1107. [30] W.L. Luyben, Simple method for tuning SISO controllers in multivariable systems, Industrial and Engineering Chemistry Process Design and Development 25 (1986) 654–660.

[31] A.G.J. MacFarlane, J.J. Belletrutti, The characteristic locus design method, Automatica 9 (1973) 575–588. [32] J.M. Maciejowski, Multivariable Feedback Design, Addison-Wesley, Reading, MA, 1989. [33] I. Masubuchi, Y. Kamitane, A. Ohara, N. Suda, H1 control for descriptor system: a matrix inequalities approach, Automatica 33 (1997) 669–673. [34] M. Morari, Robust stability of systems with integral control, IEEE Transactions on Automatic Control 30 (6) (1985) 574–577. [35] O.D.I. Nwokah, The stability of linear multivariable systems, International Journal of Control 37 (1983) 623–629. [36] Osita D.I. Nwokah, Ronald A. Perez, On multivariable stability in the gain space, in: Proceedings of the 29th conference on Decision and Control, Honolulu, Hawaii, 1990, pp. 328–333. [37] Z.J. Palmor, Y. Halevi, N. Kradney, Automatic tuning of decentralized pid controller for tito process, Automatica 31 (1993) 1001–1010. [38] D. Peaucelle, D. Arzelier, O. Bachelier, J. Bernussou, A new robust D-stability condition for real convex polytopic uncertainty, Systems Control Letters 40 (2000) 21–30. [39] D.C.W. Ramos, P.L.D. Peres, An LMI condition for the robust stability of uncertain continuous-time linear systems, IEEE Transactions on Automatic Control 47 (4) (2002) 675–678. [40] H.H. Rosenbrock, State Space and Multivariable Theory, John Wiley, New York, 1970. [41] M.G. Safonov, M. Athans, A multiloop generalization of the circle criterion for stability margin analysis, IEEE Transaction on Automatic Control 26 (2) (1981) 415–422. [42] D.E. Seborg, T.F. Edgar, D.A. Mellichamp, Process Dynamics and Control, John Wiley & Sons, New York, 1989. [43] F.G. Shinskey, Process-control Systems. Application, Design and Tuning, MaGraw-Hill, New York, 1988. [44] G.J. Silva, A. Datta, S.P. Bhattacharyya, PID Controllers for Timedelay Systems, Birkhauser, Boston, 2005. [45] S. Skogestad, I. Postlethwaite, Multivariable Feedback Control: Analysis and Design, second ed., John Wiley & Sons, 2005. [46] Alexandre Trofino-Neto, Vladimir Kucera, Stabilization via static output feedback, IEEE Transactions on Automatic Control 38 (5) (1993) 764–765. [47] J. Tsinias, N. Kalouptsidos, Output feedback stabilization, IEEE Transactions on Automatic Control 35 (8) (1990) 951–954. [48] Qing-Guo Wang, Book review of ‘structure and synthesis of PID controllers’ by A. Datta, M.-T. Ho and S.P. Bhattacharyya, Automatica 39 (4) (2003) 758–759. [49] Qing-Guo Wang, Decoupling Control, Springer-Verlag, New York, 2003. [50] Qing-Guo Wang, Book review: ‘PID controllers for time-delay systems’ by G.J. Silva, A. Datta and S.P. Bhattacharyya, SIAM Review 47 (4) (2005) 855–856. [51] Qing-Guo Wang, Biao Zhou, Tong-Heng Lee, Qiang Bi, Auto-tuning of multivariable PID controllers from decentralized relay feedback, Automatica 33 (3) (1997) 319–330. [52] Qing-Guo Wang, H.W. Fung, T.H. Lee, PID tuning for improved performance, IEEE Transactions on Control Systems Technology 7 (4) (1999) 457–465. [53] Qing-Guo Wang, H.W. Fung, Y. Zhang, PID tuning with exact gain and phase margins, ISA Transactions 38 (1999) 243–249. [54] O. Yaniv, Synthesis of uncertain MIMO feedback systems for gain and phase margins at different channel breaking points, Automatica 28 (5) (1992) 1017–1020. [55] Zhang Yu, Qing-Guo Wang, K.J. Astrom, Dominant pole placement for multi-loop control systems, Automatica 38 (7) (2002) 1213–1220. [56] F. Zheng, Q.-G. Wang, T.H. Lee, On the design of multivariable PID controllers via LMI approach, Automatica 38 (2002) 517–526.