- Email: [email protected]
The Computer Aided Design of Multivariable Non Linear Control Systems using Frequency Domain Techniques
THE COMPUTER AIDED DESIGN OF MULTIVARIABLE NON LINEAR CONTROL SYSTEMS USING FREQUENCY DOMAIN TECHNIQUES
GRAY, J . O .
and TAYLOR, P.M.
ABSTRACT
A conputer program has been developed based on s i n u s o i d a l i n p u t d e s c r i b i n g f u n c t i o n techniques f o r t h e p r e d i c t i o n o f l i m i t c y c l e s i n m u l t i v a r i a b l e non l i n e a r feedback systems. I n t e r a c t i v e d i s p l a y s a r e produced i n t h e conplex frequency domain as an a i d t o conpensator design and t h e method is a p p l i e d t o t h e design o f a c o n t r o l l e r f o r an i n d u s t r i a l a i r conpressor system. 1.
INTRODUCTION
A f e a t u r e o f many p r a c t i c a l process c o n t r o l s y s t e m i s a multi-input multi-output i n t e r a c t i v e s t r u c t u r e with one o r more gross non l i n e a r i t i e s and a p r i mary design t a s k i n such c i r c u n s t a n c e s i s t o ensure t h e avoidance o f l i m i t c y c l i n g c o n d i t i o n s . This p a p e r o u t l i n e s how such a system may be investigated using a sinusoidal input describing function philosophy and a computer program has been developed t o p r e d i c t t h e magnitude and freqaency of p o s s i b l e limit c y c l e modes and t o a i d i n t h e t h e design o f s u i t a b l e conpensation t c ensure avoidance o f c r i t i c a l r e g i o n s i n t h e complex frequency domain. T h i s program r e p r e s e n t s an e x t e n s i o n o f a p r e v i o u s l y developed design s u i t e f o r s i n g l e i n p u t s i n g l e o u t p u t non l i n e a r feedback s y s t e m ( 7 ) and can provide t h e c o n t r o l e n g i n e e r wi?h a powerful i n t e r a c t i v e d e s i p s u i t e which is capable o f d e a l i n g with a wide range of non l i n e a r f u n c t i o n s where t h e hypothesis of harmonic l i n e a r i s a t i o n can be assumed v a l i d . The b e h a v i o w o f a l a r g e conpressor f o r s u p p l y i n g a i r t o a chemical p l a n t i s s t u d i e d and a conpensator designed t o achieve c l o s e d loop s t a b i l i t y and s p e c i f i e d performance. 2.
LIMIT CYCLE PREDICTION I N MULTIVARIABLE SYSTLE1S
S i n u s o i d a l d e s c r i b i n g f u n c t i o n techniques have been used widely i n t h e s t u d y of l i m i t c y c l e behaviour i n non l i n e a r , s i n g l e loop feedback s y s t e m where t h e a p p l i c a t i o n o f conplex frequency domain ana: ys i s has represer-ted an i n t u i t i v e extension c f l i n e a r s y s t e m t h e o r y . The l i m i t a t i o n s of t h e cpproach are, o f c o u r s e , w e l l documented(1) and, a l though t h e r e c e n t emergence o f an algorithm due t o Mees(2) has done much t o enhance the v a l i r ? i t y cf t h e method, t h e assunption o f , e s s e n t i a l l y , l o ? p a s s l i n e a r system elements s t i l l remains &r i m p o r t a n t p r e r e q u i s i t e i n any s t u d y i n v o l v i n g t h e use of t h e d e s c r i b i n g f m c t i o n philosophy. Recent
i n t e r e s t h a s been shown i n t h e a p p l i c a t i o n o f s i n u s o i d a l d e s c r i b i n g f u n c t i o n methods t o m u l t i v a r i a b l e non l i n e a r feedback systems where again advantage can be taken o f i n t u i t i v e e x t e n s i o n s t o a l r e a d y developed l i n e a r s y s t e m t h e o r y . Such an e x t e n s i o n i n v o l v i n g t h e i n v e r s e Nyquist a r r a y was proposed by Rosenbrock i n 1969(3) and l a t e r i n a more d e f i n i t i v e form by Mees(4). More r e c e n t l y Ramani( 5) has proposed g r a p h i c a l methods i n v o l v i n g eigenvalue l o c i o f t h e l i n e a r system elements i n t h e s t u d y of l i m i t c y c l e o p e r a t i o n . Although t h e a p p l i c a t i o n o f t h e i n v e r s e Nyquist a r r a y does p r o v i d e a quick, systema t i c method of compensation design i n t h e complex frequency domain f o r t h e a v o i d a n ~ eof a l l regions of p o s s i b l e t h e o r e t i c a l l i m i t c y c l e b e h a v i o w , i t can be shown i n p r a c t i c e t h a t t h e method is q u i t e p e s s i m i s t i c i n its a c t u a l p r e d i c t i o n of l i m i t c y c l e s ( 6 ) . The o b j e c t o f t h i s paper is t o o u t l i n e a method f o r t h e p r e d i c t i o n o f limit c y c l e operat i o n which i s l e s s conservative than t h a t based on the i n v e r s e Nyquist a r r a y approach b u t which s t i l l maintains t h e c m c e p t u a l l y powerful Nyquist d i s p l a y s t o a i d i n compensator design. In u s i n g t h e i n v e r s e Nyquist a r r a y t h e p o s s i b i l i t y of t h e o r e t i c a l limit cycle o p e r a t i o n i n an i n t e r a c t i v e non- l i n e a r diagonally dominant mulcivariable system can be s t u d i e d g r a p h i c a l l y b ; considering t h e behaviour o f l o c i a s s o c i a t e d w i t h t h e diagonal l i n e a r and non l i n e a r elements of t h e system. Each locus has an envelope defined by a s e r i e s of d i s c s c e n t r e d on t h e locus whose r a d i i a r e r e l a t e d t o t h e magnitudes o f corresponding o f f diagonal elements. P o s s i b l e s o l u t i o n s t o t h e s i n p l f h a r ~ o n i cbalance equations can e x i s t i n any regior. of t h e corrplex plane where t h e envelopes of r e l a t e d l i n e a r and non l i n e a r l o c i i n t e r s e c t . Althcyigh :r!e applicat i o n o f s c a l i n g :echriques can be used with advanta g e ( 8 ) ( 6 ) t > a a j u s t t h e envelspt dirrer.sions of i n d i v i d ~ a ll o c i , i t i s though- ti-.: rethod o f system compensation b c e d on r7.e ,:..cl".a;.ca c f envelope o v e r l a p must t e n d t o he cori;-?i.-.t;.:t. Provide.? t h a t a s y s t e m a t i c metlir;i. c : c c . t~= t i o n is employed it i s , however, not r j ; t f i c u l t 'r c a l c u l a t i a t any freq,uency a l l p o s s i b i e s s i . ~ : l m s t t h e simple s t e a d y s t a t e harr.or.ic t;?-an;* rr.atrl. e q ~ i t i o n over the range of physica-ly e l l c i a b i e s i g n s magnitudes a: d freq-encies an< t c ;r,-~l;;ee i:. t h e complex frecuency dcmair. ir:di v i ?ual s c i u: ior. v e c t o r s which r e p r e s e n t t h e con:ribc?ic:. 3f kc::. diagonal and o f f diagcnal systerr e:emer.tz.
.
Theoretical l i m i t cycle operation can then be determined by a study of these individual s o l u t i o n vectors i n t h e v i c i n i t y of t h e (-1, jO) point. AS t h e magnitude and phase o f each p o s s i b l e s o l u t i o n vector is knm.,the a r e a o f uncertainty introduced i n the inverse Nyquist a r r a y approach by t h e a r t i f i c e o f the Gershgorin envelope is removed and a less conservative r e s u l t ensured. Such an approach is, a f course, only v i a b l e i f t h e computational e f f o r t can be constrained within reasonable born&. In t h e method adopted here t h e steady s t a t e harmoni c balance equations o f all loops a r e considered i n d i v i d u a l l y i n a s e q u e n t i a l manner such t h a t t h e number o f s o l u t i o n s e t s t o t h e o v e r a l l system harmonic balance equation diminishes r a p i d l y as t h e conputation proceeds. For convenience, the s o l u t i o n vectors of each loop can be displayed graphically as a Nyquist a r r a y i n t h e conplex frequency domain. 3.
THE PREDICTION OF LIMIT CYCLES USING A SEQUENTIAL COMPUTATIONAL- PROCEDURE
?he method is a p p l i c a b l e t o t h e system configurat i o n s shown i n Figure 1when t h e following condit i o n s a r e met: (a)
all l i n e a r elements a c t , e s s e n t i a l l y , as low pass f i l t e r s t o j u s t i f y t h e a p p l i c a t i o n o f t h e describing function philosophy;
(b)
i f a l i m i t cycle is present all loops w i l l o s c i l l a t e a t t h e same frequency. In general t h i s is a reasonable assunption p a r t i c u l a r l y i f the non l i n e a r elements a n s i m i l a r i n form a d i f the dominant l i n e a r elements have approximately t h e same frequency c h a r a c t e r i s t i c s ( 5 1.
'Ihe equations governing l i m i t cycle "peration i n a multivariable non l i n e a r system can be expressed as follows
where T is the non l i n e a r r e t u r n r a t i o system matrix with elements t i . ( a j , j w ) and I is a colunm vector o f sinusoids a t $he i n p u t t o t h e non l i n e a r e l e m n t s such t h a t a i = A i / B i . The method c o n s i s t s o f findi n g all p o s s i b l e G l u t i o n s t o equation s e t ( 1 ) over physically allowable ranges of s i g n a l magnitudes m d frequencies. I f we consider f i r s t t h e case o f Limit cycle operation where
w f can, f o r t h e range of p h y s i c a l l y allowable, d i s c r e t e values of the reference s i n u s o i d a 1 = A1 and a t a s p e c i f i c -raiue o f f r e q u e n c y n , determine those f i n i t e p s e t s o f sinusoids ( a l k , a 2 ); k = 1,2,. ,p which w i l l s a t i s f y t h e con8ition &or l i m i t cycle operation i n t h e n t h loop as given by t h e equation
..
I f we now examine t h e equation governing l i m i t cycle operation i n t h e (n-11th loop which is given by
where again al
=
A1,
We can t r y f o r s a t i s f a c t i o n o f t h i s equation using only those f i n i t e p s o l u t i o n s e t s derived from equaticm.(3) o v e r t h e range o f d i s c r e t e values o f A1A s i n p r a c t i c e h i s c r e t e values of A1 a r e used i n the determination o f t h e s e t s o f s i n u s o i d s which w i l l s a t i s f y both equations ( 3 ) and (41, t h e r e could t h e o r e t i c a l l y be s o l u t i o n s which e x i s t o w r t h e d a t a i n t e r v a l s and t h i s could give r i s e t o e r r o r s i n l i m i t cycle p r e d i c t i o n . To guard a g a i n s t t h i s and t o add a s p e c i f i e d degree o f robustness t o t h e conputation we pass from one equaticn t o t h e next, n o t cnly those s e t s o f v e c t o r s which s a t i , s f y t h e f i r s t equation e x a c t l y b u t a l s o those s e t s which y i e l d an approximate s o l u t i c n . The degree o f approximation is r e l a t e d d i r e c t l y t o t h e anplitude, frequency and phase i n t e r v a l s chosen and a l s o t o s p e c i f i e d e r r o r bomds on the system parameters. These approximate s o l u t i o n s a r e generated automatically within t h e program once t h e d a t a i n t e r v a l s and parameter e r r o r bounds are s p e c i f i e d by t h e user. I n the f i n a l g r a p h i c a l display cmly those e r r o r c i r c l e s associated with t h e system parameter e r r o r bounds need be displayed. I f q s o l u t i o n s e t s o f s i n u s o i d s a r e found (where q 6 p) which s a t i s f y both equations ( 3 ) and ( 4 ) t h e n t h e (n-21th loop equation is n e x t examined f o r harmcnic balance using these q s e t s and the process repeated u n t i l loop 1 is reached. This has a harmonic balance equation which can be p u t i n the form
which i s p l o t t e d as shown i n Figure 2 f o r a range of values o f n. I f t h e end of m e l e f t hand s i d e s o l u t i o n vectors touches t h e (-1, jO) p o i n t o r i f t h i s p o i n t l i e s within an e r r o r c i r c l e associated with one of t h e s o l u t i o n v e c t o r s , then a conditicn could e x i s t which would s a t i s f y t h e o v e r a l l steady s t a t e harmonic balance requirements o f equation (1).
It should be noted t h a t as soon as any equation is not s a t i s f i e d by a s e t of s o l u t i o n sinusoids passed from t h e previous equation then t h i s s e t is immed i a t e l y discarded. When s o l u t i o n s have been obtained f o r a l l loop equations over t h e allowable range o f A 1 and Q the i n e q u a l i t y
,...a
is n e x t considered and t h e above conputational procedure repeated in a cycLic manner f o r a range of physically r e a l i s a b l e values o f t h e reference sinus o i d a2 = A?. The d e s c r i b i n g function equation f o r loop 2 is f i n a l l y derived a s
where a = A1 1
.
i n a s i m i l a r manner t o t h a t f o r loop 1 and displayed i n the conplex frequency domain. The computatior. p r o c e e d s . u n t i 1 d i s p l a y s a r e produced f o r a l l n loops i n t h e system. Limit cycles could t h e o r e t i c a l l y occur i f any s o l u t i o n v e c t o r o f any o f t h e n loop equations touches t h e (-1, jO) p o i n t o r , f o r a more conservative p r e d i c t i o n , i f t h i s p o i n t l i e s within one of t h e e r r o r c i r c l e s a s s o c i a t e d with any solut i o n vector. A study of such a Nyquist a r r a y can quickly d e t e r
mine t h e magpitude and frequency of a l l t h e o r e t i c a l l y p o s s i b l e modes of l i m i t cycle operation and isol a t e those loops which a r e responsible f o r o v e r a l l system i n s t a b i l i t y . Of equal inportance, i t can a l s o i n d i c a t e t h e form of compensation necessary t o avoid t h e p o s s i b i l i t y of limit cycle operation i n a p a r t i c u l a r loop. I t is o f t e n convenient,when a preconpensator element is envisaged,to use a Nyquist display of t h e e l e -
ments of T. Such a display i s e a s i l y generated have been once all s o l u t i o n s e t s f o r the vector determined using t h e method o u t l i n e d above. ?he corresponding harmonic balance equation f o r loop 1 now has t h e form
. ht.
where
(j
?CZ j d j 4.
=
...n)
1..
l o c i ot ?..c rl-g,:.ae 1 e r n r . t ~ tii. (i=l,.. n ) . For any cndnge ir. loop m the diago::al locus of t h e i t h loop i s m a l t e r e d h c t the magpitude and phase o f t h e superinposed s o l u r i o n v e c t o ~ sw i l l a l t e r . Any r.ew v a i i e I ,sz- :?cme:.r? iie irithin bounds defined by a c i r c l e with r a d i u s
centred on t h e m a l t e r e d diagonal locus. This c i r c l e can e a s i l y be computed and displayed i f necessary. The corresponding l i n e a r preconpenoator
a c t s on t h e colunns of T rendering t h e Nyquist displays of T no longer very u s e f u l f o r compensator design. Such a precompensator has howev$r t h e e f f e c t o f multiplying t h e rnth row o f T by km(j w ) where again t h e value of m is chosen a s t h a t of t h e most troublesame loop. Desigp i n t h e frequency domain gan thus proceed by using displays o f t h e v e c t o r tii ( i = l , . .n) with superimposed s e t s of s o l u t i o n vectors
.
a r e t h e elements of
COMPENSATION TECHNIQUES
In Figure 1 a p o s t conpensator o f t h e form
w i l l have the e f f e c t of multiplying the rnth row o f T by Em(jw) where the value of m is chosen a s t h a t of the most t r o u b l e s o w loop. Every vector e l e m n t i n t h e rnth loop equation w i l l be s i m i l a r l y affected. As we a m dealing with a Nyquist a r r a y t h e r e q u i r e d form of (,L jw) t o s u i t a b l y compensate t h e rnth loop proceeds n a t u r a l l y from a knowledge e f c l a s s i c a l s i n g l e loop feedback theory. N l o t h e r vectors tii; ( i = l , . .n)ifm w i l l remain u n a l t e r e d a s w i l l passed t o t h e rnth equation. t h e s o l u t i o n values of However, t h e s o l u t i o n values of passed out o f t h e rnth equation w i l l be a l t e r e d which i n t u r n w i l l vary t h e vectors
.
The e f f e c t of such a compensator on t h e behaviour of t h e o t h e r loops can most e a s i l y be v i s u a l i s e d i f we consider a l l vectors
t o be superimpmed a t s p e c i f i c frequencies on t h e
5.
COWUTATIONAL ASPECTS
In t h e s e q u e n t i a l conputational procedure up t o twenty magnitude values of t h e reference s i n u s o i d and up t o twenty values o f frequency are used i n t h e determination o f s o l u t i o n s vectors f o r any loop equation. The number o f phase angles is a l s o n o r Although t h i s may maily r e s t r i c t e d t o twenty. i n i t i a l l y r e s u l t i n a l a r g e number o f s o l u t i o n s e t s most of these a r e p r o g r e s s i v e l y discarded a s the conputation proceeds and t h e r e is a r a p i d convergence on t h e f i n a l displayed s o l u t i o n s e t . The frequency and amplitude ranges of p o s s i b l e l i m i t cycle behaviour a r e automatically displayed i n t a b u l a r form along with t h e most probable l i m i t cycle a m p l i t u d e s and frequency. Once c r i t i c a l a n p l i t u d e and frequency ranges have been a s c e r t a i n e d irom t h e t a b l e and an i n i t i a l study of t h e d i s p l a y s , t h e a m p l i t u d e and frequency ranges used i n the calculaticn can be conpressed t o give a d d i t i o n a l s o l u t i o n s e t s i n t h e c r i t i c a l region o f t h e complex f r e q u e n c ~ plane i f required. The robustness of t h e s o l u t i o n i s enhanced by compensating automatically f o r the i n t e r v a l s i n t h e d i s c r e t e data used i n t h e conputation and by allowi n g f o r known e r r o r bounds on any element o r group of elements within t h e feedback system under consideration. Such bounds can be r e l a t e d d i r e c t l y t o bounds on i n d i v i d u a l s o l u t i o n vectors and i t i s thus conputationally easy t o produce corresponding e r r o r c i r c l e s on t h e conplex frequency plane. For c l a r i t y of p r e s e n t a t i o n only t h a t e r r o r c i r c l e
c l o s e s t t o t h e (-1, jO) p o i n t need be s t o r e d and displayed a t a y frequency p o i n t cn t h e c o l p l e t e s o l u t i o n w c t o r . This t e c h n i q r r ca a l s o be used t o obtain a m a s u r e o f t h e s e n s i t i v i t y o f system l i m i t c y c h behaviour t o v a r i a t i o n s i n one o r more o f t h e system's p h y s i c a l parameters. A l l design program a r e w r i t t e n in Fortran I V with a chain f a c i l i t y t o reduce core s t o m requirements t o l e s s than 10K words. Conputer questions and mswers use lrnemonic English format and a command i n t e r c e p t o r r o u t i n e provides c o q l e t e f l e x i b i l i t y i n e n t e r i n g o r l e a v i n g a n y . p r o p m block. This allow8 r a p i d change o f e i t h e r l i n e a r o r non l i n e a r d a t a and innnediate access t o t h e display p r o g r a ~ ~ t o evaluate the e f f e c t uf any d e s i w decision. Linear d a t a can e i t h e r be e n t e r e d i n pole wm o r polynomial form. Non l i n e a r d a t a can be e n t e r e d i n polynomial o r piecewise l i n e a r form o r a s a s e t o f d i s c r e t e measured d a t a points. The v i s u a l display u n i t f a c i l i t i e s include l i n e o r p o i n t format, alphzr numeric e n t r y , a x i s display and automatic s c a l i n g .
6.
THE DESIGN OF A CONTROLLER FOR AN INDUSTRIAL COMPRESSOR
The system which is s h a m diagrammatically i n Figure 3 c o n s i s t s o f a l a r g e conpressor f o r supplyi n g air t o a chemical p l a n t and has been s t u d i e d previously as a purely Linear system(9) (10). The o u t l e t p n s s u r e is c o n t r o l l e d by guide vane p o s i t i c n and t h e o u t l e t flow by means o f a blow o f f valve which' prevents surging a t i n p u t flows below 20% o f f u l l flow r a t e . When t h e blow o f f valve is closed the system c o n s i s t s o f a s i n g l e c c n t r o l loop and good r e g u l a t i o n can be obtained using a standard pneumatic c o n t r o l l e r . When t h e blow o f f valve opens however the system has two i n p u t s and e x h i b i t s l i m i t cycle behaviour when t h e loop gains are increased t o any u s e f u l value. The dynamics o f t h e l i n e a r system kre obtained(l0) with t h e blow o f f valvc open and are shown i n Figure 4, which now includes a s a t u r a t i o n eJement t o represent t h e non l i n e a r c h a r a c t e r i s t i c s o f the blow o f f valve loop. It is d e s i r e d t o design a preconpensator K(s) which w i l l guarantee avoidance of l i m i t cycle behavisur, a f a s t blow o f f valve response, low i n t e r a c t i e n and zero o f f s e t i n t h e pressure loop. Various sources of e r r o r e x i s t i n t h e s i n p l i f i e d system model shown i n Figure 4 due t o t h e following factors. (1)
The l i n e a r system elements given by t h e model only approximate t o t h e a c t u a l p l a n t behaviour.
(2)
A very s i m p l i f i e d model o f t h e blow o f f valve non l i n e a r i t y is used.
(3)
A non l i n e a r c h a r a c t e r i s t i c is known t o e x i s t
i n t h e guide vane linkage. D e t a i l s o f t h i s were n o t a v i l a b l e b u t it was known t h a t l i n e a r operation could be achieved i f required. (4)
Valve dynamics a r e not included.
To ensure a conservative design i n t h e presence o f these modelling e r r o r s and e r r o r band o f 230% was used on a l l elements o f t h e - r e t u r n r a t i o matrix T throu&out s u b s e ~ a e n tc a l c u l a t i o n s . An i n i t i a l comp u t e r i n v e s t i g a t i o n f o r l i m i t cycle operation using
a proposed c o n t r o l l e r i n t h e pressura loop with a gain o f 1 5 pmduced t h e displays shcmn i n f i g u r e 5 which include t h e appropriate e r r o r c i r c l e s . Ihe system h a s obviously a very high degme o f i n t e r a c t i o n a n d ' a l i m i t c y c l e was e d i c t e d with t h e f o l lowing moet l i k e l y paranmter values w = 1.35 Rad/sec, u = 7.04, ly = 2.75. A subsequent simw l a t i o n confirmed lidt cycle operation with pararmter values w = 1.45 Rad/sec, ul = 6.6, u2 = 1.5. To reduce i n t e r a c t i o n e f f e c t s a preconpensator o f constants w a s & s i p a d using s t a n d a r d frequency domain techniques(l0) t o mat c l o s e l y diaponalise t h e system a t a s p e c i f i c frequency and anplitude set. A frequency o f 1 Radian/sec was chosen on t h e b a s i s o f e a r l i e r work with a conpletely l i n e a r i s e d system. The displays o f the system b e h a n o w with only t h i s compensator are shown i n f i g u r e 6. It w i l l be seen t h a t no limit cycle is now predicted. A s we a r e now considering a preconr pensator element t h e inverse displays a r e now used with e r r o r c i r c l e s corresponding t o system parameter e r r o r b o m b o f t30% being determined by u s u a l mapping procedures. I n t e r a c t i o n e f f e c t s i n loop 2 a r e n e g l i g i b l e and a r e a l s o small i n loop 1 f o r low values o f u2. The gain i n loop 1 is now increased by a f a c t o r 5 and it w i l l be seen from Figure 7 t h a t t h e e f f e c t on loop 2 is small. I f t h e gain i n loop 2 is subsequently .increased by a f a c t o r 2 the system s t i l l remains f r e e from l i m i t cycle operation b u t loop i n t e r a c t i o n e f f e c t s now increase i n loop 1 as shown i n Figure 8 s o t h i s loop gain was n o t adopted. The f i n a l conpensator was o f t h e form
The multivariable c i r c l e c r i t e r i o n can now be i n voked t o ensure a s y n p t o t i c system s t a b i l i t y and the absence o f l i m i t cycles f o r a l l bomded s e t s of inp u t s and t h e appropriate displays a r e shown i n Figure 9 where, following t h e arguments o f ;ook(+l), we can use t h e s c a l e d mean Gershgorin bands of GK and the c r i t i c a l d i s c o f nqp. Examination o f t h i s p l o t shows t h a t loop 2 w i l l be s t a b l e and a l s o t h a t loop 1 w i l l be s t a b l e with the blow o f f valve loop f u l l y closed assuming t h a t the l i n e a r p a r t o f the model is s t i l l v a l i d m d e r these conditions. The response of the closed loop syster. t o a s t e p change i n pressure demand is shown i n Figure 10 and the corresponding response t o a s t e p change In a i r fiow demand is shown i n Figure 11. In both czses f a t s e t t l i n g times and low i n t e r a c t i o n e f f e c t s a r e achi$ved. Any o f f s e t i n t h e pressure loop car1 be remedied by long time constant i n t e g r a l zc'icn which should nDt a l t e r the p l o t s s i g n i f i c a n t l y . The o f f diagonal t e r n s of (i are l a r g e r than zne diagc:;ai term thus i n d i c a t i n g t h a t The icops are being e f f e c t i v e l y inter(.hangc.d. Th:s could ~ e s u i t i n i n e f f i c i e n t operation, although t h e s i t c a t i a n could be improved by allowing a g r e a t e r degree of loop i n t e r a c t i o n . I r should be noted t h a t t h e des i g n method could have been s e d t o equal e f f e c t if t h e model of t h e non l i n e a r i t y i n loop 1 had beep available f o r i n c i s i o n .
Control S y s t e m Centre Repost No. 2 6 2 , U.M.1 .S .T., September 1974. Provided t h a t t h e e s s e n t i a l p r e r e q u i s i t e s e x i s t f o r i t s a p p l i c a t i o n t h e s e q u e n t i a l conputational procedure p r e s e n t s an accurate, systematic and e f f i c i e n t method of p r e d i c t i n g l i m i t cycle operation i n a multivariable non l i n e a r feedback system. By studying the behaviour of i n d i v i d u a l s o l u t i o n vect o r s i n the region of t h e (-1, jO) point,accurate p r e d i c t i o n s can be made a s t o t h e m a e i t u d e and f r e quency o f p o s s i b l e l i m i t cycles. A s p e c i f i e d degree of robustness can be introduces i n t o t h e method by considering e r r o r bounds on t h e p h y s i c a l elements within t h e system matrix and t h e r e is a d i r e c t and obvious r e l a t i o n s h i p betweer: such bounds and t h e behaviour of s o l u t i o n vectors i n t h e conplex frequencies domain. This f a c i l i t y can a l s o b e used t o r e l a t e v a r i a t i o n s in l i m i t cycle b h a v i o u r t o variat i o n s i n one o r more system parameters t h u s allowing t h e i d e n t i f i c a t i o n of p a r t i c u l a r l y s e n s i t i v e system e l e m n t s . I f limit cycles a r e p r e d i c t e d by an examination of the s o i u t i o n v e c t o r s i n any one loop then t h e s t r u c t u r e of a s u i t a b l e preconpensator o r p o s t compensator t o a d j u s t t h e s e ve.:tor s o l u t i o n s i n t o l e s s c r i t i c a l regions of t h e complex freqwncy domain i s immediately obvious from an elemenrary knowledge of s i n g l e 1.oop feedback theory. h e influence of such a compensator on t h e o t h e r loops is constrained within e a s i l y determined bounds and, provided t h a t t h e system i s n o t grossly i n t e r a c t i v e , a viable m t h o d of system conpensation may be devised.
(7)
GRAY, J.O., SAWIDES, L.S. and TALOR, P.M., I n t e r a c t i v e graphics design using sinusoiZal describing f u i c t i o n methods, I .S.A. Transact i o n s , '101. 1 4 , No. 1, 1975. pp. 68-77.
(8)
RAMANI, N . and ATHERTON, D.F., Correspondence, Proc. IEE, Vol. 20, No. 7, i973. p . El.&.
(9)
LUNTZ, R., MUNRc3, N. a d MCLEOD, R.S., Computer aided design o f m u l t i v a r i a b i e c o n t r o l systems, IEE 4 t h WAC Convention, Manchester, 1971. pp. 55-59.
(10) ROSENBROCK, H.H., 'Ihe computer aided design of c o n t r o l s y s t e m , Academic P r e s s , 1974. (11) COOK, P. , Modified m u l t i v a r i a b l e c i r c l e t h e o r e m . Recent mathematical developments i r c o n t r o l , Academic P r e s s , 1973. pp. 367-372.
The s e q u e n t i a l computation method is thought t o be t h e only e x i s t i n g method f o r s t u d y i n g multivariahle systems with gross non l i n e a r i t i e s i n o r d e r t o a s c e r t a i n l i m i t c y c l e operation and t o derive both pre- and p o s t compensator elements which ensure closed loop s t a b i l i t y with known i n t e r a c t i o n effects.
FIGURE 1
REW ENCES
(1) ELB, A. and VANDERVELD, W.E., Multiple i n p u t describing functions and non l i n e a r s y s t e m design, McGraw H i l l , New York, 1969.
(2)
MEES , A. I . and BERGEN, A . R. , Describing F ; m c t i ~ n sRevisited, E l e c t r o n i c s Research Laboratory, University o f C a l i f o r n i a , Berkeley, 1974.
(3)
RXENBROCK, H . H . , Design o f m u l t i v a r i a b l e cont r o l s y s t e m using t h e i n v e r s e Nyquist a r r a y , Proc. IEE., Vol. 116, No. 11, 1969. pp. 1929-1536.
(4)
MEES, A. I . , Describing functions c i r c l e c r i t e r i a and multiloop feedlack systems, Proc. IEE., Vol. 120, No. 1, 1973. pp. 126-130.
(5)
(6)
,SAKANI, N. and A'TiAERTON, D.P., Frequency response nethods f o r non i i n e a r multivariable system', Proc. CanadiaCorif. Auto. Ccnt., 1973. 9.2-1 9.2-15. GRAY, J.O. and TAYLOR, P.M., & s c r i b i n g functions and t h e inverse Nyquist a r r a y ,
SOLUTION VECTORS
' I
. I
FIGUFE 2
.,
I
?5:.
: :~~..:. , 'a$~,
PLOT OF
Q
LOOP 2
GAIN X5 LOOP 1 GAIN X1
LOOP 2
+
S
LOOP 2
+
1
+ J
NYQUIST ARRAY
/
LOOP 1 M O P 1 GAIN X5 LOOP 2 GAIN X 2
LOW 2
.
FIGURE 9
;.,
9' ,
STEP CHANGE PESSURE
...
M
,:
,F
, LC
SECS
10 ,,,
VALVE POSITION
. r-o OPEN
FIGURE 10
PRESSURE
STEP CHANGE FLOW
I
> o
VALVE POSITION
OPEN
FIGURE 11
GRAY, J . O .
and TAYLOR, P.M.
ABSTRACT
A conputer program has been developed based on s i n u s o i d a l i n p u t d e s c r i b i n g f u n c t i o n techniques f o r t h e p r e d i c t i o n o f l i m i t c y c l e s i n m u l t i v a r i a b l e non l i n e a r feedback systems. I n t e r a c t i v e d i s p l a y s a r e produced i n t h e conplex frequency domain as an a i d t o conpensator design and t h e method is a p p l i e d t o t h e design o f a c o n t r o l l e r f o r an i n d u s t r i a l a i r conpressor system. 1.
INTRODUCTION
A f e a t u r e o f many p r a c t i c a l process c o n t r o l s y s t e m i s a multi-input multi-output i n t e r a c t i v e s t r u c t u r e with one o r more gross non l i n e a r i t i e s and a p r i mary design t a s k i n such c i r c u n s t a n c e s i s t o ensure t h e avoidance o f l i m i t c y c l i n g c o n d i t i o n s . This p a p e r o u t l i n e s how such a system may be investigated using a sinusoidal input describing function philosophy and a computer program has been developed t o p r e d i c t t h e magnitude and freqaency of p o s s i b l e limit c y c l e modes and t o a i d i n t h e t h e design o f s u i t a b l e conpensation t c ensure avoidance o f c r i t i c a l r e g i o n s i n t h e complex frequency domain. T h i s program r e p r e s e n t s an e x t e n s i o n o f a p r e v i o u s l y developed design s u i t e f o r s i n g l e i n p u t s i n g l e o u t p u t non l i n e a r feedback s y s t e m ( 7 ) and can provide t h e c o n t r o l e n g i n e e r wi?h a powerful i n t e r a c t i v e d e s i p s u i t e which is capable o f d e a l i n g with a wide range of non l i n e a r f u n c t i o n s where t h e hypothesis of harmonic l i n e a r i s a t i o n can be assumed v a l i d . The b e h a v i o w o f a l a r g e conpressor f o r s u p p l y i n g a i r t o a chemical p l a n t i s s t u d i e d and a conpensator designed t o achieve c l o s e d loop s t a b i l i t y and s p e c i f i e d performance. 2.
LIMIT CYCLE PREDICTION I N MULTIVARIABLE SYSTLE1S
S i n u s o i d a l d e s c r i b i n g f u n c t i o n techniques have been used widely i n t h e s t u d y of l i m i t c y c l e behaviour i n non l i n e a r , s i n g l e loop feedback s y s t e m where t h e a p p l i c a t i o n o f conplex frequency domain ana: ys i s has represer-ted an i n t u i t i v e extension c f l i n e a r s y s t e m t h e o r y . The l i m i t a t i o n s of t h e cpproach are, o f c o u r s e , w e l l documented(1) and, a l though t h e r e c e n t emergence o f an algorithm due t o Mees(2) has done much t o enhance the v a l i r ? i t y cf t h e method, t h e assunption o f , e s s e n t i a l l y , l o ? p a s s l i n e a r system elements s t i l l remains &r i m p o r t a n t p r e r e q u i s i t e i n any s t u d y i n v o l v i n g t h e use of t h e d e s c r i b i n g f m c t i o n philosophy. Recent
i n t e r e s t h a s been shown i n t h e a p p l i c a t i o n o f s i n u s o i d a l d e s c r i b i n g f u n c t i o n methods t o m u l t i v a r i a b l e non l i n e a r feedback systems where again advantage can be taken o f i n t u i t i v e e x t e n s i o n s t o a l r e a d y developed l i n e a r s y s t e m t h e o r y . Such an e x t e n s i o n i n v o l v i n g t h e i n v e r s e Nyquist a r r a y was proposed by Rosenbrock i n 1969(3) and l a t e r i n a more d e f i n i t i v e form by Mees(4). More r e c e n t l y Ramani( 5) has proposed g r a p h i c a l methods i n v o l v i n g eigenvalue l o c i o f t h e l i n e a r system elements i n t h e s t u d y of l i m i t c y c l e o p e r a t i o n . Although t h e a p p l i c a t i o n o f t h e i n v e r s e Nyquist a r r a y does p r o v i d e a quick, systema t i c method of compensation design i n t h e complex frequency domain f o r t h e a v o i d a n ~ eof a l l regions of p o s s i b l e t h e o r e t i c a l l i m i t c y c l e b e h a v i o w , i t can be shown i n p r a c t i c e t h a t t h e method is q u i t e p e s s i m i s t i c i n its a c t u a l p r e d i c t i o n of l i m i t c y c l e s ( 6 ) . The o b j e c t o f t h i s paper is t o o u t l i n e a method f o r t h e p r e d i c t i o n o f limit c y c l e operat i o n which i s l e s s conservative than t h a t based on the i n v e r s e Nyquist a r r a y approach b u t which s t i l l maintains t h e c m c e p t u a l l y powerful Nyquist d i s p l a y s t o a i d i n compensator design. In u s i n g t h e i n v e r s e Nyquist a r r a y t h e p o s s i b i l i t y of t h e o r e t i c a l limit cycle o p e r a t i o n i n an i n t e r a c t i v e non- l i n e a r diagonally dominant mulcivariable system can be s t u d i e d g r a p h i c a l l y b ; considering t h e behaviour o f l o c i a s s o c i a t e d w i t h t h e diagonal l i n e a r and non l i n e a r elements of t h e system. Each locus has an envelope defined by a s e r i e s of d i s c s c e n t r e d on t h e locus whose r a d i i a r e r e l a t e d t o t h e magnitudes o f corresponding o f f diagonal elements. P o s s i b l e s o l u t i o n s t o t h e s i n p l f h a r ~ o n i cbalance equations can e x i s t i n any regior. of t h e corrplex plane where t h e envelopes of r e l a t e d l i n e a r and non l i n e a r l o c i i n t e r s e c t . Althcyigh :r!e applicat i o n o f s c a l i n g :echriques can be used with advanta g e ( 8 ) ( 6 ) t > a a j u s t t h e envelspt dirrer.sions of i n d i v i d ~ a ll o c i , i t i s though- ti-.: rethod o f system compensation b c e d on r7.e ,:..cl".a;.ca c f envelope o v e r l a p must t e n d t o he cori;-?i.-.t;.:t. Provide.? t h a t a s y s t e m a t i c metlir;i. c : c c . t~= t i o n is employed it i s , however, not r j ; t f i c u l t 'r c a l c u l a t i a t any freq,uency a l l p o s s i b i e s s i . ~ : l m s t t h e simple s t e a d y s t a t e harr.or.ic t;?-an;* rr.atrl. e q ~ i t i o n over the range of physica-ly e l l c i a b i e s i g n s magnitudes a: d freq-encies an< t c ;r,-~l;;ee i:. t h e complex frecuency dcmair. ir:di v i ?ual s c i u: ior. v e c t o r s which r e p r e s e n t t h e con:ribc?ic:. 3f kc::. diagonal and o f f diagcnal systerr e:emer.tz.
.
Theoretical l i m i t cycle operation can then be determined by a study of these individual s o l u t i o n vectors i n t h e v i c i n i t y of t h e (-1, jO) point. AS t h e magnitude and phase o f each p o s s i b l e s o l u t i o n vector is knm.,the a r e a o f uncertainty introduced i n the inverse Nyquist a r r a y approach by t h e a r t i f i c e o f the Gershgorin envelope is removed and a less conservative r e s u l t ensured. Such an approach is, a f course, only v i a b l e i f t h e computational e f f o r t can be constrained within reasonable born&. In t h e method adopted here t h e steady s t a t e harmoni c balance equations o f all loops a r e considered i n d i v i d u a l l y i n a s e q u e n t i a l manner such t h a t t h e number o f s o l u t i o n s e t s t o t h e o v e r a l l system harmonic balance equation diminishes r a p i d l y as t h e conputation proceeds. For convenience, the s o l u t i o n vectors of each loop can be displayed graphically as a Nyquist a r r a y i n t h e conplex frequency domain. 3.
THE PREDICTION OF LIMIT CYCLES USING A SEQUENTIAL COMPUTATIONAL- PROCEDURE
?he method is a p p l i c a b l e t o t h e system configurat i o n s shown i n Figure 1when t h e following condit i o n s a r e met: (a)
all l i n e a r elements a c t , e s s e n t i a l l y , as low pass f i l t e r s t o j u s t i f y t h e a p p l i c a t i o n o f t h e describing function philosophy;
(b)
i f a l i m i t cycle is present all loops w i l l o s c i l l a t e a t t h e same frequency. In general t h i s is a reasonable assunption p a r t i c u l a r l y i f the non l i n e a r elements a n s i m i l a r i n form a d i f the dominant l i n e a r elements have approximately t h e same frequency c h a r a c t e r i s t i c s ( 5 1.
'Ihe equations governing l i m i t cycle "peration i n a multivariable non l i n e a r system can be expressed as follows
where T is the non l i n e a r r e t u r n r a t i o system matrix with elements t i . ( a j , j w ) and I is a colunm vector o f sinusoids a t $he i n p u t t o t h e non l i n e a r e l e m n t s such t h a t a i = A i / B i . The method c o n s i s t s o f findi n g all p o s s i b l e G l u t i o n s t o equation s e t ( 1 ) over physically allowable ranges of s i g n a l magnitudes m d frequencies. I f we consider f i r s t t h e case o f Limit cycle operation where
w f can, f o r t h e range of p h y s i c a l l y allowable, d i s c r e t e values of the reference s i n u s o i d a 1 = A1 and a t a s p e c i f i c -raiue o f f r e q u e n c y n , determine those f i n i t e p s e t s o f sinusoids ( a l k , a 2 ); k = 1,2,. ,p which w i l l s a t i s f y t h e con8ition &or l i m i t cycle operation i n t h e n t h loop as given by t h e equation
..
I f we now examine t h e equation governing l i m i t cycle operation i n t h e (n-11th loop which is given by
where again al
=
A1,
We can t r y f o r s a t i s f a c t i o n o f t h i s equation using only those f i n i t e p s o l u t i o n s e t s derived from equaticm.(3) o v e r t h e range o f d i s c r e t e values o f A1A s i n p r a c t i c e h i s c r e t e values of A1 a r e used i n the determination o f t h e s e t s o f s i n u s o i d s which w i l l s a t i s f y both equations ( 3 ) and (41, t h e r e could t h e o r e t i c a l l y be s o l u t i o n s which e x i s t o w r t h e d a t a i n t e r v a l s and t h i s could give r i s e t o e r r o r s i n l i m i t cycle p r e d i c t i o n . To guard a g a i n s t t h i s and t o add a s p e c i f i e d degree o f robustness t o t h e conputation we pass from one equaticn t o t h e next, n o t cnly those s e t s o f v e c t o r s which s a t i , s f y t h e f i r s t equation e x a c t l y b u t a l s o those s e t s which y i e l d an approximate s o l u t i c n . The degree o f approximation is r e l a t e d d i r e c t l y t o t h e anplitude, frequency and phase i n t e r v a l s chosen and a l s o t o s p e c i f i e d e r r o r bomds on the system parameters. These approximate s o l u t i o n s a r e generated automatically within t h e program once t h e d a t a i n t e r v a l s and parameter e r r o r bounds are s p e c i f i e d by t h e user. I n the f i n a l g r a p h i c a l display cmly those e r r o r c i r c l e s associated with t h e system parameter e r r o r bounds need be displayed. I f q s o l u t i o n s e t s o f s i n u s o i d s a r e found (where q 6 p) which s a t i s f y both equations ( 3 ) and ( 4 ) t h e n t h e (n-21th loop equation is n e x t examined f o r harmcnic balance using these q s e t s and the process repeated u n t i l loop 1 is reached. This has a harmonic balance equation which can be p u t i n the form
which i s p l o t t e d as shown i n Figure 2 f o r a range of values o f n. I f t h e end of m e l e f t hand s i d e s o l u t i o n vectors touches t h e (-1, jO) p o i n t o r i f t h i s p o i n t l i e s within an e r r o r c i r c l e associated with one of t h e s o l u t i o n v e c t o r s , then a conditicn could e x i s t which would s a t i s f y t h e o v e r a l l steady s t a t e harmonic balance requirements o f equation (1).
It should be noted t h a t as soon as any equation is not s a t i s f i e d by a s e t of s o l u t i o n sinusoids passed from t h e previous equation then t h i s s e t is immed i a t e l y discarded. When s o l u t i o n s have been obtained f o r a l l loop equations over t h e allowable range o f A 1 and Q the i n e q u a l i t y
,...a
is n e x t considered and t h e above conputational procedure repeated in a cycLic manner f o r a range of physically r e a l i s a b l e values o f t h e reference sinus o i d a2 = A?. The d e s c r i b i n g function equation f o r loop 2 is f i n a l l y derived a s
where a = A1 1
.
i n a s i m i l a r manner t o t h a t f o r loop 1 and displayed i n the conplex frequency domain. The computatior. p r o c e e d s . u n t i 1 d i s p l a y s a r e produced f o r a l l n loops i n t h e system. Limit cycles could t h e o r e t i c a l l y occur i f any s o l u t i o n v e c t o r o f any o f t h e n loop equations touches t h e (-1, jO) p o i n t o r , f o r a more conservative p r e d i c t i o n , i f t h i s p o i n t l i e s within one of t h e e r r o r c i r c l e s a s s o c i a t e d with any solut i o n vector. A study of such a Nyquist a r r a y can quickly d e t e r
mine t h e magpitude and frequency of a l l t h e o r e t i c a l l y p o s s i b l e modes of l i m i t cycle operation and isol a t e those loops which a r e responsible f o r o v e r a l l system i n s t a b i l i t y . Of equal inportance, i t can a l s o i n d i c a t e t h e form of compensation necessary t o avoid t h e p o s s i b i l i t y of limit cycle operation i n a p a r t i c u l a r loop. I t is o f t e n convenient,when a preconpensator element is envisaged,to use a Nyquist display of t h e e l e -
ments of T. Such a display i s e a s i l y generated have been once all s o l u t i o n s e t s f o r the vector determined using t h e method o u t l i n e d above. ?he corresponding harmonic balance equation f o r loop 1 now has t h e form
. ht.
where
(j
?CZ j d j 4.
=
...n)
1..
l o c i ot ?..c rl-g,:.ae 1 e r n r . t ~ tii. (i=l,.. n ) . For any cndnge ir. loop m the diago::al locus of t h e i t h loop i s m a l t e r e d h c t the magpitude and phase o f t h e superinposed s o l u r i o n v e c t o ~ sw i l l a l t e r . Any r.ew v a i i e I ,sz- :?cme:.r? iie irithin bounds defined by a c i r c l e with r a d i u s
centred on t h e m a l t e r e d diagonal locus. This c i r c l e can e a s i l y be computed and displayed i f necessary. The corresponding l i n e a r preconpenoator
a c t s on t h e colunns of T rendering t h e Nyquist displays of T no longer very u s e f u l f o r compensator design. Such a precompensator has howev$r t h e e f f e c t o f multiplying t h e rnth row o f T by km(j w ) where again t h e value of m is chosen a s t h a t of t h e most troublesame loop. Desigp i n t h e frequency domain gan thus proceed by using displays o f t h e v e c t o r tii ( i = l , . .n) with superimposed s e t s of s o l u t i o n vectors
.
a r e t h e elements of
COMPENSATION TECHNIQUES
In Figure 1 a p o s t conpensator o f t h e form
w i l l have the e f f e c t of multiplying the rnth row o f T by Em(jw) where the value of m is chosen a s t h a t of the most t r o u b l e s o w loop. Every vector e l e m n t i n t h e rnth loop equation w i l l be s i m i l a r l y affected. As we a m dealing with a Nyquist a r r a y t h e r e q u i r e d form of (,L jw) t o s u i t a b l y compensate t h e rnth loop proceeds n a t u r a l l y from a knowledge e f c l a s s i c a l s i n g l e loop feedback theory. N l o t h e r vectors tii; ( i = l , . .n)ifm w i l l remain u n a l t e r e d a s w i l l passed t o t h e rnth equation. t h e s o l u t i o n values of However, t h e s o l u t i o n values of passed out o f t h e rnth equation w i l l be a l t e r e d which i n t u r n w i l l vary t h e vectors
.
The e f f e c t of such a compensator on t h e behaviour of t h e o t h e r loops can most e a s i l y be v i s u a l i s e d i f we consider a l l vectors
t o be superimpmed a t s p e c i f i c frequencies on t h e
5.
COWUTATIONAL ASPECTS
In t h e s e q u e n t i a l conputational procedure up t o twenty magnitude values of t h e reference s i n u s o i d and up t o twenty values o f frequency are used i n t h e determination o f s o l u t i o n s vectors f o r any loop equation. The number o f phase angles is a l s o n o r Although t h i s may maily r e s t r i c t e d t o twenty. i n i t i a l l y r e s u l t i n a l a r g e number o f s o l u t i o n s e t s most of these a r e p r o g r e s s i v e l y discarded a s the conputation proceeds and t h e r e is a r a p i d convergence on t h e f i n a l displayed s o l u t i o n s e t . The frequency and amplitude ranges of p o s s i b l e l i m i t cycle behaviour a r e automatically displayed i n t a b u l a r form along with t h e most probable l i m i t cycle a m p l i t u d e s and frequency. Once c r i t i c a l a n p l i t u d e and frequency ranges have been a s c e r t a i n e d irom t h e t a b l e and an i n i t i a l study of t h e d i s p l a y s , t h e a m p l i t u d e and frequency ranges used i n the calculaticn can be conpressed t o give a d d i t i o n a l s o l u t i o n s e t s i n t h e c r i t i c a l region o f t h e complex f r e q u e n c ~ plane i f required. The robustness of t h e s o l u t i o n i s enhanced by compensating automatically f o r the i n t e r v a l s i n t h e d i s c r e t e data used i n t h e conputation and by allowi n g f o r known e r r o r bounds on any element o r group of elements within t h e feedback system under consideration. Such bounds can be r e l a t e d d i r e c t l y t o bounds on i n d i v i d u a l s o l u t i o n vectors and i t i s thus conputationally easy t o produce corresponding e r r o r c i r c l e s on t h e conplex frequency plane. For c l a r i t y of p r e s e n t a t i o n only t h a t e r r o r c i r c l e
c l o s e s t t o t h e (-1, jO) p o i n t need be s t o r e d and displayed a t a y frequency p o i n t cn t h e c o l p l e t e s o l u t i o n w c t o r . This t e c h n i q r r ca a l s o be used t o obtain a m a s u r e o f t h e s e n s i t i v i t y o f system l i m i t c y c h behaviour t o v a r i a t i o n s i n one o r more o f t h e system's p h y s i c a l parameters. A l l design program a r e w r i t t e n in Fortran I V with a chain f a c i l i t y t o reduce core s t o m requirements t o l e s s than 10K words. Conputer questions and mswers use lrnemonic English format and a command i n t e r c e p t o r r o u t i n e provides c o q l e t e f l e x i b i l i t y i n e n t e r i n g o r l e a v i n g a n y . p r o p m block. This allow8 r a p i d change o f e i t h e r l i n e a r o r non l i n e a r d a t a and innnediate access t o t h e display p r o g r a ~ ~ t o evaluate the e f f e c t uf any d e s i w decision. Linear d a t a can e i t h e r be e n t e r e d i n pole wm o r polynomial form. Non l i n e a r d a t a can be e n t e r e d i n polynomial o r piecewise l i n e a r form o r a s a s e t o f d i s c r e t e measured d a t a points. The v i s u a l display u n i t f a c i l i t i e s include l i n e o r p o i n t format, alphzr numeric e n t r y , a x i s display and automatic s c a l i n g .
6.
THE DESIGN OF A CONTROLLER FOR AN INDUSTRIAL COMPRESSOR
The system which is s h a m diagrammatically i n Figure 3 c o n s i s t s o f a l a r g e conpressor f o r supplyi n g air t o a chemical p l a n t and has been s t u d i e d previously as a purely Linear system(9) (10). The o u t l e t p n s s u r e is c o n t r o l l e d by guide vane p o s i t i c n and t h e o u t l e t flow by means o f a blow o f f valve which' prevents surging a t i n p u t flows below 20% o f f u l l flow r a t e . When t h e blow o f f valve is closed the system c o n s i s t s o f a s i n g l e c c n t r o l loop and good r e g u l a t i o n can be obtained using a standard pneumatic c o n t r o l l e r . When t h e blow o f f valve opens however the system has two i n p u t s and e x h i b i t s l i m i t cycle behaviour when t h e loop gains are increased t o any u s e f u l value. The dynamics o f t h e l i n e a r system kre obtained(l0) with t h e blow o f f valvc open and are shown i n Figure 4, which now includes a s a t u r a t i o n eJement t o represent t h e non l i n e a r c h a r a c t e r i s t i c s o f the blow o f f valve loop. It is d e s i r e d t o design a preconpensator K(s) which w i l l guarantee avoidance of l i m i t cycle behavisur, a f a s t blow o f f valve response, low i n t e r a c t i e n and zero o f f s e t i n t h e pressure loop. Various sources of e r r o r e x i s t i n t h e s i n p l i f i e d system model shown i n Figure 4 due t o t h e following factors. (1)
The l i n e a r system elements given by t h e model only approximate t o t h e a c t u a l p l a n t behaviour.
(2)
A very s i m p l i f i e d model o f t h e blow o f f valve non l i n e a r i t y is used.
(3)
A non l i n e a r c h a r a c t e r i s t i c is known t o e x i s t
i n t h e guide vane linkage. D e t a i l s o f t h i s were n o t a v i l a b l e b u t it was known t h a t l i n e a r operation could be achieved i f required. (4)
Valve dynamics a r e not included.
To ensure a conservative design i n t h e presence o f these modelling e r r o r s and e r r o r band o f 230% was used on a l l elements o f t h e - r e t u r n r a t i o matrix T throu&out s u b s e ~ a e n tc a l c u l a t i o n s . An i n i t i a l comp u t e r i n v e s t i g a t i o n f o r l i m i t cycle operation using
a proposed c o n t r o l l e r i n t h e pressura loop with a gain o f 1 5 pmduced t h e displays shcmn i n f i g u r e 5 which include t h e appropriate e r r o r c i r c l e s . Ihe system h a s obviously a very high degme o f i n t e r a c t i o n a n d ' a l i m i t c y c l e was e d i c t e d with t h e f o l lowing moet l i k e l y paranmter values w = 1.35 Rad/sec, u = 7.04, ly = 2.75. A subsequent simw l a t i o n confirmed lidt cycle operation with pararmter values w = 1.45 Rad/sec, ul = 6.6, u2 = 1.5. To reduce i n t e r a c t i o n e f f e c t s a preconpensator o f constants w a s & s i p a d using s t a n d a r d frequency domain techniques(l0) t o mat c l o s e l y diaponalise t h e system a t a s p e c i f i c frequency and anplitude set. A frequency o f 1 Radian/sec was chosen on t h e b a s i s o f e a r l i e r work with a conpletely l i n e a r i s e d system. The displays o f the system b e h a n o w with only t h i s compensator are shown i n f i g u r e 6. It w i l l be seen t h a t no limit cycle is now predicted. A s we a r e now considering a preconr pensator element t h e inverse displays a r e now used with e r r o r c i r c l e s corresponding t o system parameter e r r o r b o m b o f t30% being determined by u s u a l mapping procedures. I n t e r a c t i o n e f f e c t s i n loop 2 a r e n e g l i g i b l e and a r e a l s o small i n loop 1 f o r low values o f u2. The gain i n loop 1 is now increased by a f a c t o r 5 and it w i l l be seen from Figure 7 t h a t t h e e f f e c t on loop 2 is small. I f t h e gain i n loop 2 is subsequently .increased by a f a c t o r 2 the system s t i l l remains f r e e from l i m i t cycle operation b u t loop i n t e r a c t i o n e f f e c t s now increase i n loop 1 as shown i n Figure 8 s o t h i s loop gain was n o t adopted. The f i n a l conpensator was o f t h e form
The multivariable c i r c l e c r i t e r i o n can now be i n voked t o ensure a s y n p t o t i c system s t a b i l i t y and the absence o f l i m i t cycles f o r a l l bomded s e t s of inp u t s and t h e appropriate displays a r e shown i n Figure 9 where, following t h e arguments o f ;ook(+l), we can use t h e s c a l e d mean Gershgorin bands of GK and the c r i t i c a l d i s c o f nqp. Examination o f t h i s p l o t shows t h a t loop 2 w i l l be s t a b l e and a l s o t h a t loop 1 w i l l be s t a b l e with the blow o f f valve loop f u l l y closed assuming t h a t the l i n e a r p a r t o f the model is s t i l l v a l i d m d e r these conditions. The response of the closed loop syster. t o a s t e p change i n pressure demand is shown i n Figure 10 and the corresponding response t o a s t e p change In a i r fiow demand is shown i n Figure 11. In both czses f a t s e t t l i n g times and low i n t e r a c t i o n e f f e c t s a r e achi$ved. Any o f f s e t i n t h e pressure loop car1 be remedied by long time constant i n t e g r a l zc'icn which should nDt a l t e r the p l o t s s i g n i f i c a n t l y . The o f f diagonal t e r n s of (i are l a r g e r than zne diagc:;ai term thus i n d i c a t i n g t h a t The icops are being e f f e c t i v e l y inter(.hangc.d. Th:s could ~ e s u i t i n i n e f f i c i e n t operation, although t h e s i t c a t i a n could be improved by allowing a g r e a t e r degree of loop i n t e r a c t i o n . I r should be noted t h a t t h e des i g n method could have been s e d t o equal e f f e c t if t h e model of t h e non l i n e a r i t y i n loop 1 had beep available f o r i n c i s i o n .
Control S y s t e m Centre Repost No. 2 6 2 , U.M.1 .S .T., September 1974. Provided t h a t t h e e s s e n t i a l p r e r e q u i s i t e s e x i s t f o r i t s a p p l i c a t i o n t h e s e q u e n t i a l conputational procedure p r e s e n t s an accurate, systematic and e f f i c i e n t method of p r e d i c t i n g l i m i t cycle operation i n a multivariable non l i n e a r feedback system. By studying the behaviour of i n d i v i d u a l s o l u t i o n vect o r s i n the region of t h e (-1, jO) point,accurate p r e d i c t i o n s can be made a s t o t h e m a e i t u d e and f r e quency o f p o s s i b l e l i m i t cycles. A s p e c i f i e d degree of robustness can be introduces i n t o t h e method by considering e r r o r bounds on t h e p h y s i c a l elements within t h e system matrix and t h e r e is a d i r e c t and obvious r e l a t i o n s h i p betweer: such bounds and t h e behaviour of s o l u t i o n vectors i n t h e conplex frequencies domain. This f a c i l i t y can a l s o b e used t o r e l a t e v a r i a t i o n s in l i m i t cycle b h a v i o u r t o variat i o n s i n one o r more system parameters t h u s allowing t h e i d e n t i f i c a t i o n of p a r t i c u l a r l y s e n s i t i v e system e l e m n t s . I f limit cycles a r e p r e d i c t e d by an examination of the s o i u t i o n v e c t o r s i n any one loop then t h e s t r u c t u r e of a s u i t a b l e preconpensator o r p o s t compensator t o a d j u s t t h e s e ve.:tor s o l u t i o n s i n t o l e s s c r i t i c a l regions of t h e complex freqwncy domain i s immediately obvious from an elemenrary knowledge of s i n g l e 1.oop feedback theory. h e influence of such a compensator on t h e o t h e r loops is constrained within e a s i l y determined bounds and, provided t h a t t h e system i s n o t grossly i n t e r a c t i v e , a viable m t h o d of system conpensation may be devised.
(7)
GRAY, J.O., SAWIDES, L.S. and TALOR, P.M., I n t e r a c t i v e graphics design using sinusoiZal describing f u i c t i o n methods, I .S.A. Transact i o n s , '101. 1 4 , No. 1, 1975. pp. 68-77.
(8)
RAMANI, N . and ATHERTON, D.F., Correspondence, Proc. IEE, Vol. 20, No. 7, i973. p . El.&.
(9)
LUNTZ, R., MUNRc3, N. a d MCLEOD, R.S., Computer aided design o f m u l t i v a r i a b i e c o n t r o l systems, IEE 4 t h WAC Convention, Manchester, 1971. pp. 55-59.
(10) ROSENBROCK, H.H., 'Ihe computer aided design of c o n t r o l s y s t e m , Academic P r e s s , 1974. (11) COOK, P. , Modified m u l t i v a r i a b l e c i r c l e t h e o r e m . Recent mathematical developments i r c o n t r o l , Academic P r e s s , 1973. pp. 367-372.
The s e q u e n t i a l computation method is thought t o be t h e only e x i s t i n g method f o r s t u d y i n g multivariahle systems with gross non l i n e a r i t i e s i n o r d e r t o a s c e r t a i n l i m i t c y c l e operation and t o derive both pre- and p o s t compensator elements which ensure closed loop s t a b i l i t y with known i n t e r a c t i o n effects.
FIGURE 1
REW ENCES
(1) ELB, A. and VANDERVELD, W.E., Multiple i n p u t describing functions and non l i n e a r s y s t e m design, McGraw H i l l , New York, 1969.
(2)
MEES , A. I . and BERGEN, A . R. , Describing F ; m c t i ~ n sRevisited, E l e c t r o n i c s Research Laboratory, University o f C a l i f o r n i a , Berkeley, 1974.
(3)
RXENBROCK, H . H . , Design o f m u l t i v a r i a b l e cont r o l s y s t e m using t h e i n v e r s e Nyquist a r r a y , Proc. IEE., Vol. 116, No. 11, 1969. pp. 1929-1536.
(4)
MEES, A. I . , Describing functions c i r c l e c r i t e r i a and multiloop feedlack systems, Proc. IEE., Vol. 120, No. 1, 1973. pp. 126-130.
(5)
(6)
,SAKANI, N. and A'TiAERTON, D.P., Frequency response nethods f o r non i i n e a r multivariable system', Proc. CanadiaCorif. Auto. Ccnt., 1973. 9.2-1 9.2-15. GRAY, J.O. and TAYLOR, P.M., & s c r i b i n g functions and t h e inverse Nyquist a r r a y ,
SOLUTION VECTORS
' I
. I
FIGUFE 2
.,
I
?5:.
: :~~..:. , 'a$~,
PLOT OF
Q
LOOP 2
GAIN X5 LOOP 1 GAIN X1
LOOP 2
+
S
LOOP 2
+
1
+ J
NYQUIST ARRAY
/
LOOP 1 M O P 1 GAIN X5 LOOP 2 GAIN X 2
LOW 2
.
FIGURE 9
;.,
9' ,
STEP CHANGE PESSURE
...
M
,:
,F
, LC
SECS
10 ,,,
VALVE POSITION
. r-o OPEN
FIGURE 10
PRESSURE
STEP CHANGE FLOW
I
> o
VALVE POSITION
OPEN
FIGURE 11