A CONDITIONAL STOCHASTIC PROJECTION METHOD APPLIED TO A PARAMETRIC VIBRATIONS PROBLEM

. Parametric vibrations can be observed in cable-stayed bridges due to periodic excitations caused by a deck or a pylon. The vibrations are described by an ordinary differential equation with periodic coefficients. The paper focuses on random excitations, i.e. on the excitation amplitude and the excitation frequency which are two random variables. The excitation frequency L ω is discretized to a finite sequence of representative points, i ,L ω . Therefore, the problem is (condi-tionally) formulated and solved as a one-dimensional polynomial chaos expansion generated by the random excitation amplitude. The presented numerical analysis is focused on a real situation for which the problem of parametric resonance was observed (a cable of the Ben-Ahin bridge). The results obtained by the use of the conditional polynomial chaos ap-proximations are compared with the ones based on the Monte Carlo simulation (truly two-dimensional, not conditional one). The convergence of both methods is discussed. It is found that the conditional polynomial chaos can yield a better convergence then the Monte Carlo simulation, especially if resonant vibrations are probable.


Introduction
A dangerous phenomenon of parametric resonances can be observed among others in bridge cables during a strong wind gusts and sometimes due to a traffic load. It happens due to periodic movements of a cable support, which are usually caused by deck or rarely by pylon vibrations. There are several well-documented case studies in which parametric vibrations were observed in bridge cables, like Guazu Bridge in Argentina, Belgian bridges -Ben-Ahin and Wandre, Annacis Bridge in Canada, Faro Bridge in Denmark, Helgoland Bridge in Norway, Erasmus Bridge (Netherlands), Burlington (USA) or Second Severn Crossing (GB) Bridges. In each of the reported cases, the amplitude of vibrations increased exponentially; that is a very specific feature of parametric resonant vibrations which allows to distinguish this phenomenon from typical resonant vibrations with the amplitude increasing linearly (in absence of damping). Therefore, parametric resonant vibrations are much more dangerous than standard load-excited resonances. Moreover, parametric vibrations appear usually when an excitation frequency is near to doubled eigenfrequency -the instability zone is the largest in that case (Nayfeh, Mook 1995).
Relevant case studies of structures and the observed parametric vibrations are presented in Caetano (2007), Lilien and Pinto Da Costa (1994), Royer-Carafangi (2003). The presented paper focuses on a linear differential equation with random periodic coefficients. The Hermite-type polynomial chaos expansion is supported by the stochastic projection method (the Galerkin approach) and is implemented to a case of two input random variables. The convergence of results is presented in a context of the Monte Carlo simulation.

