Nonlinear systems driven by polynomials of filtered Poisson processes

Probabilistic Engineering Mechanics 14 (1999) 195±203

Nonlinear systems driven by polynomials of ®ltered Poisson processes Federico Waisman a, Mircea Grigoriu b,* a

b

EQE International, 44 Montgomery St. Suite 3200, San Francisco, CA 94104, USA School of Civil and Environmental Engineering, Cornell University, Ithaca, NY 14853, USA

Abstract The perturbation method is applied to determine approximately the mean, variance, skewness and kurtosis of the transient and stationary response of nonlinear systems driven by polynomials of ®ltered Poisson processes. The analysis is based on the classical perturbation method, the Itoà differentiation formula, and properties of the response of linear systems subjected to polynomials of ®ltered Poisson processes. Two examples are presented to demostrate the ef®ciency and accuracy of this approximate analysis. q 1998 Elsevier Science Ltd. All rights reserved. Keywords: Nonlinear systems; Duf®ng oscilator; Filtered Poisson processes; Polynomials; Perturbation method; Itoà formula; Moment equations; Stochastic differential equations; Probability theory; Poisson white noise

1. Introduction The perturbation method is frequently used in random vibrations to determine approximately moments of the response of weakly nonlinear systems subjected to Gaussian noise. The method approximates the solution of a nonlinear differential equation by a truncated power series with terms of decreasing order of magnitude that satisfy equations with an identical linear differential operator but different inputs. This work proposes a general method for calculating approximately the transient and stationary moments of any order for the solution of nonlinear systems driven by polynomials of ®ltered Poisson processes. The method is based on the perturbation technique, the Itoà differentiation formula, and an algorithm for ®nding exactly moments of any order of the solution of linear systems driven by polynomials of ®ltered Poisson processes [1]. The method is applied to ®nd response moments of a Duf®ng oscillator subject to (i) a Poisson white noise process, and (ii) a quadratic form of a ®ltered Poisson process. A comparison of results with the Monte Carlo simulation method is also presented. 2. Problem de®nition Let {X(t) [ R m}, t $ 0} be the solution of the nonlinear system _ 1 b‰X…t† 1 1f …X…t††Š ˆ p…Z…t†; t†; t $ 0 (1)  1 aX…t† X…t† * Corresponding author.

where f is a nonlinear differentiable function of X(t), a denotes the damping matrix, b is approximately the stiffness matrix of the system for small oscillations, 1 . 0 is a small parameter, and p(Z(t),t) is a vector of m polynomials of Z(t) with components pk …Z…t†; t† ˆ

X dk

gk ; dk …t†

n2 Y rˆ1

Zr …t†dk;r

(2)

depending on the non-negative integers dk ˆ {dk;1 ; dk;2 ; ¼; dk;n2 } and the time dependent coef®cients {g k,d}(t)}. The n2-dimensional ®ltered Poisson process Z(t) is the solution of dZ…t† ˆ b…t†Z…t†dt 1 dY…t†

(3)

in which b(t) are time dependent coef®cients and Y…t† ˆ

Zt Z 0

Y

yM…ds; dy†

(4)

is a compound Poisson process of intensity n . 0 characterized by the Poisson random measure M [2±4]. The random measure M of an arbitrary set of [0,1) £ Y is a Poisson random variable counting the number of points in this set. The mean of this Poisson random variable for (t,t 1 dt) £ (y,y 1 dy) is EM…dt; dy† ˆ n dt dF…y†

(5)

where F denotes the distribution of a Y-valued random variable. The Poisson random measure M(ds,dy) can be used to represent the homogeneous compound Poisson

0266-8920/99/$ - see front matter q 1998 Elsevier Science Ltd. All rights reserved. PII: S 0266-892 0(98)00031-9

196

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

process 8 0; > > < X Y…t† ˆ N…t† > Yk ; > :

Let N…t† ˆ 0 (6)

N…t† . 0

kˆ1

where {N(t),t $ 0} is a homogeneous Poisson counting process of intensity n . 0 and the independent Y-valued random variables {Yk} of distribution F arriving at the random times {G k} of {N(t),t $ 0}. The interarrival times {G k ± G k21} are independent exponential variables with mean 1/n . 2.1. Perturbation method The solution X(t) of Eq. (1) can be approximated by the power series in 1 X…t† . X0 …t† 1 1X 1 …t† 1 12 X 2 …t† 1 ¼

(7)

