A Hybrid Approach for the Random Dynamics of Uncertain Systems under Stochastic Loading

Abstract

This paper presents a hybrid Galerkin/perturbation approach based on Radial Basis Functions for the dynamic analysis of mechanical systems affected by randomness both in their parameters and loads. In specialized literature various procedures are nowadays available to evaluate the response statistics of such systems, but sometimes a choice has to be made between simpler methods (that could provide unreliable solutions) and more complex methods (where accurate solutions are provided by means of a heavy computational effort). The proposed method combines a Radial Basis Functions (RBF) based Galerkin method with a perturbation approach for the approximation of the system response. In order to keep the number of differential equations to be solved as low as possible, a Karhunen-Loève (KL) expansion for the excitation is used. As case study a non-linear single degree of freedom (SDOF) system with random parameters subjected to a stochastic windtype load is analyzed and discussed in detail; obtained numerical solutions are compared with the results given by Monte Carlo Simulation (MCS) to provide a validation of the proposed approach. The proposed method could be a valid alternative to the classical procedures as it is able to provide satisfactory approximations of the system response.

1. Introduction

The analysis of systems with uncertain parameters is an increasingly important subject in the treatment of a plethora of engineering problems. The dynamic response of a mechanical system with uncertain parameters possesses probabilistic features which depend on the probability distribution of the system parameters, and uncertainty in geometry or mechanical properties could affect both response and reliability of a structure [1, 2]. It is then necessary to take into account the consequence of these uncertainties in the study and design as randomness of structural parameters may lead to large and unexpected excursions of the structural response with, in some cases, drastic reductions in structural reliability. The reliability analysis for such structures strongly depends on the variation of these parameters; consequently, a probabilistic approach is necessary for adequate reliability analysis. Uncertainties must be adequately described to be included in a mathematical model that can be conveniently solved. To this aim, probability theory has been extensively adopted as a general tool able to characterize structural uncertainty by means of random fields and random variables [3].

The most common technique employed to analyse systems with uncertain parameters is the Monte Carlo simulation (MCS). This technique allows the statistical evaluation of the system response based on a certain number of analyses with different values of the random parameters. This approach could be very computationally expensive since it requires, after generating a large number of samples of uncertain parameters, the solution of the corresponding deterministic (linear or nonlinear) problem. Random numbers, representing uncertain parameters with given correlation structure, are usually generated by digital simulation [4]. The main advantage of the MCS method is that the procedures adopted in a deterministic setting can be directly applied to solve the mechanical problem at hand. In general, however, especially for nonlinear systems, the simulation procedures are quite inefficient due to the large number of samples needed to guarantee accurate statistical results. Moreover, the implementation of MCS requires the complete probabilistic description of the uncertain quantities through their probability density functions, which are often unavailable. All these remarks make evident that the Monte Carlo-based structural analysis, even if still the most used among the stochastic analysis methods for structural problems, may become very heavy from a computational point of view as the number of structure degrees of freedom (DOFs) and the number of uncertain parameters increase. In some cases, the computational effort makes it inapplicable (especially in case of nonlinearities), so the MCS method is mainly used as a test for analytical approaches.

For these reasons, alternative procedures have been proposed in the literature. Methods based on perturbation of the stochastic quantities have been presented both for linear [5, 6] and nonlinear problems and have been applied in the solutions of several engineering problems, such as frame analysis, reliability analysis, and buckling of shells with random initial imperfections [7, 8]. These approaches employ the conventional finite element method along with perturbation techniques and allow the probabilistic characterization of the response in terms of mean and covariance functions once the first- and second-order moments of the uncertainties are assigned. These methods mainly consist of a direct approach using probabilistic, instead of statistic, theory. This is usually pursued both in static and dynamic settings by using, for instance, expansion methods, where the stiffness matrix of the structural problem is split into a deterministic part (obtained with the mean value of random parameters) and into a part which accounts for the fluctuation of the random variables about its mean value. In order to evaluate the probabilistic response, Taylor expansions or Neumann expansions [9–11] are adopted to avoid inversion of matrices depending on the random parameters. Most of them, as pointed out before, are based on perturbation techniques, so that the stochastic finite element (SFE) method is often identified as the classical finite element (FE) method coupled with a perturbation approach [12, 13]. Unfortunately, some of these approaches show a drawback in that they are less and less accurate as the level of the uncertainty of the parameters increases. Moreover, even if the uncertainty is low, they ensure accuracy only for the second-order statistics of the response. Consequently, they can be applied only in case of Gaussianity of the response, which is rarely the case, even if the basic uncertain parameters are modelled as Gaussian. In fact, due to the nonlinear relationship between the system response and the system (random) parameters, the response is usually strongly non-Gaussian, even for linear systems.

To overcome these drawbacks, hybrid approaches have been proposed in the past. One of these has been proposed by Ghanem [14] by coupling the Monte Carlo simulation and the spectral stochastic finite element method (SSFEM). The author proposes to use a polynomial chaos expansion of the system response that, after substitution in the original equilibrium equation, leads to residuals that are orthogonalised with respect to the expansion basis. The original problem is then reduced to an iterative solution of a linear system of equations. Nevertheless, even if published in late 1998, no convincing application of these ideas has been published so far, and the methods seem practically limited to linear problems. A method due to Liu et al. [9] shows how an estimation of the time history of first two moments for the structural response in a linear or nonlinear system is obtained. This work has been improved by Chiostrini and Facchini [15] where a stochastic input has also been taken into account. The first two moments of the response have been evaluated taking into account a Taylor expansion of the structural response centred on the mean value of random parameters. The method is efficient when the dependence between the response and the random parameters is approximately a polynomial of the same degree as the expansion. It is then necessary to take in account more flexible expressions for the approximation of the response.