Polynomial chaos approximation -theoretical background, literature review
Generally, probabilistic analyses can be classified into two categories: methods using a statistical approach and methods based on a non-statistical approach. The statistical approach includes first of all the widely used Monte Carlo simulation (MC). Statistical methods are easy to implement but sometimes -in cases of complicated mathematical models -not efficient. Therefore, research effort has been made in developing non-statistical methods, like for instance perturbation techniques, usually addressed to "small uncertainties". Ghanem and Spanos proposed an innovative approach called polynomial chaos method (PC) and applied it successfully to various problems of mechanics (Ghamen, Spanos 1991); it is a situation where a stochastic process with finite variance can be presented in form of infinite convergent series, based on polynomials of standard Gaussian variable ξ (with zero mean value and standard deviation equal to 1).
The PC method separates stochastic components from deterministic ones for both -input data and random solutions. Random components of solutions can be calculated independently from the deterministic ones. If a random variable is a parameter taken into consideration, it can be expressed as a series of polynomials of the standard Gaussian variable (or variables) with series coefficients which are real numbers. For stochastic processes, the series of polynomials of the standard Gaussian variable are similar as for random variables, but the coefficients are time-dependent functions. The same is true for random fields but the coefficients are position-and/or timedependent. Clearly, some regularity conditions are required for both data and solutions, usually not restrictive in engineering applications.
The year 1991 -when Ghanem and Spanos applied the PC in conjunction with the finite element method (Ghamen, Spanos 1991) -was the turning point in the method development. In recent years, the PC method has been applied to a wide spectrum of stochastic problems, for instance to the theory of elasticity and plasticity (Ghamen, Spanos 1991;Anders, Hori 1999), stochastic dynamics (Ghanem, Spanos 1993;Sarkar, Ghanem 2002), wave propagation in random media (Manolis, Karakostas 2003), stochastic Stokes equations (Qu, Xu 2013) or sensitivity analysis (Sudret 2008;Crestaux et al. 2009;Blatman, Sudret 2010). In the mentioned examples, the application of the PC expansion presents the evidence of a good calculation efficiency. It refers also to relatively large random fluctuations, far beyond the scope of the simplest perturbation techniques. The success of this approach yields also from a conjunction with the Galerkin approximation, i.e. the best (least square) approximation in abstract Hilbert spaces.
Traditionally (Ghamen, Spanos 1991), the PC representations are based on the Hermite polynomials of the standard Gaussian variable ξ. However, in the case of complex nonlinear problems and non-Gaussian distributions of the input data or the system response, it is better to replace the Hermite polynomials by a more suitable function base. If the input data is very close to the Gaussian random variable then the linear approximation (with only two terms) can be accepted but the response expansion is usually much longer; and opposite, a change of the function base can simplify the response approximation though the simple Gaussian input will have a less accurate polynomial approximation. In other words, the "optimal" selection of the function base is problem-dependent. The so called generalized polynomial chaos expansions are addressed by Xiu and Karniadakis (2002a, 2002b, 2003, Le Maitre et al. (2004) and others. In recent years, the generalized polynomial chaos was successfully applied to solve the problems of diffusion and heat flow (Xiu, Karniadakis 2002a, 2003, or problems connected with interaction between structure and fluid (Xiu, Karniadakis 2002b;Lucor, Karniadakis 2004), etc.
Recently, a new insight into the background of the PC method was presented by Karniadakis (2005, 2006). The authors formulate a so called multielement generalized polynomial chaos method to deal with long-term integration and discontinuities in stochastic differential equations. On the other hand, Blatman and Sudret proposed the so called sparse polynomial chaos approach (Blatman, Sudret 2010, 2011. To decrease the computational cost of classical PC schemes in the case of large size of the input random vector, the authors presented a method that builds a sparse PC expansion.

A stochastic scalar product
Consider a sequence of uncorrelated random variables ) (γ k B , which have zero expected values and finite vari- , and assume that 1 const . It is a well-established fact that the expected value operator E[.] defines a scalar product .,. as follows: , which introduce the stochastic basis, makes it possible to expand every finite-variance random variable ) (γ X in a series (Ghamen, Spanos 1991):

The stochastic projection method
The stochastic projection method is the best mean-square approximation based on the distance defined by the scalar product. Usually, explicit form of the random variable X is unknown thus it is recommended to use its approximation in the form of the finite expansion (1), so for It is a general fact that the best mean-square m-dimensional approximation of X coincides with m-truncation of the orthogonal expansion (1). Indeed, for: the mean-square distance: can be calculated from the following expression:

The Hermite polynomials of the Gaussian random variable
The traditional Hermite polynomials Hek(x) of the real variable x are orthogonal with the standard Gaussian density ) (x ϕ taken as the weight function: The formula (4) is a simple prototype for the stochastic orthogonality condition, if the stochastic basis is introduced as ξ . This happens because -by definition -for such stochastic basis: Estimation of the approximation error (and the convergence rate) is strongly influenced by the weight function ) (x ϕ . The stochastic Hermite approximation is reliable for "small" values of ξ and is less reliable for "large" fluctuations -which are, however, less probable in case of the random variable N(0; 1) ξ ∈ .

Application of the polynomial chaos approximation
A PC-implementation is presented for a differential equation with periodic coefficients which yields from the equation of motion.

Equation of motion of a bridge cable with stochastic excitation
Transversal vibrations w(x,t) of a cable, i.e. in the y-direction in Figure 1, due to a vertical kinematic excitation ∆(t), are the physical prototype of the mathematical model. In the considered example the excitation is applied only in form of kinematic excitation with in-plane components. For that case, only in-plane cable vibrations occur. That is why only a two-dimensional dynamic problem is considered in the paper. The equation of motion (6) is derived in Brzakala (2013): where: E -the Young modulus of the cable; P -crosssectional area of the cable; m -cable mass per unit length; w c -damping coefficient for the transversal vibrations ξ ω ⋅ ⋅ ⋅ = 1 2 m c w ; 1 ω -natural eigenfrequency for the first symmetric mode of the cable in the transverse direction, for small values of the Irvine parameter it can be assumed as , like for the taut string (Irvine 1981;Caetano 2007); ξ -fraction of the critical damping; z(x) -assumed function of a sag for the inclined cable in case of the static load (dead weight of the cable as well as an initial tension force) ( ) The kinematic excitation ) cos( ) is the reason of periodic changes of coefficients. Positiondependent function z(x) and Hst are known (pre-defined). The solution w(x,t) is assumed in the separated variables, see Georgakis and Taylor (2005), as the following approximation: The second component in the formula (7) can be also assumed as a linear function, not the parabolic one (El Ouni, Ben Kahla 2012). However, Georgakis and Taylor propose the parabolic approximation which involves the sag changes due to a deck movement (Georgakis, Taylor 2005). Note that both concepts are approximate.
The single degree of freedom model is considered because for cables with a small sag, excited only parametrically, it is impossible to get complex vibrations, i.e. the ones consisting of several forms; so called combination resonance does not appear (Takahashi 1991). Selection of the excitation frequencies should be such that the parametric vibrations with the first form appear. They are close to the doubled first eigenfrequency (Nayfeh, Mook 1995). To get other forms of vibrations, higher excitation frequencies are necessary.
Probabilistic analyses of the random response of the system due to random fluctuations of the kinematic excitation (two random variables ∆ and L ω ) are the objective of the paper. The PC expansion and the MC simulation are applied to solve this problem.
Practically justified cases of small displacements are considered, so nonlinear terms in equation of motion are neglected. Due to the approximation (7) where: is a kinematic excitation due to the deck movement; Si -ten deterministic constants depending on geometrical or material parameters of the cable: cos sin / ( ) 0.167 cos / 2 cos / .

Example -deterministic and random variables
The deterministic input parameters are presented in Table 1.
Only positive values of the excitation amplitude and positive frequencies are possible, that is why the twodimensional log-normal probability distribution is assumed. It is assumed that the dimensionless random variab- are not correlated. The Y2 domain taken into consideration is truncated to finite values -from zero to the tripled expected value of the variable instead of the whole positive part of the Y2-axis. Having in mind this truncation, the value of the total probability equals 0.9999 which can be assumed as 1. The values of Y1 are not truncated. Note that the coefficients of variation µ σ / equal 0.3 thus they are not small. . Therefore, the maximal in-plane displacement of the cable (in the middle of the cable span) becomes also a stochastic process, ) , ( ) ( γ t w t w = . Deterministic Eqn (8) is linear but it is nonlinear from the stochastic point of view. Random parameters are present in the function of the kinematic excitation as well as in the solution w(t), but not in any explicit form.
The analysis uses conditional probability distributions and conditional moments assuming a sequence of deterministic values of kinematic excitation frequencies i L, ω . For each fixed i L, ω , the amplitude of kinematic excitation ∆ is the only one random variable which is considered. It means that Eqn (8) is solved repeatedly for discretized values of the excitation frequencies i L, ω . Making use of the one-dimensional solutions and their (conditional) probability distributions, random solution of the bi-variate stochastic problem can be found using the formula of total probability (or total moments), like ] Excitation frequencies i L, ω for each step are calculated as the probability mass center for each selected sub-domain (zone). There are 10 selected zones with the same probability . The idea of this partition is presented in Figure 2, in a simplified case of three zones. Ten coordinates of the right end of the selected zones (as in Fig. 2

The polynomial chaos approximation of the log-normal random variable
The input log-normal variable Y1 = ∆/l can be expressed in the form of the exponential function: where: ξ is a standard Gaussian random variable and ln ln For all fixed values of Y2 (ten zones with the same probability ), the expansion based on the Hermite polynomials is used to approximate the random variable Y1: etc. (Ghamen, Spanos 1991). The coefficients i ∆ β are calculated using the scalar product and the orthogonality of the Hermite polynomials: The integration in (13)

Solution of the one-dimensional stochastic problem
Since the random input parameter Y1 is in the form of (11), the solution is a stochastic process ) , ( ξ t w approximated as follows: Generally, in the m-dimensional subspace, the PC expansion leads to the approximate equation resulting from (8): where const L = ω .
In the approximated approach, a random residuum m ε should be as small as possible. It is assumed that the residuum is perpendicular to each base vector of the mdimensional subspace, so each scalar product ) ( He , ξ ε k m should be zero for 0, 1, ..., . k m = Note that only zero element is perpendicular to all base vectors, so the approximation becomes better if the subspace dimension m increases. The solution method justifies the term "projection" because the approximate solutions belong to subspaces spanned on the first m base vectors. The multiplication of both sides of Eqn (15) by the polynomials ) ( He ξ k and next the application of the expected value operator to both sides of (15) leads to ) 1 ( + m -dimensional system of ordinary differential equations (Ghamen, Spanos 1991). In this system, some coefficients ikmn e , ikm e appear, for example ; they are natural numbers, independent of the considered boundary problem, and they can be obtained by the use of the Gaussian integrals or the MC simulation. The search for the 1 + m unknown functional coefficients ) (t w i can use the standard Runge-Kutta method of the numerical integration available in Mathematica® software. Zero initial conditions are assumed in this situation. Finally, the solution of the random problem can be assembled by multiplication of the obtained functions ) by the orthogonal Hermite polynomials, as in (14).

Solution of the two-dimensional stochastic problem
Bearing in mind the solutions obtained for the onedimensional stochastic problem, it is possible to extend the analysis to the two-dimensional one where both ∆ and L ω are random. From a mathematical point of view, such a generalization can use a product-base of the Hermite polynomials, like ) ( He ) ( He where two generating Gaussian variables are used. Alternatively, a combination of the expansion with a discretization can be also effective. For the estimation of first two moments, the following expression yield from the total probability formula: ),