involving the unknown functions {Xi(t) [ Rm}, i $ 0, t $ 0 }. This representation of X(t) and Eq. (1) give …X 0 …t† 1 eX 1 …t† 1 e2 X 2 …t† 1 ¼† 1 a…X_ 0 …t† 1 eX_ 1 …t† 1 e2 X_ 2 …t† 1 ¼† 1 b‰…X0 …t† 1 eX1 …t† 1 e2 X 2 …t† 1 ¼† 1 ef …X 0 …t† 1 eX1 …t† 1 e2 X 2 …t† 1 ¼†Š ˆ p…Z…t†; t†

(8)

To identify the order of magnitude of 1 f(X(t)) consider the Taylor series of this function about X0(t) f …X0 …t† 1 eX 1 …t† 1 e2 X2 …t† 1 ¼† . f …X 0 …t†† 1 1

m X 2f …X 0 …t†† …eX 1 …t† 1 e2 X2 …t† 1 ¼† 2 X …t† i iˆ1 m 1 X 2f …X 0 …t†† …eX1 …t† 1 e2 X 2 …t† 1 ¼†2 1 ¼ (9) 2 i;jˆ1 2Xi …t†Xj …t†

The solution X0(t) is called the unperturbed solution of Eq. (1) and corresponds to 1 ˆ 0. Replacing f(X(t)) in Eq. (8) by its expansion given by Eq. (9) and equating terms of the same power of 1 , the in®nite set of linear ordinary stochastic differential equations X 0 …t† 1 aX_ 0 …t† 1 bX 0 …t† ˆ p…Z…t†; t† X 1 …t† 1 aX_ 1 …t† 1 bX 1 …t† ˆ 2bf …X 0 …t†† " # m X 2f …X 0 …t†† _  X 1 …t† X2 …t† 1 aX 2 …t† 1 bX 2 …t† ˆ 2b 2Xi …t† iˆ1 (10) results. These equations have the same linear differential operator L ˆ (d 2/dt 2) 1 a (d/dt) 1 b but different inputs, and can be solved sequentially.

X…p† …t† ˆ

n X kˆ0

1k Xk …t†

(11)