In the present paper, to provide an insight into the statistical response variations of nonlinear systems, a computationally efficient approach is proposed that takes into account randomness both in mechanical parameters and loads. The approximation of the system response is accomplished by means of an expansion of radial basis functions (RBFs): a novel class of approximating functions is adopted to model the dependence of the response displacements upon the uncertain variables. Such expansion can be classified among the so-called radial basis neural networks with time-dependent coefficients [16, 17]. A stochastic process which can be described by this model allows a direct evaluation of its main characteristics leading to a remarkable accuracy level for the statistics of the response and, in general, for the probability density function (PDF). First, after formulating the problem, the approach is briefly drafted with the aim of introducing some useful quantities and notations; then, it is particularized with reference to the system under examination: a strongly nonlinear SDOF system with an uncertain parameter subject to a random load. The numerical treatment of the mechanical problem is discussed in detail. In particular, to develop the approach, two methods are mainly taken into consideration in the present work: the first is described in Liu et al. [9], and subsequently enhanced by Chiostrini and Facchini [15]. It can be classified as a perturbation method and makes use of sensitivity vectors to evaluate the first two moments of the response. As already mentioned, the efficiency of this method is closely linked to the degree of the expansion and of the dependence of the response on the random parameters. This condition is difficult to check, and this fact can lead to great errors in the results. The second one makes use of a Galerkin approach with radial basis trial functions: it is more complex and leads to the solution of a great number of differential equations than the former, but the results are much more satisfactory (Facchini [18] and Betti et al. [4]). Both methods are briefly reviewed for convenience in the following sections. The proposed procedure can be seen as a compromise between such different situations: it provides a satisfactorily accurate solution by means of not so complex calculations. As shown in the numerical application, the features of the improved perturbation method are preserved even in the nonlinear dynamic analysis. The results highlight a good level of accuracy as is revealed by comparisons with the MCS results. As it provides a good compromise between computational cost and accuracy, the proposed hybrid approach could be an effective tool for the solution of a wide range of problems involving uncertain parameters.

2. Theoretical Remarks

The mechanical system considered in this study is assumed to be described by an equation of motion where a nonlinear restoring function characterizes the nonlinear structural behaviour; load is assumed as a random process:𝐌𝐱̈+𝐠(𝐱,𝐱̇,𝐛)=𝐅(𝑡).(2.1) In (2.1), 𝐱 is the vector of the system degrees of freedom (DOFs), 𝐌 is the mass matrix, and 𝐠 is the nonlinear restoring function; dot indicates differentiation with respect to time. The nonlinear function 𝐠 is assumed to be affected by randomness, that is, it changes according to the variability of some random parameters (e.g., Young’s modulus of the structural material, dimensions of the structural members, etc.). Such uncertain parameters are grouped in the vector 𝐛. The excitation 𝐅(𝑡) is a random process which can be described in terms of a parametric model by means of its mean vector 𝐅_𝑚 and some zero-mean random vectors 𝜶^(𝑘) with proper covariance structure: 𝐅(𝑡)=𝐅_𝑚+_𝑘𝜶^(𝑘)𝜒_𝑘(𝑡).(2.2) Eventually, 𝜒_𝑘(𝑡) denotes a function of time to be selected by the user.

The probabilistic properties of the response function 𝐱(𝑡) can thus be obtained by evaluating the dependence, at each instant, of the response function itself with respect to the random parameters (i.e., the analysis of the function 𝐱(𝐛,𝜶,𝑡)). This analysis, as discussed in the introduction, can be approached by several methods. Two of them are, for convenience, next reported: the perturbation method and the Galerkin approach. In particular, the Galerkin method is next reviewed considering an expansion of radial basis functions (RBFs).

2.1. Brief Review of Perturbation Method

The perturbation method approximates the dependency of the system response on the random vectors 𝜶^(𝑘) and 𝐛, under the assumption of their independency, by means of a series expansion:𝑥_𝑘𝑡,𝐛,𝜶^(ℎ)=𝑥_𝑘𝑡,𝝁_𝐛,𝝁_𝜶+𝜕𝑥^𝑘𝑡,𝝁_𝐛,𝝁_𝜶𝜕𝐛𝐛−𝝁^𝐛+𝜕𝑥^𝑘𝑡,𝝁_𝐛,𝝁_𝛼𝜕𝜶^(𝑖)+1𝜶^(𝑖)−𝝁_𝜶(𝑖)2𝐛−𝝁^𝐛^𝑡𝜕²𝑥_𝑘𝑡,𝝁_𝐛,𝝁_𝜶𝜕𝐛²+1𝐛−𝝁_𝐛2𝜶^(𝑖)−𝝁_𝜶(𝑖)^𝑡𝜕²𝑥_𝑘𝑡,𝝁_𝐛,𝝁_𝜶𝜕𝜶^(𝑖)𝜕𝜶^(𝑗)𝜶^(𝑗)−𝝁_𝜶(𝐣),(2.3) where 𝜇 is the mean of the random parameters. A repeated index implies summation over it. The response statistics can thus be obtained taking expected values of (2.3). The derivatives of the system response with respect to the random coefficients can be evaluated differentiating the equation of motion (2.1) (see e.g., Liu et al. [9] and Chiostrini and Facchini [15] for details). The method provides effective results when the dependence of the response on the random parameters is quadratic; a severe error can occur if this dependence is more complex (see [19]).

