Comparison of modal parameter estimation techniques on aircraft structural data

Mechanical Systems and Signal Processing (1990) 4(2), 157-172

COMPARISON TECHNIQUES

OF ON

MODAL

PARAMETER

AIRCRAFT

ESTIMATION

STRUCTURAL

DATA

J. E. COOPERt

Materials and Structures Department, Royal Aerospace Establishment, Farnborough, U.K. (Received 10 November 1988, accepted 27 September 1989) Five global time domain parameter estimation methods, least squares, double least squares, total least squares, correlation fit and Smith least squares were compared using the same set of impulse responses taken from the ground vibration test of an aircraft structure. The behaviour of the modal parameter estimates was investigated for changes in the solution order and number of data points analysed. The techniques produced comparable results for the dominant modes, although there were some differences in the estimates of the other modes. The correlation fit and Smith least squares approaches found equivalent estimates to the other methods while using much smaller computational model orders. 1. INTRODUCTION There have been a n u m b e r of global system identification methods developed in recent years for use with multiple response data produced from a number of input cases. Such techniques provide a single global estimate from many data sets. Comparisons between the global parameter estimation techniques are scarce. Usually, in such comparisons, the data sets used by each method are generated separately. Most popular time domain methods e.g., polyreference[1], Ibrahim time domain (ITD)[2] and the Eigensystem realisation algorithm (ERA)J3], can be shown to employ the least squares (LS) approach as the basic element in their formulations. It is well known that when noisy data is considered, overspecification of the computational model order is necessary in order to reduce the bias (statistical error) inherent in the LS estimates. Such an approach is c o m m o n in modal analysis. However, numerous techniques exist which attempt to produce unbiased estimates without having to resort to excessive overspecification. Statistical comparisons between the LS method and some of the alternative approaches have been performed using modal data [4, 5]. However, these investigations have so far been limited to the single input-single output (SISO) case using simulated data corrupted with Gaussian noise. Real data suffers not only from additive noise corruption but also from the fact that the mathematical model can never perfectly represent the real situation. In this paper, comparisons are made between the results produced by the polyreference method and versions of the polyreference technique where the LS element has been replaced by the double least squares (DLS), total least squares (TLS) and correlation fit (CF) approaches. The iterative Smith least squares (SLS) method is also included in the comparison. A multi input-multi output ( M I M O ) case is considered using test data from a Gnat aircraft which has been the subject of a previous study [6]. The primary aim is to investigate the behaviour of the parameter estimates obtained by the different techniques using the same data sets. Variations due to different data generation and collecting procedures are therefore eliminated. 157 tpresent address: Department of Aeronautical Engineering, University of Manchester, MI3 9PL UK.

158

J.E.

COOPER

2. MATHEMATICAL MODEL The methods considered in this paper all model the impulse response which can be generated either by; applying an impulse or step relaxation to the structure, by generating the impulse response function from the inverse Fourier transform of the frequency response function (FRF), or by using an averaging technique such as random decrement [7]. For a M degree of freedom (dof) system, the impulse response is modeled as a summation of M exponentially damped sinusoids, so that the output at the i-th measurement station for the j-th input case can be written as M

Yu( k At) = Y!ik = ~ e O°'g'~'[ R!i~ sin ( ¢oark A t ) q- S!i, COS ( wa,k At)] r:

(1)

I

or 2M

YUk = Y. C,j,/.Lf

(2)

r=l

where /Zr = e ( ~,0~+,,o, )a,

(3)

and for each mode, st , to and toe are critical damping, natural frequency and damped natural frequency respectively. R, S and C are amplitude terms. For the underdamped case considered here, C a n d / z occur in complex conjugate pairs. When the data from all NO measurement stations and NI input cases is considered, equation (2) becomes [8] Yt I/, • LYnolk

"'"

" "

YlNIk ] " = Yk = W D k V

(4)

yNONIk

where D = diag [/.t~,/x2,.. •, ~2m ] and Y is a (NO x NI) block data matrix consisting of all of the output values at the NO measurement stations at the kth time instant for each of the NI input cases. The use of a block matrix notation enables the SISO formulation to be expanded with ease to the M I M O case• W is a (NO x 2M) matrix containing the complex mode shapes and V is a (2M x NI) matrix containing the modal participation factors. The elements of V relate the complex amplitude terms of equation (2) so that if the first input case is taken as the reference, then V,; = C,,,/C,r,

i = 1. . . . , NO,

j = 1. . . . . NI,

r = 1. . . . . 2M.

(5)

It can be shown that equation (4) is the general solution of the block difference equation Yk+ Y k - , A , + ' "

"+ Yk-~A, = [0]¢NO×N,~

(6)

if the block polynomial equation D ' V + D" ' VA, +. . .+ DVA~ ,+VA~=[0]~2M×Nz~

(7)

is satisfied, where ( s x N 1 ) = 2 M and A~ are ( N I x N 1 ) block difference equation coefficients. Equation (7) can be rewritten into a canonical form so that

[D'-'V

D ' 2V

...

V]

A20. At

0

I.

"'"

0

.-.

=D[D~-'V

D~-2V

...V]

(8)

AIRCRAFT

STRUCTURAL

159

DATA

