Transient responses of dynamical systems with random uncertainties

Probabilistic Engineering Mechanics 16 (2001) 363±372

www.elsevier.com/locate/probengmech

Transient responses of dynamical systems with random uncertainties C. Soize* Universite de Marne-la-ValleÂe, Ing. 2000 Bat. Copernic, 5 Bd Descartes, 77454 Marne-la-Vallee Cedex 2, France

Abstract A new approach is presented for modeling random uncertainties by a nonparametric model allowing transient responses of mechanical systems submitted to impulsive loads to be predicted in the context of linear structural dynamics. The probability model is deduced from the use of the entropy optimization principle whose available information involves the algebraic properties related to the generalized mass, damping and stiffness matrices which have to be positive-de®nite symmetric matrices, and the knowledge of these matrices for the mean reduced matrix model. An explicit construction and representation of the probability model have been obtained and are very well suited to algebraic calculus and to Monte Carlo numerical simulation in order to compute the transient responses of structures submitted to impulsive loads. Finally, a simple example is presented. q 2001 Elsevier Science Ltd. All rights reserved. Keywords: Random uncertainties; Dynamical systems; Structural dynamics; Transient response; Impulsive load; Entropy optimization principle

1. Introduction This paper deals with predicting the transient responses of structures submitted to impulsive loads in linear structural dynamics. In general, this kind of prediction is relatively dif®cult because the structural models have to be adapted to large, medium and small vibrational wavelengths which correspond to the low-, medium- and high-frequency ranges. Here, we are interested in the case where the impulsive load under consideration has an energy which is almost entirely distributed over a broad low-frequency band and for which prediction of the impulsive load response can be obtained with a reduced matrix model constructed using the generalized coordinates of the mode-superposition method associated with the structural modes corresponding to the n lowest eigenfrequencies of the structure. Under the above assumptions and for a complex structure, dimension n of the reduced matrix model generally has to be high (several dozen or hundred structural modes may be necessary to predict transient responses). However, it is known that the higher the eigenfrequency of a structural mode, the lower its accuracy because the uncertainties in the model increase (in linear structural dynamics and vibrations, the effects of uncertainties on the model increase with the frequency and it should be kept in mind that the mechanical model and the ®nite element model of a complex structure * Tel.: 133-1-6095-76-61; fax: 133-1-6095-7657. E-mail address: [email protected] (C. Soize).

tend to be less reliable in predicting the higher structural modes). This is why random uncertainties in the mechanical model have to be taken into account. This is a fundamental problem in structural dynamics when the mechanical model has to be adapted to predict a transient response for which not only the low-frequency band is mainly concerned, but also the upper part of this low-frequency band and may be the medium-frequency-band have to be taken into account. Random uncertainties in ®nite element models are usually modeled by using parametric models. Concerning such a parametric approach, for general developments, we refer the reader to Refs. [1±7] and for aspects related to stochastic ®nite elements, we refer the reader to Refs. [8±13]. The structural modes corresponding to the n lowest eigenfrequencies of the complex structure are calculated using the ®nite element method (see for instance Refs. [14±18]). This paper presents a new nonparametric probabilistic model of random uncertainties for reduced matrix models of structures in order to predict transient responses due to impulsive loads. The information used does not require the description of the local parameters of the mechanical model. This nonparametric model of random uncertainties is based on a probability model introduced in Refs. [19,20] for symmetric positive-de®nite real random matrices deduced from the entropy optimization principle. The available information is only constituted of the mean value of the generalized mass, damping and stiffness matrices of the mean reduced matrix model which is deduced from the mean ®nite element model.

0266-8920/01/$ - see front matter q 2001 Elsevier Science Ltd. All rights reserved. PII: S 0266-892 0(01)00026-1

364

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

2. A probability model for symmetric positive-de®nite random matrices In this section, we recall the main results established in Refs. [19,20] concerning the construction of a probability model for random matrices with values in the set M1 n …R† of all the (n £ n) real symmetric positive-de®nite matrices using the entropy optimization principle which allows only the available information to be used. It should be noted that the results obtained and presented below differ from the known results concerning the Gaussian and circular ensembles for random matrices such as orthogonal (GOE), sympletic, unitary and antisymmetric hermitian ensembles which have been extensively studied in the literature (see for instance Ref. [21]). In addition, we complete the construction given in Refs. [19,20] in order to obtain a consistent probabilistic model which allows the convergence properties to be studied when dimension n approaches in®nity. 2.1. Probability density function on the space of positivede®nite symmetric real matrices and characteristic function MSn …R†

the mean value [A] of random matrix [A] Z ~ ˆ ‰AŠ; ‰AŠp‰AŠ …‰AŠ† dA E{‰AŠ} ˆ 1 Mn …R†

in which E denotes the mathematical expectation and where mean value [A] is given in M1 n …R†: In addition, we assume that random matrix [A] is such that E{ln…det‰AŠ†} ˆ v

Z M1 n …R†

~ ˆ 1: p‰AŠ …‰AŠ† dA

…3†