2.2. Brief Review of Galerkin-RBF Approach

The Galerkin approach describes the dependence of the response on the random parameters by means of time-dependent vectors 𝐰^(𝑘)(𝑡) and trial functions 𝜑_𝑘(𝐛,𝜶^(ℎ)):𝐱𝑡,𝐛,𝜶^(ℎ)≅𝐱̃𝑡,𝐛,𝜶^(ℎ)=^𝑁𝜑_𝑘=1𝐰^(𝑘)(𝑡)⋅𝜑_𝑘𝐛,𝜶^(ℎ),(2.4) where 𝑁_𝜑 is the number (chosen by the user) of the trial functions. A variational approach (Facchini [18] and Betti et al. [4]) can be established, to obtain the solution. Denoting by the symbol 𝐸[⋅] the expected value operator, the equation that allows to obtain the solution can be written as follows:𝐸𝐞𝜑_𝑘𝐛,𝜶^(ℎ)=0∀𝑘=1⋯𝑁_𝜑,(2.5) where the error function 𝐞 is defined as follows:𝐞(𝑡)=𝐌𝐱̈̃+𝑔𝐱̃,𝐱̇̃,𝑡−𝐅(𝑡).(2.6) In the classical Galerkin approach, the trial functions 𝜑_𝑘(𝐛,𝜶^(ℎ)) in (2.4) are fixed, and therefore, the solution is completely defined by (2.5). On the other hand, radial basis functions (RBFs) are dependent on some more parameters (the centre and the decay coefficients), which enable a further optimization of the procedure. A number of different RBFs can be found in the literature (see f.i. Broomhead and Lowe [16], Haykin [20], Poggio and Girosi [17], Gotovac and Kozulic [21], Mai-Duy and Tran-Cong [22]). In the present study, Gaussian bell-shaped RBF will be used. Their form is as follows:𝜑_𝑘𝐛,𝜶^(ℎ)=𝜙‖‖𝐯−𝐜^(𝑘)‖‖𝜎_𝑘;𝐯=𝐛𝜶⁽¹⁾⋯𝜶^(𝑛)^𝑡,(2.7) with 𝜙(𝑟)=exp−𝑟²,(2.8) where 𝐜^(𝑘) is the centre of the 𝑘th function and 𝜎_𝑘 is its decaying parameter. It is to be noted that other kinds of RB functions work similarly well, as all of them are universal approximators [21]. In other words, under general conditions, arbitrarily complex functions can be approximated to a desired tolerance with a sufficiently high number of RBF with proper parameters [18, 23]. The variability of the centres and the decay parameters of the RBFs induces a variability of the RBFs themselves, so that another set of equations can be derived from the variational approach and an optimization procedure can be performed (see Facchini [18] for details).

2.3. The Proposed Hybrid Approach

The main idea of the approach herein proposed is to combine the two previously cited methods to join their advantages. Extremely summarizing the perturbation approach, from a numerical point of view, requires both a smaller amount of memory and a computational effort. The Galerkin method is more accurate and its accuracy can be controlled by means of error functions. Moreover, the perturbation method gives satisfactory results when it is used to model the dependence of the response on the random parameters of the forcing process, while modelling the dependence of the response on the random parameters of the mechanical system generally requires the use of the Galerkin approach [24]. Consequently the hybrid approach aims to combine the advantages of the two methods modelling the two different types of dependency in the system equation. Based on this idea the forcing process and the structural response are expressed as follows:𝐹(𝑡)=^𝑁𝐹_ℎ=1𝐱𝑡,𝐛,𝜶_(ℎ)𝜶^(ℎ)𝑎_ℎ(𝑡),(2.9)≅𝐱̃𝑡,𝐛,𝜶^(ℎ)=^𝑁𝑤_𝑘=1𝐰^(𝑘)𝑡,𝜶^(ℎ)𝜑_𝑘(𝐛).(2.10) Eventually, an expansion like (2.3) is used:𝐰^(𝑘)𝑡,𝜶^(ℎ)=𝐰^(𝑘)𝑡,𝝁_𝛼+𝜕𝐰^(𝑘)𝑡,𝝁_𝛼𝜕𝜶^(ℎ)+1𝜶^(ℎ)−𝜇_𝛼(ℎ)2𝜶^(𝑚)−𝜇_𝛼(𝑚)𝜕²𝐰^(𝑘)𝑡,𝝁_𝛼𝜕𝜶^(𝑚)𝜕𝜶^(ℎ)𝜶^(ℎ)−𝜇_𝛼(ℎ).(2.11) Then, the solution of original problem of (2.1) can be solved as a set of coupled deterministic nonlinear equations, and the statistics of the response process (mean value and correlation function) can be directly evaluated by means of the following:𝐸[𝐱(𝑡)]=^𝑁𝑤_𝑘=1𝐑𝐸𝐰^(𝑘)(𝑡,𝜶)𝐸𝜑_𝑘(𝐛)_𝑋𝑋,(𝑡,𝑠)=−𝐸𝐸𝐰𝐰^(ℎ)(_(ℎ)𝑡,𝜶)𝐰^(𝑘)(𝑠,𝜶)𝐸𝜑_ℎ(𝐛)𝜑_𝑘(𝐛)(𝑡,𝜶)𝐸𝐰^(𝑘)(𝑠,𝜶)𝐸𝜑_ℎ(𝐛)𝐸𝜑_𝑘(𝐛).(2.12) Specific details of the proposed hybrid approach are reported next, where the method is applied to illustrative case studies.