Polynomial chaos expansions -results of the numerical calculations
Eqn (15) is solved 10 times for the ten deterministic values of i L, ω , calculated in Section 3.2. Next, the solution for two-dimensional stochastic problem is assembled as in (16). The solutions based on the PC expansions are compared with the results yielding from the direct MC simulation.

Results of the one-dimensional stochastic problem
The amplitude ∆ of the kinematic excitation is assumed to be the log-normal random variable, frequencies of the kinematic excitation are considered as deterministic values. Solutions -in the form of the variance -are presented in Figure 3 and Figure 4 for two selected excitation frequencies L ω . The first value 1 L 2ω ω ≈ corresponds to the resonant vibrations (Fig. 3), the second one presents stable solutions (Fig. 4). The value of 1 ω is constant, presented in Table 1. Note that the expected value equals to (4) and (14) thus it is not sensitive to the length of the polynomial chaos expansion (parameter m). For some excitation frequencies, a strong amplification of vibrations in the time domain appears, i.e. for the excitation frequencies close to the doubled natural eigenfrequency of the cable 1 ω (the instability zone). For that case the solution based on the 2-term approximation strongly differ from the one based on the 6-term approximation (Fig. 3).
However, a much better convergence is observed for a stable solution (Fig. 4), so the linear stochastic PC approximation (i ≤ 1) is sufficient but only far away from the resonant frequency (Fig. 4). For resonant vibrations, a more extended expansion is required, compared to the case of stable solutions.

