- Email: [email protected]
Canonical structures in the identification of multivariable systems
Automatwa, Vol 11, pp 361-374 Pergamon Press, 1975 Printed in Great Britain
Canonical Structures in the Identification of Multivariable Systems* Structures Canoniques dans l'Identification de Systbmes Multivariables Kanonische Strukturen bei der Identifikation mehrvariabler Systeme KaHom~,~ecKHe
cTpyKTypbI
B I ~ e H T I , IqbHKaRttI,I MHOrOCB~I3aHttI, IX C t t e T e M
ROBERTO
GUIDORZIt
The idennficatwn of multwanable systems can m no way be constdered as the simple esttmatwn of a set of parameters, the determmatwn o f a suttable model structure ts an equally fundamental step. Summarymln this paper a umtary ldentlficatmn procedure for linear muluvarmble systems is described The approach considered is based on a prehmmary canonical structure identification, i • on the d~rect determination, from inputoutput sequences, of a set of mvarmnt indexes completely describing the input-output structure of the system Reduced computation ume and storage, canonical representation of the identified model and easy estimation of the mmal state of the system are in this way achzeved The results obtained from the application of this procedure to a simulated process and to a chemical reactor are also reported
all the m o d e m available control theory regards these representatlons. (2) The class of the input slgnals should not be restricted to such partlcular inputs as impulse functions, step functions, colored or whlte nolse, smusoldal signals,P.R.B.S. etc.,since these inputs can seldom be apphed to real processes. It is on the contrary desirable to use normal operating records of the plant. T h e cholce of a state-space, z.e parametrlc, model imphes the prehnunary determmataon of the system order, Itis m fact well known [2] that parametric models can glve large errors when the order of the model does not agree wlth the order of the process. Another practlcallyrelevant reqmrement concerns the mmlmahty, or at leasta strong lim]tatlon, of the number of the slgnificantparameters of the model and this means the choice of a statable canomcal, or quas]-canomcal, form. Moreover, the s~mphclty of the hnk between the selectedcanonlcal state-space representatmn and the set of inputoutput equations used for the identification,the mmal data are input--output sequences, beawly condmons the reqmred amount of computatmn. The prewous problems have been completely solved for singleinput single-outputor multl-mput single-output systems [2]. Troubles appear when multJvarmblc, Lc multl-mput multi-output, systems are consldered smcc a umtary approach satlsfymg prcwous reqmrcments cannot be found m the lltcraturc. The purpose of thlspaper is to bridge thlsgap by exploring a class of state-space canonical models hnked by partlcularly slmple relatmns to inputoutput d~fference dcscnptmns that can be dlrectly identified from the input--output sequences. The resulting]dent]ficatmnprocedure ]s slmflarto those
1 INTRODUCTION ACCORDING to the classical d e f i n m o n gwen b y Z a d e h [I], ~dentfficatmn ~s the deterrmnaUon, on the bas~s o f i n p u t a n d o u t p u t , o f a system within a specified class o f systems, t o which the system u n d e r test ~s e q m v a l e n t Th~s d e f i n m o n ~s q m t e general a n d allows m a n y degrees o f f r e e d o m m the pracUcal f o r m u l a t m n o f the ~denUficatmn p r o b l e m since the class o f the m o d e l s to be ~dentlfied a n d the class o f i n p u t signals are n o t specified a n d also the m e a n i n g o f ' e q m v a l e n t ' ~s often referred to a s t a t a b l e o p t t m a h t y criterion F o r the pracUcal usefulness o f a n ~dentfficaUon p r o c e d u r e the following c o n s i d e r a t i o n s s h o u l d be t a k e n into account" (1) T h e o b w o u s choice for the class o f the m o d e l s ~s gwen by state-space r e p r e s e n t a t m n s since a l m o s t * Received 11 June 1973, revised 11 September 1974. revised 15 January 1975. The original version of this paper was presented at the 3rd IFAC Symposmm on Identification and System Parameter Esumatmn which was held m The Hague/Delft, The Netherlands, during June 1973 The published proceedings of this IFAC meeting may be ordered from Elsevier/North Holland P O B 211, Amsterdam, Netherlands It was recommended for publication m revised form by associate editor J Ackermann 1"Istltuto di Automatica, Umvers~ty of Bologna, Vlale Risorglmento 2, 40136 Bologna, Italy. 26
361
362
ROBERTO GUIDORZI
for single-input single-output systems and can be summarized in the following steps (1) Structural ldentaficataon, i.e. determmaUon by direct treatment of the input--output sequences of a set of indexes completely describing the structure, not only the order, of a canonical model for the process (2) Parametric ldentafication of the subsystems, one for every output component, of an inputoutput descnptaon of the system with the prewously identified structure (3) Deduction from the identified input--output representation of the state-space model
Corapartson wah premous hterature A fundamental contribution in flus field can be found m the paper by Ho and Kalman [3]; in this work a state-space model is deduced from an infinite sequence of Markov parameters Interestmg results have also been given by Tether [4] who studies the equivalence between finite sequences of Markov parameters and state-space models and by Ackermann [5, 6] who gives a canonical reahzation procedure by the use of an input--output descnptaon. A closely related approach has been described by Mayne [8]. Gopinath [9] assumes an initial knowledge of the system general inputoutput sequences; by direct elaboration of the data, a set of input-output equations is then obtmned. The procedure ~s, however, based on the knowledge of a statable data selector matrix correspondmg to a well-defined structure of the equaUons and an efficient method for the construction of this matrix is not g~ven. The approach by Gopmath is also the basis of a paper by Budin [10], who improves the computational reqtarements and points out a general procedure leading to a minimal but not canomcal reahzation. The basle difference between the approach presented in this paper and the previous ones lies in the structural identaficatlon of the system, i.e. m the direct deternunation from the input--output sequences of a set of indexes like those described in a dual context in [11, 12] which permit a decomposaaon of the system into subsystems that can be separately identified. This approach rehes heawly on the structural properties of blockcompamon state-space representations for multivariable systems and on their equivalence with input-output canonical polynormal representations Some early results on these propertaes had been obtained m the first half of 1969 [13] and have been developed in [7, 14, 15]; apphcatlons to both simulated and real processes can be found in [14, 16-18] The structural identification approach has been considered also by Tse and Wemert [19] and Chow [20].
2 ORGANIZATION OF THE PAPER
The class of linear discrete multlvariable systems considered can be represented by the usual statespace equations x(k + 1) = Fx(k) + Gu(k), y(k) = Hx(k),
(2 1a) (2 l b)
where x belongs to the n-dimensional state-space gF, u to the r-dimensional input space qz' and y to the m-dimensional output space ~ The lmtial data are quite general synchronous input-output sequences, i e normal operating records of y and u The paper is orgamzed as follows. In Section 3 some propertaes of canonical (blockcompamon) state-space representataons are recalled and the structural and parametric equivalence between this class of models and input-output difference descriptions of the type P(z) y(k) = Q(z) u(k),
(2 2)
where P(z) and Q(z) are polynomial matrices m z, where z -x is the unitary delay operator, is estabhshed In Section 4 the structural ldentlficaUon of equations (2. I) is considered, the subsequent parametric ldentafication of an input-output model of type (2.2) follows m Section 5. In Section 6 the algorithms given m Section 4 and in Section 5 are extended to the case of data corrupted by additave noise with zero mean and a maximum hkehhood estimator is described In Section 7 the derivation of a class of canomcal state-space models from the identified input-output description is discussed and m Section 8 some practical apphcations of the proposed identification procedure are reported Some concluding remarks are given in Section 9 3 CANONICAL INPUT-OUTPUT STRUCTURES OF MULTIVARIABLE SYSTEMS
Some brief recalls about state-space canomcal forms are first given because of their use m the following For a proof of the gwen results see Luenberger [21], complete and useful papers on this argument have also been written by Wang and Davison [22] and Rlssanen [23] where the meaning of the term 'canonical' ~s deeply investigated Consider the discrete, completely observable, multivanable system represented by equation (2 1) and let H =
h~T
(3 1)
hmT
Construct then the vector sequences 111,F T h1F w~hi,
,
(3 2)
hm, F T hm,
,
Canomcal structures m the ~dent~ficatmn of multwanable systems and select them m the following order h i, h~,
, h,,,, F T h 1,
, F r h,,, l zT' h i,
Remark 2 The dependent vectors g r ' ! h are a hnear combmatmn, w~th coefficients gwen by the tth slgmficant row of A, of prewously selected ones, l e
Fr~'h, = a~xahx+
The original system has been decomposed m this way into m interconnected subsystems, a direct consequence of the structure of the matrices A and C is gwen by the following lemma. Lemma 2
(3 4)
The) th subsystem m the canonacal representatmn (A,B,C) is completely observable from the j t h component of the output An input-output difference equataon of the type (2 2) wdl now be deduced from the triple (A, B, C) Remarkable features of this descnptaon are that ~t exhibits the structure of the system and that also the parameters of the triple (A, B, C) can be easily obtained from it, and vice versa. Let
the system equations became
w(k+ I) ffiAw(k) +nu(k),
(3 5)
y(k) = Cw(k),
(3 6)
and the matrices A, B and C exhibit the following structures X = T F r -~ = {A.} (3 7)
[0 ] [o o]
(t,j = 1,2,. ,m) where the (v~x v~) mamces /q, show the companion structure (3.8) and the (v~ x v~) matrices A,~ are of the type (3 9)
w(k) = (wk.l w~.~ w~,,F, y(k) = (Yk.l
o
,
a.. 1
u(k) = (Uk,~Uk,~.. Uk.,)v,
(3 8)
I e. for a simpler notation wg,1 = wx(k), Wk,~ = w~(k) and so on Consider then the j t h subsystem, because of the complete observablhty of ~ts state from t b e j t h output component and of the canomcal structure of the matrices A and C it is easy to
a~,v,
0
0 a.,~j
C = H T -1 =
0
.,
' (3.9)
0]
0
0
0
1
0
0 0
,oo
1
0
..
t 1
Y~,~)r
and
Ij,r l
A.jffi
+a~mahm+a~l,~F'rhl+
Ihm,. ,
F'r("-l) h.,},
a~,l
b~r
The number, v., of nonzero coefficients m the generic matrix A~j, because of the order followed m the selectmn of the vectors (3 2) and of the consequent structure of the matrix T, ~s at most equal to v.+ 1 f f j < l , to v~ f f j > z [7] [14]
The set of integers v1, . , v,n completely defines the structure of the couple (F, H) and Is mvarmnt with respect to a change of coordinates The proof directly follows from [11,12] by means of obwous duahty considerations As is well known, by operating the change of coordinates leading to the new state vector w = Tx, where
A,~ =
bnl
Remarl~ 1
Lemma 1
,Irr(~-l)hl[
bit (3 11)
[ bn T
retaining a vector FT'lh m (3 3) ~f and only ff ~t ~s independent from all previously selected ones Note that when F T ' ~ ~s tested all the vectors Fq'Jh, (0 ~ j ~
T T = {hvFThl,
bll
B = TG = I bit
(3 3)
363
(v1+ 1)
(vx +
+ vrn-1 + 1)
0
(3.10)
364
ROBERTO GUIDORZI
derive from equatmns (3 5) and (3 6) the following expressions
The substltutmn of (3.13) in (3 5) leads to the Input-output descnptaon
W/~Jvi+ +~t-x+l}= Yka,
{(zI-- A) V(z)} y(k) - {(zI - A) WZ(z) + B} u(k)
Wk,¢v~+ +~j_x+2~----"zy~.~-b~,+ +vj_,+l~u(k),
(3 17)
wk:,,+ +~,_,+a~ = : Y k ~ - b ~ , + +,i-~+2~u(k) -
Wk,{vx+ +i,j)=
bTx+ +,j-x+x~zu(k)
ZVj-lykj-- b~x+ a' -- b{vl+
+vj-l) u(k)-
+v.~-t+l) Z "-2
u(k).
.
(312)
It must be remembered that here and m the followmg z -t is the umtary delay operator; the notatmn z~ u(k) xs therefore eqmvalent to u(k + 2) and so on. By wntang equatmn (3 12) for j = l, ,m the whole state vector can be expressed as w(k) = V(z) y ( k ) - WZ(z) u(k),
In representatmn (3 17), however, only the vlth, (vt+va)th, ,nth equatmns are slgmficant; the remaining ones are simple ~denttt~es The s~gmficant equations m (3 17) can therefore be written in the form P(z) y(k) = Q(z) u(k), (3 18) where
[ pn(z)
,
pm (z)
where
qu(z)
qlr(z) ]
qmt(z)
qmr(z)
.
-au,~z-a~ta,
p~j(z) = - a2j,vo zV,J- 1 -
Z px-1
0
0
0
(3 14)
- a~,2 z - aoa
1
q~(z) =fl(.x+ +,,)oz~'-l+
(3.22)
+fl{.x+ +~,-,+2~,Jz
+fl(~,+ +~,_,+iI,j, ..
(3 21)
The polynomials of Q(z) are obtained by computing the right slde of expression (3.17), by simple passages it follows
z
0
(3 20)
The polynomials of P(z) can be Immediately obtained from (3.17), in fact it follows that p,~(z) = z~,-ai~,v,z ~ , - 1 -
v(z) =
(3 19)
pm, (z)
Q(z) =
(3 13)
0
px,(z) ]
P(z) =
(3 23)
z "--t
0
bT hT
0
0
W ffi
[n × r ( v u -
1)]
(3 15)
0 b,~,..+x
Z(z) ffi
ix] zI
z ~t-11
VM -- max~ (v0.
0
bnT-~.+l
(3 16)
0
0
where the coefficmnts ilia are the entries of the matrix 1] = M B =
(3 24)
Canomcal structures in the ldentificatmn of multwanable systems
365
and the matnx M is gwen by (see also [7] ) - -
a11.2
-- all,a
- -
--
a11,8
- -
all,v a
1
-- alto,2
all,,I
--
alln~pl
- - a11,~ 1
1
-- all.?t,vlm
aim~
m
0 0
°
1
0 (3.25)
M = -- amt,2
...
- - a m x , v,. x
- - arax,a
0
-- atom,2
0
-- amLp. x
-- amra, a
0
-- amm,v ."
0
1
Remark 3 It is useful to observe that (see also Remark 1) the polynomials m P(z) and Q(z) satisfy the following relations" deg {pa(z)} > deg {p~(z)} for j > i, deg {p,(z)} I>deg {p,j(z)} for j < t, deg {p,(z)} > deg (p:~(z)) for t # j , deg {p,(z)} > deg {q~(z)}. The eqmvalence between the class of canomcal state-space representatmns (3.5)--(3.6) and the input--output descnptaon (3 18) has thus been estabhshed.
Remark 6 To obtmn the matrix B from the knowledge of P(z), Q(z) it is first necessary to construct, by direct inspection, the mamces M and I]; matrix B Is then glven by B -- M -1 I] Note that M is always nonsmgular because of its structure, m fact det (M) --- 1. Example Let
Remark 4 The set of indexes vx,.., V,n can be lmmedmtely deduced by inspection, m&fferenfly from the knowledge of A or of P(z); also from the parametric standpoint A and P(z) are eqmvalent. Matrix C can be chrectly written ff va, , vm are known
A=
Remark 5 To obtain the matrix Q(z) from the knowledge of the couple (A, 13) it ~s first necessary to construct the matrix M, which can be written by direct lnspectaon of A, the elements of !] = MB are then the coefficients of the polynommls q¢t(z).
B=
B=
-- amra,v.
[Ol0i ] -
'
3
ffi
flr~
1
1
03
-
2
0
13
1
2
"
1 ,
0
1
I
1
-
0 0 0 1 0
Then vl = v~ ffi 3, v2 = 2; the coetticlents of the polynonuals q.(z) are the elements of the matrix l] gwen by
J 0
--2--10
z~+2z+ I
01
,
C--
0 1
21
i
1
1
1
0
1 --
-
1
-12
1
- 1
The equivalent input-output descriptaon is therefore -z~-2z+l
-
-1
1 0 0 0 0
0
i
- 1
-
-1
-1
][33110]o [ 8:1 8 . fl~, 831 tim 8,, f t .
02
y(k) ffi
-z-
1
z+ 4 1 u(k). z+2
]
0
1
366
ROBERT{) GUIDORZI
Remark 7
4 STRUCTURAL IDENTIFICATION
Definition 1
The integer N in (4 3) must be large enough
The structural identification of a multivanable system is the detenmnaUon of the set of integers % , v~ defining the structure of the couple (A, C) from input--output sequences without the intermedmte construction of a parametric model The structural ldentffication will be performed by taking advantage of the fact, given by Remark 4, that this mformataon concerning a canonical statespace representauon is also evident In the inputoutput difference model (3.18) Consider the sth equation in the set (2 18), i e m
P
~ps.,(z)yk,, = Y, qs,,(z)uk,,,
(4 1)
~,=I
~,~1
that, on the basis of the relations (3 21), (3.22) and (3.23), can also be written In the form 111 /,'m~ ~=1 j - I
+ ~ ~fl(,,+ +~._x+,,.,u(~+~_x),, (4.2) ~=1 1=1
vss = vs).
(In (4 2)
C o n s i d e r n o w the matrix o f
[N>n+rv M, where vM=max,(v,) ] in order to peru'at the selection of the necessary number of independent vectors, in other words, the inputoutput sequences must be of sufficient length to permit the complete structural identification of the system. Usually the number of available data ~s many times the system order so that the previous condmon is largely fulfilled
Remark 8 For the structural, and subsequent parametric, identification the input sequence must satisfy to well-known condmons [2], i e must 'excite all the modes' of the system. Also this requirement can be easily met when the length of the data is large with respect to the system order It is certainly not advisable to carry out the structural ldentificataon directly on the vectors of (4.3) because of the large amount of storage necessary, moreover, the required storage would be a function of N Since for every matrix, D, rank(D) = rank(DTD), a more useful algorithm follows
input-output data given by
llk,1
Yk+L1 Yk+~a
Yk,m Yk+l,m Yk+l,m Yk+z,m
Uk+l,1
Uk'r Ul~+l,r
YI+N,1
Yk+N,m
Mk+N,1
Uk+N,r
Yea
Yk+L1
I lym(k)lux(k) I lu,(k)
= kVl(k) y (k + l) Equation (42) shows that the dependence relations among the vectors of (4 3) are the same, taking into account also the input, as those among the vectors (3 2), as indicated by Remark 2 This property allows the determination of vx, , v,~ by selecting the vectors (4 3) according to the same selection plan considered in Section 3, the vectors of (4 3) will therefore be selected m the following order Yl(k)
Ym(k),
Ul(k)
yx(k+ 1)
]
(4 3)
Algorlthm for the structural tdenttficatwn Let L,(yj) = [y,(k)
yj(k + 1)
yj(k + t - 1)],
(4 5)
L,(uj) = [us(k)
u,(k + 1)
u;(k + t - 1)]
(4 6)
Then matrix (4 3) talong 3x vectors m the first submatnx, 82 vectors in the second, etc can be written as
R(8. 82, ,8.,+~)
um(k),
Ym(k+ l),
1
(4 4)
L,.+,(ur)} (4 7)
= (L,,(yO Define n o w the product R z R as
A vector is retained If and only if it is independent from previously selected ones, when a dependent vector y~(k+vs) ts found all the remaining vectors belonging to the same submatnx will also be dependent so that their test is unnecessary The selection ends when a dependent vector has been found m every output submatrlx, the numbers of vectors selected from these submatnces will be
S(~l,
, 8re÷r) = RT(81,
, ~m+r)R(~I,
, 8,n-r)
(4 8) S(81, ,3,,,+r) is therefore a square mamx whose dimension is given by 3t+ + 3m+r Construct then the sequence of increasing-dimension matrices S(2,1,
,1),
S(2,2,
,I), S(2,2,
, ,2),
(4 9)
Canonical structures in the ldenttficatlon of muitwanable systems and select from (4 9) nonsmgular ones When a singular matrix is found, one of the indexes vt is determined, the procedure ends, all remaining matrices are singular, when all m indexes are determaned.
Remark 9 If two adjacent matrices in sequence (4 9) are considered, all the elements of the first are present in the subsequent one that can thus be obtained by computing only a hmited number of terms Remark 10 Let S(pa,/~2, ,~m+r) be a singular matrix in (4.9) and let/~, be the index increased by one with respect to the previous (nonsmgular) matrix in the sequence. Then v, = / ~ , - 1 while the indexes v. are given by v. =/z~ (j = 1, ,m) (j:# t), these indexes estabhsh the number of nonzero elements in the s u b m a m c e s / q r 5 PARAMETRIC IDENTIFICATION The parametric ldenUficataon will be performed on the input-output descnptmn (3 18); the parameters to be ldentafied are therefore the coefficients of the polynomials p,j(z) and q~j(z). The computaUon of the parameters is earned out through m independent steps; in each step the parameters of the polynomials belongmg to one row of P(z) and to the corresponding row of Q(z) are deternuned. Consider the selection of the vectors in matrix (4.3); when a dependent vector ys(k + vs) is found, it is hnked to previously accepted vectors by equataons (4.2) or (4 1) so that it is possible to determine the coeflficients of pa(z)(l = 1, .. , m) and of qs~(z) 0 = 1, .,r) For this purpose define the vector of parameters
,=(a81.1
I
+,,).,)(5 I) and, for slmphclty of notataon, let
s8 =
,
,
(s 3)
By using the full length of the avmlable data, Y8 can be obtained by means of the least-squares estimator Ys ~- (R8r Ra) -1 Rs T Ys(k + va) = Ss -1Rs T ys(k + v,).
(5.4)
Remark 11 Note that, dunng the structural ldentaficatlon, vs is determined by recogmmng the singularity of the
367
matrix S(v,1, , v , + 1, , v ~ ) , tlus matrix contains both S 8 and the vector R,Tys(k+vs) The additional computations regard the inversion of S s and its product with the vector R, T ys(k + vs)
Remark 12 The maramal number of parameters that are identified in a single step is gwen by n + rvM This is therefore also the maximal order of the matrix to be inverted Efficwncy of the whole identtficatwn procedure It must be noted that the previous structural ~denttficataon of the system remarkably improves the algorithm efficiency since a well-defined model among many possible structures is selected. Procedures based on the a posterwri optamahty test of an assumed order]structure model reqmre the parametric identificataon of a large number of models also ff a process order Is assumed. If, for instance, a tenth order process with three outputs is considered, the possible canomcal structures are (I, 1,8), 0 , 8 , 0 , (8, 1, 1), (2, 1,7), (1,2,7), (7, 1,2), (7,2, 1), (1,7,2), (2,7, 1), (2,2,6), (2,6,2), (6,2,2), (3,2,5), (2,3,5), (5,2,3), (5,3,2), (2,5,3), (3,5,2), (3,3,4), (3,4,3), (4,3,3), (4,4,2), (2,4,4), (4,2,4), so that the a posteriori optmaality test on the system structure would reqmre the parametric estimation of 24 different models and all this assun~ng that the system order is 10. Further reasons for the numencal efficiency of the proposed algorithm are due to the use, at a generic step of the structural ldenufication, of a/l the data used at previous steps, no computations are lost, as indicated in Remark 9, and to the use, in the parametric ldentaficataon, of the same matrices constructed during the structural identification as suggested by Remark 11 The order of magnitude of the required identaficatton tame can be deduced from [18] where a test case concerning a power plant with five inputs and three outputs is solved, with the identified structure (3,3, 3), in 10 sec of computang time on the CDC 6600. The reduced computer storage necessary during the different steps is a direct consequence of the previous observataons and of the separate ldentaficatlon of the different subsystems, as lmphed by Remark 12. With reference to the work of Ackermann it can be noted that by following the procedure outhned in [5] a maxmuzataon of the lntnnsic delay of the model, vM, would be achieved, instead of a nunmuzataon as in tlus paper, and, possibly, the records of some outputs would not be used in the identification procedure. A maxamizataon of vM leads to a greater number of parameters in the matrix Q(z) and therefore to a greater storage
368
ROBERTO GUIDORZI
requirement as per Remark 12, this shortcoming, however, is shght. More important is the usefulheSS of a model with mlmmal vM in the synthesis of controllers and observers for the identified process On the other hand, the selection procedure outlined by Ackermann is advantageous when one output of the system is remarkably more free from noise than the remaining ones since the information from the noise.free sensor can be best utahzed
6 STRUCTURAL AND PARAMETRIC I D E N T I F I C A T I O N I N THE P R E S E N C E O F NOISE
In this section the algorithms previously considered are extended to the noisy case It will be assumed here that the input and output sequences are corrupted by an additive uncorrelated noise with zero mean; the noisy components of the input and output vectors will be denoted with Ykj* = Ykd+d(Yko),
(6 1)
tiled* ~- Ulco+ d(Uko )
(6 2)
In (6 1), (6 2) and m the following the starred quantataes are noisy ones. In tlus case it is well known [24, 25] that the least-squares formula (5 4) leads to a biased estimate of ~,, also for N-~ oo. The origin of the bias in (5.4) can, however, be easily detected and eliminated if the statistics of the noise are known. In fact, from the expressions (5.3) and (48) of S~, since the noise has zero mean and is not correlated with the input--output sequences of the system, it follows that [14] p h m l S s * ffi hm 1 S,+N(d~), N--~ zv
(63)
N--,~ .N
where N(ds) is the covarmnce matrix of the noise vector
I d(uk,,) d(uk+,,_l.r)}. It stands also [14] 1 , phm ~Rs T N--~o zv
,
Ys (k + vs)
1 T = h m - - Rs y~(k + v~) + n 8, N-,o: .N
(6 4) where n,= hm~ N-,Qo
IIN
Z [d(Y,.3d(Y,+v..s)],
[)sk
,
estlmatlon of the quantities S s and RsTys(k+ vs) by means of the expressions g, = S s * - NN(ds),
(6.6)
RsT~s(k+vs) = RsT*ys*(k+v~)-ns.
(6 7)
The 'compensated' least-squares estamator ?~ = ~ - x RsT ys(k + vs)
(6 8)
gives therefore a consistent estimate, "~s, of the parameters vector Ys The identification of the system structure will, slmdarly, be performed on the sequence of the matrices ~(pq, ,/*m+r) It can be noted that, when the same amount of uncorrelated zero-mean noise is added to the input-output sequences and nonoverlappmg sets of data are used, then N(ds) = a 21 where oa is the variance of the noise so that 1
,
1
.
phm_--.S~ = hm ~ S ~ - N o ~ I , N--,QON
1
N-ooo
(6 9)
1
phm--~Rs T* ys*(k + v~) = hm ~ Rs T ys(k + v~), .V-,av z v
N ~ o o -'fi/
(6 10) and the estimator (6 8) is given by ¢G = ( S s * - N o ~ I ) - l R s T * Y s * ( k +vs)
(6 11)
It is also important to note that, in this case [14, 26], an estimate of N(r~ is g~ven by the least eigenvalue of the symmetrical matrix S(tz,/z, ,/0 for/z > vM so that a consistent estimate of the structure and of the parameters of the system can be obtained startmg from a given upper bound to the order of the system. The same technique can be apphed when different amounts of noise are present on the various inputs and outputs by perfornung a previous scahng on the data; the ratio of the different noises must however be known In [14] it has been shown that estimator (6 11) is equivalent to the eigenvector estmaator proposed by Levm [27] and to the instrumental variable estamator proposed by Wong and Polak [28], the compensated estimator (6 8) gives therefore a maximum hkehhood estimate The covanance matrix of the esUmated coefficients cannot be easily obtained in the general case, when, however, the amount of adchUve noise is low it can be approximated [27] with cov(~s) = (I +'~T@s)(S~*--No~I)-1o~
(6.12)
N
Z [dO,).... ,), d(y,+,,,-la)]t
}
Z, [d(u,+~_x.~)d0,)+~..~)]
[ ..
7 F R O M I N P U T - O U T P U T TO STATE-SPACE MODELS
(6 5)
If the statistics of the noise are known it is therefore possible, for N large enough, to obtain a consistent
After the structural and parametric identification of the system carried out according to the algorithms described in Sections 4, 5 and 6, an input-output descnptaon of the type (3 18) will be
Canomcal structures m the ~dentafication of multavanable systems obtained. It is then necessary to obtain, with a m~mmal effort, a canomcal state-space representatmn. A first obwous choice regards the construetmn of a state-space model of the type (3.5)-(3 11) that can be obtmned as described in Remarks 4 and 6, ~ts mmal state can be computed by means of (3.13). This procedure has been outhned m [7] and [14], a new more efficient algonthm wdl be developed here. Tim new algorithm reqmres the development of a second class of canomcal statespace models obtmnable from P(z) and Q(z) To avoid an exeesswe length of the present paper some results are only stated here; their proof can, however, be found m [15] and [29] Let us consider the d~fference input--output descnptmn for multwanable systems g~ven by
0 ~,1 ] rvtt~
I 0 ~, =
(7.8)
)
Iv,-1 O~$gtip
(7.9)
[,,TI[-,: : -] ..
0
0
I~ =
:
A(z) y(k) = Q(z) u(k),
(7 1)
DnT
y(k) = Ky(k),
(7 2)
0
where A(z) and Q(z) are (mxm) and (mxr) matrices whose enmes Ao(z) and q , ( z ) are polynommls in z satasfymg to the following condmons:
369
(~=
..
0 0
..
:
flnl "'" 0
(7 I0)
)
fin,
0
1
. . . . . .
0
c~x
0
0
Cml . .
0
1
0
O] ..
0
Cm~ . . . . . .
1
(1) The polynommls X,(z) are momc. vI
(2) The relaUons deg {)q/(z)} > deg {)q~(z)} for i ~ j, deg{>~(z)}>deg{k~/(z)}
for ~~ j ,
deg {)~(z)} > deg {qo(z)) are vahd. K ~s a real tnangular (m × m) matrix, g~ven by 1
c~a
In Lemma 3 the word 'eqmvalent' means that the input-output model (7.I) (7.2) and the state-space model (7.7) (7.1 I) have the same external behavior, Le. thelr output is identlcal for every input sequence apphed wRh null mmal condmons The hnk between nonzero miUal condRions m (7.1)(7.2) and the mgaal state m (7.7) (7.11) is gwen by the relataon
= l~V(z)K-ly(0)-~lCZ(z)n(0),
(7 3) I
Cml
1
Write also - oq,,~z - oq,a,
)Lt.l( Z ) = __ Ot~j,v # ZvO--1 _ _ . . . __ Ot~],2 Z - - Oti.%1,
(7.12)
where V(z) and Z(z) are glven by (3.14) and (3.16) respectwely, l~l is gaven by expression (3.25) with the obwous substttution of the parameters ~ . k wlth %.k and W Is given by
1 Cm,ra_ 1
~,~(z) = z ~,- ~,,,~,z ~-I - .
tl
~(0) = ~iV(z) y(0)- ~ Z ( z ) ,(0)
1
K __ cm_ u
vI + v2
(7 4)
t)~r
...
b~?
o
..
0
...
...
0
(7.5)
q,i (z) -- fl(,,+ +~,)dz"-I +" "
-~"~(Pl+ "l')',-l+g)drZdrfl{,/.++,,-i+1)~) n = v~+v2+
+vm,
(7.6)
0
vM = m a x t ( v , ) .
It Is then possible to state the following lemma
• •
~n T
Lemma 3
An equivalent state-space descnptaon of the system (7.1)(7.2) is glven by the triple (~,B,C"~) where "~={'~o} ( i , j = l , . , m ) ,
(7.7)
0
(# x ~ , ~ ) .
(7.
!3)
370
ROBERTO GUIDORZI
The proof of Lemma 3 can be obteaned, starting from the canonical representaUon (7 7) (7 11), by means of considerations completely s~mdar to those presented in Sectaon 3, these results can, however, be found also in the work of Wolovlch [29] and in the complete version of [15] From the comparison of the Identified pair P(z), Q(z) with the pair A(z), Q(z) of equatlon (7 1) it follows that the only property of P(z) that is not shared by A(z) IS the condition deg{;~,~(z)} > deg A.(z)} f o r j < t since in P(z) the degree of pu(z) can, at most, equal the degree of p.~z) when j < t. By denoting with K -1 the real (mx m) triangular matnx of the type (7 3) whose tth row is given by the coefficients of the terms z ~, In the tth row of P(z) it is, however, possible to write the identified model in the form P(z) KK -x y(k) = Q(z) u(k),
(7 14)
i e. to write a model of the type (7 1) (7 2) by talong A(z) = P(z)K and y ( k ) = K-Xy(k). A canonical state-space model of the type (7.7) (7.11) can therefore be deduced from the identified pear P(z), Q(z) according to the following algorithm.
Algorithm (1) Write directly the matrix i~ from Q(z), where the entries of i~ are the coefficients of the polynormals in Q(z) (2) Write directly the trmngular mamx K -x from P(z), where the enmes of the ~th row of K -~ are the coefficlents of z ~, in the zth row of P(z). (3) Compute K and A(z) = P(z)K. (4) Write directly from K and A(z) the matrices and .~ usmg the fact that the nonzero columns of ~ are the columns of K and the slgmficant columns of .~ are gwen by the coefficients of the polynormals m A(z). The mltlal state of the system can then be computed, if desired, by means of the relatmn (7 12).
of operatmns actually requtred by the new algonthm depends on the number of polynomials p~(z) on the left side of the mmn diagonal of P(z) whose degree is equal to the degree of p~,(z), when the conditmn deg{p~(z)}>deg{p.(z)} for j < l is met by P(z) the new algorithm does not reqmre any computation, since K - I, and permits one to write directly the mple (,~,I~,~) from P(z), Q(z) A further feature of the new algorithm ~s that it does not reqmre the lnversmn of the (n × n) matrix 1¢¢1 but only the Inversion of an (m × m) triangular matrix and Its product for an (m ×m) polynomial matrix and can thus prove advantageous when n is large.
A numerwal example The input-output sequences of a system w~th one input and two outputs are gwen by u(0)=l,
u(1)=2, u(4) = - 5, u(7) = 50, u(lO) = 30,
u(3)=5,
u(5) = - 12, u(8) = 10,
u(6) = 15, u(9) = - 6 0 ,
u(11) = O,
u(12) = 0
y(0) = (0,0) T, y(l) = (0,0) T, y(2) = (1,0) T, y(3) = (0, 2) T, y(4) = (2, - 1)T, y(5) = (5, 4) T, y(6) = ( - 3, - 9) T, y(7) = (5, 10)T, y(8) = ( - 25, 8)T, y(9) = (56, - 7 ) T, y(10) = ( - 3 , 10)T, y(ll) = ( - 2 7 , - 4 )
T, y(12)= (40,5) T.
Step 1 Structural tdenttficatwn Takang N = 9, the sequence of matrices S(2, 1, 1),
S(2, 2, 1),
S(3, 2, 2), s(3, 3, 2),
ts constructed. The first singular m a m x found ts S(3,3,2) so that v2 = 3 - 1 = 2. Th~s m a m x is deleted from the sequence that continues as S(3, 2, 3),
Remark 13 The derived state-space model is, because of its structure, completely observable Its complete reachab~hty, controllablhty if the dynamical matrix ~s nonslngular, is however not assured since non reachable, or non controllable, initial states could be nonzero; in these cases the obt~uned model can be reduced by means of standard techmques [30]. This consideration ~s not relevant when the amount of avadable input--output data is large with respect to the system order The new algorithm presented here is more effiment than the more easdy obtainable algorithm based on Remarks 4 and 6 Moreover, the number
n(2)=4,
S(4, 2, 3),
The second singular matrix found xs S(4, 2, 3) so that v l = v 3 f = 3 Since m = 2 the structural ldenUfiCatlon is terminated.
Step 2. Parametric tdenuficatwn The input-output equations that can consequently be wntten are ( - a~l,~z 2- a2x,2z - a~l,x) Y~a + (z2 - a~ z - a~2a)Yk,z = (5~,xz + & 0 u~ 1, (2~ -- Oll,3 Z 2 - -
an,z z - alia) Yk,1-1-( --
a12,f1,7 -- a12,1)
YI,,a
= ~3,1 :" + ~,1 ~ <31,0 ,,,.1
Canonical structures m the identification of multwanable systems The associated parameters, obtained by means of the estimator (5, 4), are gwen by
Y2 = (a21,1a~l,~a21~l~,a a~.~]/~*.xfls,a)T
0-111 -211 0)T,
=(1 Y1 ----" (all,1 a11.ea11,aIa12,1al~
131,13~,x/3~)w --0 0-11-1 I I 1 0 0 :
and consequently P(z)=
, ze - 1
Q(z)=
z2+ 2 z - 1
. 1
Step 3 Constructwn of a state-space model
[1]
Following the new algorithm introduced in this paper the matrix ~ can be lmmedmtely written from Q(z) and is given by
371
given by Gaussian independent sequences of length N = 2400 with zero mean and variance ¢r~ = 0 1 In absence of additive noise the output variances %12 and cr~z~showed the same input values so that when a noise wxth variance ~rn2 is added to the input and to the output, it is possible to say that a 100(cr,Jcr) per cent noise level is present The test on the smgulanty of the matrices ~Ozl ,/~m+0 can be performed by reporting their determinant or the least e~genvalue In absence of no~se the results of the structural ~dentlfication are given in F~g 1 where, clearly, a structure given by v1 = 5, v2 = 6 is deterrmned For a 10 per cent addmve noise the results of the structural identification performed on the noisy matrices S* and on the compensated ones ~ are shown in Figs 2 and 3 The same tests with an additive no~se level of 20 per cent are illustrated m F~gs 4 and 5.
0
~=
0 1 IO s.
0 m
Matrix K -x can be directly written from P(z) and ~s gwen by
1]
c.
I0"
6
m
NoiseleSS data
10-( --
It follows immediately that K=
[1] --1
tO-S m
1
tO-t2
A = P(z)K = [ zz+z~+z-2 -2z
o:[:o 0
.~=
ko
o]
1
o
-1
0
-z+l ze+2z - 1
T 2
FIG
1
I
I
4
Structural
J
T,
I
6
,
8
identification
of a
simulated
eleventh order system m the absence of nmse
1
0
-1
0
1
1
-I
0
0
0
0
0
1
0
2
1
-2
.
The mitaal state, computed by means of (7 12), is given by ~(0) = (1 0 0 0 0) T.
|0
m
I
E
I0--
8 PRACTICAL APPLICATIONS OF THE DESCRIBED IDENTIFICATION PROCEDURE In this section two apphcatlons of the proposed identlficataon procedure to a simulated and to a real process are reported.
I
2
i
I
4
t
I
6
z
I
8
Stmulated process An eleventh order system with two inputs and two outputs has been s]mulated. The input w a s
F,G 2 Structural ]dentlficat,on of a simulated eleventh order system m the presence of a 10 por cent level of additive noise U n c o m p e n s a t e d data
372
ROBERTO Gurooazi I0z
io z
lo-
No)se = 10 %
io-=
P
2
6
4
FIG. 5. Structural zdenuficatlon of a s~mulated eleventh order system m the presence of a 20 per cent level of addmve nmse. Compensated data
-
8
FIG 3. Structural ldentflicatmn o f a simulated eleventh order system m the presence o f a 10 per cent level o f addmve norse. Compensated data
~" zc
to-&
8
I1
Io
No,se = 20 %
10 = -
I0
20 Noise,
4
6
%
FIG. 6 Comparison between the reconstrucUon errors on the first system output for least squares models (L S ) and compensated least square models (CLS)
, 2
30
8
P
FIG. 4 Structural identification of a simulated eleventh order system m the presence of a 20 per cent level of addmve noise. Uncompensated data.
The Identified parameters wlth the least squares esUmator (5 4) and with the compensated esUmator (6.11) (C.L S ) have been used for the reconstructaon of the output sequences of the system, the reconstrucUon errors for noise levels up to 30 per cent are m&cated in Ftgs. 6 and 7 for the first and second output respectively. All these expenments have been performed by using overlapping sets of data. Real process A multistage multi-output chemical reactor with agitator and three output flows of 60 l/hr has been considered. The input--output data for the
,< g
2C
§
I
1o
W
2o
I
3o
Norse, %
FIG 7 Comparison between the reconstruction errors on the second system output for least squares models (L S ) and compensated least squares models (CLS)
Canon,cal structures in the identification of multlvanable systems identtficatton have been obtained by mjectmg a tracing dye into the input flow and by measuring the dye level on the outputs wtth a color meter every 30 sec The structure determination test is shown in R g 8; also if a clear indicatton has not been obtained for the delay chain relative to the first output, a structure w~th vx = 4, v2 = 4 and va = 3 has been assumed and an eleventh order model considered A comparison between the dynamical inputoutput behavior of the Identified model and of the actual system is depicted m the Figs. 9, 10 and 11, where the pulse responses of the system and of the model are compared
373
)0(
80 >~
Model output
--
so
O
40
20
l 4
z
I e
I a
I lo
I ,z
I ~
I f" ) ~6 Ja zo
.
Mm FIG 10 Companson between the identified model output and the measured data Second output
i0oI "~\ \
•FJrst subsystem
\ \ \
10-2
o<
o
>,) _>
se=.d ,ub,y.,em
0 10-3 --
°t
80
~
-
~
Model output .
60 4O
2 l
)
2
)
)
4
i
I
i''%'
6
,
6
1
I
f
)
I
l
I
i
I
2 4 5 0 I0 IZ 14 i6 18 20 ~'9 Mm
ii
I_
-
FIG 11 Comparison between the ldenufied model output and the measured data. Third output
F FIG 8. Structural ldentlficatmn for a three-output
chemical reactor 9 CONCLUDING REMARKS
In this paper a umtary approach to the structural and parametric identdicatlon of hnear multwanable systems has been presented; this approach rehes heavily on the structural properties of blockc o m p a m o n state-space representations for mulUvariable systems. The applications of the procedure to both simulated and real processes that can be found m [14, 16-18] show good results for what concerns both the obtained models and the
IO(
80 -- Model output 8
406020
,
w '21 ) 3
I
computational efforts required The whole procedure can be easily extended, in an entlrcly obvlous way, to systems where also an algebrmcal link between the input and the output Is present.
~ d a t o
i)
4.
.5
6
7
8
9
I0 I[
_
Mm FIG 9. Comparison between the 1dent]fled model output and the measured data First output
REFERENCES [1] L. A ZAOEH From clrcuR theory to system theory. Proc IRE 50, 856-865 (1962) [2] K J. ASTROMand P EYKHOFF System identlficatmn, a survey. Proc. IFAC Symp on ldentaficanon and Process Parameter Esttmanon, Prague (1970)
374
ROBERTO GUIDORZl
[3] B L Ho and R E KALMAN Effective construction of linear state-variable models from input-output functions Regelungstechmk 14, 545-548 (1966) [4] A J TETHER Construction of minimal hnear statevariable models from finite input-output data IEEE Trans Aut Control AC-15, 427-436 (1970) [5] J E ACKERMANN Die mmimale Em-AusgangsBeschreibung von mehrgrossen Systemen und ihre Bestimmung aus Em-Ausgangs-Messungen Regelungstechmk 19, 203-206 (1971) [6] J E ACKERMANNand R S BucY Canomcal minimal realization of a matrix of impulse response sequences Inform Contr 19, 224-231 (1971) [7] C BONIVt:NTOand R GtnDoazl Canonical inputoutput descnpnon of linear multivariable systems Rtcerche dz Automattca 2, 72-83 (1971) [8] D Q MAYNE ParametnzaUon and ldentlficanon of hnear multlvanable systems Stabthty oJ Stochastic Dynamical Systems 294, Lecture Notes in MathemaUcs Springer Verlag (1972) [9] B GOPtNATR On the identification of hnear tlmemvanant systems from input-output data Bell Syst Tech J 48, 1101-1113 (1969) [10] M A BODtN Minimal realization of discrete lmear systems from input-output observations IEEE Trans Aut Control AC-16, 395--401 (1971) [11] P BRUNOVSKY A classlficaUon of linear controllable systems Kybernettka 3, 173-187 (1970) [12] V M PoPOV Invanant description of linear timemvanant controllable systems S I A M J Control 10, 252-264 (1972) [13] G MARRO et al Studio dl metodi numeno per Ihdent~ficazione dei s~stem~ dmamtc~ Tech Rept Centro Calcoh e Servomeccamsm~, Umvers~ty of Bologna (1969) [14] C BOMWNTOand R GUIDORZl Parametric Identification of Linear Multwarlable System~ Preprmts of the JACC, St Louts (1971) [15] R Gomot~.l and G MARRO Partmi identlficanon of large-scale systems Proc Allerton Conf on Circuit and System Theory, Urbana (1972) [16]A. DAL BIANCO, P GALLI and R GUIDORZt Identification of chemical reactors from inputoutput data Proc IEEE Conf on Dectswn and Control, Miami (1971)
[17] C BONIVENTOand R GUIDORZl Apphcatlon of an identification method to distillation columns P~oc Princeton Con] on Info~matron Sciences and Systems (1972) [18] R GUIDORZIand R Rossl Identification of a power plant from normal operating records Automatic Control Theory and Apphcatwns 2, 63-67 (1974) [19] E TSE and H L WEtNERT Structure Determination and Parameter ldentlficatwn for Multwartable Stochastic Linear Systems Preprmts of the JACC, Columbus (1973) [20] J C CHOW On the Structural ldentlficatwn and Parameter Esttmatwn m Linear Multwartable Dynamic Systems Preprmts of the JACC, Columbus (1973) [21] D G LUENBERGER Canonical forms for linear multlvanable systems IEEE Trans Aut Control AC-12, 290-293 (1967) [22] S H WANG and E J DAVISON Canomcalforms of hnear multwartable systems Control Systems Rep No 7203, Dept of Electrical Engineering, Umverslty of Toronto [23] J RtSSANEN Baszs of mvanants and canomcal forms for linear dynamic systems Automatwa 10, 175-182 (1974) [24]R G BAENDERand F W SMITH A note on the statistical estimation of a hyperplane IEEE Trans Aut Control AC-13, 591 (1968) [25] T KOOPMANS Identification Problems m Economw Model-Construcnon J Wiley, New York (1935) [26] C M WOODSIDE Estimation of the order of hnear systems Proc 1FAC Syrup on IdenttJicatton and Process Parameter Estimation, Prague (1970) [27] M J LEVIN Estimation of a system pulse transfer function In the presence of noise IEEE Trans Aut Control AC-9, 229-235 (1964) [28] K Y WONG and E POLAK Identification of hnear d~screte systems using the instrumental variable method IEEE Trans Aut Control AC-12, 707-718 (1967) [29] W A WOLOWCH The determination of state-space representations for hnear mult~varmble systems Automatwa 9, 97-106 (1973) [30] D Q MAYNE Computanonal procedure for the minimal realization of transfer function matrices Proc IEE 115, 1363-1368 (1968)
Canonical Structures in the Identification of Multivariable Systems* Structures Canoniques dans l'Identification de Systbmes Multivariables Kanonische Strukturen bei der Identifikation mehrvariabler Systeme KaHom~,~ecKHe
cTpyKTypbI
B I ~ e H T I , IqbHKaRttI,I MHOrOCB~I3aHttI, IX C t t e T e M
ROBERTO
GUIDORZIt
The idennficatwn of multwanable systems can m no way be constdered as the simple esttmatwn of a set of parameters, the determmatwn o f a suttable model structure ts an equally fundamental step. Summarymln this paper a umtary ldentlficatmn procedure for linear muluvarmble systems is described The approach considered is based on a prehmmary canonical structure identification, i • on the d~rect determination, from inputoutput sequences, of a set of mvarmnt indexes completely describing the input-output structure of the system Reduced computation ume and storage, canonical representation of the identified model and easy estimation of the mmal state of the system are in this way achzeved The results obtained from the application of this procedure to a simulated process and to a chemical reactor are also reported
all the m o d e m available control theory regards these representatlons. (2) The class of the input slgnals should not be restricted to such partlcular inputs as impulse functions, step functions, colored or whlte nolse, smusoldal signals,P.R.B.S. etc.,since these inputs can seldom be apphed to real processes. It is on the contrary desirable to use normal operating records of the plant. T h e cholce of a state-space, z.e parametrlc, model imphes the prehnunary determmataon of the system order, Itis m fact well known [2] that parametric models can glve large errors when the order of the model does not agree wlth the order of the process. Another practlcallyrelevant reqmrement concerns the mmlmahty, or at leasta strong lim]tatlon, of the number of the slgnificantparameters of the model and this means the choice of a statable canomcal, or quas]-canomcal, form. Moreover, the s~mphclty of the hnk between the selectedcanonlcal state-space representatmn and the set of inputoutput equations used for the identification,the mmal data are input--output sequences, beawly condmons the reqmred amount of computatmn. The prewous problems have been completely solved for singleinput single-outputor multl-mput single-output systems [2]. Troubles appear when multJvarmblc, Lc multl-mput multi-output, systems are consldered smcc a umtary approach satlsfymg prcwous reqmrcments cannot be found m the lltcraturc. The purpose of thlspaper is to bridge thlsgap by exploring a class of state-space canonical models hnked by partlcularly slmple relatmns to inputoutput d~fference dcscnptmns that can be dlrectly identified from the input--output sequences. The resulting]dent]ficatmnprocedure ]s slmflarto those
1 INTRODUCTION ACCORDING to the classical d e f i n m o n gwen b y Z a d e h [I], ~dentfficatmn ~s the deterrmnaUon, on the bas~s o f i n p u t a n d o u t p u t , o f a system within a specified class o f systems, t o which the system u n d e r test ~s e q m v a l e n t Th~s d e f i n m o n ~s q m t e general a n d allows m a n y degrees o f f r e e d o m m the pracUcal f o r m u l a t m n o f the ~denUficatmn p r o b l e m since the class o f the m o d e l s to be ~dentlfied a n d the class o f i n p u t signals are n o t specified a n d also the m e a n i n g o f ' e q m v a l e n t ' ~s often referred to a s t a t a b l e o p t t m a h t y criterion F o r the pracUcal usefulness o f a n ~dentfficaUon p r o c e d u r e the following c o n s i d e r a t i o n s s h o u l d be t a k e n into account" (1) T h e o b w o u s choice for the class o f the m o d e l s ~s gwen by state-space r e p r e s e n t a t m n s since a l m o s t * Received 11 June 1973, revised 11 September 1974. revised 15 January 1975. The original version of this paper was presented at the 3rd IFAC Symposmm on Identification and System Parameter Esumatmn which was held m The Hague/Delft, The Netherlands, during June 1973 The published proceedings of this IFAC meeting may be ordered from Elsevier/North Holland P O B 211, Amsterdam, Netherlands It was recommended for publication m revised form by associate editor J Ackermann 1"Istltuto di Automatica, Umvers~ty of Bologna, Vlale Risorglmento 2, 40136 Bologna, Italy. 26
361
362
ROBERTO GUIDORZI
for single-input single-output systems and can be summarized in the following steps (1) Structural ldentaficataon, i.e. determmaUon by direct treatment of the input--output sequences of a set of indexes completely describing the structure, not only the order, of a canonical model for the process (2) Parametric ldentafication of the subsystems, one for every output component, of an inputoutput descnptaon of the system with the prewously identified structure (3) Deduction from the identified input--output representation of the state-space model
Corapartson wah premous hterature A fundamental contribution in flus field can be found m the paper by Ho and Kalman [3]; in this work a state-space model is deduced from an infinite sequence of Markov parameters Interestmg results have also been given by Tether [4] who studies the equivalence between finite sequences of Markov parameters and state-space models and by Ackermann [5, 6] who gives a canonical reahzation procedure by the use of an input--output descnptaon. A closely related approach has been described by Mayne [8]. Gopinath [9] assumes an initial knowledge of the system general inputoutput sequences; by direct elaboration of the data, a set of input-output equations is then obtmned. The procedure ~s, however, based on the knowledge of a statable data selector matrix correspondmg to a well-defined structure of the equaUons and an efficient method for the construction of this matrix is not g~ven. The approach by Gopmath is also the basis of a paper by Budin [10], who improves the computational reqtarements and points out a general procedure leading to a minimal but not canomcal reahzation. The basle difference between the approach presented in this paper and the previous ones lies in the structural identaficatlon of the system, i.e. m the direct deternunation from the input--output sequences of a set of indexes like those described in a dual context in [11, 12] which permit a decomposaaon of the system into subsystems that can be separately identified. This approach rehes heawly on the structural properties of blockcompamon state-space representations for multivariable systems and on their equivalence with input-output canonical polynormal representations Some early results on these propertaes had been obtained m the first half of 1969 [13] and have been developed in [7, 14, 15]; apphcatlons to both simulated and real processes can be found in [14, 16-18] The structural identification approach has been considered also by Tse and Wemert [19] and Chow [20].
2 ORGANIZATION OF THE PAPER
The class of linear discrete multlvariable systems considered can be represented by the usual statespace equations x(k + 1) = Fx(k) + Gu(k), y(k) = Hx(k),
(2 1a) (2 l b)
where x belongs to the n-dimensional state-space gF, u to the r-dimensional input space qz' and y to the m-dimensional output space ~ The lmtial data are quite general synchronous input-output sequences, i e normal operating records of y and u The paper is orgamzed as follows. In Section 3 some propertaes of canonical (blockcompamon) state-space representataons are recalled and the structural and parametric equivalence between this class of models and input-output difference descriptions of the type P(z) y(k) = Q(z) u(k),
(2 2)
where P(z) and Q(z) are polynomial matrices m z, where z -x is the unitary delay operator, is estabhshed In Section 4 the structural ldentlficaUon of equations (2. I) is considered, the subsequent parametric ldentafication of an input-output model of type (2.2) follows m Section 5. In Section 6 the algorithms given m Section 4 and in Section 5 are extended to the case of data corrupted by additave noise with zero mean and a maximum hkehhood estimator is described In Section 7 the derivation of a class of canomcal state-space models from the identified input-output description is discussed and m Section 8 some practical apphcations of the proposed identification procedure are reported Some concluding remarks are given in Section 9 3 CANONICAL INPUT-OUTPUT STRUCTURES OF MULTIVARIABLE SYSTEMS
Some brief recalls about state-space canomcal forms are first given because of their use m the following For a proof of the gwen results see Luenberger [21], complete and useful papers on this argument have also been written by Wang and Davison [22] and Rlssanen [23] where the meaning of the term 'canonical' ~s deeply investigated Consider the discrete, completely observable, multivanable system represented by equation (2 1) and let H =
h~T
(3 1)
hmT
Construct then the vector sequences 111,F T h1F w~hi,
,
(3 2)
hm, F T hm,
,
Canomcal structures m the ~dent~ficatmn of multwanable systems and select them m the following order h i, h~,
, h,,,, F T h 1,
, F r h,,, l zT' h i,
Remark 2 The dependent vectors g r ' ! h are a hnear combmatmn, w~th coefficients gwen by the tth slgmficant row of A, of prewously selected ones, l e
Fr~'h, = a~xahx+
The original system has been decomposed m this way into m interconnected subsystems, a direct consequence of the structure of the matrices A and C is gwen by the following lemma. Lemma 2
(3 4)
The) th subsystem m the canonacal representatmn (A,B,C) is completely observable from the j t h component of the output An input-output difference equataon of the type (2 2) wdl now be deduced from the triple (A, B, C) Remarkable features of this descnptaon are that ~t exhibits the structure of the system and that also the parameters of the triple (A, B, C) can be easily obtained from it, and vice versa. Let
the system equations became
w(k+ I) ffiAw(k) +nu(k),
(3 5)
y(k) = Cw(k),
(3 6)
and the matrices A, B and C exhibit the following structures X = T F r -~ = {A.} (3 7)
[0 ] [o o]
(t,j = 1,2,. ,m) where the (v~x v~) mamces /q, show the companion structure (3.8) and the (v~ x v~) matrices A,~ are of the type (3 9)
w(k) = (wk.l w~.~ w~,,F, y(k) = (Yk.l
o
,
a.. 1
u(k) = (Uk,~Uk,~.. Uk.,)v,
(3 8)
I e. for a simpler notation wg,1 = wx(k), Wk,~ = w~(k) and so on Consider then the j t h subsystem, because of the complete observablhty of ~ts state from t b e j t h output component and of the canomcal structure of the matrices A and C it is easy to
a~,v,
0
0 a.,~j
C = H T -1 =
0
.,
' (3.9)
0]
0
0
0
1
0
0 0
,oo
1
0
..
t 1
Y~,~)r
and
Ij,r l
A.jffi
+a~mahm+a~l,~F'rhl+
Ihm,. ,
F'r("-l) h.,},
a~,l
b~r
The number, v., of nonzero coefficients m the generic matrix A~j, because of the order followed m the selectmn of the vectors (3 2) and of the consequent structure of the matrix T, ~s at most equal to v.+ 1 f f j < l , to v~ f f j > z [7] [14]
The set of integers v1, . , v,n completely defines the structure of the couple (F, H) and Is mvarmnt with respect to a change of coordinates The proof directly follows from [11,12] by means of obwous duahty considerations As is well known, by operating the change of coordinates leading to the new state vector w = Tx, where
A,~ =
bnl
Remarl~ 1
Lemma 1
,Irr(~-l)hl[
bit (3 11)
[ bn T
retaining a vector FT'lh m (3 3) ~f and only ff ~t ~s independent from all previously selected ones Note that when F T ' ~ ~s tested all the vectors Fq'Jh, (0 ~ j ~
T T = {hvFThl,
bll
B = TG = I bit
(3 3)
363
(v1+ 1)
(vx +
+ vrn-1 + 1)
0
(3.10)
364
ROBERTO GUIDORZI
derive from equatmns (3 5) and (3 6) the following expressions
The substltutmn of (3.13) in (3 5) leads to the Input-output descnptaon
W/~Jvi+ +~t-x+l}= Yka,
{(zI-- A) V(z)} y(k) - {(zI - A) WZ(z) + B} u(k)
Wk,¢v~+ +~j_x+2~----"zy~.~-b~,+ +vj_,+l~u(k),
(3 17)
wk:,,+ +~,_,+a~ = : Y k ~ - b ~ , + +,i-~+2~u(k) -
Wk,{vx+ +i,j)=
bTx+ +,j-x+x~zu(k)
ZVj-lykj-- b~x+ a' -- b{vl+
+vj-l) u(k)-
+v.~-t+l) Z "-2
u(k).
.
(312)
It must be remembered that here and m the followmg z -t is the umtary delay operator; the notatmn z~ u(k) xs therefore eqmvalent to u(k + 2) and so on. By wntang equatmn (3 12) for j = l, ,m the whole state vector can be expressed as w(k) = V(z) y ( k ) - WZ(z) u(k),
In representatmn (3 17), however, only the vlth, (vt+va)th, ,nth equatmns are slgmficant; the remaining ones are simple ~denttt~es The s~gmficant equations m (3 17) can therefore be written in the form P(z) y(k) = Q(z) u(k), (3 18) where
[ pn(z)
,
pm (z)
where
qu(z)
qlr(z) ]
qmt(z)
qmr(z)
.
-au,~z-a~ta,
p~j(z) = - a2j,vo zV,J- 1 -
Z px-1
0
0
0
(3 14)
- a~,2 z - aoa
1
q~(z) =fl(.x+ +,,)oz~'-l+
(3.22)
+fl{.x+ +~,-,+2~,Jz
+fl(~,+ +~,_,+iI,j, ..
(3 21)
The polynomials of Q(z) are obtained by computing the right slde of expression (3.17), by simple passages it follows
z
0
(3 20)
The polynomials of P(z) can be Immediately obtained from (3.17), in fact it follows that p,~(z) = z~,-ai~,v,z ~ , - 1 -
v(z) =
(3 19)
pm, (z)
Q(z) =
(3 13)
0
px,(z) ]
P(z) =
(3 23)
z "--t
0
bT hT
0
0
W ffi
[n × r ( v u -
1)]
(3 15)
0 b,~,..+x
Z(z) ffi
ix] zI
z ~t-11
VM -- max~ (v0.
0
bnT-~.+l
(3 16)
0
0
where the coefficmnts ilia are the entries of the matrix 1] = M B =
(3 24)
Canomcal structures in the ldentificatmn of multwanable systems
365
and the matnx M is gwen by (see also [7] ) - -
a11.2
-- all,a
- -
--
a11,8
- -
all,v a
1
-- alto,2
all,,I
--
alln~pl
- - a11,~ 1
1
-- all.?t,vlm
aim~
m
0 0
°
1
0 (3.25)
M = -- amt,2
...
- - a m x , v,. x
- - arax,a
0
-- atom,2
0
-- amLp. x
-- amra, a
0
-- amm,v ."
0
1
Remark 3 It is useful to observe that (see also Remark 1) the polynomials m P(z) and Q(z) satisfy the following relations" deg {pa(z)} > deg {p~(z)} for j > i, deg {p,(z)} I>deg {p,j(z)} for j < t, deg {p,(z)} > deg (p:~(z)) for t # j , deg {p,(z)} > deg {q~(z)}. The eqmvalence between the class of canomcal state-space representatmns (3.5)--(3.6) and the input--output descnptaon (3 18) has thus been estabhshed.
Remark 6 To obtmn the matrix B from the knowledge of P(z), Q(z) it is first necessary to construct, by direct inspection, the mamces M and I]; matrix B Is then glven by B -- M -1 I] Note that M is always nonsmgular because of its structure, m fact det (M) --- 1. Example Let
Remark 4 The set of indexes vx,.., V,n can be lmmedmtely deduced by inspection, m&fferenfly from the knowledge of A or of P(z); also from the parametric standpoint A and P(z) are eqmvalent. Matrix C can be chrectly written ff va, , vm are known
A=
Remark 5 To obtain the matrix Q(z) from the knowledge of the couple (A, 13) it ~s first necessary to construct the matrix M, which can be written by direct lnspectaon of A, the elements of !] = MB are then the coefficients of the polynommls q¢t(z).
B=
B=
-- amra,v.
[Ol0i ] -
'
3
ffi
flr~
1
1
03
-
2
0
13
1
2
"
1 ,
0
1
I
1
-
0 0 0 1 0
Then vl = v~ ffi 3, v2 = 2; the coetticlents of the polynonuals q.(z) are the elements of the matrix l] gwen by
J 0
--2--10
z~+2z+ I
01
,
C--
0 1
21
i
1
1
1
0
1 --
-
1
-12
1
- 1
The equivalent input-output descriptaon is therefore -z~-2z+l
-
-1
1 0 0 0 0
0
i
- 1
-
-1
-1
][33110]o [ 8:1 8 . fl~, 831 tim 8,, f t .
02
y(k) ffi
-z-
1
z+ 4 1 u(k). z+2
]
0
1
366
ROBERT{) GUIDORZI
Remark 7
4 STRUCTURAL IDENTIFICATION
Definition 1
The integer N in (4 3) must be large enough
The structural identification of a multivanable system is the detenmnaUon of the set of integers % , v~ defining the structure of the couple (A, C) from input--output sequences without the intermedmte construction of a parametric model The structural ldentffication will be performed by taking advantage of the fact, given by Remark 4, that this mformataon concerning a canonical statespace representauon is also evident In the inputoutput difference model (3.18) Consider the sth equation in the set (2 18), i e m
P
~ps.,(z)yk,, = Y, qs,,(z)uk,,,
(4 1)
~,=I
~,~1
that, on the basis of the relations (3 21), (3.22) and (3.23), can also be written In the form 111 /,'m~ ~=1 j - I
+ ~ ~fl(,,+ +~._x+,,.,u(~+~_x),, (4.2) ~=1 1=1
vss = vs).
(In (4 2)
C o n s i d e r n o w the matrix o f
[N>n+rv M, where vM=max,(v,) ] in order to peru'at the selection of the necessary number of independent vectors, in other words, the inputoutput sequences must be of sufficient length to permit the complete structural identification of the system. Usually the number of available data ~s many times the system order so that the previous condmon is largely fulfilled
Remark 8 For the structural, and subsequent parametric, identification the input sequence must satisfy to well-known condmons [2], i e must 'excite all the modes' of the system. Also this requirement can be easily met when the length of the data is large with respect to the system order It is certainly not advisable to carry out the structural ldentificataon directly on the vectors of (4.3) because of the large amount of storage necessary, moreover, the required storage would be a function of N Since for every matrix, D, rank(D) = rank(DTD), a more useful algorithm follows
input-output data given by
llk,1
Yk+L1 Yk+~a
Yk,m Yk+l,m Yk+l,m Yk+z,m
Uk+l,1
Uk'r Ul~+l,r
YI+N,1
Yk+N,m
Mk+N,1
Uk+N,r
Yea
Yk+L1
I lym(k)lux(k) I lu,(k)
= kVl(k) y (k + l) Equation (42) shows that the dependence relations among the vectors of (4 3) are the same, taking into account also the input, as those among the vectors (3 2), as indicated by Remark 2 This property allows the determination of vx, , v,~ by selecting the vectors (4 3) according to the same selection plan considered in Section 3, the vectors of (4 3) will therefore be selected m the following order Yl(k)
Ym(k),
Ul(k)
yx(k+ 1)
]
(4 3)
Algorlthm for the structural tdenttficatwn Let L,(yj) = [y,(k)
yj(k + 1)
yj(k + t - 1)],
(4 5)
L,(uj) = [us(k)
u,(k + 1)
u;(k + t - 1)]
(4 6)
Then matrix (4 3) talong 3x vectors m the first submatnx, 82 vectors in the second, etc can be written as
R(8. 82, ,8.,+~)
um(k),
Ym(k+ l),
1
(4 4)
L,.+,(ur)} (4 7)
= (L,,(yO Define n o w the product R z R as
A vector is retained If and only if it is independent from previously selected ones, when a dependent vector y~(k+vs) ts found all the remaining vectors belonging to the same submatnx will also be dependent so that their test is unnecessary The selection ends when a dependent vector has been found m every output submatrlx, the numbers of vectors selected from these submatnces will be
S(~l,
, 8re÷r) = RT(81,
, ~m+r)R(~I,
, 8,n-r)
(4 8) S(81, ,3,,,+r) is therefore a square mamx whose dimension is given by 3t+ + 3m+r Construct then the sequence of increasing-dimension matrices S(2,1,
,1),
S(2,2,
,I), S(2,2,
, ,2),
(4 9)
Canonical structures in the ldenttficatlon of muitwanable systems and select from (4 9) nonsmgular ones When a singular matrix is found, one of the indexes vt is determined, the procedure ends, all remaining matrices are singular, when all m indexes are determaned.
Remark 9 If two adjacent matrices in sequence (4 9) are considered, all the elements of the first are present in the subsequent one that can thus be obtained by computing only a hmited number of terms Remark 10 Let S(pa,/~2, ,~m+r) be a singular matrix in (4.9) and let/~, be the index increased by one with respect to the previous (nonsmgular) matrix in the sequence. Then v, = / ~ , - 1 while the indexes v. are given by v. =/z~ (j = 1, ,m) (j:# t), these indexes estabhsh the number of nonzero elements in the s u b m a m c e s / q r 5 PARAMETRIC IDENTIFICATION The parametric ldenUficataon will be performed on the input-output descnptmn (3 18); the parameters to be ldentafied are therefore the coefficients of the polynomials p,j(z) and q~j(z). The computaUon of the parameters is earned out through m independent steps; in each step the parameters of the polynomials belongmg to one row of P(z) and to the corresponding row of Q(z) are deternuned. Consider the selection of the vectors in matrix (4.3); when a dependent vector ys(k + vs) is found, it is hnked to previously accepted vectors by equataons (4.2) or (4 1) so that it is possible to determine the coeflficients of pa(z)(l = 1, .. , m) and of qs~(z) 0 = 1, .,r) For this purpose define the vector of parameters
,=(a81.1
I
+,,).,)(5 I) and, for slmphclty of notataon, let
s8 =
,
,
(s 3)
By using the full length of the avmlable data, Y8 can be obtained by means of the least-squares estimator Ys ~- (R8r Ra) -1 Rs T Ys(k + va) = Ss -1Rs T ys(k + v,).
(5.4)
Remark 11 Note that, dunng the structural ldentaficatlon, vs is determined by recogmmng the singularity of the
367
matrix S(v,1, , v , + 1, , v ~ ) , tlus matrix contains both S 8 and the vector R,Tys(k+vs) The additional computations regard the inversion of S s and its product with the vector R, T ys(k + vs)
Remark 12 The maramal number of parameters that are identified in a single step is gwen by n + rvM This is therefore also the maximal order of the matrix to be inverted Efficwncy of the whole identtficatwn procedure It must be noted that the previous structural ~denttficataon of the system remarkably improves the algorithm efficiency since a well-defined model among many possible structures is selected. Procedures based on the a posterwri optamahty test of an assumed order]structure model reqmre the parametric identificataon of a large number of models also ff a process order Is assumed. If, for instance, a tenth order process with three outputs is considered, the possible canomcal structures are (I, 1,8), 0 , 8 , 0 , (8, 1, 1), (2, 1,7), (1,2,7), (7, 1,2), (7,2, 1), (1,7,2), (2,7, 1), (2,2,6), (2,6,2), (6,2,2), (3,2,5), (2,3,5), (5,2,3), (5,3,2), (2,5,3), (3,5,2), (3,3,4), (3,4,3), (4,3,3), (4,4,2), (2,4,4), (4,2,4), so that the a posteriori optmaality test on the system structure would reqmre the parametric estimation of 24 different models and all this assun~ng that the system order is 10. Further reasons for the numencal efficiency of the proposed algorithm are due to the use, at a generic step of the structural ldenufication, of a/l the data used at previous steps, no computations are lost, as indicated in Remark 9, and to the use, in the parametric ldentaficataon, of the same matrices constructed during the structural identification as suggested by Remark 11 The order of magnitude of the required identaficatton tame can be deduced from [18] where a test case concerning a power plant with five inputs and three outputs is solved, with the identified structure (3,3, 3), in 10 sec of computang time on the CDC 6600. The reduced computer storage necessary during the different steps is a direct consequence of the previous observataons and of the separate ldentaficatlon of the different subsystems, as lmphed by Remark 12. With reference to the work of Ackermann it can be noted that by following the procedure outhned in [5] a maxmuzataon of the lntnnsic delay of the model, vM, would be achieved, instead of a nunmuzataon as in tlus paper, and, possibly, the records of some outputs would not be used in the identification procedure. A maxamizataon of vM leads to a greater number of parameters in the matrix Q(z) and therefore to a greater storage
368
ROBERTO GUIDORZI
requirement as per Remark 12, this shortcoming, however, is shght. More important is the usefulheSS of a model with mlmmal vM in the synthesis of controllers and observers for the identified process On the other hand, the selection procedure outlined by Ackermann is advantageous when one output of the system is remarkably more free from noise than the remaining ones since the information from the noise.free sensor can be best utahzed
6 STRUCTURAL AND PARAMETRIC I D E N T I F I C A T I O N I N THE P R E S E N C E O F NOISE
In this section the algorithms previously considered are extended to the noisy case It will be assumed here that the input and output sequences are corrupted by an additive uncorrelated noise with zero mean; the noisy components of the input and output vectors will be denoted with Ykj* = Ykd+d(Yko),
(6 1)
tiled* ~- Ulco+ d(Uko )
(6 2)
In (6 1), (6 2) and m the following the starred quantataes are noisy ones. In tlus case it is well known [24, 25] that the least-squares formula (5 4) leads to a biased estimate of ~,, also for N-~ oo. The origin of the bias in (5.4) can, however, be easily detected and eliminated if the statistics of the noise are known. In fact, from the expressions (5.3) and (48) of S~, since the noise has zero mean and is not correlated with the input--output sequences of the system, it follows that [14] p h m l S s * ffi hm 1 S,+N(d~), N--~ zv
(63)
N--,~ .N
where N(ds) is the covarmnce matrix of the noise vector
I d(uk,,) d(uk+,,_l.r)}. It stands also [14] 1 , phm ~Rs T N--~o zv
,
Ys (k + vs)
1 T = h m - - Rs y~(k + v~) + n 8, N-,o: .N
(6 4) where n,= hm~ N-,Qo
IIN
Z [d(Y,.3d(Y,+v..s)],
[)sk
,
estlmatlon of the quantities S s and RsTys(k+ vs) by means of the expressions g, = S s * - NN(ds),
(6.6)
RsT~s(k+vs) = RsT*ys*(k+v~)-ns.
(6 7)
The 'compensated' least-squares estamator ?~ = ~ - x RsT ys(k + vs)
(6 8)
gives therefore a consistent estimate, "~s, of the parameters vector Ys The identification of the system structure will, slmdarly, be performed on the sequence of the matrices ~(pq, ,/*m+r) It can be noted that, when the same amount of uncorrelated zero-mean noise is added to the input-output sequences and nonoverlappmg sets of data are used, then N(ds) = a 21 where oa is the variance of the noise so that 1
,
1
.
phm_--.S~ = hm ~ S ~ - N o ~ I , N--,QON
1
N-ooo
(6 9)
1
phm--~Rs T* ys*(k + v~) = hm ~ Rs T ys(k + v~), .V-,av z v
N ~ o o -'fi/
(6 10) and the estimator (6 8) is given by ¢G = ( S s * - N o ~ I ) - l R s T * Y s * ( k +vs)
(6 11)
It is also important to note that, in this case [14, 26], an estimate of N(r~ is g~ven by the least eigenvalue of the symmetrical matrix S(tz,/z, ,/0 for/z > vM so that a consistent estimate of the structure and of the parameters of the system can be obtained startmg from a given upper bound to the order of the system. The same technique can be apphed when different amounts of noise are present on the various inputs and outputs by perfornung a previous scahng on the data; the ratio of the different noises must however be known In [14] it has been shown that estimator (6 11) is equivalent to the eigenvector estmaator proposed by Levm [27] and to the instrumental variable estamator proposed by Wong and Polak [28], the compensated estimator (6 8) gives therefore a maximum hkehhood estimate The covanance matrix of the esUmated coefficients cannot be easily obtained in the general case, when, however, the amount of adchUve noise is low it can be approximated [27] with cov(~s) = (I +'~T@s)(S~*--No~I)-1o~
(6.12)
N
Z [dO,).... ,), d(y,+,,,-la)]t
}
Z, [d(u,+~_x.~)d0,)+~..~)]
[ ..
7 F R O M I N P U T - O U T P U T TO STATE-SPACE MODELS
(6 5)
If the statistics of the noise are known it is therefore possible, for N large enough, to obtain a consistent
After the structural and parametric identification of the system carried out according to the algorithms described in Sections 4, 5 and 6, an input-output descnptaon of the type (3 18) will be
Canomcal structures m the ~dentafication of multavanable systems obtained. It is then necessary to obtain, with a m~mmal effort, a canomcal state-space representatmn. A first obwous choice regards the construetmn of a state-space model of the type (3.5)-(3 11) that can be obtmned as described in Remarks 4 and 6, ~ts mmal state can be computed by means of (3.13). This procedure has been outhned m [7] and [14], a new more efficient algonthm wdl be developed here. Tim new algorithm reqmres the development of a second class of canomcal statespace models obtmnable from P(z) and Q(z) To avoid an exeesswe length of the present paper some results are only stated here; their proof can, however, be found m [15] and [29] Let us consider the d~fference input--output descnptmn for multwanable systems g~ven by
0 ~,1 ] rvtt~
I 0 ~, =
(7.8)
)
Iv,-1 O~$gtip
(7.9)
[,,TI[-,: : -] ..
0
0
I~ =
:
A(z) y(k) = Q(z) u(k),
(7 1)
DnT
y(k) = Ky(k),
(7 2)
0
where A(z) and Q(z) are (mxm) and (mxr) matrices whose enmes Ao(z) and q , ( z ) are polynommls in z satasfymg to the following condmons:
369
(~=
..
0 0
..
:
flnl "'" 0
(7 I0)
)
fin,
0
1
. . . . . .
0
c~x
0
0
Cml . .
0
1
0
O] ..
0
Cm~ . . . . . .
1
(1) The polynommls X,(z) are momc. vI
(2) The relaUons deg {)q/(z)} > deg {)q~(z)} for i ~ j, deg{>~(z)}>deg{k~/(z)}
for ~~ j ,
deg {)~(z)} > deg {qo(z)) are vahd. K ~s a real tnangular (m × m) matrix, g~ven by 1
c~a
In Lemma 3 the word 'eqmvalent' means that the input-output model (7.I) (7.2) and the state-space model (7.7) (7.1 I) have the same external behavior, Le. thelr output is identlcal for every input sequence apphed wRh null mmal condmons The hnk between nonzero miUal condRions m (7.1)(7.2) and the mgaal state m (7.7) (7.11) is gwen by the relataon
= l~V(z)K-ly(0)-~lCZ(z)n(0),
(7 3) I
Cml
1
Write also - oq,,~z - oq,a,
)Lt.l( Z ) = __ Ot~j,v # ZvO--1 _ _ . . . __ Ot~],2 Z - - Oti.%1,
(7.12)
where V(z) and Z(z) are glven by (3.14) and (3.16) respectwely, l~l is gaven by expression (3.25) with the obwous substttution of the parameters ~ . k wlth %.k and W Is given by
1 Cm,ra_ 1
~,~(z) = z ~,- ~,,,~,z ~-I - .
tl
~(0) = ~iV(z) y(0)- ~ Z ( z ) ,(0)
1
K __ cm_ u
vI + v2
(7 4)
t)~r
...
b~?
o
..
0
...
...
0
(7.5)
q,i (z) -- fl(,,+ +~,)dz"-I +" "
-~"~(Pl+ "l')',-l+g)drZdrfl{,/.++,,-i+1)~) n = v~+v2+
+vm,
(7.6)
0
vM = m a x t ( v , ) .
It Is then possible to state the following lemma
• •
~n T
Lemma 3
An equivalent state-space descnptaon of the system (7.1)(7.2) is glven by the triple (~,B,C"~) where "~={'~o} ( i , j = l , . , m ) ,
(7.7)
0
(# x ~ , ~ ) .
(7.
!3)
370
ROBERTO GUIDORZI
The proof of Lemma 3 can be obteaned, starting from the canonical representaUon (7 7) (7 11), by means of considerations completely s~mdar to those presented in Sectaon 3, these results can, however, be found also in the work of Wolovlch [29] and in the complete version of [15] From the comparison of the Identified pair P(z), Q(z) with the pair A(z), Q(z) of equatlon (7 1) it follows that the only property of P(z) that is not shared by A(z) IS the condition deg{;~,~(z)} > deg A.(z)} f o r j < t since in P(z) the degree of pu(z) can, at most, equal the degree of p.~z) when j < t. By denoting with K -1 the real (mx m) triangular matnx of the type (7 3) whose tth row is given by the coefficients of the terms z ~, In the tth row of P(z) it is, however, possible to write the identified model in the form P(z) KK -x y(k) = Q(z) u(k),
(7 14)
i e. to write a model of the type (7 1) (7 2) by talong A(z) = P(z)K and y ( k ) = K-Xy(k). A canonical state-space model of the type (7.7) (7.11) can therefore be deduced from the identified pear P(z), Q(z) according to the following algorithm.
Algorithm (1) Write directly the matrix i~ from Q(z), where the entries of i~ are the coefficients of the polynormals in Q(z) (2) Write directly the trmngular mamx K -x from P(z), where the enmes of the ~th row of K -~ are the coefficlents of z ~, in the zth row of P(z). (3) Compute K and A(z) = P(z)K. (4) Write directly from K and A(z) the matrices and .~ usmg the fact that the nonzero columns of ~ are the columns of K and the slgmficant columns of .~ are gwen by the coefficients of the polynormals m A(z). The mltlal state of the system can then be computed, if desired, by means of the relatmn (7 12).
of operatmns actually requtred by the new algonthm depends on the number of polynomials p~(z) on the left side of the mmn diagonal of P(z) whose degree is equal to the degree of p~,(z), when the conditmn deg{p~(z)}>deg{p.(z)} for j < l is met by P(z) the new algorithm does not reqmre any computation, since K - I, and permits one to write directly the mple (,~,I~,~) from P(z), Q(z) A further feature of the new algorithm ~s that it does not reqmre the lnversmn of the (n × n) matrix 1¢¢1 but only the Inversion of an (m × m) triangular matrix and Its product for an (m ×m) polynomial matrix and can thus prove advantageous when n is large.
A numerwal example The input-output sequences of a system w~th one input and two outputs are gwen by u(0)=l,
u(1)=2, u(4) = - 5, u(7) = 50, u(lO) = 30,
u(3)=5,
u(5) = - 12, u(8) = 10,
u(6) = 15, u(9) = - 6 0 ,
u(11) = O,
u(12) = 0
y(0) = (0,0) T, y(l) = (0,0) T, y(2) = (1,0) T, y(3) = (0, 2) T, y(4) = (2, - 1)T, y(5) = (5, 4) T, y(6) = ( - 3, - 9) T, y(7) = (5, 10)T, y(8) = ( - 25, 8)T, y(9) = (56, - 7 ) T, y(10) = ( - 3 , 10)T, y(ll) = ( - 2 7 , - 4 )
T, y(12)= (40,5) T.
Step 1 Structural tdenttficatwn Takang N = 9, the sequence of matrices S(2, 1, 1),
S(2, 2, 1),
S(3, 2, 2), s(3, 3, 2),
ts constructed. The first singular m a m x found ts S(3,3,2) so that v2 = 3 - 1 = 2. Th~s m a m x is deleted from the sequence that continues as S(3, 2, 3),
Remark 13 The derived state-space model is, because of its structure, completely observable Its complete reachab~hty, controllablhty if the dynamical matrix ~s nonslngular, is however not assured since non reachable, or non controllable, initial states could be nonzero; in these cases the obt~uned model can be reduced by means of standard techmques [30]. This consideration ~s not relevant when the amount of avadable input--output data is large with respect to the system order The new algorithm presented here is more effiment than the more easdy obtainable algorithm based on Remarks 4 and 6 Moreover, the number
n(2)=4,
S(4, 2, 3),
The second singular matrix found xs S(4, 2, 3) so that v l = v 3 f = 3 Since m = 2 the structural ldenUfiCatlon is terminated.
Step 2. Parametric tdenuficatwn The input-output equations that can consequently be wntten are ( - a~l,~z 2- a2x,2z - a~l,x) Y~a + (z2 - a~ z - a~2a)Yk,z = (5~,xz + & 0 u~ 1, (2~ -- Oll,3 Z 2 - -
an,z z - alia) Yk,1-1-( --
a12,f1,7 -- a12,1)
YI,,a
= ~3,1 :" + ~,1 ~ <31,0 ,,,.1
Canonical structures m the identification of multwanable systems The associated parameters, obtained by means of the estimator (5, 4), are gwen by
Y2 = (a21,1a~l,~a21~l~,a a~.~]/~*.xfls,a)T
0-111 -211 0)T,
=(1 Y1 ----" (all,1 a11.ea11,aIa12,1al~
131,13~,x/3~)w --0 0-11-1 I I 1 0 0 :
and consequently P(z)=
, ze - 1
Q(z)=
z2+ 2 z - 1
. 1
Step 3 Constructwn of a state-space model
[1]
Following the new algorithm introduced in this paper the matrix ~ can be lmmedmtely written from Q(z) and is given by
371
given by Gaussian independent sequences of length N = 2400 with zero mean and variance ¢r~ = 0 1 In absence of additive noise the output variances %12 and cr~z~showed the same input values so that when a noise wxth variance ~rn2 is added to the input and to the output, it is possible to say that a 100(cr,Jcr) per cent noise level is present The test on the smgulanty of the matrices ~Ozl ,/~m+0 can be performed by reporting their determinant or the least e~genvalue In absence of no~se the results of the structural ~dentlfication are given in F~g 1 where, clearly, a structure given by v1 = 5, v2 = 6 is deterrmned For a 10 per cent addmve noise the results of the structural identification performed on the noisy matrices S* and on the compensated ones ~ are shown in Figs 2 and 3 The same tests with an additive no~se level of 20 per cent are illustrated m F~gs 4 and 5.
0
~=
0 1 IO s.
0 m
Matrix K -x can be directly written from P(z) and ~s gwen by
1]
c.
I0"
6
m
NoiseleSS data
10-( --
It follows immediately that K=
[1] --1
tO-S m
1
tO-t2
A = P(z)K = [ zz+z~+z-2 -2z
o:[:o 0
.~=
ko
o]
1
o
-1
0
-z+l ze+2z - 1
T 2
FIG
1
I
I
4
Structural
J
T,
I
6
,
8
identification
of a
simulated
eleventh order system m the absence of nmse
1
0
-1
0
1
1
-I
0
0
0
0
0
1
0
2
1
-2
.
The mitaal state, computed by means of (7 12), is given by ~(0) = (1 0 0 0 0) T.
|0
m
I
E
I0--
8 PRACTICAL APPLICATIONS OF THE DESCRIBED IDENTIFICATION PROCEDURE In this section two apphcatlons of the proposed identlficataon procedure to a simulated and to a real process are reported.
I
2
i
I
4
t
I
6
z
I
8
Stmulated process An eleventh order system with two inputs and two outputs has been s]mulated. The input w a s
F,G 2 Structural ]dentlficat,on of a simulated eleventh order system m the presence of a 10 por cent level of additive noise U n c o m p e n s a t e d data
372
ROBERTO Gurooazi I0z
io z
lo-
No)se = 10 %
io-=
P
2
6
4
FIG. 5. Structural zdenuficatlon of a s~mulated eleventh order system m the presence of a 20 per cent level of addmve nmse. Compensated data
-
8
FIG 3. Structural ldentflicatmn o f a simulated eleventh order system m the presence o f a 10 per cent level o f addmve norse. Compensated data
~" zc
to-&
8
I1
Io
No,se = 20 %
10 = -
I0
20 Noise,
4
6
%
FIG. 6 Comparison between the reconstrucUon errors on the first system output for least squares models (L S ) and compensated least square models (CLS)
, 2
30
8
P
FIG. 4 Structural identification of a simulated eleventh order system m the presence of a 20 per cent level of addmve noise. Uncompensated data.
The Identified parameters wlth the least squares esUmator (5 4) and with the compensated esUmator (6.11) (C.L S ) have been used for the reconstructaon of the output sequences of the system, the reconstrucUon errors for noise levels up to 30 per cent are m&cated in Ftgs. 6 and 7 for the first and second output respectively. All these expenments have been performed by using overlapping sets of data. Real process A multistage multi-output chemical reactor with agitator and three output flows of 60 l/hr has been considered. The input--output data for the
,< g
2C
§
I
1o
W
2o
I
3o
Norse, %
FIG 7 Comparison between the reconstruction errors on the second system output for least squares models (L S ) and compensated least squares models (CLS)
Canon,cal structures in the identification of multlvanable systems identtficatton have been obtained by mjectmg a tracing dye into the input flow and by measuring the dye level on the outputs wtth a color meter every 30 sec The structure determination test is shown in R g 8; also if a clear indicatton has not been obtained for the delay chain relative to the first output, a structure w~th vx = 4, v2 = 4 and va = 3 has been assumed and an eleventh order model considered A comparison between the dynamical inputoutput behavior of the Identified model and of the actual system is depicted m the Figs. 9, 10 and 11, where the pulse responses of the system and of the model are compared
373
)0(
80 >~
Model output
--
so
O
40
20
l 4
z
I e
I a
I lo
I ,z
I ~
I f" ) ~6 Ja zo
.
Mm FIG 10 Companson between the identified model output and the measured data Second output
i0oI "~\ \
•FJrst subsystem
\ \ \
10-2
o<
o
>,) _>
se=.d ,ub,y.,em
0 10-3 --
°t
80
~
-
~
Model output .
60 4O
2 l
)
2
)
)
4
i
I
i''%'
6
,
6
1
I
f
)
I
l
I
i
I
2 4 5 0 I0 IZ 14 i6 18 20 ~'9 Mm
ii
I_
-
FIG 11 Comparison between the ldenufied model output and the measured data. Third output
F FIG 8. Structural ldentlficatmn for a three-output
chemical reactor 9 CONCLUDING REMARKS
In this paper a umtary approach to the structural and parametric identdicatlon of hnear multwanable systems has been presented; this approach rehes heavily on the structural properties of blockc o m p a m o n state-space representations for mulUvariable systems. The applications of the procedure to both simulated and real processes that can be found m [14, 16-18] show good results for what concerns both the obtained models and the
IO(
80 -- Model output 8
406020
,
w '21 ) 3
I
computational efforts required The whole procedure can be easily extended, in an entlrcly obvlous way, to systems where also an algebrmcal link between the input and the output Is present.
~ d a t o
i)
4.
.5
6
7
8
9
I0 I[
_
Mm FIG 9. Comparison between the 1dent]fled model output and the measured data First output
REFERENCES [1] L. A ZAOEH From clrcuR theory to system theory. Proc IRE 50, 856-865 (1962) [2] K J. ASTROMand P EYKHOFF System identlficatmn, a survey. Proc. IFAC Symp on ldentaficanon and Process Parameter Esttmanon, Prague (1970)
374
ROBERTO GUIDORZl
[3] B L Ho and R E KALMAN Effective construction of linear state-variable models from input-output functions Regelungstechmk 14, 545-548 (1966) [4] A J TETHER Construction of minimal hnear statevariable models from finite input-output data IEEE Trans Aut Control AC-15, 427-436 (1970) [5] J E ACKERMANN Die mmimale Em-AusgangsBeschreibung von mehrgrossen Systemen und ihre Bestimmung aus Em-Ausgangs-Messungen Regelungstechmk 19, 203-206 (1971) [6] J E ACKERMANNand R S BucY Canomcal minimal realization of a matrix of impulse response sequences Inform Contr 19, 224-231 (1971) [7] C BONIVt:NTOand R GtnDoazl Canonical inputoutput descnpnon of linear multivariable systems Rtcerche dz Automattca 2, 72-83 (1971) [8] D Q MAYNE ParametnzaUon and ldentlficanon of hnear multlvanable systems Stabthty oJ Stochastic Dynamical Systems 294, Lecture Notes in MathemaUcs Springer Verlag (1972) [9] B GOPtNATR On the identification of hnear tlmemvanant systems from input-output data Bell Syst Tech J 48, 1101-1113 (1969) [10] M A BODtN Minimal realization of discrete lmear systems from input-output observations IEEE Trans Aut Control AC-16, 395--401 (1971) [11] P BRUNOVSKY A classlficaUon of linear controllable systems Kybernettka 3, 173-187 (1970) [12] V M PoPOV Invanant description of linear timemvanant controllable systems S I A M J Control 10, 252-264 (1972) [13] G MARRO et al Studio dl metodi numeno per Ihdent~ficazione dei s~stem~ dmamtc~ Tech Rept Centro Calcoh e Servomeccamsm~, Umvers~ty of Bologna (1969) [14] C BOMWNTOand R GUIDORZl Parametric Identification of Linear Multwarlable System~ Preprmts of the JACC, St Louts (1971) [15] R Gomot~.l and G MARRO Partmi identlficanon of large-scale systems Proc Allerton Conf on Circuit and System Theory, Urbana (1972) [16]A. DAL BIANCO, P GALLI and R GUIDORZt Identification of chemical reactors from inputoutput data Proc IEEE Conf on Dectswn and Control, Miami (1971)
[17] C BONIVENTOand R GUIDORZl Apphcatlon of an identification method to distillation columns P~oc Princeton Con] on Info~matron Sciences and Systems (1972) [18] R GUIDORZIand R Rossl Identification of a power plant from normal operating records Automatic Control Theory and Apphcatwns 2, 63-67 (1974) [19] E TSE and H L WEtNERT Structure Determination and Parameter ldentlficatwn for Multwartable Stochastic Linear Systems Preprmts of the JACC, Columbus (1973) [20] J C CHOW On the Structural ldentlficatwn and Parameter Esttmatwn m Linear Multwartable Dynamic Systems Preprmts of the JACC, Columbus (1973) [21] D G LUENBERGER Canonical forms for linear multlvanable systems IEEE Trans Aut Control AC-12, 290-293 (1967) [22] S H WANG and E J DAVISON Canomcalforms of hnear multwartable systems Control Systems Rep No 7203, Dept of Electrical Engineering, Umverslty of Toronto [23] J RtSSANEN Baszs of mvanants and canomcal forms for linear dynamic systems Automatwa 10, 175-182 (1974) [24]R G BAENDERand F W SMITH A note on the statistical estimation of a hyperplane IEEE Trans Aut Control AC-13, 591 (1968) [25] T KOOPMANS Identification Problems m Economw Model-Construcnon J Wiley, New York (1935) [26] C M WOODSIDE Estimation of the order of hnear systems Proc 1FAC Syrup on IdenttJicatton and Process Parameter Estimation, Prague (1970) [27] M J LEVIN Estimation of a system pulse transfer function In the presence of noise IEEE Trans Aut Control AC-9, 229-235 (1964) [28] K Y WONG and E POLAK Identification of hnear d~screte systems using the instrumental variable method IEEE Trans Aut Control AC-12, 707-718 (1967) [29] W A WOLOWCH The determination of state-space representations for hnear mult~varmble systems Automatwa 9, 97-106 (1973) [30] D Q MAYNE Computanonal procedure for the minimal realization of transfer function matrices Proc IEE 115, 1363-1368 (1968)