3. Numerical Examples

To explain and discuss feasibility of the proposed hybrid approach, several case studies are considered. As a first step, the deterministic loading has been considered with the aim to compare the result of the approach with the existing method; next, the case with stochastic loading is reported. The cases are discussed detailing each step of the proposed approach.

3.1. Deterministic Loading

3.1.1. Duffing Oscillator

As a first application, a SDOF Duffing oscillator subjected to a time-dependent forcing process is considered. The equation of motion of the Duffing oscillator is written as follows:𝑚𝑥̈+𝑐𝑥̇+𝑘𝑥+𝛽𝑥³=𝑓(𝑡),(3.1) where the forcing process is assumed as a harmonic load:𝑓(𝑡)=𝑎⋅sin𝜔𝑡.(3.2) Uncertainties of the system are assumed to be concentrated in the stiffness 𝑘 of the system; the nonlinear restoring term in (3.1) can be written as follows:𝑔(𝑥,𝑥̇,𝑡;𝑘)=𝑐𝑥̇+𝑘𝑥+𝛽𝑥³,(3.3) and, consequently, the equation of motion in general form, accordingly to (2.1), can be rewritten as follows:𝑚𝑥̈+𝑔(𝑥,𝑥̇,𝑡;𝑘)=𝑓(𝑡).(3.4) Based on the proposed approach (2.10) can be particularised as follows: 𝑥(𝑡,𝑘)=^𝑁𝜑_𝑗=1𝑤_𝑗(𝑡)𝜑_𝑗(𝑘),(3.5) where 𝐰(𝑡) is a time-dependent vector and 𝜑_𝑗(𝑘) are assumed as Gaussian-shaped RB functions depending only on the random stiffness 𝑘 of the mechanical system. Considering (2.12), the solution of the original problem of (3.1) can be moved to the solution of a set of coupled deterministic nonlinear equations, and the statistics of the response process (mean value and correlation function) can be directly evaluated by means of the following: 𝜎𝜇_𝑥(𝑡)=𝐸[𝑥(𝑡,𝑙)]=_𝑗𝑤_𝑗(𝑡)𝐸𝜑_𝑗(𝑘)_2𝑥,𝑡₁,𝑡₂=_𝑖,𝑗𝑤_𝑖𝑡₁𝐸𝜑_𝑖(𝑘)−𝐸𝜑_𝑖(𝑘)𝜑_𝑗(𝑘)−𝐸𝜑_𝑗(𝑘)𝑤_𝑗𝑡₂.(3.6) To analyse the problem, the following parameters have been assumed: 𝑎=300 N; 𝜔=15 rad/sec; 𝐸_𝑘=800 N/cm; 𝛽=0.021/(cm²); 𝑐=2.85 Nsec/cm; 𝑚=200 Nsec²/cm. To check the efficiency of the approach, several tests with increasing number of RBFs (𝑁_𝜑) have been carried out considering 25 sec of analysis, and the results have been then compared with MCS. A set of 500 simulations (direct numerical integration of the differential equation of motion (3.1)) has been carried out by generating a vector of random values of the variable 𝑘 (assuming a coefficient of variation of 5%).

As a first case, 2-RBFs (𝑁_𝜑=2) have been considered; at each time step, a system of 𝑁_𝜑×𝑁 equations has to be solved. The proposed approximation with 𝑁_𝜑=2 matches quite well the MC simulation with respect to the mean value of displacement and velocity (Figures 1 and 2). On the contrary, if the standard deviation of both displacement and velocity are plotted together with MCS results, it is possible to see that the approximation is not able to reproduce the response process (Figures 3 and 4). The difficulty of the 2-RBFs approximation to predict the system response can be highlighted analysing the system response behaviour at several instants (Figure 5).

Figure 1 

Mean value of displacement.

Figure 2 

Mean value of velocity.

Figure 3 

Standard deviation of displacement.

Figure 4 

Standard deviation of velocity.

(a)

(b)

(c)

(d)

Figure 5 

System response (displacement 𝑥) at different instants (5, 10, 15, and 20 sec).

Based on the above results, an approximation with 5-RBFs (𝑁_𝜑=5) has been considered. In figures 6, 7, 8, and 9, the results of the approximation are again compared with MC simulation. In this case, the approximation function with 𝑁_𝜑=5 is able to match the results of the simulation also with respect to the standard deviation of both displacement and velocity. The analysis of the system response behaviour at several instants (Figure 10) shows the ability of the approximation to match the actual system response.

Figure 6 

Mean value of displacement.

Figure 7 

Mean value of velocity.

Figure 8 

Standard deviation of displacement.

Figure 9 

Standard deviation of velocity.

(a)

(b)

(c)

(d)

Figure 10 

System response (displacement 𝑥) at different instants (5, 10, 15, and 20 sec).

3.1.2. Mass Pendulum

As a second case study, a mass pendulum system is considered (Figure 11) with an aim to offer an effective comparison with the existing methods to solve dynamic stochastic problems. The strongly nonlinear SDOF system is a mass pendulum system under a time-space-dependent load acting on the system mass. The parameters characterizing the system are the truss length and the mass. Naming 𝜃(𝑡) the Lagrangian parameter, the equation of motion can be written as follows:𝑚𝑙𝜃̈+𝑚𝑔sin(𝜃)=𝐹(𝑡)cos(𝜃),(3.7) where 𝑙 is the pendulum length, 𝑚 is the system mass, 𝜃 is the rotation, 𝑔 the gravity acceleration, and 𝐹(𝑡) the external stochastic forcing load (“sin” and “cos” denote the standard trigonometric functions); dot superscripts indicate differentiation with respect to time.

