- Email: [email protected]

A N E W A P P R O A C H FOR T H E C A L C U L A T I O N OF RESPONSE SPECTRAL DENSITY OF A LINEAR S T A T I O N A R Y R A N D O M MULTIDEGREE OF FREEDOM SYSTEM A. M. SHARAN, S. SANKAR AND T. S. SA~KArt Department o[ Mechanical Engineering, Concordia University,Montreal, Canada H3 G 1M8 (Received 12 May 1981, and in revised[orm 8 December 1981) A new approach for the calculation of response spectral density for a linear stationary random multidegree of freedom system is presented. The method is based on modifying the stochastic dynamic equations of the system by using a set of auxiliary variables. The response spectral density matrix obtained by using this new approach contains the spectral densities and the cross-spectral densities of the system generalized displacements and velocities. The new method requires significantly less computation time as compared to the conventional method for calculating response spectral densities. Two numerical examples are presented to compare quantitatively the computation time.

1. INTRODUCTION In random vibration of systems, the calculation of the response spectral density is of prime interest due to its explicit relationship to the mean square value of system amplitude. Several approaches have been used in the past to calculate the response of linear stationary single- and multidegree of freedom systems [1--4]. The most commonly used approaches are the impulse response method and complex frequency response method. For an n degree of freedom mechanical system, the governing dynamical equation can be written as [m] (nxn)

i(t) (n•

+

[c] (n•

:k(t) + [ k ] (nxl) (n•

x(t) (n•

=

f(t) , (nxl)

(1)

where [m], [c] and [k] are the mass, damping and stiffness matrices respectively and fit) is the vector of random excitation forces. If [~bf(to)] is the n • n excitation spectral density matrix of the excitation force f(t), then the response spectral density is given by [5] [~b~(co)] = [H(jco)][~br(co)][H*(j~o)] T,

(2)

where [H(j~o)] is a (n • n) matrix of the transfer function of the system in equation (1), [H*(j~o)] T is the transpose of the complex conjugate of [H(jto)] and [~bx(co)] is an n x n response spectral matrix consisting of diagonal elements representing the spectral density of generalized displacement and off-diagonal elements representing the cross-spectral density of generalized displacements. [H(jco)] can be found by calculating the Fourier transform of equation (1) and can be represented as [H(jco)] = [-coZEm] + jo~[c] + [k]] -1. (3) 513 0022-460x/82/160513 + 07 $03.00/0 9 1982 AcademicPress Inc. (London)Limited

514

A . M . SHARAN, S. S A N K A R A N D T. S. S A N K A R

The calculation of the response spectral density [r (o9)] in general is carried out over all frequencies of interest and it involves matrix addition, multiplication and inversion. In particular, during the caleulation of [r the matrix inversion is usually a time consuming operation in a computer. Further, if the degrees of freedom or the frequencies of interest are large in number, then the computer CPU time for the evaluation of [r (~o)] also becomes large. In this paper, a new method for the calculation of the response spectral density [r (w)] is introduced which eliminates any matrix inversion. Further, this method provides, in addition to [r the response spectral densities [r [r and [r in a single calculation. 2. MATHEMATICAL FORMULATION By using a new set of variables, Y(t) and Y(t), equation (1) can be rewritten as [6]

(4)

[M]Y(t) + [K]Y(t) = F(t),

:rro

r-,l I

COl1

[M] L[m] [c]J' (2n x2n) I'x("!t

Y(t) = , (2n x l ) Lx(/)J

[K] =t LuJ [k]J' (2n x2n) t 'x!t!t

~r(t) = , (2n x l ) L~(/)J

F(t) = { ; i ; } . (2n x1) t

Let iV] be the model matrix (2n x 2n) consisting of the eigenvectors of equation (4). Defining u V]Z(t) and ~'(t)= [V]7,(t), and premultiplying each term in equation (4) by t V]T gives [M*]:Z(t) + [K*]Z(t) = [ v]TF(t),

(5)

where

[M*]=[vjT[M][V]

and

[/,2*]=Iv]TIE]IV]

(5a)

are diagonal matrices of size (2n x2n). Assuming the initial conditions are zero and taking the Laplace transformation of equation (5) gives Z(s) = [s[M*] + [K*]]-I[ V]TF(s).

(6)

Premultiplying by [ V] converts equation (6) to

Y(s)=[A(s)]F(s),

[A(s)]=[V][G*(s)][V] T,

[G*(s)]=[s[M*]+[K*]] -1.

(7a, b) The matrix [G*(s)] is a diagonal matrix with its diagonal elements being the reciprocals of the diagonal elements of the matrix [s[M*] + [K*]]. From equation (7) it is evident that [A(s)] represents the transfer function of the system given in equation (4). Then the expression for the response spectral density of the system (governed by equation

(4)) is [r = [A(joJ)] (2n • 2n) (2n • 2n) where

[r (2n x 2n)

[A*(jo:)] T, (2n • 2n)

(8)

R A N D O M M U L T I D E G R E E OF FREEDOM SYSTEM

515

is the excitation spectral matrix corresponding to the force vector F(t), [A*(joJ)] r is the transpose of the complex conjugate of A(j~o), and [~bv(~o)]is a 2n x 2n response spectral matrix. If [$v(ta)] is partitioned, then [r162

LCx~(~o) Cx(o~)J'

(9)

From equation (9) it is evident that the calculation of [4,Y(~o)]provides response spectral matrices [~bx(~o)], [~b~(~o)], [qSx~(~o)]and [qS~(~o)]. The procedure outlined in this section not only provides four times the information as compared to the solution obtained by using equation (2) but also the necessity of inverting the matrix corresponding to the system transfer function is avoided. 3. STEP-BY-STEP PROCEDURE OF THE NEW METHOD 1. Formulate the modified governing equation of the system as shown in equation (4). 2. Obtain the modal matrix IV], and the corresponding natural frequencies of the system represented by equation (4). 3. Calculate the diagonal matrices [M*] and [K*] given by equations (5a). 4. Calculate the matrices [A(s)] and [G*(s)] given by equations (7a) and (7b), respectively. 5. Calculate the spectral density of the response matrix [$v(~o)] by using equation (8). The required computations for each of the steps, starting from step 2, were done by using a package written in FORTRAN language. The modal matrix [V] and the corresponding natural frequencies were obtained by using a subroutine EIGCC from the IMS library of subroutines. 4. NUMERICAL EXAMPLES 4.1.

THE

INERTIA

M A T R I X D I A G O N A L , AND THE

STIFFNESS

AND THE DAMPING

MATRICES, N O N - D I A G O N A L

The first system chosen to provide a numerical example is shown in Figure 1, and the details of the various parameters are given in Appendix A. In practice, for mechanical systems, the matrix [$s(oJ)] can be determined experimentally. The necessary procedure for its determination can be found in reference [7]. XI

----I~ e, __

X3

C6 E

, ,,.,,

r-1~'2

---~a

Figure 1. Three degrees of freedom system.

Numerical computations were carried out with the parameter set as in Appendix A to obtain, (a) the matrix [~bx(o:)] by using equation (2), and (b) the matrix [~by(~o)]by using equation (8), respectively. The calculations were repeated with the maximum frequenb~, limit set of 2000, 3000 and 4000 Hz. In all cases the lowest frequency was selected as 5 Hz with the frequency step also chosen as 5 Hz. The computations were

516

A. M. SHARAN. S. SANKAR AND T. S. SANKAR

carried out on a CDC Cyber 174 digital computer. For a pre-selected maximum frequency limit, the matrix [&~(~o)] calculated by using the conventional approach (equation (2)) was identical to that obtained by using the method proposed in this paper (equation (8)). The computation times for each of the three sets mentioned above are given in Table 1. TABLE 1

Computation time for the response power spectral density matrix of the three degrees of freedom mechanical system Computation time (s)

Matrix

[~b~(aJ)] (3 x 3)

[r

Equation(s) used

Case 1: max. frequency limit 2000 Hz

Case 2: max. frequency limit 3000 Hz

Case 3: max. frequency limit 4000 Hz

(2)

8'821

13"345

19.960

(8) and (9)

3.501

5"195

7.005

(8)

11"488

16'991

22"550

(3 • 3) ICy(ca)] (6 X 6) t ~r

is obtained by using equation (8) and calculating only the eIements a,+1.,+1 to a2,,2,.

4.2. THE INERTIA MATRIX, NON-DIAGONAL, AND THE STIFFNESS AND THE DAMPING MATRICES, DIAGONAL

The second system considered in numerical work is shown in Figure 2. This system represents an active suspension system of a tractor cab [8]. This problem involves the solution of a four degrees of freedom system. The four suspension units under each corner of the cab base plate isolate the vibrations in the vertical, pitch and roll modes. The fifth suspension unit parallel to the cab base plate provides isolation in the side-to-side

r

J cob-/

suspension~ 17 RoseJ "

'~x 2 x=

I

X~IX3

Figure 2. Geometric representation of the four degrees of freedom cab system.

RANDOM MULTIDEGREE

517

OF F R E E D O M S Y S T E M

or the lateral direction. The forces present are along each of the co-ordinates xt, x2, x3 and x4. The system is subjected to base excitation. This system can also be represented by the matrix differential equation (1) [8]. The details of the various parameters of this system are given in Appendix B. The response spectral density matrices [q~x(ca)] and [4~v(w)], mentioned in section 4.1, were computed in a similar manner. Once again, there was no complex matrix inversion involved in the computation of the matrix [4~, (~o)]. The details of the results of these computations are given in Table 2. It should be noted that the elements of [4~x(oo)] calculated by using the conventional approach (equation (2)) were identical to those obtained by using the proposed method. TABLE

2

Computation time for the response power spectral density matrix of the cab suspension system Computation time (s)

Matrix [4,.(o~)]t (4 x 4) [r (4 • 4) [g,v(~o)] (8 x 8)

Equation(s) used

Case i: max. frequency limit 2000 Hz

Case 2: max. frequency limit 3000 Hz

Case 3: max. frequency limit 4000 I-Iz

(2)

13'577

20"359

26.876

(8) and (9)

5.517

8.024

10.712

(8)

23.958

35"370

47"338

The results of computation indicate that [4~x(~o)]tand [~b~(~o)]r were identical. 5. RESULTS AND DISCUSSION From Tables 1 and 2 it can be seen that there is a significant reduction in the computation time in obtaining the matrix [~b.(oJ)] if the new approach is used instead of equation (2). Therefore, one can infer that the reduction in computer time realized in using the new approach would become large as the size of the matrix [~bx(o~)] increases. As a point of clarification it should be pointed out that the matrix [~b.(~)] in equation (9) has been computed without actually computing the matrix [r v(~O)]. This was possible by partitioning the matrices involved in the calculation. In this new method, when [r is calculated, the natural frequencies and the corresponding mode shapes of the system are automatically obtained as a by-product. In some vibrating systems such as the one in section 4.2, the m o d e shapes and natural frequencies characterize the rigid body motion of the vehicle under consideration and thus provide very valuable information in the overall suspension design of the vehicle. In these situations, using the method proposed in this paper, one can obtain the vehicle response to random road surface undulations, [~b~(o~)], and vehicle natural frequencies and mode shapes in one single calculation. Further, it is possible by using the m e t h o d proposed to obtain the spectral density [&~], and the cross-spectral densities [4~x~] and [4~] without any added computational difficulty.

518

A.M.

SHARAN,

S. S A N K A R

AND

T. S. S A N K A R

6. CONCLUSIONS In this paper a new method for the calculation of response spectral density of a linear stationary random multidegree of freedom system has been presented. T h e m e t h o d is based on modifying the governing equations of motion of the system with a set of auxiliary variables. As a first step in the calculation of [q~ (~o)], the natural frequencies and the modal matrix of the modified system were obtained, which in turn were used to obtain the transfer function of the modified system. The procedure used in the new m e t h o d is such as to make any matrix inversion corresponding to the system transfer function completely unnecessary. Further, this new method requires significantly less computation time for the calculation of [$~ (~o)] than the conventional method. Two numerical examples are given to compare quantitatively the computation time. T h e proposed m e t h o d can be used, if needed, to calculate the response spectral density of the generalized velocity [$~ (~o)] and the cross-spectral densities [q~x~(~o)] and [ ~ (to)] without any a d d e d computational difficulty.

REFERENCES 1. S. H. CRANDALL and W. D. MARK 1963 Random Vibration in Mechanical Systems. New York: Academic Press Inc. 2. J. D. ROBSON 1963 A n Introduction to Random Vibration. Edinburgh University Press. 3. G. E. UHLENBECK and M. C. WANG 1945 Reviews of Modern Physics 17, 323-342 On the theory of Brownian motion II. 4. T. K. CAUGIqEY and H. J. STUMPF 1961 Journal of Applied Mechanics 28, 563-566. Transient response of a dynamic system under random excitation. 5. I. M. YANG 1970 Ph.D Thesis, California Institute of Technology, Pasadena, California, U.S.A. Stationary random response of multidegree-of-freedom systems. See p. 12. 6. L. MEIROVITCH 1969 Analytical Methods in Vibrations. Toronto: Collier-MacMillan Canada Limited. See p. 410. 7. J. S. BENDAT and A. G. PIERSOL 1971 Random Data: Analysis and Measurement Procedure. Toronto: John Wiley & Sons, Inc. See pp. 271-284. 8. D. G. ROLEY 1975 Ph.D. Thesis, University of California, Davis, California, U.S.A., Tractor cab suspension modeling. See p. 16.

APPENDIX A: VARIOUS MATRICES OF THE MECHANICAL SYSTEM 1. Inertia matrix (kg):

0 ~ 1 0 45.454

-45.454 [m] =

0.

0

2. Stiffness matrix (N/m):

90.909

1 !]

[k] = 176'458 x 106 I - !

2

--

.

-1 . Damping matrix (N s/m): F 3.3644 [c]--- 1750.1 I - 0 . 3 7 4 7

L-O.O581

-0,3747 3.3337 -0.4079

-0.0581] -0'4079 / . 6.3788.J

519

R A N D O M MULTIDEGREE OF FREEDOM SYSTEM

.

Excitation spectral density matrix (NZ/Hz): - 1 0 . 0 Io0go) + 100

0 - 9 , 5 log aJ +

[#)f(co)] =

o

0 - 9 . 5 log ~o+ 10

100

0

APPENDIX B: VARIOUS MATRICES OF THE TRACTOR CAB SUSPENSION SYSTEM [8] 1. Inertia matrix (kg): "681,82 0 [m] = 0 0

0 736.67 0 0

0 0 395.9 519,5

0 0 519.5 681.82

2. Stiffness matrix (N/m):

r 26844

0 5606.37 0 0

0 0 5606.34 0

0 0 0 9483.s

0 426,63 0

i]

3. Damping matrix (N s/m):

[c] =

-2042. $ 0 0 0

426.63 0 0

~

~

"

729.51

4. Excitation spectral density matrix (N2/Hz) (assumed):

[~f(~o)]

=

- 1 0 . 0 log co + 100

0

0

0 0 0

-9"5 log co + 100 0 0

0 -9"0 log co + 100 0

0 ]

0 0 -8.5 log~ +100