Results of the two-dimensional stochastic problem
The first two moments of the two-dimensional stochastic problem for the random variables L ω and ∆ are calculated with the use of (16) according to the pre-selected ten values of  The results of the numerical calculations yielding from the PC approximation are checked using the MC simulation for both situations presented in the paragraph 4, i.e. for one and two log-normal input variables, respectively. The first six Hermite polynomials are used (m = 5).

Validation of the one-dimensional stochastic problem
For the first case, excitation frequency is far away from the resonant zone, so the stable vibrations are analyzed for 1 L 31 . 1 ω ω = (Fig. 6). The second case refers to the resonant vibrations, excitation frequency is close to the doubled natural frequency of the system, 1 L 11 . 2 ω ω = (Fig. 7).  The variance of the solutions obtained by the use of the PC approximation and the MC simulations is presented in Figure 8. in Figure 8 reveal the same tendency. For both, the PC and the MC results, the variance strongly increases in time. Phase shifting between the solutions may be caused by the different approaches used for describing the random character of the variable L ω , i.e. by a discretization technique (a continuous spectrum of the excitation frequencies is considered in the MC simulation).
In the presented example, the PC expansion is more efficient (faster) than the MC simulation. A good convergence of the variance is reached for 20000 MC simulations, compared to only 6 polynomials in the PC expansion. For the solutions outside the resonant zone, a better convergence of the PC approximation is still more evident.

Summary and general conclusions
The parametric vibrations belong to one of the less recognized parts of the traditional structural dynamics. However, the parametric resonance is often reported as a potentially significant danger, mostly in bridge cables. The projection method based on the finite dimensional polynomial subspace can be an effective alternative to popular perturbation methods, the method of moments or to the Monte Carlo simulation.
The presented examples proved that the polynomial chaos expansion is more efficient (or simply faster) than the Monte Carlo simulation, especially if the resonant vibrations appear. Furthermore, the advantage of the polynomial chaos expansion is obvious if the random variables of interest (including the response variables) are Gaussian or close to the Gaussian distribution. Indeed, a very small number of expansion terms is required, say m ∼ 1÷2. On the other hand, the Monte Carlo simulation is less insensitive to probability distributions of the input variables as well as the increasing number of the input random variables. Furthermore, the numerical implementation of the Monte Carlo simulation is also much easier.
The used conditional probabilities and conditional moments can reduce the dimension of the model (number of considered random variables) with no significant loose of accuracy. However, the expansion method is also effective for multidimensional spaces with polynomial product bases. Note also that it is much easier to obtain results in the form of the probability histograms, by making use of the pre-calculated deterministic wi(t) and the Monte Carlo simulations -applied this time only to ) ( He ξ i .