Figure 11 

The mechanical nonlinear system.

It is assumed that 𝑙, the pendulum length, is an uncertain parameter; consequently, the rotation 𝜃 is depending on the length 𝑙, on the time 𝑡, and on the forcing process 𝐹(𝑡). As a first application term 𝐹(𝑡) in (3.7) has been assumed as a harmonic loading: 𝐹(𝑡)=𝐹₀sin(𝜔𝑡); 𝑚𝑙𝜃̈+𝑚𝑔sin(𝜃)=𝐹₀sin(𝜔𝑡)cos(𝜃).(3.8)

The stochastic dynamic problem described by (3.8) has been solved firstly by means of the perturbative approach, and results have been compared both with the proposed hybrid one (particularised to the case of deterministic load) and with the Monte Carlo simulation (MCS) assuming a set of 500 simulations.

Based on the perturbative approach, the system response is approximated by means of a series expansion with respect to the random parameters; the stochastic dynamic problem has been solved assuming following second-order Taylor series expansion of the solution about the mean values of the uncertain parameter:𝜃(𝑡,𝑙)≅𝜃𝑡,𝜇_𝑙+𝜕𝜃𝑡,𝜇^𝑙𝜕𝑙𝑙−𝜇^𝑙+12𝑙−𝜇^𝑙^𝑡𝜕²𝜃𝑡,𝜇_𝑙𝜕𝑙²𝑙−𝜇_𝑙.(3.9) Substituting (3.9) into (3.8) and collecting terms of the same order will yield a set of deterministic equations. In case of linear problems, this equation set has identical homogeneous parts subjected to different forcing terms, and great efficiency is achieved by the method as the system can be solved sequentially. In case of nonlinear problems, as the present case, this advantage is generally lost due to the resulting coupling of nonlinear terms; therefore, this set of equations must be solved simultaneously.

Figures 12 and 13, show the system response in terms of mean value and standard deviation obtained. Perturbative approach solutions are compared with the solution obtained by MCS. Analysing Figures 12 and 13 it is possible to observe that the solution obtained by the perturbative approach is able to match the actual system solution only up to about 20 s; afterwards, the approximation becomes unacceptable.

Figure 12 

Mean value of 𝜃: comparison between MCS and perturbative approach results (0–50 sec).

Figure 13 

Standard deviation of 𝜃: comparison between MCS and perturbative approach results (0–50 sec).

This is due to the fact that the accuracy of the results relies on the accuracy of the quadratic approximation adopted for the dependence of the response on the random parameter. Since the degree of fluctuation within this range is generally not known in advance (especially when system nonlinearity and time factors take effect) the methods can produce misleading results. In particular, the failure of the perturbative approach in the present case can be better explained analysing the dependency of the system response (rotation 𝜃) on the pendulum length 𝑙. Figure 14 shows that the pendulum response depends approximately in a linear manner on the length of the pendulum up to about 5 s. From 10 s onwards, the quadratic approximation becomes unacceptable, and this causes the highlighted errors. It is then clear that approximation (3.9) could provide acceptable results only when the dependency of the response process by the random parameters of the system is roughly quadratic, which results into a restriction on the typology of nonlinearities that affect the system behaviour or into the variability of uncertainties.

(a)

(b)

(c)

(d)

Figure 14 

System response (rotation 𝜃) at different instants (5, 10, 25, and 50 sec).

The same problem is next solved with the proposed approach (it is worthwhile noting that, in this case of deterministic loading, the hybrid approach coincides with the Galerkin/RBF); based on the assumption on deterministic loading, (2.10) is rewritten as follows: 𝜃(𝑡,𝑙)=^𝑁𝜑_𝑘=1𝑤_𝑘(𝑡)𝜑_𝑘(𝑙),(3.10) where 𝐰(𝑡) is a time-dependent vector, 𝜑_𝑘(𝑙) are functions depending only on the random parameters 𝑙 of the mechanical system. Functions 𝜑_𝑘(𝑙) are assumed to be Gaussian-shaped RBF [21]:𝜑_𝑘(𝑟)=exp−𝑟²2𝜎².(3.11) With these assumptions, (2.12) becomes 𝜎𝜇_𝜃(𝑡)=𝐸[𝜃(𝑡,𝑙)]=_𝑘𝑤_𝑘(𝑡)𝐸𝜑_𝑘(𝑙)_2𝜃,𝑡₁,𝑡₂=_ℎ,𝑘𝑤_ℎ𝑡₁𝐸𝜑_ℎ(𝑙)−𝐸𝜑_ℎ(𝑙)𝜑_𝑘(𝑙)−𝐸𝜑_𝑘(𝑙)𝑤_𝑘𝑡₂.(3.12) Figures 15 and 16 report the mean value and standard deviation of the rotation 𝜃, respectively, comparing the solution of the proposed approach with MCS; the proposed approximation is able to follow the system response even in the range where the perturbative approach fails. Furthermore, if the response of the pendulum is analysed at several time instants (Figure 17), the proposed method shows its capacity to approximate with confidence the actual system behaviour.

Figure 15 

Mean value of 𝜃: comparison between MCS and hybrid approach results (0–50 sec).

Figure 16 

Standard deviation of 𝜃: comparison between MCS and hybrid approach results (0–50 sec).

(a)

(b)

(c)

(d)

Figure 17 

System response (rotation 𝜃) at different instants (10, 20, 40, and 50 sec): comparison between MCS and hybrid approach (HA) approximation.