which can be seen to be the transposed eigensolution of a block sparse Hessenberg matrix. The roles of input and output can be reversed in the above formulation using the principle of reciprocity. 3. METHODS CONSIDERED 1N THE COMPARISON Only a short description of the methods will be given here. Full details can be found in the references. 3.1. POLYREFERENCE The polyreference method [ 1] is the most commonly used global time domain approach. A very lucid account of its derivation is given by Deblauwe et al. [8]. The technique consists of three stages: (i) the block difference equation (6) is extended columnwise so that +

tv~J

=

i

Y~.:,

Y,

lY;-,

Y;-2

...

liAl

I"2

.....

-A2

-S

or

[ Y] = ~ o

(9)

and the LS method is used to estimate the block A~ coefficients by minimising the sum of the squares of the errors on the [ Y] vector so that

__O= ( ~ r ~ ) - l ~ r [

Y].

(10)

(ii) the sparse block Hessenburg matrix in equation (8) is formed and the eigensolution found. The eigenvalues are the roots of the block polynomial in equation (7) and lead to the system frequencies and dampings via equation (3). The modal participation factors are found from the eigenvectors. (iii) the mode shapes are found by expanding equation (4) columnwise so that

IYol[ vTI'T L~)~j

LV~D~J

or

[fq~=[v]~w T

(11)

and this gives the solution to the mode shapes as W ~ = {[ V][ V]~}-'[ V][ ~']L

(12)

A common alternative approach to this third stage is to perform the calculation in the frequency domain [9]. It can be noted that the eigenvalues tzi of equation (8) appear in the calculation of all the modal parameters, thus any error in their calculation will result in erroneous estimates. The techniques described in sections 3.2-3.4 are used to replace the LS element of the first stage of the polyreference formulation in order to find better estimates of the Ai parameters and hence the system eigenvalues.

160

J.E. COOPER

3.2. DOUBL E LEAST SQUARES The DLS approach used here employs the same philosophy of the DLS technique applied by Ibrahim to the ITD [10] and Pappa and Juang to the ERA [1 l], but is applied to the LS technique [4]. The LS method described by equation (10) usually results in a positive error upon the damping estimates. However, an alternative estimate of the difference equation parameters 0__ = ( c ~ r O ) - ' q~ T[

y]

(13)

where

[Y,+, Y~ ... ~__ .Y,+2

Y,+,

...

YN-,

...

Y. ] Y:

(14)

Y~

generally gives a negative bias on the damping estimates. In the DLS method, the average of the two estimates given by equations (10) and (13) is found. It is assumed that the bias on the two separate solutions will cancel out. When applied to the LS approach, the DLS method can be formulated to use only slightly more computation than that of the LS technique. 3.3. TOTAL LEAST SQUARES In the TLS technique [12], the squares of the errors on both the [ Y] vector and the columns of the • matrix are minimised. The method can be formulated either by: (i) finding the NI eigenvectors corresponding to the smallest N | eigenvalues of the symmetric matrix [4]

OTO H0 =

[r]r~

OT[ Y] ] [y]T[y]j

(15)

so that

[[+'<'> +TIll lr-o]=[-o],,. v]To [Y]T[Y]JL

<16

(ii) by finding the N | right hand singular vectors corresponding to the smallest NI singular values from the singular value decomposition (SVD) of the matrix

H,=[OI[Y]]

(17)

as described by Zhang et al. [ 13]. The correspondence between the two approaches is easy to show as H~H,

= Ho.

(18)

Theoretically, when the correct model order is used on uncorrupted data, the smallest NI eigenvalues (or singular values) will be zero. However, when an overspecified model is used on real noisy data there will be a number of "zero" eigenvalues (singular values) that will have some non-zero value. The corresponding eigenvectors (right hand singular vectors) will each lead to estimates of the Ai parameters [4]. It is the SVD approach that has been used in this paper.

AIRCRAFT 3.4•

CORRELATION

STRUCTURAL

DATA

161

FIT

The CF method [4, 14] is based upon an analysis of the approximate data correlations used in the curve fit. Unlike the LS method, an arbitrary number of correlations can be used and those correlations that are prone to be corrupted by noise may be omitted from the curve fit. Defining the approximate autocorrelation for the j-th lag as 1

12j = N - j

N-j

k=l ~ r~Yk+j

(19)

then equation (6) can be rewritten as =0

O,+I2,_,A,+...+O,_~A~

(20)

and by extending equation (20) columnwise, hut assigning the variables NA to determine the number of correlations used, and IP to set the number of correlations omitted from the fit, then the CF equivalent of equation (9) is Oip+s

OIp+s-I

OIp+~-2

~"~lP+s+l

~"~lP+s

"(~IP+s-I

f~IP+s+NA-I

~IP+s+NA-2

f~IP+s+NA

...

IIAI]

f2/ff+I [

-A2

• "" ~I;+NAI

-As

or

_~= q,19

(21)

The CF estimate of the block A~ coefficients follow as 19 = [ ~ r ~ ] - , ~ r_r" 3.5.

SMITH

(22)

LEAST SQUARES

Unlike the other methods considered in this paper, the SLS technique [15] uses a two stage iterative procedure based upon equation (1), with a constant offset term R0o , included. Due to possible convergence problems, overspecification of the solution model order cannot be used. Generally, the technique has been neglected by the modal analysis community on the grounds that it is computationally expensive• In studies using simulated data sets [5] the technique was found to give extremely good parameter estimates• The extension to the MIMO case follows directly from the multiple response case discussed by Smith [15]. The first stage takes assumed frequency and damping values and minimises the squares of the differences between the measured and estimated output values. The minimisation leads to the set of simultaneous linear equations "Rj

N

E

GoGo

GoGt

•""

GoG2M

G1Go

GIGt

"""

GIG2~

Rijl

YijkGo

] (23)

R~jM

k=l

k=l

-~2MGo

G2MGI

"""

G2MG2M

S o,

"

LYi)kG2M.l

