On the Identification of Subsystems of Decomposed Large Scale Systems

Copyright © IFAC Identification and Sy,,,'m Estimation 19H5, York, lIK, I'll':;

I'aramct~r

ON THE IDENTIFICATION OF SUBSYSTEMS OF DECOMPOSED LARGE SCALE SYSTEMS T. Fresewinkel and H. Unbehauen DeIHLr/nll'n/

or t;/ff/riml EnginefTing, Automatif Control Laboratory,

Ruhr-University /Joch1l./II, /)-4630 Bochum J, Federal Republic of Germany

Abstract. For an efficient identification of large scale systems it is advisable to decompose the entire system into separately treatable subsystems of lower complexity. In most applications these subsystems are of a certain structure. In this paper a method is proposed for identifying the parameters of these subsystems by decomposing them into interconnected transfer functions which are linear in their parameters. The identification algorithm operating recursively/iteratively is derived by partitioning and optimizing a cost function with respect to the parameter vectors of the transfer functions. The advantage of this method is that the identification procedure for the subsystems is feasible by applying linear filters to each transfer function without calculating a common denominator, such that the number of parameters to be identified, is minimal and the system structure is saved during the identification procedure. Keywords. Large scale systems, identification, decomposition, recursive functions, iterative methods, invariant inbedding, multivariable systems.

methods (Fresewinkel 1984a) results in the conclusion that for identification purposes it is very helpful first to study the system structure very intensively. For an efficient identification it is proposed by Fresewinkel (1984b) to introduce a special multistrata hierarchical decomposition of the system due to Mesarovic (1970). In this multilevel system representation, the subsystems of the second level of description have the advantage that they are identifiable one by one. These subsystems, in most practical applications, are of a certain structure (Fresewinkel 1984b).

INTRODUCTION A good system model is a prerequisite for nearly all applications in system and control theory. The problem of estimating such a model for large scale systems in particular, has not yet been solved satisfactorily. A close look at the existing identification algorithms for large scale systems suggests the presence of a common methodology based on the following four steps (Fresewinkel 1984a): step step step step

1: 2: 3: 4:

A schematic interaction flow diagram An appropriate decomposition Verification of the identification conditions Parameter identification

The algorithm proposed in this paper is applicable to these special structures; but it is not restricted to these subsystems alone. In the following the subsystems of the mentioned second level of description are referred to as the system. This system is itself an arbitrary combination of several subsystems, which should be linear in their parameters and whose interactions are partially unknown. The parameters of the subsystems are identified by a two level algorithm such that the number of parameters to be estimated, is minimal.

The choice of the algorithm in the fourth step is determined by the type of decomposition. In general two types of algorithms for identification purpose of complex systems exist. They are: (i)

sequential algorithms The aim of these algorithms is to find an optimal global model of the complex system. This requires that all interconnections between the subsystems are measurable and that the subsystems are ordered in an hierarchical manner without any feedback. (Bubnicki 1975, 1980, Hasiewicz 1978, 1979). In general these subsystems are still very complex and there is no efficient algorithm for parameter identification until now.

FORMULATION OF THE IDENTIFICATION PROBLEM Suppose that the system outputs are given by

(1) where

E(k) = E(k-l) + y(k-l)

(ii) parallel algorithms The decomposition into parallel order subsystems forms the basis for most of the decentralized state estimators in nonlinear state space descriptions (Arafeh 1974, Hassan 1977, El-Sherief 1982, Chemouil 1981). Necessary conditions for estimation and identification are that all the subsystem state vectors are locally observable. These algorithms in genera:l involve a lot of computations.

(2)

is the system para:meter vector which is allowed to vary stochastically. The noise is described by ~(k) and v(k) with mean values and covariances E{~(k)} ~

Q

E{£(k)£T(k-i)} ~{y(k)} =

E{y(k)yT(k-i)} = 6(k-i)

The investigation of the different identification

267

~(k)6(k-i)

Q

(3)

~(k)6(k-i)

is the Kronecker Delta Function.

268

T. Fresewinkel and H. Unhehauen