be the solution of order p of Eq. (1) depending on the ®rst p solutions of Eq. (10). Denote by Pp(t) ˆ {Xp(t) 0 , Ç 1 …t† 0 , X0 …t† 0 , XÇ0 …t† 0 } 0 the …2p 1 1†mÇ }p(t) 0 ,¼, X1 …t† 0 , X X dimensional vector of the solutions Xk(t), k ˆ 0,¼,p, and their derivatives. In this paper the prime denotes matrix transposition. De®ne the right-hand side of Eq. (10) as p*(Z(t), Pp21 …t†† ˆ fp…Z…t†; t† 0 ; …2bf f …X0 …t††† 0 ¼; g 0 . The n ˆ (2p 1 1)m 1 n2-dimensional augmented state vector {Pp(t) 0 ,Z(t) 0 } 0 is the solution of the stochastic differential equation (Eqs. (3) and (10)) " # " # a…t†Pp …t† 1 p p …Z…t†; Pp21 …t†† Pp …t† ˆ dt d Z…t† b…t†Z…t† " # Z 0 1 M…dt; dy† (12) y where a(t) is a (n1,n1) matrix, n1 ˆ (2p 1 1)m, and M is the Poisson random measure de®ned in Eqs. (4) and (5). The objective is to ®nd moments of any order for the approximate solution X (p)(t). If f(X(t)) is a polynomial function, Eq. (10) is a linear system of differential equations driven by polynomials of ®ltered Poisson processes. The Itoà differentiation formula can be applied to derive equations for the moments of {Pp(t) 0 ,Z(t) 0 } 0 based on an algorithm developed previously [1]. This algorithm ®nds exactly transient and stationary moments of any order of linear systems driven by polynomials of ®ltered Poisson processes. The convergence of the series in Eq. (7) is essential for using the perturbation method. For the deterministic case, the convergence of the series in Eq. (7) to X(t) has been proven by Poincare provided that the perturbation parameter 1 is small [5±7]. When the perturbation method is applied for a particular ®nite value 1 , its result is expected to agree well with numerical simulations provided that the system is perturbed in the neighborhood of a stable equilibrium point of the unperturbed system, 1 ˆ 0 [8]. The dynamic stability of nonlinear systems cannot be obtained using this perturbation method. Another approach, also known as the Melnikov method, is commonly used in obtaining the dynamic stability of deterministic nonlinear systems [6, 7, 9]. Convergence in the stochastic case can be achieved by requiring, for each, that the series in Eq. (7) converges to X(t) almost surely, i.e. the set N of realizations of Z for which the time series does not converge has probability zero. 3. Method of analysis The Itoà formula provides a useful mathematical framework to obtain differential equations for functions of the solution of a stochastic differential equation. In this work the Itoà differentiation formula is used to obtain moment

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

197

equations of any order for the response of linear systems subjected to polynomials of ®ltered Poisson and Gaussian processes. The moment equations are exact and form a closed set of linear differential equations that can be solved by any conventional numerical method. Moreover, stationary moments can be obtained by solving linear algebraic equations. Differential equations for the characteristic function, Lyapunov exponents and the Fokker±Planck equation are examples of the wide use of this mathematical tool.

with the differential form

3.1. The Itoà differentiation formula

Any vector Q (t) satisfying a stochastic differential equation of the type of Eqs. (13) and (14), and the Itoà differentiation formulas of Eqs. (15) and (16) can be used to develop differential equations for differentiable functions of this vector. This observation is used to develop differential equations for the moments of the augmented state vector {Pp(t) 0 ,Z(t) 0 } 0 in Eq. (12).

Consider the solution {Q (t) [ R n,t $ 0} of the stochastic differential equation dQ…t† ˆ m…Q…t†; t†dt 1

Z Y

c…Q…t†; t; y†M…dt; dy†; t $ 0 (13)

where m…u; t†; and c(u ,t,y) are (n,1) and (n,q) matrices, respectively, and M is the Poisson random measure of Eqs. (4) and (5). It is usually convenient to consider the integral form of Eq. (13)

Q…t† ˆ Q…0† 1 1

Zt 0

Zt Z 0

Y

m…Q…s†; s†ds c…Q…s†; s; y†M…ds; dy†; t $ 0

1

0

Y

1

Z Y

‰g…Q…t† 1 c…Q…t†; t; y†; t†

2 g…Q…t†; t†ŠM…dt; dy†:

(16)

3.2. Moment equations Let Q (t) ˆ {Pp(t) 0 ,Z(t) 0 } 0 be the solution of Eq. (12) and consider the differentiable function n

depending on the components of this process and positive integers {qi}. The Itoà differentiation formula gives (16) dg…Q…t†† ˆ

(n "n 1 1 X X kˆ1

lˆ1

# ak;l …t†Ql …t† 1 pk p …Z…t†; Pp21 …t††

 qk Q21 k g…Q…t†† n2 X

n X

1

kˆn1 1 1 lˆ1

1

2

Z Y n Y iˆ1

) bk2n1 ;l …t†Qn1 1l …t†qk Q21 k g…Q…t††

dt

2 n1 n Y Y 4 Qi …t†qi …Qi …t† 1 yi2n1 †qi iˆ1

iˆn1 1 1

#

Qi …t†qi M…dt; dy†

(18)

Let f(X(t)) be a polynomial function. Thus, it can be also written as f …X…t†† ˆ

X ek

(15)

(17)

iˆ1

(14)

‰g…Q…s† 1 c…Q…s†; s; y†; s†

2 g…Q…s†; s†ŠM…ds; dy†

dg…Q…t†; t† ˆ

) n X 2g…Q…t†; t† 2g…Q…t†; t† 1 dt mk …Q…t†; t† 2t 2Qk …t† kˆ1

g…Q…t†† ˆ P Qi …t†qi

The ®rst integral can be de®ned in the Riemann sense for every sample of Q(t) because the process has continuous samples with probability one and the function m(x,t) is a smooth function. The second integral is de®ned in the Itoà sense [2±4, 10]. Existence and uniqueness of the solution of Eqs. (13) and (14) is guaranteed if these three conditions are satis®ed [2, 4, 11]: (i) c(x,t,y) is a continuous function of x and Y and left continuous with respect to t, (ii) the functions m(x,t) and c(x,t,y) are Lipschitz continuous, and (iii) the function c(x,t,y) is bounded in the time-space rectangle of [0,1) £ Y for every x. The conditions (i) and (ii) guarantee the Rt R existence with probability one of the integral 0 Y c…X…s†; s; y†M…ds; dy†. The last condition assures the uniqueness of the solution X(t). Let g(Q (t),t) be a differentiable function and Q (t) de®ned by Eqs. (13) and (14). The process g(V (t),t) is the solution of the stochastic integral equation [2±4, 10, 11] Zt  2g…Q…s†; s† g…Q…t†; t† ˆ g…Q…0†; 0† 1 2s 0 ) n X 2g…Q…s†; s† ds mk …Q…s†; s† 1 2Qk…s† kˆ1 Zt Z