3.2. Stochastic Loading

3.2.1. Mass Pendulum

To detail the proposed approach with respect to the case of stochastic loading, the term 𝐹(𝑡) in (3.7) is assumed to be a stochastic load. Under very general conditions, the Karhunen-Loève (KL) expansion for the forcing process [15, 24] can be used:𝐹(𝑡)=^𝑁𝛼_ℎ=1𝛼_ℎ⋅𝑓_ℎ(𝑡),(3.13) where 𝛼_ℎ are the KL coefficients of the forcing process. The perturbation method is adopted to approximate the dependency of the system response on the random vector collecting the terms of KL expansion. Taking into account (2.9), the following series expansion is used for the generalized displacement:𝜃(𝑡,𝑙,𝜶)≅^𝑁𝜑_𝑘=1𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙),(3.14) where vector 𝜶 groups 𝛼_ℎ terms and 𝑤_𝑘 functions in (3.14) are expressed by the perturbative approach:𝑤_𝑘(𝑡,𝜶)=𝑤_𝑘𝑡,𝝁_𝜶+𝜕𝑤^𝑘𝑡,𝝁_𝜶𝜕𝛼_𝑖𝛼_𝑖+1−𝜇_𝛼𝑖2𝛼^𝑖−𝜇_𝛼𝑖𝜕²𝑤_𝑘𝑡,𝝁_𝜶𝜕𝛼_𝑗𝜕𝛼_𝑖𝛼^𝑗−𝜇_𝛼𝑗.(3.15) It is noteworthy to observe that also the possibility to use the KL expansion for the response process too would minimize the number of equations to be solved as this expansion is optimal [14]. Clearly, this is not possible since the covariance function of the system response is not known a priori.

Inserting (3.14) in (3.7), one obtains𝑚𝑙=𝐹(^𝑁𝜑_𝑘=1𝑡𝑤)̈c_𝑘o(s𝑡,𝜶)⋅𝜑_𝑘(_{𝑁𝜑𝑘=1}𝑙)+𝑚𝑔sin^𝑁𝜑_𝑘=1𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙)𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙).(3.16) and multiplying (3.16) by 𝜑_ℎ(𝑙) and rearranging, one obtains^𝑁𝜑_𝑘=1𝑚𝑙𝑤̈_𝑘(𝑡,𝜶)=⋅𝐹𝜑_𝑘(𝑡()𝑙)c⋅o𝜑s_ℎ(𝑙)+𝑚_{𝑁𝜑𝑘=1}𝑔sin^𝑁𝜑_𝑘=1𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙).(3.17) Denoting the expected value operator by the symbol 𝐸[⋅], evaluating the expected value of both members of (3.17), the following system of equations is obtained:𝑚𝐸𝑙𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)⋅𝑤̈_𝑘(𝑡,𝜶)+𝑚𝑔=𝐸𝐹s(i𝑡n)𝐸𝑤c_𝑘o(s𝑡,𝜶𝑤)_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)(𝑡,𝜶)⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙),(3.18) and assuming the following notation:𝑔𝑀_ℎ𝑘=𝑚𝐸𝑙𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)_ℎ,=𝑚𝑔𝐸sin𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙),𝐹_ℎ(𝑡,𝐰)=^𝑁𝛼_𝑗=1𝛼_𝑗⋅𝑓_𝑗(𝑡)⋅𝐸cos𝑤_𝑘(𝑡,𝜶)⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)=^𝑁𝛼_𝑗=1𝛼_𝑗⋅𝑓_𝑗(𝑡)⋅𝑇_ℎ(𝐰),(3.19) the original problem of (3.7) can be written in vectorial form as follows:𝐌𝐰̈+𝐠(𝐰)=𝐅(𝑡,𝐰),(3.20) where the effects of the load randomness have been evaluated by means of the perturbation method, and the effects of the randomness of the uncertain parameter 𝑙 are taken into account by means of the Galerkin approach. Using (3.15), evaluating (3.20) with respect to the mean value of the KL coefficients 𝛼_ℎ of the forcing process, and collecting terms, it is possible to obtain the solving systems of equations.