If the system is a combination of transfer functions or in general a composition of N linear subsystems, then each subsystem is characterized by its parameter vector ~i' As the interconnections between the subsystems are arbitrary, the parameter vectors are normally combined in a nonlinear way for calculating the system outputs according to eq. (I). By uSing the maximum a-posteriori (MAP) identification method, i t can be shown (Sage 1971) that the most probable parameter vector which causes y(k) k~k2kf to appear is the vector which minimizes the cost function

aH(k) ap. (k-I)

1. (k-l)

-1

-1

- y. (k)

(Bb)

-1

(Bc) aH (k) ay. (k)

pO(k) - p. (k-I) - v. (k-I)

0

-1

-1.

-1

T aE. (k)

N

~~ 0

1

~ E~ (k)E. (k)+y. (k)

I j;1

ar~1 (k)

(Bd)

-1.

apo(k) -J

-J

(Be)

-1

-1

jfi

Let us define the matrices

M. (k)

Minimizing eq. (4) with respect to p(k) in general results in a highly non linear probl~m. To circumvent this, the costfunction has to be decomposed into N subfunctions each describing the part of the whole estimation problem which is related to one subsystem. By introducing a vector po of new variables, so called coordinating variables, it is possible tg fix all parameters ~ in eq. (4) to some value ~ except the subsystem parameter vector ~i' By definition each subsystem is linear in its parameter vector p .. So minimization of the cost func1 tion, which is now a function of the interaction variables ~o, becomes a linear problem with respect to ~i:

1

-