(

hk ; ek …t†

m Y rˆ1

Xr …t†ek;r

(19)

depending on the non-negative integers ek ˆ {ek,1,ek,2,¼,ek,m} and the time dependent coef®cients {h k,e(t)}. The right-hand

198

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

side of Eq. (10) then becomes 2

1

3

n2 Y

n1 X X kˆ1 ek

Sdk gk ; dk …t† Zr …t†dk;r 7 6 7 6 rˆ1 7 6 7 6 7 6 m Y X p 6 ek;r 7 p …Z…t†; Pp21 …t†† ˆ 6 2b 7 h ; e …t† X …t† e k k k 0;r 7 6 7 6 rˆ1 7 6 5 4 .. .

qk hk ; ek m…q1 ; ¼; qk 2 1; ¼; q2n1 ; q2n1 11

1 ek;n1 11 ¼; q3n1 1 ek;n ; q3n1 11 ; ¼; qn ; t† n X

1¼ 1

n X

kˆn1 1 1 lˆn1 1 1

qk bk2n1 ;l2n1 m

 …q1 ; ¼; qn1 ; ¼; qk 2 1; ¼; qn ; t† (20) (

where Xk,r(t) is the rth element of Xk(t). By averaging Eq. (18), it is found

1n

n1 X n1 X

ˆ

kˆ1 lˆ1

1

£ my …qn1 11 2 rn1 11 ; ¼; qn 2 rn ; t† 2 m…q1 ; ¼; qn ; t†

n1 X X

qk gk ; dk m…q1 ; ¼; qk 2 1; ¼; qn1 ; qn1 1 1

1 dk;n1 11 ¼; q2n1 1 dk;n ; q2n1 11 ; ¼; qn ; t† n1 X X

qk hk ; ek m…q1 ; ¼; qk 2 1; ¼; q2n1 ; q2n1 11

kˆ1 ek

1 ek;n1 11 ¼; q3n1 1 ek;n ; q3n1 11 ; ¼; qn ; t† 1 ¼ 1

n X

n X

kˆn1 1 1 lˆn1 1 1

(

1 n dt E

n1 Y iˆ1

) qk bk2n1 ;l2n1 m…q1 ; ¼; qn1 ; ¼; qk21 ; ¼; qn ; t† dt

Qi …t†qi

Z

i 2ri  Qi …t†ri yqi2n dF…y† 2 1

n Y

qi X

Y iˆn 1 1 r ˆ0 1 i n Z Y Y iˆ1

qi ! ri !…qi 2 ri †! )

Qi …t†qi dF…y†

…21†

Q where m (q1,q2,¼,qn;t) ˆ E niˆ1 ui …t†qi denotes moments of the state vector. The average of the last term of this equation can be expressed in terms of moments of Q (t) because it involves powers of the components of the state vector. Dividing by dt, the ordinary differential equation

m_ …q1 ; q2 ¼; qn ; t† ˆ

n1 X n1 X kˆ1 lˆ1

1

r[q iˆn1 1 1

qi ! m…q1 ; ¼; qn1 ; rn1 11 ; ¼; rn ; t† ri !…qi 2 ri †! (22)

qk ak;l m…q1 ; ¼; ql 1 1; ¼; qk 2 1; ¼; qn1 ; ¼; qn ; t†

kˆ1 dk

1

n Y

)

d m…q1 ; ¼; qn ; t† (

X

qk ak;l m…q1 ; ¼; ql 1 1; ¼; qk21 ; ¼; qn1 ; ¼; qn ; t†

n1 X X kˆ1 dk

qk gk ; dk m…q1 ; ¼; qk 2 1; ¼; qn1 ; qn1 11

1 dk;n1 11 ¼; q2n1 1 dk;n ; q2n1 11 ; ¼; qn ; t†

is obtained for the moments of Q (t), where r [ q means that ru ˆ {0,1,¼,qu} for u ˆ n1 1 1,¼,n. An algorithm had been developed [1] for solving exactly Eq. (22) based on three observations:P (i) the moments n m (0,¼,0,rn111,¼,rn) of any order iˆn1 11 ri can be calculated exactly because they correspond to the ®ltered Poisson process Z(t) [12], (ii) the moments m (q1,¼,qn;t) vanish if one or more arguments qi are strictly smaller than zero, P n and (iii) the set of moment equations of order qi generated by Eq. (22) involve the moments s ˆ iˆ1 m…q1 ; ¼; ql 1 1; ¼; qk 2 1; ¼; qn1 ; ¼; qn ; t†; m…q1 ; ¼; qk 2 1; ¼; qn1 ; qn111 1 dk,n111 ¼, qn 1 dk,n;t), m…q1 ; ¼; qn1 ; ¼; qk 2 1; ¼; P qn ; t†; and m (q1,¼,qn1,rn111,¼,rn) of order s, s 2 1 1 niˆn1 11 dk;i ; s 2 1, and # s, respectively. The set of moment equations ofPorder s involves unknown moments of order s 2 1 1 niˆn1 11 dk;i . However, these unknown moments can be determined prior to solving the set of moment equations of order s from the set of P moment 1 qi 0 ˆ equations for m (q1 0 ,¼,qn1 0 ,rn111,¼,rn;t), where niˆ1 s 2 1 and {rn111,¼,rn} are arbitrary non-negative integers. These moments can be obtained exactly for every s by starting from s ˆ 1 and increasing this parameter up to the required value because the moments m (q1 0 ,¼,qn1 0 ,rn111,¼,rn;t) are known for s ˆ 1. Therefore, the moment equations given by Eq. (22) are closed and can be solved exactly. There is no need to use closure techniques for solution. Initial conditions are needed for solving Eq. (22) and these conditions need to be speci®ed. If the coef®cients of the differential equation de®ning Q (t) are time invariant and Q (t) becomes stationary as t ! 1; the moments m (q1,¼,qn;t) ˆ m (q1,¼,qn) do not depend on time so that m_ (q1,q2,¼,qn;t) ˆ 0 in Eq. (22). In this case, the moments of Q (t) are the solution of algebraic equations. Moments of any order of the approximate solution of order p of the nonlinear stochastic differential equation in Eq. (1) can be found exactly using Eq. (22),

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

as long as the system nonlinearity is of polynomial form. The algorithm to obtain moments equations for nonlinear systems with polynomial nonlinearities is presented in the next section for a Duf®ng oscillator [6].

4. Numerical examples Let {X(t) [ R,t $ 0} be the solution of the stochastic differential equation _ 1 v2 X…t†…1 1 1X…t†2 † ˆ Z…t†d;  1 2zvX…t† X…t†

t$0 (23)

where v is the linearized frequency, z is the damping ratio, 1 . 0 is a small parameter controlling the cubic nonlinearity, and Z(t) is a ®ltered Poisson process de®ned by dZ…t† ˆ bZ…t†dt 1 dY…t†

(24)

where

X 0 …t† 1 aX_ 0 …t† 1 bX0 …t† ˆ Z…t†d X 1 …t† 1 aX_ 1 …t† 1 bX1 …t† ˆ 2v X0 …t† 2

7 6 X_ 2 …t† 7 6 7 6 7 6 7 6 2 2 3 6 2v X2 …t† 2 2zvX_ 2 …t† 2 3v X0 …t† X1 …t† 7 7 6 7 6 7 6 7 6 _ …t† X 7 6 1 7 6 7 6 7 m…Q…t†† ˆ 6 6 2v2 X …t† 2 2zvX_ …t† 2 v2 X …t†3 7 7 6 1 1 0 7 6 7 6 7 6 7 6 _ 0 …t† X 7 6 7 6 7 6 7 6 2 _ 0 …t† 1 Z…t†d 7 6 X v X …t† 2 2 zv 2 0 7 6 5 4 bZ…t† (27) 4.1. Moment equations Consider the continuous and differentiable function g…Q…t†† ˆ

7 Y

Qi …t†qi

(28)

of the solution Q (t) of Eq. (26). Moments of any order of the solution Q (t) of Eq. (26) satisfy the differential equation (22)

m_ …q1 ; q2 ; q3 ; q4 ; q5 ; q6 ; q7 ; t† ˆq1 m…q1 2 1; q2 1 1; q3 ; q4 ; q5 ; q6 ; q7 ; t† 2 ‰2zv…q2 1 q4 1 q6 † 2 bq7 Šm…q1 ; q2 ; q3 ; q4 ; q5 ; q6 ; q7 ; t† 2 3v2 q2 m…q1 ; q2 2 1; q3 1 1; q4 ; q5 1 d; q6 ; q7 ; t† 1 q3 …q1 ; q2 ; q3

3

2 1; q41 1; q5 ; q6 ; q7 ; t† 2 v2 q4 m…q1 ; q2 ; q311 ; q4

X 2 …t† 1 aX_ 2 …t† 1 bX2 …t† ˆ 22v2 X0 …t†2 X1 …t†

2 1; q5 ; q6 ; q7 ; t† 2 v2 q4 m…q1 ; q2 ; q3 ; q4 2 1; q5

X 3 …t† 1 X_ 3 …t† 1 bX3 …t† ˆ 23v2 …X0 …t†X1 …t†2 1 X0 …t†2 X2 …t†† .. .

1 3; q6 ; q7 ; t† 2 q5 m…q1 ; q2 ; q3 ; q4 ; q521 ; q611 ; q7 ; t† 2 v2 q6 …q1 ; q2 ; q3 ; q4 ; q5 1 1; q6 2 1; q7 ; t†

(25) and form a linear system of differential equations driven by polynomials of ®ltered Poisson processes. P2 k 1 Xk(t) Consider the second order solution X (2)(t) ˆ kˆ0 of the Duf®ng oscillator in Eq. (14) and the augmented state vector Q…t† ˆ {P2 …t† 0 ; Z…t†} 0 where P2 …t† ˆ _ 2 …t†; X1 …t†; X _ 1 …t†; X0 …t†; X _ 0 …t†} 0 . The set of {X2 …t†; X equations in Eq. (25) can be written in state space form as dQ…t† ˆ m…Q…t††dt 1

3

2

iˆ1

where b , 0 is a coef®cient, and dY(t) is de®ned in Eqs. (4)±(6). The system in Eq. (23) is called a Duf®ng oscillator [5±7]. The case of z . 0, v . 0 is considered. The nonlinear function f in Eq. (1) for the Duf®ng oscillator in Eq. (23) is f(X(t)) ˆ X(t) 3, and the approximation of 1 f(X(t)) in Eq. (9) becomes an equality. The set of linear differential equations resulting from the application of the perturbation method is

199

" # Z 0 y

M…dt; dy†

(26)

1 q6 m…q1 ; q2 ; q3 ; q4 ; q5 ; q6 2 1; q7 1 d; t† 1n

q7 X kˆ1

q7 ! …q ; q ; q ; q ; q ; q ; q 2 k; t†EY1k k!…q7 2 k†! 1 2 3 4 5 6 7 (29)

Q7 qi where m (q1,q2,q3,q4,q5,q6,qP 7;t) ˆ E iˆ1 Qi …t† denotes the 7 moment of order s ˆ iˆ1qi of Q (t). The moments m (0,0,0,0,0,0,k;t) can be obtained from the cumulants xk …t† ˆ n

Zt 0

E‰Y1 exp…b…t 2 s††Šk ds ˆ

nEY1k ‰exp…2kbt† 2 1Š kb (30)

200

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

Fig. 1. Transient response moments of example 1: (a) variance; (b) coef®cient of kurtosis.

of Z(t) or from the characteristic function of this process [12]

w…u; t† ˆ Eexp…iuZ…t†† ˆ exp{ 2 nt‰1 2

1 Zt Eexp‰…iuY1 exp…b…t 2 s††ŠdsŠ} t 0 (31)

and the equality i 2kd kw (u;t)/du k at u ˆ 0. Approximate momemts of the solution X(t) of Eq. (23) can be obtained from the moments of the state vector Q (t). The ®rst two central moments of the second order perturbation approximation of X(t) and their relationship with the moments of the state vector Q (t) are listed below. E‰X …2† …t†Š ˆ m…2† …t† ˆ E‰X0 …t† 1 X1 …t† 1 12 X2 …t†Š ˆ E‰X0 …t†Š 1 E‰X1 …t†Š 1 12 E‰X2 …t†Š ˆ m…0; 0; 0; 0; 1; 0; 0; t† 1 1m…0; 0; 1; 0; 0; 0; 0; t† 1 12 m…1; 0; 0; 0; 0; 0; 0; t†

(32)

v ˆ 1.0, z ˆ 0.05, d ˆ 1, b ˆ 0, n ˆ 1, and Gaussian random variables {Yk} with mean m y ˆ 0 and variance s2y ˆ Y12 2 …EY1 †2 ˆ 0.04. Consider the second order approximation X …2† …t† ˆ X0 …t† 1 1X1 …t† 1 12 X2 …t† of X…t†. Eqs. (25) and (28) become 3 2 7 6 X_ 2 …t† 7 6 7 6 7 6 7 6 2 2 3 6 2v X2 …t† 2 2zvX_ 2 …t† 2 3v X0 …t† X1 …t† 7 7 6 7 6 7 6 7 6 _ X 1 …t† 7 6 7dt 6 dQ…t† ˆ 6 7 7 6 6 2v2 X …t† 2 2zvX_ …t† 2 v2 X …t†3 7 7 6 1 1 0 7 6 7 6 7 6 7 6 _ …t† X 7 6 0 7 6 5 4 2v2 X0 …t† 2 2zvX_ 0 …t† 1 Z…t† " # Z 0 1 M…dt; dy† y (34) and

E‰…X …2† …t† 2 m…2† …t††2 Š ˆ E‰…X0 …t† 1 X1 …t†

m_ …q1 ; q2 ; q3 ; q4 ; q5 ; q6 ; t† ˆq1 m…q1 2 1; q2 1 1; q3 ; q4 ; q5 ; q6 ; t†

1 12 X2 …t† 2 m…2† …t††2 Š

2 2zv…q2 1 q4 1 q6 †…q1 ; q2 ; q3 ; q4 ; q5 ; q6 ; t†

ˆ m…2† …t†2 2 2m…2† …t†m…0; 0; 0; 0; 1; 0; 0; t†

2 3v2 q2 m…q1 ; q2 2 1; q3 1 1; q4 ; q5 1 2; q6 ; t†

1 m…0; 0; 0; 0; 2; 0; 0; t† 1 21‰m…0; 0; 1; 0; 1; 0; 0; t† 2 m…2† …t†…0; 0; 1; 0; 0; 0; 0; t†Š 1 12 ‰m…0; 0; 2; 0; 0; 0; 0; t† 2 2m…2† …t†m…1; 0; 0; 0; 0; 0; 0; t† 1 2m…1; 0; 0; 0; 1; 0; 0; t†Š 1 213 m…1; 0; 1; 0; 0; 0; 0; t† 1 14 m…2; 0; 0; 0; 0; 0; 0; t† (33) 4.2. Example 1 Let X(t) be the response of a Duf®ng oscilator of Eq. (23) driven by a Poisson white noise process [12, 13]. The solution X(t) is de®ned by Eqs. (4), (23) and (24) with

1 q3 m…q1 ; q2 ; q3 2 1; q4 1 1; q5 ; q6 ; t† 2 v2 q2 m…q1 1 1; q2 2 1; q3 ; q4 ; q5 ; q6 ; t† 2 v2 q4 m…q1 ; q2 ; q3 1 1; q4 2 1; q5 ; q6 ; t† 2 v2 q4 m…q1 ; q2 ; q3 ; q4 2 1; q5 1 3; q6 ; t† 2 v2 q6 m…q1 ; q2 ; q3 ; q4 ; q5 1 1; q6 2 1; t† 1 q5 m…q1 ; q2 ; q3 ; q4 ; q5 2 1; q6 1 1; t† 1n

q6 X kˆ1

q6 ! m…q1 ; q2 ; q3 ; q4 ; q5 ; q6 2 k; t†EY1k k!…q6 2 k†! (35)

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

201

Fig. 2. Stationary response moments of example 1: (a) variance; (b) coef®cient of kurtosis.

respectively. Stationary and transient moments of any order of X (2)(t) can be found exactly using Eq. (35). The same methodology can be used to obtained moments for the ®rst order approximation X (1)(t) of X(t). Fig. 1 shows as solid and dotted lines the transient variance and coef®cient of kurtosis of the ®rst and second order perturbation solutions of X(t), respectively. These moments correspond to zero initial conditions Q (0) ˆ 0 and 1 ˆ 0.025. Fig. 1 also shows as broken line the time history of the variance and coef®cient of kurtosis of X(t) in Eqs. (23) and (24) obtained via the Monte Carlo simulation method. The estimated moments of X(t) are based on 20 000 samples of the process. Fig. 2 shows as solid, dotted and broken lines the stationary variance and coef®cient of kurtosis of the ®rst and second order perturbation solutions of X(t), and by the Monte Carlo simulation of Eqs. (23) and (24), respectively. These stationary values are plotted for different values of the small parameter 1 that controls the cubic nonlinearity of the

oscillator. Results show that the accuracy of the approximate moments increase with the order p of the perturbation approximation. This observation agrees with PoincareÂ's monotonic convergence theorem of the deterministic pertubation method [5±7]. It is also noted that for small values of 1 the ®rst order perturbation solution provides satisfactory ®ndings. The second order perturbation is satisfactory for larger values of 1 but diverges from the Monte Carlo solution as 1 increases. 4.3. Example 2 Let X(t) be the response of a Duf®ng oscillator of Eq. (23) driven by a quadratic form of a ®ltered Poisson process. The solution X(t) is de®ned by Eqs. (4), (23) and (24) with v ˆ 1.0, z ˆ 0.05, d ˆ 2, b ˆ 3.0, n ˆ 1.0, and Gaussian random variables {Yk} with mean m y ˆ 0 and variance s y2 ˆ EY12 2 (EY1) 2 ˆ 0.4. Stationary and transient moments of any order of X (2)(t)

Fig. 3. Transient response moments of example 2: (a) mean; (b) variance; (c) coef®cient of skewness; (d) coef®cient of kurtosis.

202

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203

Fig. 4. Stationary response moments of example 2: (a) mean; (b) variance; (c) coef®cient of skewness; (d) coef®cient of kurtosis.

can be found exactly using Eq. (35). The same methodology can be used to obtain moments for the ®rst order approximation X (1)(t) of X(t). Fig. 3 shows as solid and dotted lines the transient mean, variance, coef®cient of skewness, and coef®cient of kurtosis of the ®rst and second order perturbation solutions of X(t), respectively. These moments correspond to zero initial conditions Q (0) ˆ 0 and 1 ˆ 0.025. Fig. 3 also shows as a broken line the time history of the transient mean, variance, coef®cient of skewness, and coef®cient of kurtosis of X(t) in Eqs. (23) and (24) obtained via the Monte Carlo simulation method. The estimated moments of X(t) are based on 20 000 samples of the process. Fig. 4 shows as solid, dotted and broken lines the stationary mean, variance, coef®cient of skewness, and coef®cient of kurtosis of the ®rst and second order perturbation solutions of X(t), and by the Monte Carlo simulation of Eqs. (23) and (24), respectively. These stationary values are plotted for different values of the small parameter 1 . Results show that the accuracy of the approximate moments increases with the order p of the perturbation approximation. It is also noted that for small values of 1 the ®rst order perturbation solution provides a satisfactory ®nding. The second order perturbation is satisfactory for larger values of 1 but diverges from the Monte Carlo solution as 1 increases. 5. Conclusions A method was developed for calculating approximately moments of any order of the state of nonlinear systems driven by polynomials of ®ltered Poisson processes. The method is based on the classical perturbation method, a

generalized version of the Itoà differentiation formula, and an algorithm for ®nding exactly moments of any order of the solution of linear systems driven by polynomials of ®ltered Poisson processes. Two numerical examples were presented to demonstrate the use and accuracy of the proposed method. The examples show that the accuracy of the approximate moments increase with the order p of the perturbation approximation. The ®rst order perturbation solution provides a good approximation for the moments of the solution of systems with small nonlinearities.

References [1] Grigoriu M, Waisman F. Linear systems with polynomials of ®ltered poisson processes. Probabilistic Engineering Mechanics 1997;12(2):97±103. [2] Gihman II, Skorohod AV. Stochastic differential equations. New York: Springer Verlag, 1972. [3] Ikeda N, Watanabe S. Stochastic differential equations and diffusion processes. New York: North-Holland, 1981. [4] Snyder DL. Random point processes. New York: Wiley, 1975. [5] Poincare H. Les meÂthodes nouvelles de la meÂcanique ceÂleste, vol I. Paris, France: Librairie Albert Blanchard, 1987. [6] Stoker JJ. Nonlinear vibrations. New York: Interscience, 1950. [7] Verhulst F. Nonlinear differential equations and dynamical systems. New York: Springer Verlag, 1985. [8] Waisman, F. Applications of non-Gaussian models to the solution of structural engineering problems, PhD Thesis, Cornell University, 1998. [9] Guckenheimer J, Holmes P. Nonlinear oscillations, dynamical systems and bifurcations of vector ®elds. New York: Springer Verlag, 1983. [10] Samorodnitsky G, Taqqu MS. Stable non-Gaussian random

F. Waisman, M. Grigoriu / Probabilistic Engineering Mechanics 14 (1999) 195±203 processes. Stochastic models with in®nite variance. London: Chapman and Hall, 1994. [11] Snyder DL, Miller M. Random point processes in time and space. New York: Springer Verlag, 1991. [12] Grigoriu M. Applied non-Gaussian processes: examples, theory,

203

simulation, linear random vibration, and MATLAB solutions. Englewoods Cliffs, NJ: Prentice Hall, 1995. [13] Stratonovich RL. Topics in the theory of random noise, vol 1. New York: Gordon and Breach, 1981.