The first term in the expansion is the response calculated for mean value of the KL coefficients, so that𝑀𝐌𝐰̈+𝐠𝐰−𝐅(𝑡,𝐰)=𝟎,_ℎ𝑘𝑤̈_𝑘+𝑚𝑔𝑆_ℎ−^𝑁𝛼_𝑙=1𝜇_𝛼𝑙⋅𝑓_𝑙(𝑡)⋅𝐶_ℎ=0,(3.21) with𝐶𝑆_ℎ=𝐸sin𝑤_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)_ℎ,=𝐸cos𝑤_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙),(3.22) where the underbar denotes evaluation of corresponding quantities with respect to the mean value of the KL coefficients [9, 15]. The first derivatives—or, in other words, the sensitivity vectors—can be derived quite easily by differentiating with respect to the generic coefficient 𝛼_𝑖:𝐌𝜕𝐰̈𝜕𝛼_𝑖+𝜕𝐠𝐰𝜕𝐰𝜕𝜕𝛼𝐰_𝑖−𝜕𝐅𝑡,𝐰𝜕𝛼_𝑖𝑀=𝟎,_ℎ𝑘𝑤̈_𝑘,𝑖+𝑚𝑔𝐶_ℎ𝑚𝑤_𝑚,𝑖+⎡⎢⎢⎣^𝑁𝛼_𝑙=1𝜇_𝛼𝑙⋅𝑓_𝑙(𝑡)⎤⎥⎥⎦⋅𝑆_ℎ𝑘𝑤_𝑘,𝑖−𝑓_𝑖(𝑡)⋅𝐶_ℎ=0,(3.23) where𝐶𝑆_ℎ𝑚=𝐸sin𝑤_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)⋅𝜑_𝑚(𝑙)_ℎ𝑚,=𝐸cos𝑤_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)⋅𝜑_𝑚(𝑙).(3.24) And eventually a similar relation holds for the second derivatives:𝐌𝜕²𝐰̈𝜕𝛼_𝑖𝜕𝛼_𝑗+𝜕𝐠𝐰𝜕𝐰𝜕²𝐰𝜕𝛼_𝑖𝜕𝛼_𝑗+𝜕𝐰𝜕𝛼_𝑖𝜕²𝐠𝐰𝜕𝐰²𝜕𝐰𝜕𝛼_𝑗−𝜕²𝐅𝑡,𝐰𝜕𝛼_𝑖𝜕𝛼_𝑗𝑀=𝟎,_ℎ𝑘𝑤̈_{𝑘,𝑖𝑗}+𝑚g𝐶_ℎ𝑚𝑤_{𝑚,𝑖𝑗}+∑^𝑁𝛼_𝑙=1𝜇_𝛼𝑙⋅𝑓_𝑙(𝑡)⋅𝑆_ℎ𝑘𝑤_{𝑘,𝑖𝑗}−𝑚𝑔𝑆_ℎ𝑚𝑙𝑤_𝑙,𝑗𝑤_𝑚,𝑖+𝑓_𝑖(𝑡)𝑆_ℎ𝑘𝑤_𝑘,𝑗−𝑓_𝑗(𝑡)𝑆_ℎ𝑘𝑤_𝑘,𝑖+∑^𝑁𝛼_𝑙=1𝜇_𝛼𝑙⋅𝑓_𝑙(𝑡)⋅𝐶_ℎ𝑘𝑘𝑤_𝑘,𝑖𝑤_𝑘,𝑗=0,(3.25) with𝐶𝑆_ℎ𝑚𝑛=𝐸sin𝑤_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)⋅𝜑_𝑚(𝑙)⋅𝜑_𝑛(𝑙)_ℎ𝑚𝑛,=𝐸cos𝑤_𝑘⋅𝜑_𝑘(𝑙)⋅𝜑_ℎ(𝑙)⋅𝜑_𝑚(𝑙)⋅𝜑_𝑛(𝑙).(3.26) Equations (3.21), (3.23), and (3.25) constitute a system of 𝑁_𝜑+𝑁_𝜑×𝑁_𝛼+𝑁_𝜑×𝑁^2𝛼 equations (the degree of approximation by the hybrid Galerkin/perturbation approach) that allows to evaluate the solution of the original stochastic problem described by (2.1) by solving a set of coupled deterministic nonlinear equations. This system can be solved by means of the standard techniques available for deterministic equations. After the coefficients in the hybrid Galerkin/perturbation expansion have been calculated, they allow obtaining the realization of the response process by direct evaluation of (3.14)-(3.15). Mean value and standard deviation of system response process, taking into account the property of the KL, are calculated according to (2.12) as follows: 𝜇_𝑥(𝑡)==𝐸[𝑥(𝑡,𝑏)]=∑^𝑁𝜑_𝑘=1𝐸𝑤_𝑘(𝑡,𝛼)⋅𝐸𝜑_𝑘(𝑏)_𝑁𝜑_𝑘∑₌₁𝑤^𝑘𝑡,𝜇_𝛼+12∑^𝑁𝛼_𝑖,𝑗𝐸𝛼_𝑖−𝜇_𝛼𝑖𝛼_𝑗−𝜇_𝛼𝑗𝜕²𝑤_𝑘𝜕𝛼_𝑖𝜕𝛼_𝑗⋅𝐸𝜑^𝑘𝜎(𝑏)_2𝑥,=^𝑁𝜑_ℎ=1𝐸𝑤^2ℎ⋅𝐸𝜑^2ℎ−𝐸𝑤_ℎ²⋅𝐸𝜑_ℎ².(3.27) Based on the results discussed in previous paragraphs, the case of stochastic loading is analysed. Preliminarily a set of 500 simulations has been carried out by generating a vector of random values of the variable 𝑙 with assigned mean and variance. A Gaussian distribution has been assumed, considering a mean value 𝜇_𝑙=5 m and a coefficient of variation for the uncertain parameter of 5% and 10%. The response of the system has been obtained by direct numerical integration of the differential (3.7) governing the problem. The equations of motions were integrated in the case of a wind-pressure-type random excitation, with a Ficthl-Mc Vehil [25] spectral density function defined by𝑓𝑆_𝐹𝐹(𝑓)𝜎^2𝐹=𝐴𝑓1+𝐵𝑓^𝛽^𝛼𝐵𝐴==0.0216,𝐵𝛼_1/𝛽=1.3608,𝛽=,1.7034,Γ(𝛼)Γ(𝛼−1/𝛽)Γ(1+1/𝛽)=0.0879.(3.28) The excitation was described by 200 Karhunen-Loève eigenfunctions (Figure 18 reports a typical realization). Integration lasted for 50 seconds.

Figure 18 

A realization of the excitation process employing 200 terms in Karhunen-Loève expansion.

