Strojniški vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 © 2016 Journal of Mechanical Engineering. All rights reserved. D0l:10.5545/sv-jme.2016.3752 Original Scientific Paper Received for review: 2016-03-03 Received revised form: 2016-05-17 Accepted for publication: 2016-05-23 Nonlinear MDOF system Survival Probability Determination Subject to Evolutionary Stochastic Excitation Ioannis P. Mitseas14* - Ioannis A. Kougioumtzoglou2 - Pol D. Spanos3 - Michael Beer145 1 Leibniz University Hannover, Faculty of Civil Engineering and Geodetic Science, Germany 2 Columbia University, Department of Civil Engineering and Engineering Mechanics, USA 3 Rice University, Department of Mechanical Engineering, USA 4 University of Liverpool, Institute for Risk and Uncertainty and School of Engineering, United Kingdom 5 Tongji University, School of Civil Engineering and Shanghai Institute of Disaster Prevention and Relief, China An approximate technique for assessing the reliability of nonlinear multi-degree-of-freedom (MDOF) systems subject to a non-stationary stochastic excitation vector is developed. The proposed technique can be construed as a two-stage approach. First, relying on statistical linearization and utilizing a dimension reduction approach the nonlinear n-degree-of-freedom system is decoupled and cast into (n) effective single-degree-of-freedom (SDOF) linear time-variant (LTV) oscillators. Second, utilizing the effective SDOF LTV oscillator time-varying stiffness and damping elements in conjunction with a stochastic averaging treatment of the problem, the MDOF system survival probability and firstpassage PDF are determined. Overall, the developed technique appears to be efficient and versatile since it can handle readily at a low computational cost, a wide range of nonlinear/hysteretic behaviors as well as various stochastic excitation forms, even of the fully non-stationary in time and frequency kind. A 3-DOF system exhibiting hysteresis following the Bouc-Wen model is included in the numerical examples section. Comparisons with pertinent Monte Carlo simulations demonstrate the accuracy of the technique. Keywords: first-passage problem, nonlinear stochastic dynamics, evolutionary stochastic processes, nonlinear/hysteretic systems, survival probability Highlights • A stochastic dynamics methodology for determining the survival probability of nonlinear MDOF systems. • Approximate analytical expressions provided for estimating the time-varying survival probability. • Survival probabilities and the associated first-passage PDFs are determined at a low computational cost. • The developed technique is characterized by enhanced versatility as it can handle readily a wide range of nonlinear behaviors as well as various stochastic excitations with arbitrary EPS forms. 0 INTRODUCTION Excitations acting upon dynamical systems such as wind, wave, and seismic loads commonly exhibit evolutionary features. In this setting, not only the intensity of the excitation but also its frequency content exhibit strong variability. This fact necessitates the representation of this class of structural loads by non-stationary stochastic processes. Further, structural systems under severe excitations can exhibit significant nonlinear behavior of the hysteretic kind. Thus, of particular interest to the structural dynamics community is the development of techniques for determining the response and assessing the reliability of nonlinear/hysteretic systems subject to evolutionary stochastic excitations (e.g., [1] to [3]). Further, in engineering dynamics, the evaluation of the probability that the system response stays within prescribed limits for a specified time interval is advantageous for reliability based system design applications. In this regard, the first-passage problem, that is, the determination of the above time-variant probability known as survival probability, has been a persistent challenge in the field of stochastic dynamics for many decades. Monte Carlo simulation techniques are among the most potent tools for assessing the reliability of a system (e.g. [4]). Nevertheless, there are cases where the computational cost of these techniques can be prohibitive, especially when large-scale complex systems are considered; thus, rendering the development of alternative efficient approximate analytical/numerical techniques for addressing the first-passage problem necessary. Indicatively, one of the early approaches, restricted to linear systems, relies on the knowledge of the mean up-crossing rates and on Poisson distribution based approximations (e.g., [5] to [7]). Further attempts to address the firstpassage problem range from analytical ones (e.g., [8]) to numerical ones (e.g., [9]). Furthermore, techniques based on the concepts of the numerical path integral (e.g., [10] to [13]), of the probability density evolution (e.g., [3]), or of stochastic averaging/linearization (e.g., [14]) constitute some of the more recent approaches. 440 *Corr. Author's Address: Institute for Risk and Reliability, Faculty of Civil Eng. and Geodetic Science, Callinstr. 34, 30167, Hannover, Germany, mitseas@irz.uni-hannover.de Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 In this paper, an approximate analytical technique for determining the survival probability and firstpassage probability density function (PDF) of nonlinear multi-degree-of-freedom (MDOF) systems subject to an evolutionary stochastic excitation vector is developed. Specifically, first relying on a statistical linearization based dimension reduction approach the original MDOF system is decoupled and cast into (n) effective single-degree-of-freedom (SDOF) linear time-variant (LTV) systems corresponding to each and every degree of freedom of the original MDOF system. Second, a stochastic averaging based approximate technique is utilized to derive the nonlinear MDOF system survival probability and first-passage PDF at a low computational cost. The remainder of this paper is organized as follows: In section 1.1 the statistical linearization technique for nonlinear MDOF systems is presented. Next, in section 1.2 a stochastic averaging/statistical linearization treatment of the problem, through a system dimension reduction approach is briefly delineated. In section 1.3, it is shown that the nonlinear MDOF system non-stationary marginal, transition and the joint response amplitude probability density functions (PDFs) can be approximated by closed-form expressions. Further, section 2 provides analytical closed-form expressions for the time-dependent survival probability of the nonlinear MDOF structural system as well as for the corresponding first-passage PDF. In section 3, illustrative examples comprising a 3-DOF system exhibiting Bouc-Wen hysteresis and subject to evolutionary stochastic excitations are considered. Pertinent MCS data demonstrate the reliability of the proposed technique. Finally, section 4 provides with concluding remarks. 1 MDOF SYSTEM DIMENSION REDUCTION In this section, the basic elements of an approximate dimension reduction/decoupling technique developed by some of the authors for determining the non-stationary response amplitude PDF of nonlinear MDOF systems subject to evolutionary stochastic excitation are reviewed for completeness; see [15] and [16] for a more detailed presentation. 1.1 Statistical Linearization Treatment Consider an n-degree-of-freedom nonlinear system governed by the equation: where y denotes the response acceleration vector, y is the response velocity vector, y is the response displacement vector, defined in relative coordinates; M, C and K denote the (n * n) mass, damping and stiffness matrices, respectively; g(y,y) is an arbitrary nonlinear (n * 1) vector function of the variables y and y . F(t)T = (/1(t), /2(t), ...,fn(t)) is a (n * 1) zero mean, non-stationary stochastic excitation vector process defined as F(t) = ya (t) where Yt = (7i, Y2, •••, Yn) is an arbitrary (n * 1) vector of constant weighting coefficients, and a (t) is a non-stationary process with an evolutionary power spectrum (EPS) Sa (m, i). In this regard, F(t) possesses the EPS matrix: SF (o, t) Y Ss (o, t) 0 0 y2 (o, t) 0 0 YS, (o, t) .(2) Further, the non-stationary stochastic excitation process is regarded to be a filtered stationary stochastic process according to the concept proposed by Priestley [17]; see also [18]. Thus, the excitation EPS matrix of Eq. (2) takes the form: SF (m, t) = A(m,t^SF (m)A(m,t (3) where the superscripts (T) and (*) denote matrix transposition and complex conjugation, respectively; A(m, t) is the modulating matrix which serves as a time-variant filter; and SF (co) is the power spectrum matrix corresponding to the stationary stochastic vector process F (t). Note that both separable and non-separable EPS can be defined considering Eq. (3). In this manner, excitations exhibiting variability in both the intensity and the frequency content can be considered. Focusing next on the frequency domain, the response determination problem is defined as seeking the corresponding system response EPS matrix of the form: S, (o,t ) = 'S™ (o t) »(ot) (ot) t) t) ^ (ot) .(ot) t) ^ (o t) .(4) My + Cy + Ky + g (y, y ) = F (t). (1) According to the statistical linearization method (e.g., [1] to [3]), a linearized version of Eq. (1) takes the form: Nonlinear MDOF system Survival Probability Determination Subject to Evolutionary Stochastic Excitation 441 Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 My + (C + Ceq )) + (K + Keq )y = F (t). (5) where H(m) is the frequency response function (FRF) matrix defined as: Next, adopting the standard assumption that the response processes are Gaussian, the time-dependent elements of the equivalent linear matrices Ceq and Keq are given by the expressions: and «=* It | (6) (7) Further, for a linear MDOF system subject to evolutionary stochastic excitation a matrix input-output spectral relationship of the form: Sy (1) = Hgen (( t )Sf ((a) Hen ((, t), (8) can be derived (e.g., [1] and [3]), where Hgen (®, t) = jh (t-t)A (m,t)e-'^dt. (9) 0 In Eq. (9) denotes the impulse response function matrix. Furthermore, the time dependent cross-variance of the response can be evaluated by the expression: K> E[yy ]=isy,y] (®,t(!°) — K> It can be readily seen that Eqs. (6) to (10) constitute a coupled nonlinear system of algebraic equations to be solved numerically for the system response covariance matrix. Note in passing that instead of the frequency domain Wiener-Khinchin relationship of Eq. (8), a state-variable formulation can be adopted yielding a system of differential equations of the Lyapunov kind (e.g., [1] and [19]) for the system response covariance matrix. Nevertheless, although a pre-filtering treatment can be applied for considering non-stationary stochastic excitation processes of the separable kind (e.g., [1]), excitations possessing a non-separable EPS (e.g. realistic cases of earthquake excitations) cannot be accounted for, at least in a straightforward manner. Next, omitting the convolution of the impulse response function matrix with the modulating matrix can lead to substantial reduction of computational effort, especially for the case of MDOF systems (e.g., [16] and [20]). In this manner, Eq. (9) takes the form: H(«) = ((M + m(C + C„) + (( + Keq))-1. (12) Consequently, taking into account Eqs. (3) and (8), Eq. (11) becomes: s; (ft),t) = H(a)SF (m,t)H(a). (13) Note that the Eq. (13) can be regarded as a quasi-stationary approximate relationship which, in general, yields satisfactory accuracy in cases of relatively stiff systems (e.g., [20] to [22]). Note in passing that the spectral input-output relationship of Eq. (13) is exact for the case of stationary processes (e.g., [1] to [3]). Further, adopting the aforementioned quasi-stationary approach, it can be readily seen that for the ith degree of freedom, using Eqs. (2), (10) and (13) yields: E [ y2 (t )]=J (|Ha (®)f Y +... — IX> + |H n (m) y2JSa (m,t)dm, (14) E [ y,2 (t )]=]a2 (| Ha (a) yY +... — K> + |H n (ffl)|2 y2 )(a, t)da, (15) and Hgen (, t ) = H (CO) A ((, t), (11) Eqs. (14) and (15) hold true in the approximate quasi-stationary sense delineated earlier. Clearly, Eq. (13) constitutes an approximate formula for determining the MDOF system response EPS matrix at a low computational cost; thus, circumventing computationally intensive Monte Carlo simulations. 1.2 Effective SDOF Linear Time-Variant System Following next the system dimension reduction approach developed in [16], an auxiliary effective SDOF LTV system corresponding to the ith degree of freedom can be defined as: ^ + Peq,i (t) y + (t) y, = a, (t) (16) where the time-varying equivalent stiffness and damping elements of the effective LTV system can be determined by equating the variances of the response displacement and velocity expressed utilizing the quasi-stationary FRF of Eq. (16) with the 442 Mitseas, I.P. - Kougioumtzoglou, I.A. - Spanos, P.D. - Beer, M. Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 corresponding ones determined via Eqs. (14) and (15); this yields: * [ y (t )]= to = i (m2^ (t)-m2)2 + (PtqJ (t)m)2 xyi2Sa (m,t)dm,(17) and (t W)2 + (P^ (tt xy2 Sa (to, t )dm, (18) Clearly, Eqs. (17) and (18) in conjunction with Eqs. (14) and (15) constitute a nonlinear system of two algebraic equations to be solved for the evaluation of the LTV system time-varying equivalent stiffness mlq,i (t) and damping i (t) coefficients. Note that determining the time-varying natural frequency meq,i (t) is especially important for a number of reasons such as tracking and avoiding moving resonance phenomena (e.g., [23]), determining peak system response estimates based on design spectrum compatible excitation power spectra (e.g., [24]), or developing efficient approximate techniques for determining nonlinear system survival probability and first-passage PDF (e.g., [14]). Next, a stochastic averaging technique (e.g., [15] and [16]) is applied for casting the second-order stochastic differential equation (SDE) of Eq. (1) into a first-order SDE governing the evolution in time of the response amplitude ai(t). In this regard, and based primarily on the assumption of light damping, it can be argued that the response yi (t) of the effective LTV system of Eq. (16) exhibits a pseudo-harmonic behavior described by the equations: y, (t) = at (t)cos((t)t + ft (t)), (19) and y, (t) = -(°eqJ (t(t)x sm((t)t+% (t)). (20) In Eq. (19) the response amplitude ai (t) is a slowly varying function with respect to time defined as: ( ■ V «2 (t) = y2 (t) y (t) (t). m \eq, (21) whereas (t) stands for the phase of the response yi (t). Further, relying on a combination of deterministic and stochastic averaging (e.g., [16]) a first-order SDE governing each and every degree-of-freedom response amplitude process ai (t) takes the form: ■f\ ^a t\ t \ nSF ((t)'t) a (t) = -1 Peqj (t )a (t) + 0 - t2 ) = 0. (28) Relying on the assumption that the equivalent damping and stiffness coefficients follow a slowly varying with respect to time behavior, the following approximations over a small time interval [t^, tj are introduced; i.e., Pev(thj)=^(ty-i) and Meg, i(tj=i(tij-i) for 16 [t^,ty]. Next, based on the slowly varying with time behavior of the EPS, Sf,, (a>,t) is also treated as a constant over the interval [ty-i, tiJ]. Further, based on the above assumptions, introducing the variable Tij=tiJ- - tiJ--1, and applying a first-order Taylor expansion around the point tij = 0, Eqs. (27) and (28) become (see [14] for a detailed derivation): , ( tu ) = nSF ( ( j ) j ) . , (29) meq,i ,j-1 ) and hi (j-l. j ) = -Peqj (j K (30) Furthermore, considering Eqs. (25) and (29) and applying a first-order Taylor expansion for the response variance ci(t) around the point t=tij-1 yields: c, (-„ t,, j ) = C ( j)-c, (j )(l-Peq, (j ).) .(31) 444 Relying next on the Markovian assumption for the process ai, the joint-response amplitude PDF P(ai,j-i, tj-i ; ai,j, tj) is given by: p (ai ,j-i, h-i; ai, j > tj) =p (ai ,j-i, O-i) p(ai,j, h1 ai, j-i, O-i) .(32) Utilizing Eqs. (24) and (26), Eq. (32) becomes: t \ a, j_iai j p {( tj-i; a J , j = (. ) ({ . )x C (j-i )c, ((-1, tjj ( -afjci (j-1 ) - atj-iC (-1' tJ ) exp (( -i )c ((tj ) (tj)(-1) 1J (auh(tj) (j-1 )c ( tj ) 01 c ( tj ) . Further, setting (33) 2 ri.j = ^ (1 -K, (( ,-,),, ). Eq. (31) yields: C (( j-„ t,, j ) = Ci ((,, )(l- rl ). (34) (35) Next, considering Eqs. (29) and (30) and Eqs. (34) and (35), the joint response amplitude PDF p(ai]-1, tij-1 ; aij,tj of Eq. (33) is given in the form: P («y-P ti,j-{;ai, j, t, j ) = ( 2 a, x exp f x * (, j - ) ( j )(( - ^ 2j ) ( j-( ) + aa j-iC (( ^ ) 2c- (j )c (, j )(( - r 2j ) A Io x h (tj ) (t,, j ) (1 - rf2, ) (36) 2 NONLINEAR MDOF SYSTEM RELIABILITY ASSESSMENT In this section the approximate analytical technique developed by some of the authors in [14] for nonlinear SDOF survival probability determination is generalized herein to account for MDOF systems by utilizing the dimension reduction/decoupling technique outlined in section 1. In this regard, the survival probability Pf is defined as the probability that the system response 54 Mitseas, I.P. - Kougioumtzoglou, I.A. - Spanos, P.D. - Beer, M. a Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 amplitude a, stays below a prescribed barrier B over the time interval [0, T], given that a(t=0) < B. Further, the first-passage PDF and the survival probability pi (T) are related according to the expression: dPf (T) pB (t ) = - dT (37) Next, adopting the discretization scheme employed in [9] yields intervals of the form: \}i,j-i, fi, j ], j = 1,2,..., m, tio = 0, Teq,i (-1 ) km = T, h, j - ti,i-1 =■ 2 (38) where the response amplitude ai is assumed to be constant over [t^, tiJ] due to its slowly varying in time behavior. In Eq. (38) Teqi represents the LTV system equivalent natural period given by: Teqi (t )=- 2n ,, (t )" (39) Note in passing that a smaller time interval can be chosen if higher accuracy is required. In this regard, the survival probability Pf is assumed to have a constant value over the same time interval as well. Obviously, the survival probability is given by p- (t)=n [i - ^ ], (40) j=i where Ffj is defined as the probability that the response amplitude ai will exceed the prescribed barrier B over the time interval [tij_j, tj, given that no crossings have occurred prior to time tij_1. Next, invoking the Markovian property of the response amplitude ai, one gets: fb = Prob[a j'u B n a (j )< B] = j, (41) Prob[at [tiJ_1 ) < B] H j-i where n denotes the intersection symbol. Utilizing Eq. (24) Hfj_1 can be determined analytically in a straightforward manner; that is, HU-1 = \P (-i> ti,j-l )i,j-1 =1 - exp B \ 2C ( j ). ,(42) whereas HB_x j is defined as a double integral of the form: œ B Hu-1,j = \dai,j x \p(ai,j-i ' ti,j-1;ai,j 'j ). (43) Further, taking into account Eq. (36) and expanding the Bessel function I0(x) in the form (e.g., [27]): ^ (x /2 ) 1o (x ) = X k\r(K+1) (44) analytical treatment of the involved integrals is possible yielding: K= 4o+£4,, (45) where A ,0 = exP - Bl (,- )(1 - r 2- ) 1 - exp -B2 2c,- (— )( - r 2- ) x(l - r—, (46) and Ain =- (((j) j)(1 -^r with Ln = 4"ci (,1)"+1 C ( j )(l - r) X , (47) r 1 + n B '2c, (j-i)( - r2 ) -^[1 + n,0] ( (j)(1 - r2, ))n r[l + n] / \-n B2 ( X v v c ( j )(i - r 2, ) r, [i+n]-r, 1 + n, B2 2c, (, )(1 - £)_ X B 2 ,(48) In Eq. (48) r [y,z] represents the incomplete Gamma function defined as r [y,z] = j'r-1e-dt. z Concisely, the developed technique comprises the following steps: i. Determination of the MDOF system non-stationary response covariance matrix (Eqs. (10) and (13)) via a statistical linearization treatment of the problem. ii. Determination of the equivalent linear time-varying elements fieqi(t) and meqi(t) by solving the system of algebraic equations (Eqs. (17) and (18)). Nonlinear MDOF system Survival Probability Determination Subject to Evolutionary Stochastic Excitation 445 k=0 n=l 2n L X + X 0 B 0 Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 111. IV. Determination of c(t) via numerically integrating the first-order ODE Eq. (25). Determination of the equivalent natural period Teq i(t) (Eq. (39)) and discretization of the time domain via Eq. (38). v. Determination of the parameters Hi H j-i,j via Eqs. (42) and (43). j-i and vi. Determination of the survival probability Pf (T) via Eq. (40) and of the corresponding firstpassage PDF pi (T) via Eq. (37). 3 NUMERICAL APPLICATIONS In this section, a nonlinear three-degree-of-freedom system following the Bouc-Wen hysteretic model (e.g., [28] and [29]) subject to evolutionary stochastic excitation is considered to demonstrate the reliability of the technique. The survival probabilities and the first-passage PDFs obtained via the developed approximate technique are compared with survival probability and first-passage PDF estimates obtained via pertinent Monte Carlo simulations (10,000 realizations). The Monte Carlo simulations were conducted by utilizing a spectral representation methodology; additional details can be found in [30]. Further, a standard fourth-order Runge-Kutta numerical integration scheme is employed for solving the nonlinear system differential equation of motion (Eq. (1)), whereas the barrier level B is expressed as a fraction X of the maximum over time and over DOF value of the non-stationary response displacement standard deviation, i.e. B = Xmax(ai (t)) with .__i and t ui (t) = yjci (t). Considering displacements defined in relative coordinates, the 3-DOF nonlinear system is governed by Eq. (1) where yT =(y y2 ys zi z2 zs ) where M = Mn = M11 M12 M21 M22. (49) (50) Further, K = Kii Ki K„ K,. where K„ = akj -ak2 0 ak -ak 0 ak (53) (54) K12 = (1 - a) -(1 - a)k. 0 0 0 (1 - a)k2 -(1 - a)k3 0 (1 - a)k3 and K21 K 22 03,3' (55) (56) In Eqs. (54) and (55) stands for the rigidity ratio which can be viewed as a form of post-yield to pre-yield stiffness ratio (a = 1 corresponds to the linear system). Further, the damping matrix of the structural system C is assumed to be proportional to the stiffness matrix; that is, C = C C ^11 12 C C 21 22 where Cn = cKn, C = C = 0 12 '-'21 03,3 ' and C = ^22 1 0 0 0 1 0 0 0 1 (57) (58) (59) (60) In Eq. (58) c is taken equal to 0.2 x 102. For the specific example yt = mt, and the loading vector becomes F(t)T = (( (') f (t) f (t) 0 0 0). (61) Further, ml 0 0 m2 m2 0 , (51) g (y, y )T = m3 m3 m3 = (000 - gl (jp Z1 )-g 2 (j^ Z2 )-g 3 (j^ Z3 and M12 = M21 = M22 = 03,3 . (52) In the Bouc-Wen model the additional state zt is associated with the displacement y via the equation: = gi (, ), (63) where 446 Mitseas, IP. - Kougioumtzoglou, I.A. - Spanos, P.D. - Beer, M. Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 gi(,zi) = -r\.^iIzilziI"1 -Py>Zf + Ay'i. (64) The parameters y, ft, A and n are capable of representing a wide range of hysteresis loops (e.g., [28] and [29]). In this example the values a=0.15, ft = y = 0.5, n = 1 and A = 1 are considered. The equivalent linear matrices take the form (e.g., [1] to [3]): Ceq = C C eqll eq12 C C eq21 eq22 where and C = C = C = 0 eqll ^eq12 eq 22 "3,3' (65) (66) Further, where and Ceql 0 0 " Ceq 21 0 C eq 2 0 . (67) 0 0 c , eq3 Keq = Keq11 _ Keq 21 K eq12 Keq 22 (68) Keq11 Keq12 Keq 21 = O33 (69) keql 0 0 Keq22 0 keq 2 0 . (70) 0 0 k 3 eq3 The elements c and k in Eqs. (67) and (70) are given by the expressions: E(y*), "¿2 and E fa - A, (71) (72) respectively. 3.1 A 3-DOF Hysteretic System under Evolutionary Stochastic Excitation of the Separable Form In this example, the excitation EPS Sa (to, i) takes the form S„ () represents the widely used in engineering applications Clough-Penzien power spectrum (e.g., [31]) and w (t) denotes a time-modulating envelope function given by: w(t) = k( -e-2), (74) where b1 = 0.1 and b2 = 0.3; and k is a normalization constant so that w (t)max = 1. The Clough-Penzien spectrum is given by: SCP (®) = S0 wg 4+4(^g v -2-X (2 -®2) + 2®gV x_mf)4_, (l -(/ mf )2) + 4£f 2 (/ mf )2 (75) where S0 is the amplitude of the excitation spectrum, modeled as a white noise process. The parameters values used are ;g=0.7, mg= 2 rad s-1, 0.6, mf= 12.5 rad s-1. The total duration of the excitation is 20 seconds. Further, the hysteretic 3-DOF system has the properties mx = 2.0615x105 kg, m2 = 2.0559x105 kg, m3 = 2.0261x105 kg, k = 3.9668x10« Nm-1 k2 = 3.5007x10« Nm-1 and k3 = 2.6927x10« Nm-1. In Fig. 1 the EPS of Sa (a,i) is plotted for S0 = 20 m2s 3. Time [s] (73) Fig. 1. Separable excitation evolutionary power spectrum In Figs. 2 and 3 the equivalent time-varying natural frequency meq i(t) and fteqJ(t) the damping element corresponding to each DOF are plotted, respectively. Note that the hysteretic/degrading behavior of the system is captured by the decreasing with time trend of the stiffness element, as well as the increasing with time trend of the damping element. Further, in Figs. 4 and 5 the survival probabilities Pf (T) and indicatively the corresponding firstpassage PDFs pi (T) for the first DOF of the hysteretic MDOF system are plotted for various barrier levels, respectively. The value N = 30 is chosen regarding the number of terms to be included in Eq. (45). Comparisons between the analytical approximate technique and MCS data (10,000 realizations) Nonlinear MDOF system Survival Probability Determination Subject to Evolutionary Stochastic Excitation 447 Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 Fig. 2. Equivalent natural frequency mt Fig. 4. Survival probability for various values of the parameter X for the first DOF; comparisons with MCS (10,000 realizations) Fig. 3. Equivalent damping coefficient fteqj(t) demonstrate a satisfactory degree of agreement. Note that the irregular/non-smooth shape of the approximate technique based first-passage PDFs is due to the differentiation of the survival probability (Eq. (37)). In this regard, the survival probability Eq. (40) is assumed to have constant values over the time intervals [t^, ty], resulting in a non-smooth representation. Obviously, the level of non-smoothness increases when differentiation takes place. Furthermore, in Figs. 6 and 7 the survival probabilities Pf (T) corresponding to the second and third DOF of the system are plotted for various barrier levels. Comparisons with MCS demonstrate a satisfactory degree of accuracy for these cases as well. Fig. 5. First-passage PDF for various values of the parameter X for the first DOF; comparisons with MCS (10,000 realizations) Fig. 6. Survival probability for various values of the parameter X for the second DOF; comparisons with MCS (10,000 realizations) 3.2 A 3-DOF Hysteretic System under Evolutionary Stochastic Excitation of the Non-Separable Form The excitation EPS Sa (c, i) is assumed to have the non-separable form: (œ, t ) = 50ebt btt2e (76) with S0 = 10 m2s-3 and b = 0.5. This spectrum comprises some characteristics of particular interest, such as decreasing of the dominant frequency with respect to 8max(n, ft)) - 3 DOF - Analytical Smaxfn^t)) - 3rd DOF - MCS (10,000) .2max(l(t) and damping fieq>l(t) elements corresponding to each DOF are plotted, respectively. Underlying the analytical approximate approach is the attempt to capture the time evolution as well as the essential characteristics of the frequency content of the nonlinear system response. Note that the ability of the technique to provide with time-varying natural frequencies weqi(t) can be of particular importance if seen in conjunction with recent theoretical developments regarding the concept of the mean instantaneous frequency (MIF) (e.g., [33] to [35]). In this regard, weqJ(t) together with the MIF of the excitation can be potentially employed for evaluating the effects of temporal non-stationarity in the frequency content of the excitation on the system response as well as for tracking moving resonance phenomena (e.g., [23] and [36]). Further, in Figs. 11, 12 and 13 the survival probabilities Pf (T) for every DOF of the hysteretic MDOF system are plotted for various barrier levels, respectively; comparisons with MCS (10,000 realizations) demonstrate a satisfactory degree of accuracy. Fig. 11. Survival probability for various values of the parameter X for the first DOF; comparisons with MCS (10,000 realizations) Fig. 12. Survival probability for various values of the parameter X for the second DOF; comparisons with MCS (10,000 realizations) Fig. 13. Survival probability for various values of the parameter X for the third DOF; comparisons with MCS (10,000 realizations) Nonlinear MDOF system Survival Probability Determination Subject to Evolutionary Stochastic Excitation 449 Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 4 CONCLUDING REMARKS An approximate analytical technique for determining the time-varying survival probability and associated first-passage PDF of nonlinear/hysteretic MDOF systems subject to evolutionary stochastic excitation has been developed. Specifically, based on an efficient dimension reduction approach and relying on the concepts of stochastic averaging and statistical linearization, the original nonlinear w-degree-of-freedom system has been decoupled and cast into (w) effective single-degree-of-freedom (SDOF) linear time-variant (LTV) oscillators corresponding to each and every DOF. In this regard, time-varying effective stiffness a>2q ¡ (t) and damping (t) elements corresponding to each and every DOF have been defined and computed, while the non-stationary marginal, transition and joint response amplitude PDFs have been efficiently determined in closed-form expressions. Finally, the MDOF system survival probability and first-passage PDF have been determined approximately in a computationally efficient manner. Overall, the developed technique exhibits enhanced versatility since it can handle readily a wide range of nonlinear behaviors as well as various stochastic excitations with arbitrary non-separable EPS forms that exhibit strong variability in both the intensity and the frequency content. A 3-DOF system exhibiting hysteresis following the Bouc-Wen model has been included in the numerical examples section. Comparisons with pertinent Monte Carlo simulations have demonstrated the reliability of the technique. 5 ACKNOWLEDGEMENTS The first author acknowledges the financial support of the State Scholarships Foundation - IKY in Greece. 8 REFERENCES [1] Roberts J.B., Spanos P.D. (2003). Random Vibration and Statistical Linearization. Dover Publications, New York. [2] Soong, T.T., Grigoriu, M. (1993). Random Vibrations of Mechanical and Structural Systems. Prentice-Hall, New Jersey. [3] Li, J., Chen, J. (2009). Stochastic Dynamics of Structures, Dover Publications, New York, D0l:10.1002/9780470824269. [4] Schueller, G.I., Pradlwarter, H.J., Koutsourelakis, P.S. (2004). A critical appraisal of reliability estimation procedures for high dimensions. Probabilistic Engineering Mechanics, vol. 19, no. 4, p. 463-474, D0I:10.1016/j.probengmech.2004.05.004. [5] Corotis, R., Vanmarcke, E.H., Cornell, C.A. (1972). First passage of nonstationary random processes. Journal of the Engineering Mechanics Division, vol. 98, no. 2, p. 401-414. [6] Vanmarcke, E.H. (1975). On the distribution of the firstpassage time for normal stationary random processes. ASME Journal of Applied Mechanics, vol. 42, no. 1, p. 215-220, D0I:10.1115/1.3423521. [7] Barbato, M., Conte, J.P. (2001). Structural reliability applications of non-stationary spectral characteristics. ASCE Journal of Engineering Mechanics, vol. 137, no. 5, p. 371-382, D0I:10.1061/(ASCE)EM.1943-7889.0000238. [8] Kovaleva, A. (2009). An exact solution of the first-exit time problem for a class of structural systems. Probabilistic Engineering Mechanics, vol. 24, no. 3, p. 463-466, D0I:10.1016/j.probengmech.2009.01.002. [9] Solomos, G.P., Spanos, P.D. (1983). Structural reliability under evolutionary seismic excitation. International Journal of Soil Dynamics and Earthquake Engineering, vol. 2, no. 2, p. 110116, D0I:10.1016/0261-7277(83)90007-4. [10] lourtchenko, D.V., Mo, E., Naess, A. (2006). Response probability density functions of strongly non-linear systems by the path integration method. International Journal of NonLinear Mechanics, vol. 41, no. 5, p. 693-705, D0I:10.1016/j. ijnonlinmec.2006.04.002. [11] Naess, A., lourtchenko, D., Batsevych, 0. (2011). Reliability of systems with randomly varying parameters by the path integration method. Probabilistic Engineering Mechanics, vol. 26, no. 1, p. 5-9, D0I:10.1016/j.probengmech.2010.05.005. [12] Kougioumtzoglou, I.A., Spanos, P.D. (2013). Response and first-passage statistics of nonlinear oscillators via a numerical path integral approach. ASCE Journal of Engineering Mechanics, vol. 139, no. 9, p. 1207-1217, D0I:10.1061/ (ASCE)EM.1943-7889.0000564. [13] Kougioumtzoglou, I.A., Spanos, P.D. (2014). Stochastic response analysis of the softening Duffing oscillator and ship capsizing probability determination via a numerical path integral approach. Probabilistic Engineering Mechanics, vol. 35, p. 67-74, D0I:10.1016/j.probengmech.2013.06.001. [14] Spanos, P.D, Kougioumtzoglou, I.A. (2014). Survival probability determination of nonlinear oscillators subject to evolutionary stochastic excitation. Journal of Applied Mechanics, vol. 81, no. 5, 051016, D0I:10.1115/1.4026182. [15] Spanos, P.D., Lutes, L.D. (1980). Probability of response to evolutionary process. Journal of Engineering Mechanics, vol. 106, no. 2, p. 213-224. [16] Kougioumtzoglou, I.A., Spanos, P.D. (2013). Nonlinear MD0F system stochastic response determination via a dimension reduction approach. Computers and Structures, vol. 126, p. 135-148, D0I:10.1016/j.compstruc.2012.11.020. [17] Priestley, M.B. (1965). Evolutionary spectra and non-stationary processes. Journal of the Royal Statistical Society, vol. 27, no. 2, p. 204-237. [18] Dahlhaus, R. (1997). Fitting time series models to non-stationary processes. The Annals of Statistics, vol. 25, no. 1, p. 1-37, D0I:10.1214/aos/1034276620. [19] Gajic, Z., Qureshi, M. (1995). Lyapunov Matrix Equation in System Stability and Control. Academic Press, New York. [20] Mitseas, I.P., Kougioumtzoglou, I.A., Beer, M. (2016). An approximate stochastic dynamics approach for nonlinear structural system performance-based multi-objective optimum 450 Mitseas, I.P. - Kougioumtzoglou, I.A. - Spanos, P.D. - Beer, M. Strojniski vestnik - Journal of Mechanical Engineering 62(2016)7-8, 440-451 design. Structural Safety, vol. 60, p. 67-76, D0l:10.1016/j. strusafe.2016.01.003. [21] Hammond, J.K. (1973). Evolutionary spectra in random vibration. Journal of the Royal Statistical Society, vol. 35, no. 2, p. 167-188. [22] Jangid, R.S., Datta, T.K. (1999). Evaluation of the methods for response analysis under non-stationary excitation. Shock and Vibration, vol. 6, no. 5-6, p. 285-297, D0I:10.1155/1999/312010. [23] Tubaldi, E., Kougioumtzoglou, I .A. (2015). Nonstationary stochastic response of structural systems equipped with nonlinear viscous dampers under seismic excitation. Earthquake Engineering & Structural Dynamics, vol. 44, no. 1, p. 121-138, D0I:10.1002/eqe.2462. [24] Giaralis, A., Spanos, P.D. (2010). Effective linear damping and stiffness coefficients of nonlinear systems for design spectrum based analysis. Soil Dynamics and Earthquake Engineering, vol. 30, no. 9, p. 798-810, D0I:10.1016/j.soildyn.2010.01.012. [25] Spanos, P.D., Solomos, G.P. (1983). Markov approximation to transient vibration. Journal of Engineering Mechanics, vol. 109, no. 4, p. 1134-1150, D0I:10.1061/(ASCE)0733-9399(1983)109:4(1134). [26] Spanos, P.D. (1978). Non-stationary random vibration of a linear structure. International Journal of Solids and Structures, vol 14, no. 10, p. 861-867, D0I:10.1016/0020-7683(78)90076-8. [27] Abramowitz, M., Stegun, I.A. (1970). Handbook of Mathematical Functions. Dover Publications, New York. [28] Wen, Y.K. (1980). Equivalent linearization for hysteretic systems under random excitation. Journal of Applied Mechanics, vol. 47, no. 1, p. 150-154, D0I:10.1115/1.3153594. [29] Ikhouane, F., Rodellar, J. (2007). Systems with Hysteresis: Analysis, Identification and Control using the Bouc-Wen Model. John Wiley and Sons, Chichester, D0l:10.1002/9780470513200. [30] Shinozuka, M., Deodatis, G. (1991). Simulation of stochastic processes by spectral representation. Applied Mechanics Reviews, vol. 44, no. 4, p. 191-204, D0I:10.1115/1.3119501. [31] Clough, R.W., Penzien, J. (1993). Dynamics of Structures, McGraw-Hill, New York. [32] Liu, S.C. (1970). Evolutionary power spectral density of strong-motion earthquakes. Bulletin of the Seismological Society of America, vol. 60, no. 3, p. 891-900. [33] Qian, S. (2002). Introduction to Time-Frequency and Wavelet Transforms. Prentice Hall, New Jersey. [34] Kijewski-Correa, T., Kareem, A. (2006). Efficacy of Hilbert and wavelet transforms for time-frequency analysis. Journal of Engineering Mechanics, vol. 132, no. 10, p. 1037-1049, D0I:10.1061/(ASCE)0733-9399(2006)132:10(1037). [35] Spanos, P.D., Giaralis, A., Politis, N.P., Roesset, J.M. (2007). Numerical treatment of seismic accelerograms and of inelasticseismic structural responses using harmonic wavelets. Computer-Aided Civil and Infrastructure Engineering, vol. 22, no. 4, p. 254-264, D0I:10.1111/j.1467-8667.2007.00483.x. [36] Beck, J.L., Papadimitriou, C. (1993). Moving resonance in nonlinear response to fully non-stationary stochastic ground motion. Probabilistic Engineering Mechanics, vol. 8, no. 3-4, p. 157-167, D0I:10.1016/0266-8920(93)90011-J. Nonlinear MDOF system Survival Probability Determination Subject to Evolutionary Stochastic Excitation 451