2(~i (ko)-~i (k o »

k 1 + -2

I

f

k~k

+1

~i

(k o )

{ET(k)E~I(k)E. -1

-1

dET(k)

_ ---=2-

Because of eq.

-lJ

T

-~

p. (k-I)

1. (k)

-I 1. (k-I)+y. (k)+M. (k)E. (k)E. (k)

-1

M. (k)

-1

(5)

-1

-1

(6)

to the minimization problem of eq. (5) it can be shown (Schoeffler, 1971) that minimizing eq. (4) is equivalent to minimizing eq. (5) if eq. (6) is satisfied. From these considerations the identification problem is described by the Hamiltonian 1 T -I .IN(2I T ~i (k-I)~i~k-I)~i (k-I)~i (k)~i (k)~i (k)

1~1

T + 1 (k) (p. (k-I)+v. (k_I»)+yT(k) -1

-1

-1

-1

- v. (k-I»)

(p~(k)-p. -1

-1

-1

~i

0

(k ) (~i (ko ) -~i (k o ) ) o

(7c)

T

-1

d~i (k) -1 E. (k) E. (k) ~i (k-I)~i(k-I)lav (k-I) -1 -1 -i + 1. (k)-y. (k) (Ba) -1

(lOc)

p~(k) -1

~

P. (k)

(IOd)

-1

To obtain an optimal recursive solution, we first assume that the vector rO(k) in eq. (10) is known. Then y(k) is known, too. So by invariant inbedding we can derive the following recursive equations to solve (IOa,b) as shown in the appendix.

P.

(k)

-1

-I

p.

-1

~P.

-1

(k-I)+p. (k) [y. (k)+M. -1

-1

-1

(k)E~I(k)e. -1

-1

(klk-I)] (lla)

(k)

( llc)

The optimal trajectory of the parameter vector is then defined by the following derivatives of H(k):

o

-1 E. (k) E. (k) -J -J

with

(7b)

O.

dH(k) av. (k-I)

(lOb)

-1

(k-I)

with the transversality conditions -1

-1

(7a)

-1

1. (k )

I

j~1

-1.

(loa)

These equations together with the transversality conditions (7b, 7c) describe a nonlinear Two Point Boundary Value Problem (TPBVP). A pure recursive solution of a similar problem in a suboptimal way was proposed by Arafeh (1974).

If we add the constraint

~

-1

-1.

y. (k)

with

H(k)

- v. (k-I)1. (k-ll -1 -1

p. (k)

jfi

1 (k-l)V.- (k-l)v. (k-l)} -1

j;I, ... ,N

Thus the optimization problem results in the following system of equations:

N

(k)

"

is valid.

o + v

(Bd) the identity

M. (k) ~ MO. (k) ,

-1

-1

-1

(9b)

ap~(k) -1

-1

T -I

(9a)

a~i (k-I)

-1

-1

~

-1

If we iterate now to solve eq. (IOd) in each sample, we get an optimal recursive/iterative solution which can be realized by the two level algorithm shown in Fig. I. On the first level, which represents the subsystem level, all subsystem parameter vectors are estimated by N parallel filters described by eq. (11). In each sample and each iteration the new estimate for the parameter vectors ~ . (k) is supplied to the 2. level which represents tfie coordinator. On this level the coordinating variables p~(k) are updated for the next iteration step 1+1, by (Pearson, 1971) :

e:g.

(12 )

269

Identification of Large Scale Systems coordination level

T

!!1 U ll (s)

•

subsystem level ~

Two level identification algorithm

with these new variables, the vector y . (k) puted. These two vectors are then supplied the first level to calculate the parameter ~. (k) for the next iteration step. l. 1+1 So the flowchart of the algorithm shown in is valid.

I

is comback to vector

U

LN

(5)

L

T

+

!!L

Fig. 2 + ~

Special system structure in large scale systems The outputs of the subsystems i=l, ... ,L, which are transfer functions, can be calculated by -1 Bi (q ) -d i zi (k) -------1- q u (k), i=l, ... ,L (13) i Ai (q )

Initial ize : £j(k 0)' ~j(ko)

i

= 1 ...

N

-1 -1 where Bi (q ) and Ai (q ) are the following polynomials of order m and n respectively : i i m. l. -1 -j b q (14a) Bi (q ) ij j=O n. l. -j -1 a q (14b) 1 + Ai (q ) ij j=l

L

L

level 2

q-1 is the backward shift operator and d. is the time delay of subsystem i. The parameterl.vector of subsystem i is defined by level 1 (15)

The gain vectors ~i represent the subsystems j=L+1, .. ,N with outputs no

T

'ESj(k)~j

Zj (k)

(16)

with 'E: (k) j

= [u

T ~j = [h i1

i1

(k) ••• u

iNi

(k)

1

••• hiNil

Using eqs. (13,16) the system output in Fig. 3 is determined by the equation -1

L

~

Bi(q) -d. T ----1- q l.[u. (k)+m .(k)p.l l.0 -SJ -J i=l Ai (q )

L

Flowchart of the two level algorithm COMPUTATION FOR A SPECIAL SUBSYSTEM STRUCTURE

In Fresewinkel (19B4a) it has been found out that a lot of identifiable subsystems as parts of large scale systems, are of a certain structure. One of these structures is shown in Fig. 3. In this section a specific algorithm is derived from the general algorithm of the second section for application to the system of Fig. 3. The system consists of L transfer functions. The input of each transfer function is the sum of measurable signals of the process, which are amplified by an unknown constant gain vector. The system output is disturbed by some noise n(k) . So the system can be decomposed into N=2L subsystems which are linear in their parameters.

(17)

with j=i+L. Expanding eq. (17) by the product of the denominator polynomials, we get L

1

L

IT A (q- )YM(k) \J=1 \J

L

i=l

L

IT A (q \J=1 \J

-1

)u. (k)

l.

\Jfi with

u. (k)

(lB)

l.

Supposing that the system noise n(k) L

E(k) =

eq.

is given by

-1

IT A\J(q \J=1

(17) results in

)n (k),

(19)

270

T. Fresewinkel and H. Unbehauen L

-I

A (q

Il

E (k)

~

~=I

L

)y(k)-

I

-I

B. (q

i=1

)q

-d i

1

L

Il A (q ~ =I ~

~+i

-I

with )u. (k) 1

(20)

T o o m . (k)=[-!. (k-I) ... -2. (k-n.) u . (k-d .) ... u . (k-ml.-d .)] -z~ -~ -1 1 l l l l and

If the system parameters are constant or slowly time varying, we get the derivatives of eq. (20) according to eqs. 69,10), calculated with the coordinating vector ~ (k), as M. (k) -1

m -a .L

m. (k)

(k)-I

= [U

(k)

~i (k)

il

for i=l, .. . ,L,

(k) ••• u

iNi

(k)]

j=i+L.

REMARK

(21 )

f

-1

~:j

To accelerate the convergence of the algorithm, we multiply each vector Y (k) by an additional weighti ing factor p(k). Simulation results show that a good choice for p(k ) is about p(k)~I/N for k>50 and p( k)=o for 0
with

m . (k)

-al

EXAMPLE As an example, we consider a system according to Fig. 3 of two transfer functions with six inputs. The transfer functions are given by G I (s)

= (1+2s) (1+1.5s)

G (s)

=

2

~i(k)

(I+s) (1+2.5s)

and the gain vectors are hT -I T

-

and -d

°

-1 m.(k)=B (q)q -J i

L i

0

~2

u il (k)

-I

. .

Il A ( q )

r

~+i

(22)

•

~=I ~

u

iNi

-

[-0.9

-1.4]

Simulation conditions

output noise weighting factors p . (k) J

for i=l, ... ,L, j=i+L. The superscript '0' indicates that the coefficients of the polynomials and gain vectors are taken from the parameter vectors p. instead of p .. For the identification error of-the subsystems; we get L

o -I T IT A~(q )Y(k)-~i(k)ei(k-I)

ei(klk-I) =

0.9]

The input signals are uncorrelated discrete signals. The simulation conditions are listed in table I. TABLE I

(k)

-

[0.6

~=I

5% 0
E. J V (k) -J initial values p. (0) -J P. (0) -J sampling time T = 0.4s.

-]

0 10

4

,

j=l, ... ,4

The identification error is shown in Fig. 4.

~+i .03

(23a)

and e. (k\k-I) J

25

50

75

[s] 100 t

Fig. 4. Estimation error of the example o

-I

-B. (q)q

-d i L

1

0

IT A (q

~=I

~

-I

(k)

)u

The identified transfer functions are

10

~fi

1.02 (1+2.ls) (1+1.45s) (23b)

G (s) 2

for i=l, ... ,L,

o¥

An on-line computation of the subsystem outputs is possible by using the equation

1

1.04 (1+0.98s) (1+2.67s)

and the identified gain vectors are

j=i+L.

To calculate ev(k), the vector p (k-I) is substituted by pO(k), v=I, ... ,N. So the ~vposteriori errors e (k) all subsystems are identical. v

2. (k)

=

(2 4)

[0.609 [-0.902

0.903] -1.401]

The parameters converge very fast in the neighbourhood of the true values, as can be seen in Fig. 5. The number of iterations, which are necessary in each sample, decrease from 5 iterations, when

Identification of Large Scale Systems starting the algorithm, to 1 or 2 iterations at the end (Fig. 6).

271

with the split boundary conditions A (k )

-

oot: o03pt

0

7 I

c.

I 0'2

l"

-2.0

O~O~

~(ko)

I

b ll

o.oo:-=::

A, band c are known and the initial vector ~nd-the terminal vector

(A.4 )

I

I

(A.3)

b'2

are unknown . Assume that the problem is already solved and let us step from k=k to k=k +l . Then eq. (A.4) bef f comes

7 /

I 0 22

(.Q

21

-2.0

~3F

=

l

I

I

bn

b 2,

0.00

2

o ~(~,kf)

W: ::;=::

+

00 1.

Z

I

I

hll

h12

with the first partial difference

0.00

Om[ -2.0 0 ~

oa(b)

~

1

2

50

25

'

~

/22

7

75

(A.S)

T t::.c

o~okf

~(~+t::.~)-~(~)

t::.b

i

Eqs. (A.l) and (A.2), when substituted in (A.S), give [s]

100

Curves of the identified discrete parameters

B (A.6)

2·

~

As there is no analytic solution for the problem defined by eq. (A.6) we approximate a linear solution by

. - ... " .. . .... . .... _ ....... .. .. .

[s) 100 t 75 50 25 Number of iterations at each sample

(A.7) With eq.

(A.7), eq.

(A.6) becomes

CONCLUSIONS The identification algorithm proposed in this paper is a new method to determine the parameters of subsystems of a large scale system on-line. It has the advantage that we can identify the parameters of transfer functions within an arbitrary system structure without transforming the system into one of the classical multi-input multi-output representa-tions. Therefore it is possible to save the system structure during the identification procedure and to estimate certain signals which are outputs of the subsystems.

Because ;:>(k ) is unknown, we expand eq. f Taylor series about ~(kf) for small c :

3!(~(kf)-~(kf)~,~,kf)

(A . B) in a

I

ac-

c

c=O

-

ACKNOWLEDGEMENT

(A.9)

The authors are very grateful to the DFG (German Research Association) for a support under Un2S/20.

APPENDIX

If we substitute the equations (10a,b) of our original problem into eq. (A.9), we get two equations in the zero and first powers of ~: ~(kf+l)=~(kf)+~(ktl)(r(kf+l)+~(kf+l)~

The method of invariant inbedding is a very powerful tool for solving nonlinear TPBV problems (Sage 1971, 1974, Lee 196B). Suppose that there are in general two nonlinear equations ;:>(k)

!(;:>(k-l),~

~ (k )

2(;:>(k-l),~(k-l),k-l)

(k-l) ,k-l)

(A . l) ko+l~k~kf

(A.2)

~(kf+llkf)) ~

-1

(k + l ) = (~(kf)+~(kf)) f

-1

-1

(k +l) f (A. loa)

- - - - [~(kf+l) 32(k ) f (A . lOb)

T. Fresewinkel and H . Unbehauen

272

Substitution of k to k-1 in eq . (A.10) results in f eq. (11).

REFERENCES Arafeh, S.A., A.P. Sage ( 1974). Multi lev el discrete time system identification in large scale sys tems. Int. J. of Syst. Science,S, 753 - 791. Bubnicki, z. (1975). Global and l oca l identification of complex systems with cascade structure. Systems Science, 1, 55 - 65. Bubnicki, Z. ( 1980). Identification of control plants. Elsevier, Amsterdam. Chemouil, P., M.R. Katebi, D. Sastry, M. Singh (1981). Maximum a posteriori parameter estimation in large scale systems. 8th IFAC World Congress, Kyoto, XII 21 -30. El Sherief, H., M.A. Maud (1982). Microprocessor based parameter and state estimator and its application to load forecasting of interconnected power systems. 6th IFAC-Symposium on Identification and System Parameter Estimation , 1520-1525. Fresewinkel, T., H. Unbehauen (1984a) . The identification of complex systems - problems and solution approaches. Regelungstechnik 32 , 3-7. Fresewinkel, T. (1984b). Modelling of comp lex systems by sequential identification of subsystems . Regelungstechnik 32, 51-55. Hasiewicz, Z. (1978). Progressive identification of complex systems I. Statement of the problem. System Science, 4, 301-313. Hasiewicz, z. (1979). Progressive identification of complex systems 11. Validity of the method . System Science,S, 63-79. Hassan, M., M. Singh (1977). A two-level costate prediction algorithm for nonlinear systems. Automatica, 6 , 629-634. Lee, E.S. (1968). Quasilinearization ana invariant imbedding. Academic Press, New York. Mesarovic, M.D., D. Macko, Y. Takahara (1970). Theory of hierarchical multi level systems. Academic Press, New York. Pearson, J .D . (1971). Dynamic decomposition techniques. In D.A. Wismer: Optimization methods for large scale systems. McGraw Hill, 121-190. Sage, A.P., J . L. Melsa (1971). System identification. Academic Press, New York. Schoeffler, J.D. ( 1971). Static multilevel systems. In D.A. Wismer: Optimization methods for large scale systems. McGraw Hill, 1-46 .