.SoM.

where Go = 1 with Gi = e -¢~°'k'a'

sin (toaik At)

i = 1, 2 , . . . , M

162

J. E, C O O P E R

and

Gi+M = e %,o,ka, COS (walk At)

i = 1, 2 . . . . , M

and these equations are solved separately for each combination of input case and measurement station. The second stage of the iteration finds updates of the frequency and damping values by linearising the model using Taylor's Theorem about the estimates so that

N

[Ei"Eu'

Ei"EiP-

"'"

Ei"Eu2M

Eij2Eijl

Eij2Eij2

"""

Eij2Eff2M

...

Eij2MEij2~

[Eu2MEu, Eij2MEij2

[Ei"] d~2wa

N

Eli2

LEu2Mj AO,)dM

or

[EE]ij[A]=[E]! i i = l , . . . , N O j = l

(24)

.... ,NI

where E o, = k A t e-¢%ka'[R u, sin (tOd,k At) + Sur cos (wd,k At)],

r=l,M

with

Eij~+M = k At e-¢~ka'[ Rijr COS (tOdrk At)

-- S!jr

sin (tOdrk At)],

r-l,M

and

Fqk= Y o k - R u o - { i (

E u r } / k At.

The update terms in vector [A] are found by solving

[El,,

[EE]I,

[EE]INI

[a] =

[E],NI

[EE]21

[E]2,

[EE] NONt

[ E ] ~o~,.

in a LS sense. The updated frequency and damping values are now used as the starting estimates for stage 1. Steps one and two are repeated until convergence is achieved. Assuming that the constant offset term, Ru0, is zero, then the R and S amplitude terms can be converted to the complex C amplitude terms through manipulation of equations (1) and (2). 4. RESULTS The impulse responses analysed in this study were generated from the FRFs obtained using burst random excitation on a Gnat aircraft [6]. The impulse responses used were transformed directly from the frequency domain to the time domain without any windowing. The comparison was carried out using four impulse responses taken from two measurement stations for two excitation cases. The FRF using 512 points from one of the responses is shown in Fig. 1. The mode at 3.2 Hz is a support mode with a very low damping ratio.

",(l!l!qms ' 0 '(%~) ad~qs apom putt %1 ,(ouanbaaj u! £~!l!qms 'r7 '%; ~u!dtut~p pul~ %I £auonbaaJ u! ,(~!l!qms '41, '%1 £auonba~j u! £1!l!q~ls ' x 'uo!1eS!l!qg]s ou '+ "0~E = N "jopDN SU.t,(a~A 'S"I 'laid '(~![!q~s "E aan~!:-I

~<>c<~x xe)((+

¢ +

X X X X X + X~exxx x x o x "¢0

+X

+)eRe xx.x +X.K X<.¢<, x +

xxx X++

xe+x+e xeo XO

*x x

e

+

0

e~l* x.ee

i;'i e~

•

ox+e

o

o

x.~+x .e X

e

x

e x x

x+ o

•

~XX

xe<

eel+

o)e

+ x

X

+¢))e+~

xQ~

xx

"~ +

+

.

x .

.

x+

.

+

x

.

.

+

+

+x

x+

x

x

+

QeO

e~K

-l< x x ~

x

+ + +

o +

o+ 0 o

O+

IL

I~

~I =I

Z

•

.:

.0£

=~ x

.0~

+ x

o +

"÷ +

x

÷ x

OHo+ o+

o o x

x x

+

x+

.. ~x

o

o o x

+

+ +

+ ÷ +

o o

x x +

+

+

++

+ "+

+

•

Cl+O.HI-

+)e

++e~o 2 ..... :x:°o

+

X+ +~.

*

+

÷

÷x~

x)e +

x+ ++

o

xx

~

x+~ + + x

++ + +

+~-

x+ x

+

+

x

x

+

.

~

+ +

x

÷

.,~ao+~( ,lO¢-K3e

XOI

•

-I~+ xx x ex x +x ............ ; +÷ 'x . " . . . . . . . . . . . .

+x

+

+

OL

+

+ +

+

+

++

0 (ZH) ~ouenbe.~d

uo!lpu!qtuoa tuntu!ldo aql au!maalap ol q~poaddp up pug o1Su!paaaoad s! ~laOhA"dI pup V N j o suo!lpu!quaoo luaaoB! p aoj pau!plqo sllnsaa oql u! uoD~!apA atuos aq u~o oaaq.L "sapotu mols,(s oql qs!n~u!:ls!p dla q ol lp!ayouoq s! uo!lpog!oadsaaAO autos qSnoq] 'anb!uqoal d;D aql liutsn uaqA~ lgau0q OU Jo s! Iopotu iPuo!l~lndtuoo pa!j!oadsaaAO ,OOA P Jo aSh ottL "pa!.II]A S! uo!lPlno[pa aql u! pasn slu!od p]pp jo aaqtunu oql pup ',lop 0E lp 'J,up~suoo ~,do~I s! aapao lapotu iL,uo!] -mndtuoa aql om!l s!ql £1uo 'L-~ slt!d u! UA~oqs OaP SlOld £l!l!qpls apl!tu!S "S~Old aql tuoaJ pall!tuo oap pup snounds aq ol poaop!suoo aap si~u!dtupp aA!lpSau ql!A~, sapotu pO1Ptu!lSTI "adpqs apotu pup ltu!dtupp '£ouonbaaj u! (uo.ual.ua pou!tuaa]ap-aad amos u!q:l!m s! alptu!lsa aalotupapd u! ai~upqo aql "o'!) pas!i!qms oap sapom aPlna!l.tpd uaqA~ A~oqs SlOId or,LL "aspo snotAaad aql moaj osoql ql!m poapdtuoo oap saalatu~attd paletu!lsa j o 1as qopa j o salptu!~sa •aaaSap alqpoo!lou p ol £aPA £[[njadoq ii!A~ sopotu snounds aql ]nq sapotu tua:lS£S uo ]oaBa alll!l OAleL1 plnoqs SUO!lPUPA qon S "sa]Ptu!lso ~alatupapd pau!plqo aq] uo (s!s,(IPUP oq] aoj pasn slu!od plt~p j o aaqtunu ao aapao uo!:mlos '-~'o) 01qP!aPA pougap aosn auo ~Ul/LIPA Jo loaBa aql A~oqs 01 posn oap S~old qon S "aap~o uo!lnlos ~utsttoaout up aoj UA~OqS 0ap sanb!uqoal S'II pup S'IU 'S'I aql £q pau!mqo saletu!lsa jo [91] SlOld ,(l!l!qttl s 't,-i; s$!~I uI •sosuodsaa aslndtu ! ~ql jo auo jo :I~1:-I "l aan$!d ( Z H ) , ' I o u e n b e J -i

06

OZ

08

09

Og

O~

0£

O~

O~

0

90

O.L ~

0"~

~9|

VJ.VO 1V~ITIIDfi~I.LS a._4V~ID~IIV

164

J. E. C O O P E R Frequency

0 O"

10

20 +

+

+

+

+

×

x

o

+o

x o

o o

+o + x

.ID x~

o

o o

o

ee

o o+x++xo + o x + X.K

~0"

+

o

v

Z

+

+o +o

o++ +'N

+

401+.,, +

5

+

I~+

+

lU+

+II-

+

++

+ ÷

÷ x+

~-xxe+x~

•

•

+

•

x+ +x +

+

x++ x+ +

+ +x

+

~-~ +

.~< +e+

+

o+x X+

+

~+ +e
+ +

4<+

+

+

+

++ ~+ x

o )1 ~

+•

~÷

+x

e + + ~ e o

x+

+

XK++

O+

~+

x +...xxxo++

xt

+

X

+x x ~++

+

:,,,xx~<

x

+

o X

x X÷

x

xo~

x

x ~ +

•+

+~

x

~( x

~o

x

i~x

•

Xl(+X

Xl< x + x

xx

[] + ~

x

x+x

x++x-

m ~o

•1<

e

x

x
+xx++

x- )el

x

~

x X

+ + x

x+x+

~O+O

~x

+

)~ +

x+xx xxxx xx

*

Xxx

+x

~x

x~ xx+

xx

×x ee

+

x+xl-+

xxdee<+

x* x~

+

+ + ÷ x

x~o

0

x+ +

+

+

+~+

+

+

+ +

x+e(

÷ ×+

x+

XX

+× +~

+ +

+x

+

+

+

O.~e~X> +Ix

+

×+ ×

++

+

+

++

+~+ +

x+~+x~ o-ex+e>l< o

+

+

+

x+ ++x+

•

+

70

+

*+

+

++

+

+

÷

+

.<.+p

÷.BI+

÷

÷+

+

o o~

60

+ +

÷

+

x

÷

50

+

+

+

÷

40

+ +

10"

(Hz)

30

+e

+~+~

+

.1~

x3

e<

~

x

•

~

xl

:,1-~

e

+ x

Figure 3. Stability plot, DLS, varying NCdof, N = 250. Key as Fig. 2.

Frequency

0

20

30

40

i

I

L

i

0

+ +

x

÷

o

+ +

x

+

e

x-

o

x

x

+

++

o

+

x

x+

+e

4.0

o

o

x

o ~

o

o

$

++$+

o

o

x++x+

o

o

x+O~K

o

o

XXX

,m

~"

o

X',~X~+

:iz,

o

.~

ex

+o

+o

cx

+0+ f,,,~o,-~< "~

0

>r,

+ .I<

x +

x

+

+

x

+

+

÷

++

x

+

+

~0+

+

,,K+

x +

÷

x

O+x

+

,,
x

+

÷

+ +

+ +

+

+ x

+ o -I<+

+ +

+

x

x

x

x

*

x ~

x

x- x

x

~,l~-X+

÷ +

+ x +++

+ x+o

x

x + x

~ ~

÷×x

x.- + x

x

~

x

xxx+xx

>~x x x~++x

+

++xo

x~~ x

x

xxx

xxK

x+

+

x~+.lK

sx+

x

x x+

+ x x x ~

~I< x

e<

+ x

x xe

$~ix

+x~+~

+

X

x-+~+-K

xl-+

++ +

.$<+ x . o K +

x

+x

+

xx

x-~ ~ x

x+ +

x +

+

+

x +

+ +

+x x

+ +

+

x x

O+-K"

+o

+ +

x

$X$ +~-

+

x

+

x

o

x

+

+ + x +

~ o + +

+

~

÷

xl< x

~.XOXl(

~,

+

4+

)0+

,,

+

+

)o

+ +

+

~<

• x

+ x

x

.B

+

+

++

o

~l 5

÷ ÷

x

$

407

+ x

+

+

+

o

30Z

+ +

+

x

70

i

+

~ o

+

60

~

+

+ ÷

x

10-

50

+

+÷ +

20-

(Hz)

10

:,.(

x- + x xxx+x x

x~x

Figure 4. Stability plot, TLS, varying NCdof, N = 250. Key as Fig. 2.

of the two variables. A stability plot for the CF method with a 20 d o f model varying the number o f data points in the analysis can be seen in Fig. 8. The singular values o f the square matrix generated in the calculation o f the A~ parameters can be used to give an indication of the number o f modes present in the data. Such singular values are a by-product o f the TLS method. (To provide a proper comparison between the singular values obtained by the various techniques, the singular values of the TLS method have to be squared and those of the CF method square rooted.) Plots showing the ratios o f consecutive singular values obtained by the LS, DLS, TLS and CF techniques appear in Figs 9-12. Theoretically, a sudden drop in the singular values (i.e., a peak o f the plotted ratio) indicates the rank of the matrix. However, it can be seen that numerous peaks occur in all o f the plots. An accumulative sum o f the singular values is also shown, normalised so that the total sum o f the singular values is unity. The

165

AIRCRAFT STRUCTURAL DATA Frequency (Hz) 0

10

20

30

I

I

I

0

100

200

+

40

50

I

70I

60 I

I

+

+

+

+

÷

++

~-

+

+

+

+

+

+

+÷

o

[]

o

o

o$

xo

¢~

x

o

*

*

x

x

~z~

D

o

o

o

oo

4-~

oi-$

o

o

o

x

x

ox

o

o

o

o

<.o

o*

o~*

o

o

o

o

a

>o<

o

o

¢

o

o0

oao

0+o

o

o

¢

*

o

o

<.

o

o

o

o

o0~

oxo

o

o

*

•

o

0

0

0

0

0

000

OXO

0

0

$

O

¢

)O<

0

0

0

.0¢0

000

0

0

0

0

0

40

i.-

n 500 Z

400-

O

500-

Figure 5. Stability plot, LS, varying N, NCdof= 20. Key as Fig. 2.

Frequency (Hz) 0 0

100.

200'

10

20

30

40

50

60

70

I

I

I

I

I

I

I

++

+

+

+

+

+ +

++

++

+

+

+

+

+

+

++

D+

cl

¢

x

o

#t3

xx

.x

x

o

o

x

x

x

xx

o

o

o

o*

~

ox

o

o

t3

•

$

*

*G

o

o

o

o

oo

<~

oo

o

o

o

o

o

*

*x

o

o

o

o

*x

~<.

oo

o

o

o

o

o

o

~<

~

<~

O0

~IO

<>0

0

•

X

X

•

X

40

o+

o+ t~

~. 3 0 0

D~

Z

400

500

O

i

Figure 6. Stability plot, DLS, varying N, NCdof= 20. Key as Fig. 2.

rank o f the system can be determined by taking the peak that corresponds to an accumulative sum o f almost unity. Accurate starting estimates of frequency and damping are needed in order to reduce significantly the amount o f computation that is required to use the SLS method. In the author's experience, there is only a danger o f divergence in the computational solution when there is overspecification of the solution order, or when there are no modes near to the initial frequency estimates which have been given. Using the FRF and estimates from the other techniques, it is possible to start with a few modes and build up the solution order mode by mode. Estimates found for one computational order can be used as starting estimates for the next. Such an approach greatly reduces the amount of computation needed if an arbitrary set of initial frequencies and dampings were taken.

166

J. E. C O O P E R F r e q u e n c y (Hz}

0

O#. . . . .

,ooi

10

20

30

40

50

70

60

+

200-

Q_ 3 0 0 ~

o

o

.

o

.

x~.

o

o

2

o

o

o

o

xo

+

*

o

[]

o

o

o

*.

×o

x

o

o

o

o

o

o

o

x~.

~

o

o

e

o

o

ot,

xo

x

o

x

x

+

+

<

•

•

×

o

~.

x

~,

x

+

.

.

×

o

*

+

,,c

+

×

e

x

x

400 ~ 50O

Figure 7. Stability plot, TLS, Varying N, NCdof= 20. Key as Fig. 2.

Frequency (Hz) 0

100-

200-

10

20

30

+

+

+

40

+

+

+

+

+

+

o

o

o

o

oo

o$

[]


o

$

o o

$0

o

o

o

x

oo

xx

x

o

o

o

×

~o

xx

o

~

50

+

60

70

÷

+

+÷

+

+

~÷*

o

x

~

~

*

x

oxx

~.

x

*

~+

~

x

m

x

o

+

x

~o

x

x

x

x

+

O3 I-" a_ 3007' +

400-

500-

Figure 8. Stability plot, CF(NA=35 IP=6), varying N, NCdof= 20. Key as Fig. 2.

It was not p o s s i b l e to find m o r e t h a n nine m o d e s with the SLS technique. At h i g h e r o r d e r s a large a m o u n t o f c o m p u t a t i o n was n e e d e d , a n d the extra m o d e s f o u n d h a d very high d a m p i n g values which were a s s u m e d to be s p u r i o u s . T a b l e 1 s h o w s the e s t i m a t e s o f the m o d a l p a r a m e t e r s for w h a t are c o n s i d e r e d to be valid system m o d e s , p r o d u c e d by the v a r i o u s m e t h o d s which have b e e n studied.

5. DISCUSSION O F RESULTS T h e F R F d i a g r a m suggests that there are g o i n g to be a b o u t 11 m o d e s with a s c a l e d a m p l i t u d e o f o v e r 0.1. H o w e v e r , there were differences b e t w e e n e a c h o f t h e F R F plots, a l t h o u g h t h e y were all d o m i n a t e d by the m o d e s at 11.5 a n d 15.7 Hz.

AIRCRAFT STRUCTURAL DATA

1

"0t

167

================================================ o° ~

0.8

4 ~2

.

~ o.6

3~

~ o.4

2§ "6 o

0.2

0"0 0

I

I

I

8

16

24

I

I

I

I

I

32 4 0 48 56 64 Singulor volue number

I

I

I

72

80

88

Figure 9. Singular value plot, LS, N C d o f = 40, N = 250.

1 0 :::::::::::::::::::::::::::::::::::::::::::::::::::::::::::

5

0.8 ;=

0.6

c o

2 o

0"4

o

0

8

16

7;

24

32 4 0 48 56 64 Singular volue number

i

72

80

88

Figure 10. Singular value plot, DLS, N C d o f = 40, N = 250.

1.0

5 o

0.8

4 "5

.~ o 6

3g

~ 04

2"6

8 0.2 0.0 0

I

I

I

8

16

24

I

I

I

I

I

32 40 48 56 64 Singular value number

I

I

I

72

80

88

Figure 11. Singular value plot, TLS, N C d o f = 4 0 , N = 250.

168

J.E. COOPER

O~¢OQ~O~000~000000~0000000000~O~O O

10

5

47

08

= 0 ~o " i 0-4-o 02 00 0

6

o=

i 2 " 86 o

,io 4

8

12

16

20

24

28

32

Singulaval r uenumber

56

40

44

Figure 12. Singular value plot, C F ( N A = 3 5 ] P = 6 ) , N C d o f = 4 0 , N =250.

The stability plots (Figs 2-8) show that when there is a large number of modes this type of plot can be very difficult to interpret. One feature of all three plots is that pole splitting occurs. It is believed that this behaviour is due partly to the separate impulse responses having slightly different resonance frequencies, so that when the model order is high enough the techniques tend to identify the parameters of each of the responses separately. Pole splitting is also a characteristic phenomenon when higher solution orders are used. More modes stabilise when the LS technique is used than for the other methods. C o m p a r e d with the other techniques, many more modes are found between 40 and 50 Hz using the LS approach. The highest modal density is in the 18-23 Hz region and when a large model order is used, all the techniques produce estimates for five or six modes in this region. In this region the results from the LS approach stabilises slightly better than the other approaches. Five modes are found in the 18-23 Hz region using the DLS and CF techniques. It is not straightforward to determine the number of modes in the system. The dominant modes are all found easily enough but it is not possible to determine which of the remaining modes are real or spurious. It is possible that the stability criteria are too stringent. More stabilised modes would be found if more lenient criteria were used. The stability plots should be used in conjunction with techniques such as modal confidence factors [17, 18] which indicate the degree of confidence in the estimate of each mode. There is considerably less variation in the stability plots when the number of data points used is varied. The LS and DLS plots are very insensitive to the number of points used in the calculation and this remains so when a larger number of modes is considered. It should be emphasised that just because a mode has stabilised it does not mean that it is a valid system mode nor that the estimated parameters are unbiased. If a smaller solution model is used, all of the LS parameters stabilise throughout the plot. The singular value ratio plots (Figs 9-12) all show a peak at the thirteenth point corresponding to an accumulative sum of nearly unity, and therefore indicate that there are six d o f in the system. Such plots would seem only useful in determining the number of dominant modes in the system rather than showing the total number of modes that are thought to occur. In Table 1 it can be seen that the estimates obtained by the LS, DLS and TLS approaches vary somewhat when the solution order is increased. More modes are found by the LS

11.54 0.42

11.53 0-45

11.53 0.43

11.53 0.47

11-53 0.50

11.53 0.44

11.53 0.45

3"21 0-24

3-21 0.24

3"22 0.42

3"20 0" 17

3"18 0.36

3.21 0-38

3.21 0.24

DLS (40)

TLS (40)

LS (20)

DLS (20)

TLS (20)

C F (20)

SLS

15.67 0"90

15'67 0"91

15.66 0.75

15.67 0"90

15.67 0"91

15"67 0-90

15"68 0.88

15-67 0-90

18.56 4"31

18"38 4-20

18"67 2-09

18"59 3"25

18.45 4.70

18"48 3-76

18"42 3.94

18.48 3.77

20.31 4-97

19"51 2.01

19-60 1.70

19"53 1 "86

21.47 1.92

21.53 1"68

21.45 3"06

21"26 1 "96

21"20 2.20

21.46 1"88

21.44 2-03

21-46 3-74 21"83 5.46

23.39 2.63

Figures in parentheses indicate number of modes in computational model.

11-53 0.46

3.21 0.24

LS (40)

28-47 2.01 28.59 2"20

23"81 2.11 23-74 2.75

38"32 1-39

36"38 2.84

37.54 2-34

28-87 1"56

23.88 2.70 23 "98 5.69

37"97 1"80

41 "42 3-01

41-83 1.44

41.50 1"82

40.48 0.73

37.78 1.49

28.51 4-15

40.48 0-98

37-78 1"44

23.82 3-51

30.92 3"56

30.84 3.60

37.79 1'89

28.90 2"84

29-01 0-89

28.34 1"85

28"42 1 "88

23 "97 0.94

24.02 1.77

24.01 1.93

Frequencies and dampings obtained by methods

TABLE 1

47-92 1.18

45.33 2.60

64.17 (Hz) 3.56 (%)

64-40 (Hz) 2.82 (%)

45.27 2.32 45-75 3.27

(Hz) (%)

65.49 (Hz) 1.03 (%)

(Hz) (%)

(Hz) (%)

65.21 (Hz) 1-76 (%)

(Hz) (%)

45.94 5.42 47.96 0.68

48.15 0.98

45.69 3.23

43.52 2.85

48.08 0-76

48"14 0"79

45-68 1.70

45"68 1"80

43.44 1.79

43.41 1-77

~3 -]

r-

,--t

), "11 .-]

7=

70

170

J.E. COOPI:R

technique. The CF method finds nearly as many modes as the LS method using half the size of computational model. The SLS technique produces estimates of what can be considered to be the nine most dominant modes. The results from all the methods agree within a small amount of scatter on the frequencies and dampings of the dominant modes. As it is real data that is being used, it is impossible to determine which of the methods has found the modal parameters closest to truth. The LS, DLS, TLS and CF methods are all non-iterative so no interaction from the user is required once the user defined parameters (order of the solution and the data points that are considered) have been set. The CF method also requires the NA and IP parameters to be defined, so extra analysis cases have to be performed in order to find the best combination of the parameters. When the TLS technique is used only one eigenvector corresponding to any of the zero eigenvalues needs to be chosen. It is easiest to choose the eigenvector corresponding to the smallest eigenvalue. The SLS method requires a lot of user interaction in the procedure described above, if the amount of computation required is to be significantly reduced. Initial estimates of frequency and damping have to be supplied as well as a convergence criterion. The time required for convergence is related to the accuracy of the initial estimates. Whereas it is relatively easy to implement the non-iterative techniques, due to the overspecification that is required the interpretation of the results can be difficult. However, although more user interaction is required with the SLS approach, as overspecification of the solution order is not needed, the results that are obtained are easy to interpret. Exact timings of the non-iterative methods discussed in this paper are described elsewhere [4]. A direct comparison of the computation times required by the methods is made more difficult as a greater solution order is required when some techniques are used in order to obtain the same degree of accuracy. Also, as the SLS method is iterative, the computation time that is required depends upon the convergence criterion and the starting estimates. The amount of computation for the CF technique depends upon the NA and IP parameters as well. Assuming that the LS, DLS and TLS methods require a solution order about twice that required by the CF and SLS techniques in order to get the same degree of accuracy, then the computational requirement is roughly in the ratio LS : DLS : TLS : CF: SLS ~ 1 : 1 : 3 : 1 : 30. As the order of the system increases, then a much greater computational requirement is required for the TLS and SLS methods than the other techniques.

6. CONCLUSIONS The results that were obtained show that the modal parameters obtained from the same set of test data can vary considerably when the number of data points and solution order are varied. Although the techniques that were used all produced comparable estimates of the dominant modes, there were a number of differences in the estimates of the other modes. The estimates found with the LS approach were less sensitive to changes than the other methods, although some doubts can be expressed as to whether all of the modes that stabilise are system modes. Stability plots should not be the sole approach used to distinguish real from spurious modes. The C F technique, using a much smaller solution order than the other non-iterative techniques, produced comparable estimates to the LS method and avoided the pole splitting problems that the LS, DLS and TLS techniques displayed. There did not appear to be any advantage in using the TLS and DLS approaches rather than the simpler LS

AIRCRAFT STRUCTURAL DATA

171

m e t h o d . T h e SLS m e t h o d p r o d u c e d estimates o f the nine m o s t d o m i n a n t m o d e s w i t h o u t n e e d i n g to use o v e r s p e c i f i c a t i o n b u t f a i l e d to e s t i m a t e a n y o t h e r m o d e s . H o w e v e r , a m u c h g r e a t e r a m o u n t o f c o m p u t a t i o n is r e q u i r e d to use the SLS t e c h n i q u e t h a n the o t h e r methods. T h e use o f the s i n g u l a r v a l u e d e c o m p o s i t i o n to d e t e r m i n e the o r d e r o f the system w o u l d a p p e a r to d e t e c t o n l y the m o s t d o m i n a n t m o d e s in the system.

ACKNOWLEDGMENT O C o n t r o l l e r , H e r M a j e s t y ' s S t a t i o n e r y Office, L o n d o n , 1990.

REFERENCES 1. H. VOLD, J. A. KUNDRAT, G. T. ROCKLIN and R. RUSSEL 1982 SAE Paper 820194. A multi-input estimation algorithm for mini-computers. 2. R. S. PAPPA and S. R. IBRAHIM 1981 Shock and Vibration Bulletin 51(3), 43-72. A parametric study of the Ibrahim time domain modal identification algorithm. 3. J. n . JUANG and R. S. PAPPA 1985 Journal of Guidance, Control and Dynamics 8, 620-627. An Eigensystem realization algorithm for modal parameter identification and model reduction. 4. J. E. COOPER 1988 Proceedings of the 6th International Modal Analysis Conference, Kissimmee, Florida, 1589-1595. Comparison of some time domain system identification techniques using approximate data correlations. Accepted for publication in the International Journal of Analytical

and Experimental Modal Analysis. 5. J. E. COOPER and J. R. WRIGHT 1986 Proceedings of the 4th International Modal Analysis Conference, Los Angeles, California, 831-836. Comparison of some time domain system identification methods for free response data. 6. F. LEMBREGTS, J. LIPKENS, J. LEUR1DAN and H.VAN DER AUWERAER 1988 Proceedings of the 6th International Modal Analysis Conference, Kissimmee, Florida, 1706-1713. Comparison of stepped sine and broad band excitation to an aircraft modal test. 7. H. A. COLE 1973 NASA-CR-2205, On-line failure detection and damping measurement of aerospace structures by random decrement signatures. 8. F. DEBLAUWE, D. L. BROWN and R. J. ALLEMANG 1987 Proceedings of the 5th International Modal Analysis Conference, London, UK, 832-845. The polyreference time domain technique. 9. S. M. CROWLEY, O. L. BROWN and R. J. ALLEMANG 1985 Proceedings of the 3rd International

Modal Analysis Conference, Orlando, Florida, 80-87. The extraction of valid residual terms using the polyreference technique. 10. S. R. IBRAHIM 1986 AIAA Journal 24, 499-503. Double least squares approach for use in structural modal identification. 11. R. S. PAPPA and J. N. JUANG 1985 Journal of Astronautical Sciences 33, 95-118. Galileo spacecraft modal identification using an eigensystem realization algorithm. 12. G. H. GOLUB and C. F. VAN LOAN 1983 Matrix Computations, pp 420-425. Oxford: North Oxford Academic Press. 13. L. ZHANG, Y. YAO and M. Lu 1987 Mechanical Systems and Signal Processing !, 399-414. An improved time domain polyreference method for modal identification. 14. J. R. WRIGHT 1974 Journal of Aircraft I1, 774-777. Flutter test analysis in the time domain using a recursive system representation. 15. W. R. SMITH 1981 AIAA Paper 810530. Least squares time domain method for simultaneous identification of vibration parameters from multiple free response records. 16. F. LEMBREGTS, R. SNOEYS and J. LEURIDAN 1987 International Journal of Analytical and Experimental Modal Analysis 2, 19-31. Application and evaluation of multiple input modal parameter estimation. 17. S. R. IBRAHIM 1978 Journal of Spacecraft and Rockets 15, 313-316. Modal confidence factor in vibration testing. 18. H. VOLD and J. CROWLEY 1985 Proceedings of the 3rd International Modal Analysis Conference, Orlando, Florida, 305-310. A modal confidence factor for the polyreference method.

172

J. E. C O O P E R

APPENDIX: NOMENCLATURE block difference equation coefficient complex amplitude term CF correlation fit D diagonal matrix of system eigenvalues dof degrees of freedom (number of modes) DLS double least squares term used in second stage of SLS technique Eit,. [E] vector of E,, terms [EE] matrix of E~, terms ERA Eigensystem realisation algorithm F,, term used in second stage of SLS technique frequency response function ERF term used in first stage of SLS technique Gi Ho, H~ data matrices for the TLS technique variable allowing corrupted correlations to be omitted from CF technique IP ITD lbrahim time domain method LS least squares M number of modes MIMO multi-input-multi-output N number of data points considered NA variable controlling number of correlations used by CF method NCdof number of modes in computational solution NI number of input cases number of measurement stations NO r vector of correlation blocks constant offset term Riio R q, real amplitude term (r = 1 , . . . , M) real amplitude term (r = 1 , . . . , M) Sip single input-single output SISO Smith least squares SLS SVD singular value decomposition TLS total least squares V modal participation factor matrix [ v] vector used in estimation of mode shape matrix W complex mode shape matrix output at kth time instant, ith measurement station, jth input case Y(ik r~ block output matrix at kth time instant { Y] { ~'] vector of block response matrices At sampling interval {a] vector of updated frequency and damping terms critical damping, jth mode at matrix of block response matrices time shifted block response matrices 0 vector of difference equation coefficients natural frequency, jth mode o) i damped natural frequency, jth mode (odj system eigenvalues output autocorrelation at jth lag A diagonal matrix of NI smallest eigenvalues

ai CHr