In order to estimate the quality of the approximation offered by the proposed hybrid Galerkin/perturbation approach, several cases are investigated where different values of the radial basis functions numbers and KL terms in (3.14) and (3.15) have been taken into account. As for the simulation task, a time-length of 50 seconds of the analysis has been considered. The outputs of the outlined procedure, the sensitivity coefficients, are obtained by integration of (3.21)–(3.23)–(3.25). Indicating with 𝑁_𝜑 the number of RBFs adopted and with 𝑁_𝛼 the number of KL terms, at each time steps a system of 𝑁_𝜑+𝑁_𝜑×𝑁_𝛼+𝑁_𝜑×𝑁^2𝛼 equations has to be solved. After integration, upper and lower bounds were obtained for both rotation and angular velocity of the pendulum by means of [9]:ULopwpeerrbboouunndd∶∶siunpf((𝜃𝜃))≈≈𝜇𝜇_𝜃𝜃+3𝜎_𝜃,−3𝜎_𝜃.(3.29)

As a first test, a case where wind load is approximated by 200 terms in KL summation (3.13) assuming a coefficient of variation of 5%, has been considered. Several tests with different number of RBFs have been taken into account.

Figures 19 and 20 report, respectively, upper and lower bounds for pendulum rotation and for angular velocity assuming 5 RBF. Results of the proposed approach compared with Monte Carlo simulation show that the level of approximation by radial basis functions is not satisfactory. Increasing the number 𝑁_𝜑 of radial basis functions up to 10, no relevant differences appear with respect to the case with 5 RBFs. A better approximation of the PDF is, of course, obtained, but the approximation is able to reproduce the simulation only for a limited time range. A good estimation of the system response both in terms of rotation and angular velocity is possible to obtain taking into account 20 RBFs (Figures 21 and 22).

Figure 19 

Upper and lower bounds for pendulum rotation (𝑁𝜑=5; 𝑁𝛼=200).

Figure 20 

Upper and lower bounds for angular velocity (𝑁𝜑=5; 𝑁𝛼=200).

Figure 21 

Upper and lower bounds for pendulum rotation (𝑁𝜑=20; 𝑁𝛼=200).

Figure 22 

Upper and lower bounds for angular velocity (𝑁𝜑=20; 𝑁𝛼=200).

As a second case, the system response with a coefficient of variation of 10% has been analysed (the numbers of 𝑁_𝛼 terms have been considered again 200). Investigations have been made with respect to the number 𝑁_𝜑 of radial basis functions that varies from 5 to 40. Taking into account a small number of RBFs (𝑁_𝜑=5), the results of the proposed approach are not satisfactory as shown in Figures 23 and 24 where upper and lower bounds for pendulum rotation and for angular velocity are reported compared with MCS results. Increasing the number of RBF up to 10, no significant improvement of the system response is obtained (Figures 25 and 26). It is possible to obtain a good estimation of the system response both in terms of rotation and angular velocity if the number of RBF is increased to 40. Figure 27 reports the upper and lower bounds for pendulum rotation, and Figure 28 reports the upper and lower bounds for angular velocity for this case.

Figure 23 

Upper and lower bounds for pendulum rotation (𝑁𝜑=5; 𝑁𝛼=200).

Figure 24 

Upper and lower bounds for angular velocity (𝑁𝜑=5; 𝑁𝛼=200).

Figure 25 

Upper and lower bounds for pendulum rotation (𝑁𝜑=10; 𝑁𝛼=200).

Figure 26 

Upper and lower bounds for angular velocity (𝑁𝜑=10; 𝑁𝛼=200).

Figure 27 

Upper and lower bounds for pendulum rotation (𝑁𝜑=40; 𝑁𝛼=200).

Figure 28 

Upper and lower bounds for angular velocity (𝑁𝜑=40; 𝑁𝛼=200).

As expected, as the coefficient of variation of uncertain parameters increases, an increased number of RBFs need to be taken into account in the response approximation to match the actual system behaviour, hence increasing the overall computational cost of the proposed procedure. Nevertheless, the hybrid approach seems to be able to reproduce with sufficient care main statistics of the system response. Further developments of the current version of the proposed hybrid approach will be aimed at a proper definition of an error function to address the number of terms to be considered in the RBF approximation. It seems, in fact, that the approach is able to match the response process of arbitrarily complex nonlinear systems provided that a sufficient number of terms is included in the approximation. This is, in fact, a general property of neural networks.

4. Concluding Remarks

As randomness of the response process of a mechanical system can be due either to randomness in structural parameters and/or in the forcing process (f.i. wind or earthquake), the paper presents a combined hybrid Galerkin/perturbation approach based on radial basis functions for the dynamic analysis of nonlinear mechanical systems affected by randomness both in system parameters and loads. The paper suggests that the two aspects could be separated, treating the effects of randomness of mechanical parameters with a Galerkin-radial basis functions approach while the effects of the load randomness are investigated by means of a perturbation technique. This is due to the characteristics of the two methods; on the one hand, the Galerkin—RBF approach enables to evaluate with very good approximation the dependence of the structural response on the random parameters, but it increases dramatically the degrees of freedom of the examined problem. On the other hand, the perturbation method enables to evaluate a second-order expansion of the dependence of the structural response on the random parameters, thus leading to big errors in the case of strong deviations from the parabolic law. The resulting approach is quite general, and the proposed expansion for the solution contains all the necessary probabilistic information that allows to characterize the response process. As a test, the proposed procedure has been applied to a strongly nonlinear single degree of freedom problem, and the results have been compared with Monte Carlo simulation. The comparisons demonstrate that the proposed method, although limited for the time being, is a promising candidate to approach uncertain dynamic problems under stochastic loading even if it deserves further investigation to assess its efficiency with respect to other applications.