For all [Q ] in M nS …R†; the characteristic function of S random matrix [A] with values in M1 n …R† , Mn …R† is de®ned by F ‰AŠ …‰QŠ† ˆ E{exp…iR‰QŠ; ‰AŠS†} in which R‰QŠ; ‰AŠS ˆ tr{‰QŠ‰AŠ T } ˆ tr{‰QŠ‰AŠ where tr denotes the trace of matrices and where [A] T denotes the transpose of matrix [A]. We then have

F ‰AŠ …‰QŠ† ˆ

Z M n1 …R†

~ exp…iR‰QŠ; ‰AŠS†p ‰AŠ …‰AŠ† dA:

…4†

2.2. Available information for construction of the probability model We are interested in the construction of the probability distribution of a second-order random variable [A] with values in M1 n …R† for which the available information is

…7†

Mn …R†

Z

This probability density function is such that

…6†

in which g $ 1 is a positive integer and where i‰AŠiF ˆ …tr{‰AŠ‰AŠT }†1=2 is the Frobenius norm of matrix [A] in Mn(R) (the set of all the (n £ n) real matrices). Consequently, from Eqs. (3), (5) and (6), we deduce that the constraints imposed for the construction of the probability model of random matrix [A] with values in M1 n …R† are Z ~ ˆ 1; p‰AŠ …‰AŠ† dA …8† 1

~ P‰AŠ ˆ p‰AŠ …‰AŠ† dA

1#i#j#n

uvu , 11:

E{i‰AŠ21 igF } , 11;

Z

is de®ned by a probability density function ‰AŠ 7 ! p‰AŠ …‰AŠ† 1 from M1 n …R† into R ˆ ‰0; 11Š with respect to the ~ on MSn …R† de®ned [19,20] by measure (volume element) dA Y ~ ˆ 2n…n21†=4 dA d‰AŠij : …2†

with

We proved in Refs. [19,20] that the constraint de®ned by Eq. (6) allows us to obtain the existence of moments related to the inverse random matrix [A] 21

Let be the set of all the (n £ n) real symmetric matrices. Let [A] be a random matrix with values in S M1 n …R† , Mn …R† whose probability distribution …1†

…5†

M1 n …R†

M1 n …R†

~ ˆ ‰AŠ [ M1 ‰AŠp‰AŠ …‰AŠ† dA n …R†; ~ ˆ v; ln…det‰AŠ†p‰AŠ …‰AŠ† dA

…9† …10†

in which uvu , 11. 2.3. Probability model using the maximum entropy principle The measure of entropy [22] and the maximum entropy principle [23±25] are introduced to construct the probability model of random matrix [A] with values in M1 n …R† based only on the use of the available information de®ned by Eqs. (8)±(10). Let l A be the real parameter such that 1 2 l A is the Lagrange multiplier corresponding to the constraints de®ned by Eq. (10). It can then be proved [19,20] that, for l A . 0 and ‰QŠ [ M nS …R†; probability density function p[A]([A]) and characteristic function F [A]([F ]) of positivede®nite random matrix [A] are written as p ‰AŠ …‰AŠ† ˆ 1M1n …R† …‰AŠ† £ cA £ …det‰AŠ†lA 21   …n 2 1 1 2lA † tr{‰AŠ21 ‰AŠ} ; £ exp 2 2  F ‰AŠ …‰QŠ† ˆ det…‰I n Š 2

…11†

 2…n2112l A †=2 2i ‰AŠ‰QŠ† …n 2 1 1 2lA † …12†

in which det is the determinant of the matrices, [In] is the (n £ n) unity matrix and where 1M1n …R† …‰AŠ† is equal to 1 if 1 ‰AŠ [ M1 n …R† and is equal to zero if ‰AŠ Ó Mn …R†: When l A is an integer, the probability distribution de®ned by

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

Eq. (11) or (12) coincides with a Wishart distribution [26]. If l A is not an integer, then the probability distribution de®ned by Eq. (11) or (12) is not a Wishart distribution. In Eq. (11), positive constant cA is written as   n 2 1 1 2lA n…n2112lA †=2 …2p†2n…n21†=4 29 cA ˆ 8 ;  = n
…13†

where GR(z) is the gamma function de®ned for Re z . 0 by 1 z21 2t G…z† ˆ 1 t e dt: The range of parameter l A satisfy0 ing Eq. (7) in which g $ 1 is a positive integer has to be determined. For g $ 1, it can be proved that

lA . g 1 1 )

Efk‰AŠ21 gF

, 11;

g $ 1:

…14†

lA . 0 ) E{k‰AŠkhF } , 11;

; h $ 0:

…15†

Eq. (15) means that for l A . 0, all the moments of random matrix [A] exist (h is any positive integer). The A covariance Cjk;j 0 k 0 ˆ E{…‰AŠjk 2 ‰AŠjk †…‰AŠj 0 k 0 2 ‰AŠj 0 k 0 †} of random variables [A]jk and [A]j 0 k 0 is written as A 21 0 k 0 ˆ …n 2 1 1 2lA † Cjk;j {‰AŠj 0 k ‰AŠjk 0 1 ‰AŠjj 0 ‰AŠkk 0 }; …16† A of random variable [A]jk is such and the variance VjkA ˆ Cjk;jk A that Vjk ˆ …n 2 1 1 2lA †21 {‰AŠ2jk 1 ‰AŠjj ‰AŠkk }: Since [A] is a positive-de®nite real matrix, there is an upper triangular matrix [LA] in Mn(R) such that

‰AŠ ˆ ‰LA ŠT ‰LA Š;

…17†

which corresponds to the Cholesky factorization of matrix [A]. Considering Eq. (17), random matrix [A] can be written as ‰AŠ ˆ ‰LA ŠT ‰GA Š‰LA Š

…18†

in which matrix [GA] is a random variable with values in M1 n …R†: From Eqs. (5) and (18), we deduce that the mean value [GA] of random matrix [GA] is such that ‰GA Š ˆ E{‰GA Š} ˆ ‰In Š:

…19†

The probability density function p[GA]([G]) with respect to ~ on MSn …R† of random matrix [GA] with values in measure dG 1 Mn …R† is given by Eqs. (11) and (13) in which [A] has to be replaced by [In]. We then have p‰GA Š …‰GŠ†

in which positive constant CGA is such that   n 2 1 1 2lA n…n2112lA †=2 …2p†2n…n21†=4 2 8 9 CGA ˆ :   n
ˆ 1M1n …R† …‰GŠ† £ CGA £ …det‰GŠ†lA 21   …n 2 1 1 2lA † tr‰GŠ ; £ exp 2 2

…20†

…21†

From Eqs. (16) and (18), we deduce that the covariance G Cjk;j 0 k 0 of random variables [GA]jk and [GA]j 0 k 0 , de®ned G by Cjk;j is 0 k 0 ˆ E{…‰GA Šjk 2 ‰GA Šjk †…‰GA Šj 0 k 0 2 ‰GA Šj 0 k 0 †}; written as G 21 Cjk;j {…‰GA Šj 0 k ‰GA Šjk 0 † 1 …‰GA Šjj 0 ‰GA Škk 0 }: 0 k 0 ˆ …n 2 1 1 2lA †

…22† Since ‰GA Š ˆ ‰In Š; the variance variable [GA]jk is such that VjkG ˆ …n 2 1 1 2lA †21 …1 1 djk †;

In addition, we have

365

VjkG

ˆ

G Cjk;jk

of random …23†

in which djk ˆ 0 if j ± k and djj ˆ 1: Let dA . 0 be de®ned by ( )1=2 2 E{i‰GA Š 2 ‰GA ŠiF } dA ˆ : …24† i‰GA Ši2F P P Eq. (23) yields E{i‰GA Š 2 ‰GA Ši2F ˆ j k VjkG ˆ 2 21 n…n 1 1†…n 2 1 1 2lA † and since i‰GA ŠiF ˆ i‰In Ši2F ˆ n; we deduce that 1=2  n11 dA ˆ ; …25† n 2 1 1 2l A and consequently

lA ˆ `A …n†;

…26†

in which n 7 ! `A …n† is the mapping de®ned on the set N p of all positive integers such that `A …n† ˆ

1 2 d2A 1 1 d2A n 1 : 2d2A 2d2A

…27†

From Eqs. (14) and (25), we deduce that parameter d A has to be such that s n11 , 1; g $ 1; ; n $ 1: …28† 0 , dA , n 1 1 1 2g Eq. (28) shows that g has to be chosen as small as possible in order to increase the domain of possible values for d A. From convergence considerations when n ! 11 and from Eq. (28), it can be deduced that g ˆ 2 is an optimal value (see Ref. [27]). Let n0 $ 1 be a ®xed integer. Taking the value g ˆ 2; we then deduce that, if parameter d A satis®es s n0 1 1 ; …29† 0 , dA , n0 1 5 then, ; n $ n0 ; we have lA ˆ `A …n† . g 1 1 ˆ 3 and consequently, Eq. (14) holds. These equations will be

366

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

used as follows. The lower bound n0 of positive integer n is ®xed. Then, the dispersion of the probability model is ®xed by giving parameter d A, independent of n, a value such that Eq. (29) is satis®ed. For each value of integer n $ n0 ; parameter lA ˆ `A …n† is then calculated by using Eq. (27). Consequently, lA ˆ `A …n† appears as a function of n.

3. For j ˆ j 0 ; positive-valued p random variable [LA]jj can be written as ‰LA Šjj ˆ vYj in which v is given by Eq. (34) and where Yj is a positive-valued gamma random variable whose probability density function with respect to dy is given by

G j …y† ˆ

2.4. Monte Carlo simulation of random matrix [A] when l A is an integer When lA ˆ `A …n† is a positive integer, we introduce the positive integer mA such that mA …n† ˆ n 2 1 1 2`A …n†:

…30†

Substituting Eq. (27) in the right-hand side of Eq. (30) yields mA …n† ˆ …n 1 1†=d2A :

…31†

Since mA(n) is a positive integer, it can be veri®ed that the probability distribution de®ned by Eq. (11) or (12) is a Wishart distribution [26] and that random matrix [A] can be written [19,20] as ‰AŠ ˆ

A …n† 1 mX …‰LA ŠT Xj †…‰LA ŠT Xj †T ; mA …n† jˆ1

…32†

in which [LA] is the upper triangular matrix de®ned by Eq. (17) and where X1 ; ¼; XmA …n† are independent random vectors, each vector Xj being an R n-valued second-order Gaussian random variable, centered and whose covariance matrix is ‰CXj Š ˆ E{Xj XTj } ˆ ‰In Š: Consequently, Eq. (32) gives an ef®cient procedure for the Monte Carlo simulation of random matrix [A]. 2.5. Monte Carlo simulation of random matrix [A] when l A is a not an integer Let us now assume that lA ˆ `A …n†; given by Eq. (27), is a positive real number (the particular case for which l A is a positive integer is presented above in Section 2.4). Since [GA] de®ned by Eq. (18), is a random matrix with values in M1 n …R†; the Cholesky factorization allows us to write ‰GA Š ˆ ‰LA ŠT ‰LA Š

a:s:;

…33†

in which [LA] is an upper triangular random matrix with values in Mn …R†: The following results, which allow a procedure for the Monte Carlo simulation of random matrix [A] to be de®ned, are proved [19,20] 1. Random variables {‰LA Šjj 0 ; j # j 0 } are independent. 2. For j , j 0 ; real-valued random variable [LA]jj 0 can be written as ‰LA Šjj 0 ˆ 221=2 ‰L~ A Šjj 0 in which ‰L~ A Šjj 0 is a real-valued Gaussian random variable with zero mean and variance given by v ˆ 2…n 2 1 1 2`a …n††21 :

…34†

1‰0;11‰…y†  y……n2j12`A …n††=2†21 e2y : …35† n 2 j 1 2` A …n† G 2 

4. We have ‰GA Š ˆ ‰LA ŠT ‰LA Š and ‰AŠ ˆ ‰LA ŠT ‰GA Š‰LA Š:

2.6. Probability model of a set of positive-de®nite symmetric real random matrices Let us consider n random matrices ‰A1 Š; ¼; ‰An Š with values in M1 n …R† such that for each j in {1; ¼; n}; the probability density function of random matrix [Aj] satis®es Eqs. (8)±(10). This means that only the mean values of the random matrices are known. Applying the maximum entropy principle, it can be proved that the probability density function …‰A1 Š; ¼; ‰An Š† 7 ! ¼ £ M1 p‰A1 Š;¼;‰An Š …‰A1 Š; ¼; ‰An Š† from M1 n …R† £ n …R† into 1 ~ 1£ R with respect to the measure (volume element) dA ~ n on MSn …R† £ ¼ £ MSn …R† is written as ¼ £ dA p‰A1 Š;¼;‰An Š …‰A1 Š; ¼; ‰An Š† ˆ p‰A1 Š …‰A1 Š† £ ¼ £ p‰An Š …‰An Š†; …36† which means that ‰A1 Š; ¼; ‰An Š are independent random matrices. 3. Transient responses of structural dynamical systems with random uncertainties 3.1. Introduction of the mean ®nite element model for transient responses to impulsive loads Let us consider the linear transient reponse of a damped ®xed structure around a static equilibrium con®guration considered as a natural state without prestresses, submitted to an impulsive load. We introduce the ®nite element model considered as the `mean ®nite element model' of this mechanical system. The transient response {y…t†; t $ 0} of this mean ®nite element model is the solution of the following time evolution problem:  ‡ ‰DŠy…t† _ ‡ ‰KŠy…t† ˆ x…t†; ‰MŠy…t†

t $ 0;

…37†

with the initial conditions y…0† ˆ y0 ;

_ ˆ y1 ; y…0†

…38†

in which y ˆ …y1 ; ¼; ym † is the transient response vector of the m DOFs (displacements and/or rotations), x ˆ …x1 ; ¼; xm † is the impulsive load vector of the m inputs (forces and/or moments), y0 and y1 are the given initial

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

_ conditions for displacement and velocity vectors y and y; respectively. The mass, damping and stiffness matrices ‰MŠ; ‰DŠ and ‰KŠ are positive-de®nite symmetric …m £ m† real matrices (the structure being assumed to be ®xed, there are no rigid body displacements). Note that underlined quantities refer to the `mean ®nite element model'.

in which yS ˆ ‰KŠ21 xmax with xmax ˆ maxt$0 x…t†; corresponds to the maximum of the quasi-static response constructed with the mean ®nite elemnt model. The dynamic magni®cation factor b n of the mean model (associated with the elastic energy) is de®ned by b n ˆ max r n …t†: t$0

3.2. Introduction of the mean reduced matrix model The mean reduced matrix model is constructed using the mode-superposition method. We consider the generalized eigenvalue problem ‰KŠw ˆ l‰MŠw associated with the mean ®nite element model. Since [K] is a positive-de®nite matrix, we have 0 , l1 # l2 # ¼ and the associated eigenvectors {w1 ; w2 ; ¼} are such that k‰MŠwa ; wbpl ˆ ma dab and k‰KŠwa ; wb l ˆ ma v2a dab ; in which va ˆ l a is the eigenfrequency of structural mode wa whose normalization is de®ned by the generalized mass ma and where ky; xl ˆ y1 x1 1 ¼ 1 ym xm : The mean reduced matrix model related to structural modes {w1 ; ¼; wn } with n p m is then written as yn …t† ˆ ‰FŠq…t† ˆ

n X

aˆ1

q a …t†wa ;

…39†

in which ‰FŠ is the …m £ n† real matrix whose columns are constituted of structural modes {w 1 ; ¼; wn } and where q…t† ˆ …q1 …t†; ¼; qn …t†† is an n real vector of the generalized coordinates which is the unique solution of the time evolution mean reduced matrix problem  1 ‰DŠq…t† _ 1 ‰KŠq…t† ˆ F…t†; ‰MŠq…t†

t $ 0;

…40†

with the initial conditions q…0† ˆ q0 ;

_ q…0† ˆ q1 ;

…41†

in which the generalized force F(t) is an n real vector such that F…t† ˆ ‰FŠ T x…t†:

…42†

The generalized mass, damping and stiffness matrices ‰MŠ; ‰DŠ and ‰KŠ are positive-de®nite symmetric …n £ n† real matrices such that ‰MŠab ˆ ma dab ; ‰DŠab ˆ k‰DŠwb ; wa l and ‰KŠab ˆ ma v2a dab : The initial conditions q0 and q1 are such that q0 ˆ ‰MŠ21 ‰FŠ T ‰MŠy0 and q1 ˆ ‰MŠ21 ‰FŠ T ‰MŠy1 : It is assumed that n is such that yn …t†; given by the mean reduced matrix model de®ned by Eqs. (39)±(42), is a good approximation of y…t† given by the mean ®nite element model de®ned by Eqs. (37)±(38). The response ratio r n …t† at time t, associated with the elastic energy of transient response yn …t† of the mean model and calculated with the mean reduced matrix model, is de®ned by r n …t† ˆ k‰KŠy…t†; y…t†l1=2 k‰KŠyS ; yS l21=2 ˆ k‰KŠq…t†; q…t†l1=2 k‰KŠyS ; yS l21=2 ;

…43†

367

…44†

For n ˆ m; r m …t† is the response ratio associated with transient response {y…t†; t $ 0} of the mean ®nite element model (see Eqs. (37) and (38) and is written as r m …t† ˆ k‰KŠy…t†; y…t†l1=2 k‰KŠyS ; yS l21=2 :

…45†

The convergence with respect to n can be analyzed in studying response ratio r n …t† at time t and dynamic magni®cation factor bn : 3.3. Construction of a nonparametric model of random uncertainties for the reduced matrix model In this section, we introduce the principle of construction of a nonparametric model of random uncertainties, the available information being constituted of the mean reduced matrix model of the structure. It should be noted that the mean ®nite element model de®ned by Eqs. (37) and (38) is not able to predict the transient response due to impulsive load whose energy is distributed over a very broad frequency band, i.e. over the low-, medium- and highfrequency ranges (for instance, if there is energy in the medium-frequency range, more advanced probabilistic mechanical models such as the fuzzy structure theory have to be used to take into account the role played by the structural complexity [14]); the most that this kind of deterministic mean ®nite element model is able to predict is the transient response due to impulsive loads whose energy is mainly distributed over a broad low-frequency range for which the mean reduced matrix model de®ned by Eqs. (39)±(42) is suitable and allows the transient response to be predicted with good accuracy. This means that the mean ®nite element model does not constitute available information for constructing the nonparametric model of random uncertainties. However, the mean reduced matrix model de®ned by Eqs. (39)±(42) (with n not too large) does constitute the available information for constructing the transient response of the mean model, then constructing the probability model of random uncertainties. This probabilistic model is a nonparametric model of random uncertainties because the sources of random uncertainties in the mechanical model which are due to uncertain mechanical parameters such as geometrical parameters, boundary conditions, junction stiffness, mass density, Young's modulus, etc., are not directly modeled by random variables or stochastic ®elds. These random uncertain geometrical and mechanical parameters mean that the generalized mass, damping and stiffness matrices of the reduced matrix model are random matrices. The nonparametric model of random uncertainties that is proposed consists in introducing a direct construction

368

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

of a probabilistic model of these random generalized matrices. This random reduced matrix model associated with the mean reduced matrix model is then written as Yn …t† ˆ ‰FŠQ…t† ˆ

n X

aˆ1

Q a …t†wa ;

…46†

in which Q…t† ˆ …Q1 …t†; ¼; Qn …t†† and where {Q…t†; t $ 0} is an R n-valued stochatic process such that  1 ‰DŠQ…t† _ 1 ‰KŠQ…t† ˆ F…t†; ‰MŠQ…t† Q…0† ˆ q0 ;

t $ 0;

_ Q…0† ˆ q1 ;

…47† …48†

in which generalized force F(t) is the R n-valued vector de®ned by Eq. (42), where initial conditions q0 and q1 have been previously de®ned and where [M], [D] and [K] are the random generalized mass, damping and stiffness matrices with values in space M1 n …R†: The basic available information is the mean reduced matrix model which is constituted of mean generalized mass, damping and stiffness matrices ‰MŠ; ‰DŠ and ‰KŠ de®ned in Section 3.2 and which belong to M1 n …R†: Random generalized mass, damping and stiffness matrices [M], [D] and [K] are second-order random variables with values in M1 n …R† such that E{‰MŠ} ˆ ‰MŠ;

E{‰DŠ} ˆ ‰DŠ;

E{‰KŠ} ˆ ‰KŠ: …49†

In addition, in order to obtain a consistent probabilistic model and in particular, to obtain convergence properties of stochastic transient response {Yn …t†; t $ 0} when dimension n approaches in®nity, we need to introduce information relative to the existence of moments of random variables [M] 21, [D] 21 and [K] 21 (such as second-order moments). It should be noted that since random matrices [M], [D] and [K] are almost surely positive de®nite, the inverse matrices exist almost surely, but the existence of second-order moments does not follow. We therefore introduce the following constraints: E{i‰MŠ21 i2F } , 11; 2

E{i‰KŠ21 iF } , 11:

E{i‰DŠi21 i2F } , 11;

…50†

We then have to construct a probability model for symmetric positive-de®nite real random matrices [M], [D] and [K] with the available information de®ned by Eqs. (49) and (50). This construction is performed using the results presented in Section 2. 3.4. Nonparametric probability model of the reduced matrix model In this section we complete the construction of the probability model introduced in Section 3.3 using the developments of Section 2. Let n0 $ 1 be a ®xed integer and n $ n0 : We apply the results of Section 2 to the set of positive-de®nite symmetric real random matrices {[M], [D], [K]} de®ned in Section 3.3, for which the available

information is described by Eqs. (49) and (50). As indicated in Section 2.3, the levels of dispersion of random matrices [M], [D] and [K] are controlled by parameters d M, d D and d K, respectively, which are independent of n and are chosen such that (see Eq. (29)) s n0 1 1 : …51† 0 , dM ; dD ; dK , n0 1 5 Parameters l M, l D and l K are de®ned by Eq. (26)

lM ˆ `M …n†;

lD ˆ `D …n†;

lK ˆ `K …n†;

…52†

in which `M …n†; `D …n† and `K …n† are given by Eq. (27) `M …n† ˆ

1 2 d2M 1 1 d2M n 1 ; 2d2M 2d2M

…53†

`D …n† ˆ

1 2 d2D 1 1 d2D n 1 ; 2d2D 2d2D

…54†

`K …n† ˆ

1 2 d2K 1 1 d2K n1 : 2 2dK 2d2K

…55†

From Section 2.6, we deduce that random matrices [M], [D] and [K] are independent random variables with values in M1 n …R† and their probability density functions p[M]([M]), p[D]([D]) and p[K]([K]) with respect to the measures (volume ~ dD ~ and dK ~ on MSn …R† are given by Eqs. (11) elements) dM; and (13), and their characteristic functions by Eq. (12). 3.5. Construction of the stochastic transient response For ®xed positive integer n $ n0 ; we have to construct stochastic processes {Yn …t†; t $ 0} de®ned by Eqs. (46)± (48), stochastic process {Rn …t†; t $ 0} de®ned (see Eq. (43)) by Rn …t† ˆ k‰KŠYn …t†; Yn …t†l1=2 k‰KŠyS ; yS l21=2 ˆ k‰KŠQ…t†; Q…t†l1=2 k‰KŠyS ; yS l21=2 ;

…56†

and ®nally, random variable Bn de®ned (see Eq. (44)) by Bn ˆ max Rn …t†:

…57†

t$0

Below, we present a formulation which is adapted to Monte Carlo numerical simulation. For given matrices [M], [D], [K] in M1 n …R† let t 7 ! qF;q0 ;q1 (t;[M],[D],[K]) be the solution from R 1 into R n of the deterministic second-order differential equation ‰MŠq F;q0 ;q1 …t† 1 ‰DŠq_ F;q0 ;q1 …t† 1 ‰KŠqF;q0 ;q1 …t† ˆ F…t†; t $ 0;

…58†

with the initial conditions qF;q0 ;q1 …0† ˆ q0 ;

q_ F;q0 ;q1 …0† ˆ q1 :

…59†

We deduce that stochastic process {Q…t†; t $ 0} which is the solution of the stochastic dynamical problem de®ned by

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

369

Eqs. (47) and (48), can be written as Q…t† ˆ qF;q0 ;q1 …t; ‰MŠ; ‰DŠ; ‰KŠ†:

…60†

It should be noted that if q0 ˆ q1 ˆ 0; then Q(t) can usually be written as Q…t† ˆ

Zt 0

‰h…t 2 t†ŠF…t† dt;

…61†

in which t 7 ! ‰h…t†Š is the matrix-valued impulse response function of the linear ®lter associated with second-order differential Eq. (58). If q 7 ! ‰u…q†Š is a mapping from R n into the set Mn1 ;n2 …R† of all the …n1 £ n2 † real matrices, we have E{‰u…Q…t††Š} ˆ

Z

Z M1 n …R†

Z M1 n …R†

M1 n …R†

Fig. 1. Geometry of the mean structure.

‰u…qF;q0 ;q1 …t; ‰MŠ; ‰DŠ; ‰KŠ††Š ‰KŠ ˆ

~ dD ~ dK: ~ £ p‰MŠ …‰MŠ† £ p‰DŠ …‰DŠ† £ p‰KŠ …‰KŠ† dM …62† For instance, Rn(t) de®ned by Eq. (56) can be written as Rn …t† ˆ ‰u…Q…t††Š with n1 ˆ n2 ˆ 1: Calculation of the stochastic transient response of the dynamical system with random uncertainties requires the numerical construction of mapping t 7 ! qF;q0 ;q1 …t; ‰MŠ; ‰DŠ; ‰KŠ†† as the solution of the deterministic Eqs. (58) and (59). Since matrices [M], [D] and [K] are full matrices (not diagonal) as samplings of random matrices [M], [D] and [K], any integral representation of Q(t) (for instance such as Eq. (61) when q0 ˆ q1 ˆ 0) is not really ef®cient but second-order differential Eq. (58) is solved directly using an unconditionally stable implicit step-bystep integration method (such as the Newmark integration scheme [16]) with initial conditions de®ned by Eq. (59). In addition, we have to calculate multiple integrals in a higher dimension (see Eq. (62)) for which a well suited method consists in using a Monte Carlo calculation with or without variance reduction procedures [28±31]. This method is very ef®cient if there is a Monte Carlo simulation procedure for random matrices [M], [D] and [K] which is the case of the method presented in Sections 2.4 and 2.5. It should be noted that for many applications, integer n is suf®ciently high that l M, l D and l K can be considered as positive integers without introducing any signi®cant limitation in the model. Applying Eqs. (31) and (32) to random matrices [M], [D] and [K] yields ‰MŠ ˆ

‰DŠ ˆ

M …n† 1 mX …‰LM ŠT Xj †…‰LM ŠT Xj †T ; mM …n† jˆ1

D …n† 1 mX …‰LD ŠT Yj †…‰LD ŠT Yj †T ; mD …n† jˆ1

…63†

…64†

K …n† 1 mX …‰LK ŠT Zj †…‰LK ŠT Zj †T ; mK …n† jˆ1

in which

! n11 mM …n† ˆ Fix ; d2M ! n11 mK …n† ˆ Fix ; d2K

! n11 mD …n† ˆ Fix ; d2D

…65†

…66†

where Fix(x) is equal to x when x is an integer and Fix(x) rounds down x 1 1 to the nearest integer when x is not an integer. In Eqs. (63)±(65), ‰LM Š; ‰LD Š and ‰LK Š are upper triangular matrices in Mn …R† corresponding to the Cholesky factorization of symmetric positive-de®nite matrices ‰MŠ; ‰DŠ and ‰KŠ ‰MŠ ˆ ‰LM ŠT ‰LM Š; ‰KŠ ˆ ‰LK ŠT ‰LK Š:

‰DŠ ˆ ‰LD ŠT ‰LD Š;

…67†

The set of all the components of vectors X1 ; ¼; XmM …n† ; Y1 ; ¼; YmD …n† and Z1 ; ¼; ZmK …n† with values in Rn is constituted of mM …n† £ n 1 mD …n† £ n 1 mK …n† £ n independent random variables, each of which is a real-valued secondorder normalized Gaussian random variable (zero mean value and unit variance). 4. Simple example The mean structure is constituted of a rectangular homogeneous and isotropic plate located in the plane (Ox, Oy) of a Cartesian coordinate system (Oxyz), in bending mode (the outplane displacement is z), with constant thickness 4 £ 10 24 m, width 0.40 m, length 0.50 m, mass density 7800 kg/m 3, Young's modulus 2.1 £ 10 11 N/m 2 and Poisson ratio 0.29. This plate is simply supported on three edges and

370

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

Fig. 3. Graph of the modulus of the Fourier transform of impulse function. Fig. 2. Graph of impulse function t 7 ! p…t†:

free on the fourth edge corresponding to y ˆ 0 (see Fig. 1). The mean ®nite element model is constituted of a regular rectangular mesh with a constant step size of 0.01 m in x and y (41 nodes in the width, 51 nodes in the length). Consequently, all the ®nite elements are the same and each one is a four-nodes square plate element. There are 2000 ®nite elements and m ˆ 6009 degrees of freedom (z-translations and x- and y-rotations). The eigenfrequencies va ˆ 2pf a of the mean ®nite element model are f 1 ˆ 1:94; f 2 ˆ 10:28; f 3 ˆ 15:47;¼,f 25 ˆ 167:44;¼,f 80 ˆ 527:29 Hz: The initial conditions de®ned by Eq. (38) are such that y0 ˆ y1 ˆ 0: External impulsive load vector x(t) de®ned in Section 3.1 is written as x…t† ˆ p…t†g: Spatial part g ˆ …g1 ; ¼; gm † [ Rm is independent of time t and is such that gj ˆ 0 for all j in {1,¼,m} except for the nine DOFs corresponding to the nodes whose (x,y) coordinates are (0.30,0.25), (0.30,0.26), (0.30,0.27), (0.31,0.25), (0.31,0.26), (0.31,0.27), (0.32,0.25), (0.32,0.26) and (0.32,0.27), for which gj ˆ 1 (see Fig. 1). Impulse function t 7 ! p…t† is the rectangular impulse function de®ned by p…t† ˆ pmax 1‰0;t1 Š …t† with pmax ˆ 1 and t1 ˆ 0:4 s; in which 1‰0;t1 Š …t† ˆ 1 if t [ ‰0; t1 Š and 1‰0;t1 Š …t† ˆ 0 if t Ó ‰0; t1 Š: Fig. 2 shows the graph of impulse function t 7 ! p…t† and Fig. 3 shows the graph of the modulus of its Fourier transform. It can be seen in Fig. 3 that the main part of the energy of impulse function p is distributed over the [215,15] Hz frequency band in which there are three structural modes of the mean model. Consequently, only the ®rst structural modes signi®cantly contribute in the dynamical response and the structural modes whose eigenfrequencies are greater that 168 Hz (n . 25) contribute mainly in the quasi-static response and not in the dynamical response. The mean damping matrix is ‰DŠ ˆ 2jV ref ‰MŠ with V ref ˆ 2p £ 10 and j ˆ 0:001: The transient response of the mean ®nite element model is calculated by solving the time evolution problem de®ned by Eqs. (37) and (38) with y…0† ˆ y_ …0† ˆ 0 using an unconditionaly stable implicit step-bystep integration method (Newmark integration scheme) with a time-step size Dt ˆ 1=1300 s: The maximum 0:5k‰KŠy S ; yS l of the quasi-static elastic energy of the

mean ®nite element model is equal to 2.552. The dynamic magni®cation factor bm ˆ maxt$0 r m …t†; in which r m …t† is calculated by Eq. (45), is equal to 1.595. For ®xed positive integer n, the mean reduced matrix problem de®ned by Eqs. _ (40) and (41) with q…0† ˆ q…0† ˆ 0 is solved with the same Newmark integration scheme with the same time step size Dt ˆ 1=1300 s: Fig. 4 shows the convergence of the dynamic magni®cation factor bn of the mean model, de®ned by Eq. (44), as dimension n of the mean reduced matrix model increases. From Fig. 4, it can be deduced that the transient response of the mean model is reasonably converged when n ˆ 80 for which b n ˆ 1:487: Fig. 5 shows the graph of function t 7 ! r n …t† for n ˆ 80: Concerning the structure with random uncertainties, we choose n0 ˆ 1 which allows the convergence analysis with respect to dimension n of the reduced matrix model with random uncertainties to be performed for n $ n0 ˆ 1 The dispersions of the generalized mass, damping and stiffness random matrices of the reduced matrix model with random uncertainties, are controlled by parameters d M, d D and d K

Fig. 4. Graph of the convergence of dynamic magni®cation factor bn with respect to dimension n of the mean reduced matrix model.

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

Fig. 5. Graph of function t 7 ! r n …t† for n ˆ 80 corresponding to the response ratio for the mean reduced matrix model.

introduced in Section 3.4, which have to verify the constraints de®ned by Eq. (51) 0 , dM ; dD ; dK , 0:57:

…68†

The numerical simulations presented below correspond to the values

dM ˆ 0:1;

dD ˆ 0:1;

dK ˆ 0:1;

…69†

which verify Eq. (68). We are interested in the random response ratio Rn(t) de®ned by Eq. (56) and the random dynamic magni®cation factor Bn de®ned by Eq. (57). The transient response of the structure with random uncertainties is calculated using the Monte Carlo numerical simulation method. For given generalized mass, damping and stiffness matrices, the time evolution problem de®ned by Eqs. (58)

Fig. 6. Graph of the convergence of E{Bn} (mathematical expectation of the random dynamic magni®cation factor) with respect to dimension n of the random reduced matrix model.

371

Fig. 7. Graph of the convergence of Bn,max(u) maximum of the random dynamic magni®cation factor) with respect to dimension n of the random reduced matrix model.

and (59) is solved with the same Newmark integration scheme and with the same time step size Dt ˆ 1=1300 s: The Monte Carlo numerical simulation is carried out with nS ˆ 1000 samples, denoted as u1 ; ¼; unS ; for which the samples t 7 ! Rn …t; u1 †; ¼; t 7 ! Rn …t; unS †; are numerically calculated. For t ®xed, the mean value of random variable Rn(t) is estimated by E{Rn …t†} .

nS 1 X R …t; uj †: nS jˆ1 n

…70†

The samples of random variable Bn are such that Bn …uj † ˆ max Rn …t; uj †: t$0

…71†

Fig. 8. Transient responses t 7 ! rn …t† (lower irregular thin solid line), t 7 ! E{Rn …t†} (lower smooth thick solid line) and t 7 ! Rn;max …t; u† (upper irregular thin solid line) for n ˆ 80 and nS ˆ 1000:

372

C. Soize / Probabilistic Engineering Mechanics 16 (2001) 363±372

The mean value of random variable Bn is estimated by E{Bn } .

nS 1 X B …u †: nS jˆ1 n j

…72†

Finally, we introduce the function t 7 ! Rn;max …t; u† and the real number Bn,max(u ) de®ned by Rn;max …t; u† ˆ max Rn …t; uj †;

…73†

Bn;max …u† ˆ max Bn …uj †;

…74†

jˆ1;¼;ns

jˆ1;¼;nS

in which u ˆ …u1 ; ¼; unS †: Fig. 6 shows function n 7 ! E{Bn } calculated by Eq. (72) and Fig. 7 shows the function n 7 ! Bn;max …u† calculated by Eq. (74). Figs. 6 and 7 show that a reasonable convergence is obtained for n ˆ 80 and for n ˆ 80; the value of Bn,max(u) is 1.676. This value has to be compared to the value for the mean model, which is 1.487. Fig. 8 is relative to n ˆ 80 and shows three curves: the lower irregular thin solid line corresponds to the graph of function t 7 ! r n …t†; the lower smooth thick solid line to the graph of function t 7 ! E{Rn …t†} calculated by Eq. (70) and the upper irregular thin solid line to the graph of function t 7 ! Rn;max …t; u† de®ned by Eq. (73). 5. Conclusions We have presented a new approach allowing the random uncertainties to be modeled by a nonparametric model for prediction of transient responses to impulsive loads in linear structural dynamics. The information used does not require the description of the local parameters of the mechanical model. The probability model is deduced from the use of the entropy optimization principle whose available information is constituted of the fundamental algebraic properties related to the generalized mass, damping and stiffness matrices which have to be positive-de®nite symmetric matrices, and the knowledge of these matrices for the mean reduced matrix model. An explicit construction and representation of the probability model have been obtained and are very well suited to algebraic calculus and to Monte Carlo numerical simulation in order to compute the transient responses of structures submitted to impulsive loads. The numerical analysis related to the convergence of the stochastic solution with respect to the dimension of the random reduced matrix model has been carried out. References [1] Chen PC, Soroka WW. Multi-degree dynamic response of a system with statistical properties. J Sound Vibrat 1973;37(4):547±56.

[2] Soong TT. Random differential equations in science and engineering. New York: Academic Press, 1973. [3] Haug EJ, Choi KK, Komkov V. Design sensitivity analysis of structural systems. San Diego: Academic Press, 1986. [4] Ibrahim RA. Structural dynamics with parameter uncertainties. Appl Mech Rev 1987;40(3):309±28. [5] Shinozuka M. Structural response variability. ASCE Journal of Engineering Mechanics 1987;113(6):825±42. [6] Spanos PD, Zeldin BA. Galerkin sampling method for stochastic mechanics problems. ASCE J Engng Mech 1994;120(5):1091±106. [7] Lin YK, Cai GQ. Probabilistic structural dynamics. New York: McGraw-Hill, 1995. [8] Vanmarcke E, Grigoriu M. Stochastic ®nite element analysis of simple beams. ASCE J Engng Mech 1983;109(5):1203±14. [9] Shinozuka M, Deodatis G. Response variability of stochastic ®nite element systems. ASCE J Engng Mech 1988;114(3):499±519. [10] Spanos PD, Ghanem RG. Stochastic ®nite element expansion for random media. ASCE J Engng Mech 1989;115(5):1035±53. [11] Ghanem RG, Spanos PD. Stochastic ®nite elements: a spectral approach. New York: Springer, 1991. [12] Kleiber M, Tran DH, Hien TD. The stochastic ®nite element method. Chichester: Wiley, 1992. [13] Ditlevsen O, Tarp-Johansen NJ. Choice of input ®elds in stochastic ®nite elements. Probab Engng Mech 1998;14(1±2):63±72. [14] Ohayon R, Soize C. Structural acoustics and vibration. San Diego: Academic Press, 1998. [15] Clough RW, Penzien J. Dynamics of structures. New York: McGrawHill, 1975. [16] Bathe KJ, Wilson EL. Numerical methods in ®nite element analysis. New York: Prentice Hall, 1976. [17] Zienkiewicz OC, Taylor RL. 4. The ®nite element method, vol. 1. New York: McGraw-Hill, 1989 see also vol. 2, 1991. [18] Argyris J, Mlejnek HP. Dynamics of structures. Amsterdam: NorthHolland, 1991. [19] Soize C. ISSN 1159-09747. A nonparametric model of random uncertainties in linear structural dynamics. Prog Stochast Struct Dynam 1999;152:109±38. [20] Soize C. A nonparametric model of random uncertainties for reduced matrix models in structural dynamics. Probab Engng Mech 2000;15(3):277±94. [21] Mehta ML. Random matrices, revised and enlarged second edition. London: Academic Press, 1991. [22] Shannon CE. A mathematical theory of communication. Bell Syst Tech J 1948;27:379±423 p. 623±59. [23] Jaynes ET. Information theory and statistical mechanics. Phys Rev 1957;106(4):620±30. [24] Jaynes ET. Information theory and statistical mechanics. Phys Rev 1957;108(2):171±90. [25] Kapur JN, Kesavan HK. Entropy optimization principles with applications. San Diego: Academic Press, 1992. [26] Anderson TW. Introduction to multivariate statistical analysis. Chichester: Wiley, 1958. [27] Soize C. Maximum entropy approach for modeling random uncertainties in transient elastodynamics. J Acoust Soc Am 2001;109(5). [28] Cochran WG. Sampling techniques. New York: Wiley, 1977. [29] Rubinstein RY. Simulation and the Monte Carlo method. New York: Wiley, 1981. [30] Kalos MH, Whitlock PA. Monte Carlo methods, Basics, vol. 1. New York: Wiley, 1986. [31] Bratley P, Fox BL, Schrage EL. A guide to simulation. 2nd ed.. New York: Springer, 1987.