225
Ind. Eng. Chem. Res. 1995,34, 225-236
Geometric Interpretation and Measures of Dynamic Interactions in Multivariable Control Systems Shyh-Hong Hwangt Department of Chemical Engineering, National Cheng-Kung University, Tainan, Taiwan 70101, Republic of China
A technique is developed to evaluate the interactions of a multivariable control system. The technique can be employed t o select the best pairing configuration of the control system, but more than that, it is capable of providing guidelines for controller design. The technique arises from a geometric analysis of the stability boundaries in the proportional gain space. As a result, the relative gain array (RGA) is given a geometric interpretation and demonstrated to be a convenient steady-state interaction measure. A convenient dynamic counterpart is then quantified to predict dynamic interactions based on simple response experiments. The validity of the proposed dynamic interaction measures is verified by numerous literature examples.
Introduction
dn
Interactions arising in a multivariable environment produce several undesirable effects for feedback design. Cross-coupling between loops prevents the control engineer from designing each loop independently. Changes in controller parameters of one loop affect the performance of another, sometimes to the extent of destabilizing the entire system. Where interactions are significant, the experienced control engineer may possibly arrive at a successful feedback design via adroit tuning of the control loops but not without a certain loss in performance. Consequently,proper pairing of the input/ output variables to minimize interactions is of paramount importance. To achieve this, the relative gain array (RGA) (Bristol, 1966) has long been a popular variable pairing technique because of its simplicity. The RGA is based on a steady-state analysis, but it is often used to infer dynamic properties of the multiloop control system qualitatively (Bristol, 1977). However, this qualitative information is not suficient for feedback design and sometimes gives an incorrect estimate of dynamic interactions (Witcher and McAvoy, 1977; McAvoy, 1983; Friedly, 1984). The shortcoming of the RGA led to the development of many dynamic interaction measures; many of them have been discussed in detail by Jensen et al. (1986)and Grosdidier and Morari (1986). These measures usually require detailed dynamic models to analyze dynamic interactions on a quantitative basis. Unfortunately, they are seldom intended t o interpret clearly the way in which interactions affect the dynamics of the closed-loop system to facilitate feedback design. Often, it is not worthwhile to develop accurate dynamic models for the process industry. Hence, interaction measures derived from detailed dynamic models would have limited applications. In practice, fast engineering judgment in pairing the inputioutput variables of a plant should be based on easily obtainable information. For example, steady-state gains, ultimate gains and frequencies, and closed-loop peak response data can be measured quite easily even under noisy conditions and for nonlinear systems (Yuwana and Seborg, 1982; McDermott et al., 1984; Astrom and Hagglund, 1988). While these pieces of information are insufficient to construct the full model, they have been traditionally applied to measure steady-state interact E-mail
address:
[email protected].
Figure 1. Feedback structure of the multiloop control system.
tions (RGA),determine PID controller settings (ZieglerNichols rules), or construct a simple first-order plus dead-time model (Yuwana and Seborg, 1982). It is desirable t o seek dynamic interaction measures which are merely based on such information. In this study, a technique is developed to evaluate the interactions of a multiloop control system. The technique can be employed to select the best pairing configuration of the control system, but more than that, it is capable of providing guidelines for controller design. The technique arises from a geometric analysis and requires only the aforementioned information to offer quantitative measures of dynamic interactions. The RGA, which estimates steady-state interactions, is also given a clear geometric interpretation. Derivation of dynamic interaction measures is based on an analysis of stability boundaries in the proportional gain space. In particular, it is demonstrated that interactions are related t o the variations of stability boundaries at the intersections with the gain axes of all loops.
Closed-Loop Interactions Consider an N-input N-output multivariable process, y(s) = G,(s) u(s)
+4s)
(1)
where y(s) is the output vector, u(s)the input vector, and d(s)the disturbance vector. The closed-loop system consists of N single-input'single-output (SISO) loops which stipulate a diagonal controller G,(s),
[G,(s)lG = g,i(SM,j.
(2)
as depicted in Figure 1. The matrix P is a permutation matrix which determines the pairing configuration of the closed-loop system. For instance, if P is an identity matrix, the inputioutput pairing is (ui- yi). The closedloop system is described by
0888-5885/95/2634-0225$09.00/0 0 1995 American Chemical Society
226 Ind. Eng. Chem. Res., Vol. 34, No. 1, 1995
+ GpPG,)-lGpPG,y, + (I + GPPGc)-'d = (I + GG,)-lGG,y, + (I + GGJ-ld (3)
y = (I
where y,(s) is the reference vector and G(s) [=GpP(s)l represents the plant with a pairing configuration in question. The transfer function matrix relating the disturbances to the outputs,
Td(S) = (I
+ GG,)-l
(4)
Nyquist Locus
Stability Boundaries on k-axis
mRe A
I
Im
_- 1
is sometimes referred to as the sensitivity matrix and it governs disturbance rejection properties of the closedloop system. Alternatively, the complementary sensitivity matrix,
T(s) = (I
+ GG,)-lGG,
(5)
measures the sensitivity of the references to the outputs. One possible approach for analyzing closed-loop interactions would be to define a measure for the diagonal dominance of Td and T for all frequencies. However, such techniques usually require accurate dynamic models and are too complicated to reveal the way in which interactions affect the closed-loop dynamics in a simple manner. The present analysis views interactions from another perspective. An important aspect is to obtain closedloop information about a system from steady-state gains, ultimate gains and frequencies, and peak response data. The idea is first illustrated with a singleloop system. The closed-loop poles are defined by
1+gg,(s)= 0
(6)
Assume that g(s)is a stable transfer function which can be destabilized by increasing the proportional gain as stipulated by the Ziegler-Nichols tuning method. The steady-state gain g(0) is a real number, and this indicates that, under proportional control, there exists a proportional gain ks,
k,
= -l/g(O)
(7)
such that zero is a closed-loop pole at this setting. This corresponds to a simple bifurcation point and has definite nonlinear implications (Chang and Chen, 1984; Hwang and Chang, 1986). Also from eq 6, if Im(g(jw)} vanishes at a particular nonzero frequency, WH, a pair of purely imaginary closed-loop poles appear, under P control, at
kH
-l/g(jW,)
(8)
This is a Hopf bifurcation point and corresponds to the ultimate gain of the Ziegler-Nichols rules with the ultimate frequency W H . It is clear that a unique ks exists on the k-axis in the gain space but there may be multiple or no k~ for a given system. This can be conveniently seen from the open-loop Nyquist locus of g(jo). The intercept of the locus with the real axis corresponds to ( - l / k s ) or (-l/ KH). As is seen from Figure 2, the stability boundaries bounding the origin of the k-axis correspond to the two outermost intercepts of the Nyquist locus. As is wellknown, the Nyquist locus of a high-order or time-delay system will usually intercept the real axis with a nonzero o. Such a system then possesses at least one ultimate gain k ~ .
Figure 2. Open-loop Nyquist loci and their equivalent stability diagrams under proportional control.
The multiloop case is now treated. The closed-loop poles are determined by
det[I
+ GG,(s)l = 0
(9)
For purely proportional controllers, G, = K,
one can plot the bifurcation loci I' (which consist of all simple and Hopf bifurcation points, KS and KH) in the N-dimensional proportional gain space. These loci would form stability boundaries that demarcate the stable and unstable regions in the gain space and are defined by det[I
+ G(jw)Kl= 0
(11)
To distinguish between the two components of I' arising from ks and k ~ one , can define
(12) Along Ts, the leading (dominant) pole is zero (s = 0) and it may form static (simple) stability boundaries and rH may form oscillatory (Hopf) stability boundaries associated with a pair of purely imaginary poles (s = hjjw~). This implies that Ts contains important information about steady-state interactions, while rH contains that about dynamic interactions. Two examples of the 2 x 2 case are depicted in Figure 3. Figure 3a corresponds to a totally decoupled (912 = g 2 1 = 0) or a partially decoupled (gl2 or g 2 1 = 0) system. The simple stability boundaries (Ts) intersect the gain axes at k$ and kE2, which represent the inverse steady-state gains of the diagonal elements of G(s), -gii-l(O) (see eq 7). The Hopf stability boundaries (rH) intersect the gain axes at
Ind. Eng. Chem. Res., Vol. 34,No. 1, 1995 227 k2
4
Double
Zero
(b) Interacting Loops
(a) Noninteracting Loops
Figure 3. Two examples of 2 x 2 systems illustrating the evolution of stability boundaries under interactions.
kyl and kf,, each of which is the ultimate gain of a freely standing loop (with the other loop open) with the ultimate frequency UI: as given by (see eq 8)
(13) Note that for the decoupled case the simple and Hopf stability boundaries are straight lines perpendicular t o the gain axes. This implies that the inverse steady-state gain, ultimate gain, and ultimate frequency of each individual loop are invariant with changes in the controller gains of the other loop. Thus, as far as stability is concerned, the system is regarded as being free of interactions. When g12 and g21 f 0, the stability boundaries are no longer straight lines as seen in Figure 3b. Given the nominal settings & I and & p , one can draw two straight lines parallel t o the gain axes through the nominal settings. The line parallel to the ki axis intersects the simple stability boundary at k: and the Hopf stability boundary at kH with the ultimate frequency of. It appears that k,3! , k:, and W: with loopj closed can be quite different from their values for the freely standing S H H loop, i.e., kiL, kii, and wii. Consequently, the interaction present in each loop can be uncovered to a certain extent by examining variations in these quantities along the stability boundaries. It is postulated that one may quantify the steady-state interaction of loop i at the nominal gain setting of loop j, kj, via the relative difference between ks and k: and the dynamic interaction via the relative difference between k: and kf andl or that between UI? and a:. For instance, one may quantify the steady-state and dynamic intercctions of loop i at the nominal gain setting of loopj, kj, as the ratios:
The above statements for steady-state and dynamic interaction estimates are not very convenient to use since the nominal gain settings are not necessarily given before the interaction problem is analyzed. Instead, we shall consider two extreme cases, namely kj approaching 0 and for general estimates of interactions. The geometric approach is intended to decipher the way in which interactions affect ks, k:, and m y of the interacting loops and develop subsequently measures of steady-state and dynamic interactions based on easily obtainable response data. It is adopted as a convention that gii(0) is positive for all i. Loops with negative steady-state gains can be transformed to positive ones by redefining input or output variables. With this convention, k: is always negative. Hence, negative feedback in loop i corresponds to a positive ki and positive feedback corresponds t o a negative ki. Insofar as positive feedback is usually avoided, the most important region in the gain space is the first quadrant where all gains are positive and all loops are under negative feedback. QO,
Steady-State Interaction Measures and the Topology of rs
To conveniently illustrate the topological features of Ts and its relationship with steady-state interactions, a 2 x 2 system is first examined. The steady-state version of eq 11 for 2 x 2 systems is
[1 + k1g11(0)1[1 + h&22(0)1 - k1kg,,(O) g21(0) = 0 (14) Introducing scaled gains, each normalized with respect to lk3,
El
= k1/1kSl(
(154
62
=k d l k 3
(15b)
steady-state interaction estimate of loop i =
kf at kj 12: with loopj open (=k:) dynamic interaction estimate of loop i =
kr at
k? with loopj open (=kt)
the diagonal elements of the 2 x 2 relative gain array A become
228 Ind. Eng. Chem. Res., Vol. 34, No. 1, 1995 rS
Unstable
I
Unstable
Stable
Unstable
)
Unstable
Unstable
Unstable
Stable
Unstable
't
Unstable
Unstable
proportional control systems. Figure 4d shows that a negative A stipulates the appearance of a simple curve in the 'negative feedback region (the first quadrant where kl, k z > 0). This greatly narrows the stable region of this practically important quadrant and deteriorates control quality substantially since a slow real pole is introduced by Ts to this quadrant. Furthermore, with the introduction of integral action, a system with a negative A will become unstable as shown by Grosdidier et al. (1985). It is proposed to estimate steady-state interactions by means of variations in the scaled inverse steady-state gain of loop i, 17, along rs. Since &: is generally dependent on the nominal gain of loop j, a precise measure of the steady-state interaction that loop j imposes on loop i can be defined only for a particular gain setting of loop j. For simplicity, one merely considers two extreme cases _in the negative feedback region of loopj (kj > 01, i.e., kj approaching 0 and -: 1. kj 0, loopj is marginally open: It is assumed that the effect of the controller of the marginally open loopj on the inverse steady-state gain of loop i, &;, can be estimated by the slope of Ts at the intersection with the ki axis. The slope, which is zero for perpendicular intersection, is calculated from eq 17 as
-
1
(c) x > I (d) k < O Figure 4. Topology of static stability boundaries as a function of
1.
The off-diagonal terms of A are (1 - A). Upon rearrangement, eq 14 reduces to
( E , + A)(&
+ A) = A(A - 1)
[;](Ei=-l,&j=o)
= *-1 - 1
A
, i , j = 1, 2 0' f
Z)
(17)
which defines Ts in the scaled gain space of & I and &Z for a given A. It is a canonical equation for a hyperbola, and four possible topologies are depicted in Figure 4. When A = 1, implying no steady-state interactions according to the RGA, the simple stability boundaries (rs)are two straight lines, orthogonal to the gain axes as already indicated in Figure 3a where g 1 2 andlor gzl = 0. These two lines intersect at a double-zero singularity; viz., two zero poles exist a t this point. For 0 < A < 1, two simple curves are split from the intersection; a Hopf curve must appear t o bound the stable region as depicted in Figure 4b. This resultant Hopf curve ends at each simple curve t o form a new double-zero singularity. This split of the double-zerosingularity as A is decreased from 1 can be demonstrated rigorously by modern dynamic singularity theory (Guckenheimer and Holmes, 1983). The additional rH that is created is a particular "unfolding" of a double-zero singularity. However, one can argue its existence on more intuitive and physical grounds. As one changes A from 1 to a smaller value, the stability of the system should only change by an equally small extent; viz., the stability boundaries should only shift slightly. Hence, the stability regions in the first and third quadrants should retain their stability assignments. Consequently, one must cross a stability boundary as one passes diagonally in the gain space from the third quadrant to the first. This stability boundary cannot be Ts since all simple loci are defined by eq 17 and shown in Figure 4b. Hence, a rH must appear to demarcate the stability region. For 1 > 1, one can show that a simple curve serves as the stability boundary while the other appears in the unstable region as illustrated in Figure 4c. There is no change in the stability of the system (it remains unstable) as one crosses the second curve (hatched in Figure 4c). The case of 1 < 0, corresponding t o negative steady-state interactions, exhibits a peculiar feature for
where ~ ( 0 is ) the Rijnsdorp interaction quotient, ~ ( s ) , evaluated at s = 0 for steady state. When ~ ( 0= ) 0, the simple curve Ts is perpendicular to the ki axis and the inverse steady-state gain of loop i is scarcely affected by clcsing loop j. 2. kj -, loopj is under perfect control: From eq 17, the large gain asymptote of Ts is found to be
-
&!=-A
for
5-m
(19)
A convenient steady-state interaction measure can then be defined as
k: with loopj perfect (=-A) =A k: with loopj open (=-1)
(20)
It is interesting to note that this definition leads to the diagonal elements of the 2 x 2 RGA, 2. When 2 = 1, the system is free of steady-state interactions for large gains. When A < 0, the asymptote will appear in the negative feedback region (5: > 0)as indicated by eq 19, causing a catastrophic effect on the system dynamics as explained previously. Equations 18 and 20 show that the two commonly ) A, used steady-state interaction measures, ~ ( 0 and represen5 respectively the slope of Ts corresponding to loop i at kj = 0 and its asymptote for large kj, a geometric interpretation that is different from the original motivations for their construction. Since the two measures are not independent ( 4 0 ) = (A - 1)12), either one can be employed to estimate steady-state interactions. The present analysis then provides another justification of the RGA. It should be emphasized that the two measures, eqs 18 and 20, are symmetric in that the steady-state effect
Ind. Eng. Chem. Res., Vol. 34, No. 1, 1995 229 of controller 1on loop 2 is identical to that of controller 2 on loop 1. Moreover, since Ts does not exist in the first quadrant (negative feedback region) of the gain space except for the prohibitive case of 1 0, the steadystate interaction measures based on rs provide only indirect information regarding the closed-loop dynamics. The definition of eq 20 can be extended readily to N x N systems. The steady-state interaction measure of loop i is calculated from eq 11with all kj (j t i) 00 as
-
Ai = [G(0)Iii[G(O)Tlii-l
(21)
which is again a diagonal element of the N x N RGA.
Dynamic Interaction Measures For a complete dynamic interaction analysis, the Hopf stability boundaries ( r H ) must be overlaid on the simple stability boundaries as shown in Figure 3. Unfortunately, the topology of the Hopf loci under interactions is not as simply classifiable as that of the simple loci. As a result, the construction of the complete Hopf loci for dynamic interaction analysis is more difficult, if not impossible, even with full knowledge of the transfer function matrix. On the contrary, the slopes of rH at the gain axes and the associated frequency changes can be accurately estimated from easily obtainable response data as will be expounded later. They reveal important features of the Hopf stability boundaries and represent convenient estimates of dynamic interactions. Again, one starts with a 2 x 2 system. The Hopf loci are defined by
[1 + k,g,,0'w)l[l
+ kg220'dl k,kgg,,O'w) gp10'w) = 0 (22)
Scaling the gains with respect t o the ultimate gains H
kii 9
k1 = k,/kyl k, = k4kf2
(234 (23b)
one can define, analogous to the analysis of steady-state interactions, the following dimensionless version of the slope along r H at the ki axis:
s, = (&:/&j)(ki=l;kj=o);
i , j = 1,2 0'f i)
where
si,
The sign and magnitude of each may provide certain information regarding integrity and stability robustness of the multiloop control system. The associated rpte of change of the ultimate frequency of loop i with kj can be obtained as
_wyl
where df = oFlwf. Note khat as and are different, the two measures S1p and So1 are seldom or never identical, and neither are Q1z and Qp1. This asymmetric property of dynamic interactions is a distinct difference from the steady-state case. Geometrically, this implies that rH, doesAnot possess reflective symmetry across the line k l = kz in the gain space as Ts. In other words, the relative effects of one loop on the other may not be the same. It is possible that one loop is decoupled from the other but not vice versa. This unique frequency dependence has been observed by McAvoy (1983) and O'Reilly and Leithead (1993). The dynamic interaction present in each loop is quantified in this work as
Di = Wlls12I + w21Q121
(26a)
= w11~2,I + w21Q211
(26b)
D2
where the measures &jand S2g represent, respectively, the gain and frequency variations caused by interactions. Choice of the weighting functions w1 and wz reflects the relative effects of the two on the coupled loops. If such information is lacking, it is suggested to use w1 = wp = 0.5 to place equal penalties on the two measures. (Niederlinski (1971) used the product of the fundamental critical loop gain and the corresponding critical frequency as a measure of control susceptibility for a multiloop system by assuming implicitly that the loops are affected equally by the gain and frequency.) Often, the overall magnitude of dynamic interactions is of primary concern; then it is convenient to define the dynamic interaction index (DIX) as
DIX=D, + D 2
(27)
Upon the use of eqs 24-27, the problem of measuring dynamic interactions reduces to the inspection of the slopes and associated frequency changes of the Hopf stability boundaries at the gain axes. This consideration becomes more persuasive by introducing the idea of dominant poles. In general, the closed-loop dynamics are largely determined by sets of dominant poles (Hwang and Chang, 1987;Astrom and Hagglund, 1988). These poles destabilize at ultimate controller gains which form the Hopf stability boundaries in the gain space. Thus, it is obvious that the examination of the Hopf stability boundaries is able to provide certain information about the dominant poles and then accounts directly for dynamic properties of the system. For example, the Hopf stability boundaries of a totally decoupled or a partially decoupled system (glz andlor gal = 0) consist of two disjoint, orthogonal lines as depicted in Figure 3a. This implies that the dominant poles of one loop tend to be independent ofJhe controller setting of the other loop. The measures Sc and SZv are Eimply zeros, indicating this noninteracting feature. As Sg and SZu increase in magnitude, dynamic interactions, causing variations in the dominant poles, are then known to increase accordingly. The relationship of
230 Ind. Eng. Chem. Res., Vol. 34,No. 1,1995 these measures with the Hopf stability boundaries will be thoroughly investigated in simulation work. Since the exact transfer function matrix G(s) is generally unknown, the measures Sij and S2, are preferably derived from response experiments. To this end, each element of G(s)is approximated with a secondorder plus dead-time transfer function,
which can be constructed directly from easily obtainable response data as will be elaborated in a later section. Substitution of eq 28 in eqs 24 and 25 yields
This implies that if of the first loop is much larger than that of the second loop, the dominant poles and dynamic performance of the first loop could be relatively independent of the controller setting of the second loop (D1 is much smaller than D2) but not vice versa. Hence, the first loop with a sufficiently higher uFl. may be regarded as a decoupled loop dynamically. This unique ability of distinguishing the relative dynamic interaction effects between two loops is indeed an advantage over the other measures, which consider the system interactions as a whole. The above dynamic interaction results can be e?tended readily to an N x N system. The measures SQ and become
S , = (dy/&j)(,&i= l;kj = 0; all other loops open); i,j = 1,N (i f i) (30a) all other loops open); . . S,J = 1,N (i t i) (30b)
S2, = (adf/&j)(ki=l;kj=O; and
The dynamic interaction index (DM)is obtained as N
D E=x D i
(314
i=l
where
Construction of the Contrived Model Based on Response Experiments
q i= d i p : Iii = - 1.5
+ 0.5 sin(2diio;) +
+ 0.5 c0s(2diiw3 -I-
2 sin(dii&
kfk;
2 cos(d,,o:)
k:ki
To obtain dynamic interaction measures, each element of G(s)is approximated with the contrived model, (294 kpe-ds
l/;eir
g(s)= as2
(290
where
H
H
y,(oii) H = sin(d,ou)(oii/w,) H
(29i)
Note that kH and WH represent the ultimate gain and frequency of&), wkch are calculated from the identified second-order plus dead-time model. Equations 29 are also applicable to open-loop processes modeled as first-order plus dead-time transfer functions since they are expressed by eq 28 with a, = 0. By examining eqs 29, the asymmetric frequency effects of dynamic interactions are again evident as opposed to steady-state interactions. Note that in general Ag(w) decreases with an increase in frequency.
+ bs + 1
(32)
Hwang and Shiu (1993)have proposed an identification scheme to obtain the parameters in the above model. They employed a 3/3 Pad6 approximation for the timedelay term. Here, the use of a 2/1 Pad6 approximation for the time-delay term is found appropriate. The procedure is as follows. A response experimept is performed under proportional feedback (K = K ) by applying a step signal of amplitude A in set point t o obtain an underdamped response as seen in the paper of Yuwana and Seborg (1982). The required data measured from the test response curve are the first peak value Cpl,the first minimum value Cml,the second peak value Cp2, and the corresponding times tp1, t m l , and tp2. In most situations, the test response at k = is dominated by a pair of complex conjugate poles, b &j&, and hence the following estimates can be obtained:
Ind. Eng. Chem. Res., Vol. 34, No. 1, 1995 231 Table 1. Ultimate Gains and Frequencies of All Examples
(33d) where C, represents the steady state of the response. Since the dominant poles b k j h must satisfy the closed-loop characteristic equation, the model parameters a and b can be expressed as
1 a=
+ &kpe-dd[cos(dh)+ S sin(dh)/hl P+h2
b=
[kkpe-dosin(d&)- 2aShl h
(34b)
where
+ d ( P + h2)/3] h2+ (q - S)[S + d ( P + h2)/3] Or26 - q
+ = - d a ( B + h2)
k:: kY2 k;l e 2 5.29 28.7 3.19 9.75 364 139 66.6 1553 9 -4.56 4.5 9 469 200 143 338 7.13 7.44 -3036 1.37 1.39 -993 -0.06 -0.054 12.4 74.6 -2.12 -1.78 -1219 23.3 -3.15 -13.8 1.07 1.21
(
)
&
0 : :
6
1
4 2
1.66 2.55 11.2 1.5
5.33 0.935 4.56 0.173 0.0706 0.654 4.77 7.07 11.2 1 1.01 0.55 0.686 0.512 1.64 0.352 0.627 1.39 0.225 0.211 1.7 2.22 0.028 0.743 0.684 1.09 1.78 0.228 0.154 1.04
(
1
(34a)
Note that only one parameter, the dead time d, remains to be determined. It can be calculated by fitting the f i s t peak time to the analytically derived value as
tan 4 =
example 1 2 3 4 5
frequencies for dynamic interactions. However, the resultant measures do not require complete dynamic models. Another advantage of the measures is that they are dimensionless and invariant under scaling. In contrast, the interpretation of interactions from the Gershgorin bands would vary with rescaling of input or output variables. For instance, by rescaling an interacting input variable, the previously-regardednoninteracting loop may become significantly interacting. This is certainly an undesirable feature. Most importantly, the present measures provide a clear and quantitative description of steady-state and dynamic interactions in terms of the concept of stability boundaries, which cannot be deciphered from the DNA method with Gershgorin bands.
(35b)
Pairing Strategy and Examples
3(1 kk,)
(35c)
Note that eq 35a is essentially a nonlinear equation of d, whose solution can be easily found by an iteration procedure using (3t,1 - tp2)/2as an initial guess for d.
Comparison with the Direct Nyquist Array and Gershgorin Bands The direct Nyquist array (DNA) method of Rosenbrock (1974) was recommended by Jenson et al. (1986) to measure system interactions. The DNA includes a set of Nyquist plots for all elements of GO'@). Offdiagonal elements are compared with the corresponding diagonal elements to obtain various interaction measures. This method, of course, requires complete knowledge of the dynamic model. A disadvantage of this method is that it cannot easily compare the magnitudes of all elements of G(s)a t the same frequency. This can be resolved by superimposing Gershgorin bands on the Nyquist plots of diagonal elements giiGo). The Gershgorin bands on each diagonal element are formed by drawing circles of radii X&s,,igyl with centers on the gii locus for several frequencies. The radius is evaluated at the corresponding frequency along the gii locus and equal to the sum of the magnitudes of all off-diagonal elements in the same column. Interactions are interpreted as being significant at a certain frequency range when the Gershgorip bands encircle the origins at that frequency range. Note that this is essentially the 1-norm adopted by Grosdidier and Morari (1986) for estimating the diagonal dominance of a control system. The proposed interaction approach utilizes the motivation of the direct Nyquist array with Gershgorin bands for frequencies of practical interest, viz., the zero frequency €or steady-state interactions and ultimate
The selection of the best pairing configuration for N singe input-single output loops is based on the dynamic interaction index (DIX). The primary objective is to minimize dynamic interactions so that each individual loop can respond as a freely standing loop. In this respect, the DM should be as small as possible. In addition, one should always avoid pairings with a negative Ai. If any Ai is negative, the closed-loop system lacks either stability or integrity (unstable with integral action or some loops open (Niederlinski, 1971; Grosdidier et al., 1985)). Also, it follows from a previous section that a negative Ai may introduce a static stability boundary into the stable negative feedback region. This greatly narrows the available stable region for controller tuning and a dominant pole emerging from a zero pole of rs would slow down the closed-loop dynamics considerably. The validity of the proposed approach for measuring dynamic interactions is verified on several examples. The Hopf stability boundaries are plotted for all 2 x 2 examples to illustrate the specific topological features of rH, which may be quite different from those of the simple stability boundaries. The ultimate gains and frequencies of all freely gtanding loops are listed in Table 1. The values of Ai, Sg, Qy, Di,and D E of all 2 x 2 examples are derived via eqs 16, 26, 27, and 29 together with simulated response experiments and are compared to the exact values for the last two examples. The results are listed in Table 2. It is found that the proposed dynamic interaction measures based on response experiments are particularly suited for systems with appreciable dynamic interactions. They may give less accurate estimates for systems with extremely small dynamic interactions as is evidenced by the first pairing configuration of example 3 in Table 2. This shortcoming, however, does not diminish the applicability of the proposed measures. Note that a system with
232 Ind. Eng. Chem. Res., Vol. 34, No. 1, 1995 Table 2. Steady-State and Dynamic Interaction Measures of 2 /?
example 1 estimates example 2 estimates example 3 (1-1/2-2) pairing estimates exact values ( 1-2/2- 1) pairing estimates exact values example 4 estimates exact values
SlZ
QlZ
Di
D2
DM
-0.133
-0.147
-0.0124
0.627
0.0797
0.706
3.42
0.0
0.0
-0.113
-0.103
0.0
0.108
0.108
0.45 0.45
0.055 0.112
0.184 0.089
0.120 0.101
0.120 0.101
0.240 0.202
0.055 0.112
0.184 0.089
0.55 0.55
-6.19 -6.90
-2.21 -1.86
-0.514 -1.01
-0.68 -0.758
4.20 4.38
0.597 0.884
4.80 5.26
3.12 3.60
-0.226 -0.127
-0.356 -0.176
-4.96 -4.06
-2.13 -1.59
0.291 0.152
3.55 2.83
3.84 2.98
1.3e-0.3S1
G(s)=
Example 2 (Tyreus, 1979):
+ l)e-'.ls + u3
[0.1153(10s (4s
0.0887e-12.6S [(43s 1)(22s 1)
+
QZl
-1.12
r2.2e-ls
+
S2l
1.63
extremely small dynamic interactions is simply considered as noninteracting loops, which does not require quantitative accuracy for controller design. On the other hand, the estimates are always accurate enough for pairing structure selection as will be elaborated in a later section. Table 3 presents the interaction measures of two 3 x 3 examples. For all examples, both weighting functions w1 and w2 are chosen to be 0.5 to penalize equally the gain and frequency variations caused by interactions. The static stability loci Ts are also presented, but as expounded previously, their description is exact and require no further verification. The individual examples are discussed in detail below. Example 1 (Vinante and Luyben, 1972):
G(s)=
x 2 Examples
1
0.2429e-2s
(33s
+ 112
0.2429e-0.17s (44s 1)(20s 1)J
+
+
Two "real" distillation systems are used as opposed to dreaming up unrealistic transfer function matrices. The first distillation column model, reported by Vinante and Luyben (19721, describes how the temperatures on the 17th and 4th trays from the bottom of the column are affected by manipulating the reflux flow rate and the steam flow rate to the reboiler. The second model was obtained experimentally from an industrial 30-tray distillation column by Tyreus (1979). These two litera-
L I Figure 5. Stability boundaries of example 1.
I
Table 3. Steady-State and Dynamic Interaction Measures of 3 x 3 Example I, Di D2 03 DM example 5 (2.12 1.85 1.67) 0.511 0.352 0.0482 0.911 example 6 (1-1/2-2/3-3)
(1.09 0.1 0.1)
(1-1/2-3/3-2)
(1.06 0.867 0.974) 0.00422 0.0241 4.82
pairing pairing
0.0769 0.279 0.675 1.03 4.85
ture examples illustrate typical behavior of Hopf stability boundaries for 2 x 2 systems under proportional control. Each system has two possible pairing configurations; the first (1-1/2-21 pairs on the diagonal elements of G(s)whereas the second (1-2/2-11 pairs on the off-diagonal elements. For both examples, the second pairing has a negative 13. and should be avoided. The stability boundaries of two pairings of example 1 are plotted in Figure 5. As seen in Figure 5b, the negative ;1 of the second pairing introduces a simple stability boundary in the first quadrant, causing undesirable dynamic characteristics. The first pairing has moderate steady-state and dynamic interactions as indicated by 13. = 1.63 (compared to 1)and DIX = 0.706 (compared to 0). Owing t o two widely different ultimate frequencies = 1.66, = 4.561, the Hopf stability boundaries consist of two disjoint Hopf curves. This implies that in the stable region bounded by the Hopf stability boundaries and the gain axes, there exist two pairs of dominant complex poles that govern the closedloop dynamics; each pair is destabilized when crossing its corresponding Hopf stability boundaries. Both pairs of complex poles must be taken into account for controller design of the interacting loops. Loop 1 with the lower ultimate frequency is more significantly affected by closing loop 2 as is evidenced by the larger measure D1 (0.627) and the moderately varying Hopf curve corresponding to loop 1. On the other hand, the second loop with the higher ultimate frequency is relatively
(~7~
Ind. Eng. Chem. Res., Vol. 34, No. 1, 1995 233 2.0
1
I
I
(a)UMS\epChangein y,, 1.5
1.0
0.5
0.0
-0.5
0
5
10
15
20
1
I
I
Figure 7. Stability boundaries of example 2.
Time
2.0
1.5 1.o
0.5
0.0
-0.5
0 0
5
10
15
50
x,
100
150
200
Time
Time Figure 6. Set-point responses of example 1 under PI control (k = 2.38,TI = 3.16,for loop 1; k = 4.39,ZI = 1.15,for loop 2;- - -, freely standing loop; -, interacting loop).
insensitive to the first loop as indicated by the much smaller measure D2 (0.0797) and the slowly varying Hopf curve corresponding to loop 2. The asymmetric property and amount of dynamic interactions are easily uncovered by the present approach and verified in Figure 6. Figure 6 compares the set-point responses of the two coupled loops with those of the freely standing loops. Each loop is under PI control with the settings proposed by the original authors. For the response of y1 subject to a unit step change inylr, i.e., the set-point response of loop 1, the coupled loop behaves quite differently from the freely standing single loop (loop 2 open) as seen in Figure 6a. For the response of y2 subject to a unit step change in yzr, i.e., the set-point response of loop 2, the coupled loop behaves rather similarly to the single loop a t least during the initial stage as seen in Figure 6b. The cross responses 011 to ~ 2 r y2 , to ylr) also confirm this moderate amount of dynamic interactions. The RGA (A = 3.42) of example 2 indicates that the first pairing suffers from severe steady-state interactions, implying that the design of the control system will encounter a serious interaction problem. The present dynamic interaction measures (Dl = 0, DZ= 0.108, D M = 0.1081, however, reveal that the system possesses minor dynamic interactions and can be designed in a simple manner. The stability boundaries are shown in Figure 7. The Hopf curve corresponding t o loop 1is a straight line perpendicular to the 121-axis as clearly predicted by the zero D1 measure, implying a noninteracting loop. The results are again verified with PI
0
50
100
1%
200
Time
Figure 8. Set-point responses of example 2 under PI control (k = 15.7,TI = 4.35,for loop 1; k = 17.7,TI = 42,for loop 2; - - -, freely standing loop; -, interacting loop).
controllers in Figure 8. Loop 1 can be tuned independently and the set-point response of y1 is scarcely affected by closing the other loop as illustrated in Figure 8a. Loop 2 is subsequently tuned with loop 1 closed at the predetermined settings. The y2 response of the coupled loop subject to a step change inyzr is somewhat different from that of the single loop as implied by the DZ measure. Also as anticipated, the cross responses are both insignificant. These results clearly demonstrate that the present dynamic interaction measures are more suited t o designing interacting loops than the RGA, which raises unnecessary fear for controller design.
234 Ind. Eng. Chem. Res., Vol. 34, No. 1,1995 Im
Im
....
Figure 10. Direct Nyquist array with Gershgorin bands of example 3.
Figure 9. Stability boundaries of example 3.
Example 3 (Niederlinski, 1971):
G(s)= 1
[(O.ls
ko.1,
+ 1)(0.2s+ 112 1 + 1)2(0.2s+ 112
-1.2
1 J
+ 1)(0.2s+ 1)2(0.5s+ 1) 1 (0.1s + ”2s + 112 (0.1s
Two pairings of this example have comparable %’s (0.45 and 0.55). However, the first pairing has a much smaller DM than the second pairing (0.24 compared t o 4.8) as shown in Table 2, indicating that while steadystate interactions are rather equivalent for both pairings, the first one has far less dynamic interactions and should be selected for controller design. This result is evident in Figure 9, in which the first pairing resembles dynamically noninteracting loops. Note that the first pairing is indeed the better pairing chosen by the original author according to his detailed analysis. Example 3 is different from the previous ones in that, for both pairings, the two ultimate frequencies wfl and mFz are the same. As a result, the same Hopf curve intersects both gain axes to bound the stable region. However, in the first pairing, another Hopf curve slightly intrudes in the bounded region as illustrated in Figure 9a. This intruding Hopf curve has frequencies (around 9) less than those associated with the previous Hopf curve (around 11)and slows down the closed-loop response somewhat. However, this phenomenon of the intruding Hopf curve is not observed for the second pairing, which possesses two widely different ultimate frequencies. Nevertheless, the first pairing is far superior t o the second pairing dynamically. The conclusion obtained from % and DM is consistent with the DNA with Gershgorin bands which assumes that the detailed dynamic model is known as depicted
c
1
+ l)(s + l)(O.ls+ u4 0.85 .(167s + 1)(0.5s + lI4(s+ 1) (167s
0.85
(83s
+ l)(s+ 1)2 1
(167s
+ l)(s+ 1)’
This system has to be paired on the diagonal elements since the other pairing has a negative A. However, the diagonal pairing yields excessively large % (3.12) and D E (3.841, implying extremely severe steady-state and dynamic interactions. As a result, a catastrophic intruding Hopf curve appears in the stable region as shown in Figure 11. Different from the previous case, this Hopf curve has frequencies at about 0.01. The dominant poles associated with this Hopf curve certainly deteriorate the system dynamics considerably (as compared to the two ultimate frequencies 1.5 and 1). Note that this situation occurrs only when both loops are closed. Each single loop does not “feel”the existence of the intruding Hopf curve. A multivariable controller (Ray, 1981) may be required t o achieve good perfor-
Ind. Eng. Chem. Res., Vol. 34,No. 1, 1995 235 mance for this extremely ill-behaved system, This difficulty,predicted readily by the present analysis, was already observed by the original author but without a simple explanation for its occurrence. Upon examining a considerable number of 2 x 2 examples, certain rules of thumb can be summarized below: (i)When dynamic interactions are negligible ( D E>1), one must take its possible existence into account and consider employing a multivariable controller. Example 5 (Ogunnaike and Ray, 1979):
G(s)=
+
+
+
8.64s 1 9.06s 1 6.7s 1 1.11e-6.5s 2.36e-% -0.01 5s 1 7.09s 1 3.25s 1 -34.68e-'.% -46.2e-9.4s 0.87(11.61~4-l)e-s 8.15s 1 10.9s 1 (3.89s 1)(18.8s 1)
+
+
+
+
+
+
+
1
Example 6 (Tyreus, 1982):
I
1.986e-0.71s -5.24e-60S 400s 1 14.29s 1 66.67s 1 -0 .0204e-0.59S 0.33e-0.68" -2.38e-0.42S G(s)= (7.14s 1)2 (2.38s 1)' (1.43s U2
+
+
+
+
+
+
The two 3 x 3 distillation systems are to demonstrate the validity of the approach for N x N systems. The Di and DM measures are shown in Table 3. The first example has only one feasible pairing (u1 - y1, u2 - yz, u3 - y3) with a positive Ai. The D M does not imply significant dynamic interactions. The second example has two feasible pairing configurations. The RGA indicates that the first pairing (u1 - y1, uz - yz, u3 373) has serious steady-state interactions while the second pairing (u1 - y1, u~ - y3, u3 - y2) is almost steady-state decoupled. On the contrary, the DM measures ascertain clearly that the first pairing has
much less dynamic interactions than the second pairing and should be selected for subsequent controller design.
Conclusions The proposed approach for measuring dynamic interactions accounts directly for dynamic properties of the multiloop control system. In contrast to the RGA, which can give misleading indications about the amount of interactions present in dynamic systems, the proposed D M measure always yields quantitatively reliable interaction information regarding the system. Furthermore, the approach gives a clear geometric interpretation of steady-state and dynamic interactions, which leads to useful guidelines for controller design. The approach is presented for continuous-time systems. However, it can be extended readily to sampleddata systems for which discrete-time proportional controllers are used. The Hopf stability boundaries differ little for small sample times, whereas the stability regions may decrease considerably for larger sample times. The dynamic interaction measures for multiloop control systems require merely simple transient response experiments to construct second-order plus deadtime models. By virtue of eqs 29, the measures have the advantage of being dimensionless like the RGA. This indicates that they are invariant under scaling and represent meaningful quantitative estimation. Although the proposed interaction measures are derived from proportional controllers, their application is not restricted to proportional control systems. The major reason is that proportional feedback always constitutes a large portion of any control. Hence, systems with large interactions under proportional control must possess large interactions under control of any other kind. This is exactly true for steady-state interaction since the RGA, which is obtainable from proportional control systems according to the present approach, is surely not restricted to any specific controllers in applications. For convenience,the Di and D M measures are usually sufEcient to provide fair quantitative estimates for dynamic Lnteractions. On the other hand, the sign of the slope S i may provide certain information regarding integrity and stability robustness- of the interacting system. For example, a positive S i implies a convex stable re@on. If the controller setting can be chosen such that ki is larger than unity, the interacting system will be unstable if that Controller or sensor fails, viz., it lacks integrity. However, if S i is negative, the stable region in the first quadrant is narrow and all settings are close to the stability boundary and the system stability can be destroyed easily by a perturbation that shifts the stability boundaries, viz., it lacks stability robustness with respect to system perturbations.
Acknowledgment The author wishes to acknowledge helpful discussions with Professor Hsueh-Chia Chang at the University of Notre Dame. This work is supported by the National Science Council of ROC under Grant No. NSC-83-0402E006-018.
Nomenclature C,l = first minimum value of test response C,1 = first peak value of test response
238 Ind. Eng. Chem. Res., Vol. 34,No. 1, 1995
C,, = second peak value of test response C, = steady-state value of test response d = disturbance vector G = transfer function matrix of the plant with a pairing configuration in question G, = process transfer function matrix G, = controller transfer function matrix go = elements of G gp = elements of G , I = identity matrix K = proportional controller matrix k = proportional controller gain kH = ultimate gain of freely standing loop i = inverse steady-state gain of freely standing loop i 4 = nominal proportional controller setting k = normalized proportional controller gains defined in eqs
ki
15
& = normalized proportional controller gains defined in eqs 23 P = permutation matrix So = dimensionless slope along Hopf locus at ki-axis s = Laplace transform variable T = complementary sensitivity matrix Td = sensitivity matrix t,l = first minimum time of test response tpl = first peak time of test response t,, = second peak time of test response u = input vector y = output vector yr = reference vector Greek symbols r = bifurcation locus
A = relative gain array K
= Rijnsdorp interaction quotient
d = relative static gain 9.. = change of ultimate frequency of loop i with kj
= ultimate frequency of freely standing loop i Indices H = Hopf bifurcation P = process
S = simple bifurcation
Literature Cited Astrsm, K J.; Htigglund, T. Automatic Tuning of PZD Controllers; Instrument Society of America: Research Triangle Park, NC, 1988. Bristol, E. H. On a New Measure of Interaction for Multivariable Process Control. ZEEE Trans. Autom. Control 1966,11,133. Bristol, E. H.RGA Dynamic Effects of Interaction. ZEEE Conf. Decision Control. 1977,16th. Chang, H. C.; Chen, L. H. Bifurcation Characteristics of Nonlinear Systems under Conventional PID Control. Chem. Eng. Sei. 1984,40,2191. Friedly, J. C. Use of the Bristol Array in Designing Noninteracting Control Loops A Limitation and Extension. Znd. Eng. Chem. Process Des. Dev. 1984,23,469.
-
Grosdidier, P.; Morari, M. Interaction Measures for Systems Under Decentralized Control. Automatica 1986,22,309. Grosdidier, P.; Morari, M.; mot, B. R. Closed-LoopProperties from Steady State Gain Information. I d . Eng. Chem. Fundam. 1985, 24,221. Guckenheimer, J.; Holmes, P. J. Nonlinear Oscillations, Dynamical Systems and Bifurcation of Vector Fields; Springer-Verlag: New York, 1983. Hwang, S. H.; Chang, H. C. Process Dynamic Models for Heterogeneous Chemical Reactors-An Application of Dynamic Singularity Theory. Chem. Eng. Sei. 1986,41,953. Hwang, S.H.; Chang, H. C. A Theoretical Examination of ClosedLoop Properties and Tuning Methods of Single-Loop PI Controllers. Chem. Eng. Sei. 1987,42,2395. Hwang, S. H.; Shiu, S. J. A New Auto-Tuning Method with Specifications on Dominant Pole Placement. Int. J . Control 1994,60, 265. Jensen, N.; Fisher, D. G.; Shah, S. L. Interaction Analysis in Multivariable Control systems. MChE J . 1986,32,959. McAvoy, T.J . Interaction Analysis: Principles and Applications; Monograph Series 6; Instrument Society of America: Research Triangle Park, NC, 1983. McDermott, P.; Bonvin, D.; Mellichamp, D.; Chang, H. C. Eigenvalue Spectra and Modal Contributions for Counterflow Reaction Models. Chem. Eng. Sci. 1984,31,263. Niederlinski, A. A Heuristic Approach to the Design of Linear Multivariable Interacting Control Systems. Automatica. 1971, 7, 691. Ogunnaike, B. A.; Ray, W. H. Multivariable Controller Design for Linear Systems Having Multiple Time Delays. AZChE J . 1979, 25, 1043. OReilly, J.; Leithead, W. E. Multivariable Control by Individual Channel Design. Int. J . Control 1991,54, 1. Peltomaa, A.; Koivo, H. N. Tuning of a Multivariable Discrete Time PI Controller for Unknown Systems. Int. J . Control 1983,38, 735. Ray, W. H. Advanced Process Control; McGraw-Hill: New York, 1981. Rosenbrock, H. H. Computer-Aided Control System Design; Academic Press: London, 1974. Tyreus, B. D. Multivariable Control System Design for an Industrial Distillation Column. Ind. Eng. Chem. Process Des. Dev. 1979,18,177. Tyreus, B. D. Paper presented a t the Lehigh University Distillation Control Short Course, Bethlehem, PA, 1982. Vinante, C. D.; Luyben, W. L. Experimental Studies of Distillation Decoupling. K e n . Teollisuus. 1972,29,499. Witcher, M. F.; McAvoy, T. F. Interacting Control Systems: Steady State and Dynamic Measurement of Interaction. ZSA Trans. 1977,16,35. Yuwana, M.; Seborg, D. E. A New Method for On-Line Controller Tuning. AlChE J . 1982,28,434.
Received for review December 1, 1993 Revised manuscript received August 30, 1994 Accepted September 14,1994@
IE9306077 @
Abstract published in Advance ACS Abstracts, December
1, 1994.
