Time response of structure with interval and random parameters using a new hybrid uncertain analysis method

Accepted Manuscript

Time Response of Structure with Interval and Random Parameters Using a New Hybrid Uncertain Analysis Method Xingxing Feng , Jinglai Wu , Yunqing Zhang PII: DOI: Reference:

S0307-904X(18)30359-7 https://doi.org/10.1016/j.apm.2018.07.043 APM 12396

To appear in:

Applied Mathematical Modelling

Received date: Revised date: Accepted date:

18 November 2017 25 May 2018 25 July 2018

Please cite this article as: Xingxing Feng , Jinglai Wu , Yunqing Zhang , Time Response of Structure with Interval and Random Parameters Using a New Hybrid Uncertain Analysis Method, Applied Mathematical Modelling (2018), doi: https://doi.org/10.1016/j.apm.2018.07.043

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

Time Response of Structure with Interval and Random Parameters Using a New Hybrid Uncertain Analysis Method 1

Xingxing Feng, 2Jinglai Wu, 3*Yunqing Zhang

1,2,3

State Key Laboratory of Digital Manufacturing Equipment and Technology, School

of Mechanical Science and Engineering, Huazhong University of Science and

*Corresponding author: [email protected]

CR IP T

Technology, Wuhan, Hubei 430074, China 1. Time responses of engineering structures with interval and/or random parameters are investigated systematically.

2. Polynomial-chaos-Legendre-metamodel method is presented for structure with

AN US

hybrid interval and random parameters.

3. Legendre metamodel method is presented for structure with interval parameters. 4. Polynomial chaos theory is applied for structure with random parameters.

M

Abstract

Practical structures often operate with some degree of uncertainties, and the

ED

uncertainties are often modelled as random parameters or interval parameters. For realistic predictions of the structures behaviour and performance, structure models

PT

should account for these uncertainties. This paper deals with time responses of

CE

engineering structures in the presence of random and/or interval uncertainties. Three uncertain structure models are introduced. The first one is random uncertain structure

AC

model with only random variables. The generalized polynomial chaos (PC) theory is applied to solve the random uncertain structure model. The second one is interval uncertain structure model with only interval variables. The Legendre metamodel (LM) method is presented to solve the interval uncertain structure model. The LM is based on Legendre polynomial expansion. The third one is hybrid uncertain structure model with both random and interval variables. The polynomial-chaos-Legendre-metamodel

1

ACCEPTED MANUSCRIPT (PCLM) method is presented to solve the hybrid uncertain structure model. The PCLM is a combination of PC and LM. Three engineering examples are employed to demonstrate the effectiveness of the proposed methods. The uncertainties resulting from geometrical size, material properties or external loads are studied.

CR IP T

Keywords: structure; time response; uncertainty; polynomial chaos; finite element analysis

1. Introduction

Uncertainty and imprecision in structural parameters and in environmental conditions

AN US

and loads are challenging phenomena in engineering analyses [1]. For instance, variations of inertial, stiffness and damping properties due to uncertainties may significantly affect the dynamics characteristics of structures. Therefore, it is of

M

significant importance to estimate the effects of these uncertainties on structural dynamics. The most common approaches to structures with uncertainties includes

ED

probabilistic approaches and non-probabilistic approaches. The probabilistic approaches demand the probability density functions of the structural parameters. The common

PT

probabilistic approaches include generalized polynomial chaos method [2,3] and Monte

CE

Carlo, etc. Unfortunately, the knowledge of probability density function is not always available. As alternative tools, the non-probabilistic approaches, such as fuzzy finite

AC

element analysis [4,5],

interval model [6] and convex model [7,8], have been

extensively used for handling uncertainties arising in engineering problems. The interval model is widely used when only the bounds of uncertain parameters

are required. Moore [9] and Alefeld [10] done the pioneering work in the field of interval analysis. When the uncertain parameters are modelled as interval parameters, the mass, stiffness and damping matrices of the structure are also interval matrices; the

2

ACCEPTED MANUSCRIPT time responses and eigenvalues (natural frequencies) are solved by so-called interval analysis methods. Both eigenvalue problem and time response problem of structures have attracted much research attention in the last decades. The evaluation of natural frequencies of structure with interval uncertainty is crucial in vibration analysis. In order to evaluate the natural frequencies, Qiu et al. [11] introduced

CR IP T

the eigenvalue inclusion principle; Chen et al. [12] applied matrix perturbation theory; Wang et al. [13] presented a modified interval perturbation finite element method ; Gao [14] presented the interval factor method. In addition, Li et al. [15] divided the structural eigenvalue problem with interval parameters into a series of QB (quadratic

AN US

programming with box constraints) problems by the use of information of the secondorder partial derivatives of eigenvalues. Sofi et al. [16] presented an efficient procedure to seek the bounds of the eigenvalues, where interval uncertainties were handled

M

following the improved interval analysis via extra unitary interval [17,18]. The derivatives of mass and stiffness matrix with respect to interval variables are

ED

required in the aforementioned perturbation based methods. Although they may be unavailable analytically for some complicated problems, they can be solved by using

PT

centre difference method [19]. The drawback of perturbation based methods is that they

CE

are limited to solve the small uncertainty extent problems. Otherwise they may produce a large error when the uncertain range of interval variables is large.

AC

Time response of structure may seriously suffers from uncertainties associated

with material properties, geometrical properties and external loads, etc. Time response of structure with interval uncertainty is an important problem and has been widely investigated. Qiu et al. [20,21] presented first-order perturbation methods to estimate the range of the dynamic response of structures. Xia et al. [22] presented a perturbation method based on the vertex solution theorem for the first-order derivation of the

3

ACCEPTED MANUSCRIPT dynamic response from its central value and avoided interval extension/overestimation problems. Muscolino et al. [23] presented a procedure for estimating the lower and upper bounds, whose key idea is to adopt a first-order approximation of the random response derived by properly improving the ordinary interval analysis, based on the philosophy of affine arithmetic. The perturbation methods require derivatives of the

CR IP T

system responses w.r.t. uncertain variables. Wu et al. [6] presented the Chebyshev interval method for multibody dynamics to achieve shaper and tighter bounds for meaningful solutions. In comparison with perturbation methods, Chebyshev interval method can reduce the interval overestimation and it doesn’t require derivatives of

AN US

system responses w.r.t. uncertain variables. Xia et al. [24] presented a time-variant interval process model for time-variant uncertainties, in which the Chebyshev polynomial expansion is used. Muhanna et al. [25] presented an interval-based finite

M

element formulation for uncertainty in solid and structural mechanics, which avoids most sources of overestimation and computes a very sharp hull for the solution set of

ED

interval linear equations with wide interval quantities. Practical structures may contains random uncertainty and interval uncertainty

PT

simultaneously, so hybrid uncertain analysis methods are developed for hybrid

CE

uncertain structures. Xia et al. [26,27] presented several methods such as RM-CVISPM, HPMCM and HPVM to obtain the sound pressure of structural-acoustic system with

AC

interval random variables. Chen et al. [28] presented a hybrid perturbation method for the prediction of exterior acoustic field with interval and random variables. Feng et al. [29] presented the unified perturbation mathematical programming (UPMP) approach for static analysis of hybrid uncertain engineering structure. Gao et al. [30] presented a mix perturbation Monte-Carlo method for static analysis of bar structures with a mixture of random and interval parameters. These methods [24-28] devoted to interval

4

ACCEPTED MANUSCRIPT and random analysis of uncertain structure are based on perturbation theory. The perturbation methods require derivatives information of system and are limited to problems with small uncertainty level. Yin et al. [31] presented the Gegenbauer series expansion method (GSEM) for response prediction of the structure-acoustic system with interval and/or random variables. GSEM is a kind of orthogonal polynomial

CR IP T

approximation method for hybrid uncertain problems. The PCCI method may be the first orthogonal polynomial approximation method for hybrid uncertain problems, and it is originally applied in vehicle dynamics [32].

Although many methods are presented for engineering structures with hybrid

AN US

uncertainties, most of them focus on static analysis [29][30][33], natural frequency analysis [34] or acoustic field analysis of structural-acoustic system [24-26,29]. The orthogonal polynomial approximation methods have shown good accuracy and

M

efficiency in uncertainty analysis, and the application of orthogonal polynomial approximation method for time response analysis of hybrid uncertain structures is

ED

gaining popularity. Therefore, the major goal of this study is to develop a methodology based on orthogonal polynomial approximation for engineering structures with hybrid

PT

interval and/or random variables. The dynamic behaviour and performance of uncertain

CE

structures are investigated systematically. The uncertainties resulting from geometrical size, material properties or external loads are studied systematically.

AC

The paper is organized as follows. Firstly, three types of uncertain structure models

are introduced. Secondly, the methods to solve the uncertain models are introduced and some new methods are presented. Thirdly, the methods are applied to predict the dynamic responses of engineering structure. Finally, main conclusions are drawn.

2. Problem formulation It is known the equation of motion for a deterministic structure system is formulated by. 5

ACCEPTED MANUSCRIPT Mu(t )  Cu(t )  Ku(t )  P(t )  u(0)  u0 , u(0)  u0

(1)

where M, C, K are the constant mass, damping and stiffness matrices, respectively; P  t  is the external load vector; u 0 and u 0 are the initial conditions.

Practical structures often operate with some degree of uncertainty, which can result

CR IP T

from poorly known or variable parameters (e.g., variation in Young’s modulus and density of material) or from uncertain inputs (e.g. variation in external force). The uncertainties in the practical structures are often modelled as random parameters or

AN US

interval parameters. For realistic predictions of the structures behaviour and performance, structure models should account for these random and interval uncertainties. For this reason, a hybrid uncertain structure model is introduced:

,n ]T is the n-dimensional random variable and   =[-1,1]m

ED

In Eq. (2), =[1 ,2 ,

(2)

M

 M  ,   u(,   , t )  C  ,   u(,   , t )  K  ,   u(,   , t )  P(,   , t )   u(,   , 0)  u0  ,   , u(,   , 0)  u 0  ,  

is the m-dimensional interval variable (interval variables will be introduced in Section

PT

3.2). These uncertainties may result from inertial, stiffness and damping properties of

CE

the structural system or from the external force, so the mass matrix M(,  ) , the damping matrix C(,  ) , the stiffness matrix K(,  ) and the external force

AC

P  ,   ,t  are functions of random variable =[1 ,2 ,

  =[-1,1]m .

,n ]T and interval variable

In addition, u0  ,   and u0  ,   can also be considered as the

uncertain initial conditions, but in this study the initial conditions are deterministic. Equation (2) is a hybrid uncertain structure model in which the random variables and interval variables exist simultaneously. In some cases only interval uncertainties exist in the structure, so the interval uncertain structure model is introduced: 6

ACCEPTED MANUSCRIPT  M   u(  , t )  C   u(  , t )  K   u(  , t )  P(  , t )   u(  , 0)  u0   , u(  , 0)  u 0  

(3)

In other cases only random variables exist in the structure, so the random uncertain structure model is introduced:

M    u(, t )  C    u(, t )  K    u(, t )  P(, t )  u(, 0)  u0    , u(, 0)  u0   

CR IP T

(4)

It is noted that Eq. (3) and (4) are special cases of Eq. (2). The major goal of this study is to establish a framework of analysis of uncertain structures with random and/or

AN US

interval uncertainties. In the next section, uncertain analysis methods to solve Eq. (2-4) will be introduced and some new methods are presented.

3. Analysis methods for uncertainty

M

3.1. Polynomial chaos method for random uncertainty

The polynomial chaos (PC) method [35][36] is applied to solve Eq. (4) with random

ED

parameters. The basic idea of polynomial chaos approach is that random process of interest can be approximated by sums of orthogonal polynomials chaos of independent

PT

variables. A random process Y   , viewed as a function of the random event  , can be

AC

CE

expanded in terms of the orthogonal polynomial chaos as 

Y     yii     

(5)

i 0

where y j represents the deterministic coefficients to be estimated;  j    is the generalized Askey-Wiener polynomial chaos of degree j ;   1 , 2 ,

, n  is the T

multi-dimensional random variable. For different types of random variables, the corresponding basis functions are different. Polynomial basis and corresponding random variables are summarized in Table 1 [35][36]. 7

ACCEPTED MANUSCRIPT Table 1 Random variables and corresponding polynomial chaos Polynomial chaos

Uniform

Legendre

Beta

Jacobi

Gaussian

Hermite

Gamma

Laguerre

CR IP T

Random variable

In practice, a truncated expansion is used with a finite number of terms s 1

Y     yii      i 0

(6)

AN US

This means a finite number of s   n  p !/  n ! p ! terms are considered; n is the number of the random variable and p is the maximum order of the polynomial basis. Several methods can be used to estimate the unknown coefficients in the

M

expansion, such as the collocation method [35]. In the collocation method, the model outputs at some selected collocation points are used to congress the coefficients. The

ED

collocation points are selected from the roots of the polynomial, which is one order

PT

higher than the polynomial chaos expansion. For one-dimensional problem, the available interpolation points are the roots of (p+1)th order polynomial. For n-

CE

dimensional problems, the interpolation points are the tensor product of each dimensional interpolation points, so the total number of available interpolation points is n

AC

 p  1

. It is not necessary to evaluate all the interpolation points, meanwhile the

number of sampling points must be higher than the number of coefficients to be estimated; selecting the number of points equalling twice the number of coefficients is recommended for obtaining robust estimates [37]. A basic rule for selecting interpolation points is summarized:

8

ACCEPTED MANUSCRIPT If 2  n  p !/ n! p !   p 1 , then 2  n  p !/  n! p ! interpolation points are n

recommended for stable calculation [37], and if 2  n  p !/ n ! p !   p  1 , then all the n

 p  1

n

interpolation points should be selected.

Take the problem with n  1 uniformly distributed random variable and p  4 for example, the available  p  1  5 interpolation points -0.9061,-0.5384, 0, 0.5384 and

CR IP T

n

0.5384 are the roots of the (p+1)th order Legendre polynomial, while the number of coefficients is s   n  p !/ n! p !  5 . Since 2  n  p !/ n ! p !   p  1 , N=5 interpolation n

AN US

points are selected.

Take the problem with n  2 uniformly distributed random variables and p  4 for example, the available  p  1  25 interpolation points the tensor products of the roots n

M

of one-dimensional polynomials, i.e. (-0.9061, -0.9061), (-0.9061, -0.5384), (-0.9061, 0), (-0.9061, 0.5384), (-0.9061, 0.9061), (-0.5384, -0.9061),…,(0.9061, 0.5384) ,

ED

(0.9061, 0.9061), while the number of coefficients is s   n  p !/ n ! p !  15 . Since 2  n  p !/ n ! p !   p  1 , N=25 interpolation points are selected.

PT

n

CE

Take the problem with n  4 uniformly distributed random variables and p  4 for example, the available  p  1  625 interpolation points the tensor products of the n

AC

roots of the one-dimensional polynomials. The number of coefficients is s   n  p !/ n ! p !  70 . Since 2  n  p !/ n ! p !   p  1 , N  2  n  p !/ n! p !  140 n

interpolation points are selected. Once the collection points are selected, the least square method can be used to produce the unknown coefficients in Eq. (6) [32] :

9

ACCEPTED MANUSCRIPT



y  T  T 

where y   y0 , y1 ,

T



1

T 

T

 0  1   Y, T      0   N 

s 1  1  

  s 1   N  

(7)

, ys 1  is the vector of coefficients to be estimated; s is the number T

of coefficients; N denotes the number of interpolation points; 1 interpolation points; Y  Y  1  ,

 N denote the

, Y   N  denotes the vector of model output at the

CR IP T

T

interpolation points; T    is the transform matrix. After deriving the coefficients, the mean and variance of Y can be easily obtained [2,32,36]: s 1

AN US

  y0 ,  2   yi2 i2

(8)

i 1

3.2. Legendre metamodel method for interval uncertainty

M

The Legendre metamodel (LM) method is presented to solve Eq. (3) with interval parameters. For each interval parameter, there are two numbers that represent the lower

ED

and upper bounds of the parameter. A real interval  x  in real set R is defined as

x  R

(9)

PT

 x    x, x    x  x  x

where x is the lower bound and x is the upper bound. Any interval  x    a, b can be

AC

CE

transformed to the expression of    -1,1 , e.g.

 x 

ab ba    2 2

(10)

Thus, only the interval    -1,1 is considered. Consider a function f   in which    -1,1 is an interval variable, the basic task of interval analysis is to obtain the bounds of f   . Many methods are presented

10

ACCEPTED MANUSCRIPT to solve this problem. In this paper, the one-dimensional interval function f   is expanded by Legendre polynomials in consideration of interval variable    -1,1 :

 f    

p

f 0   fi Li  

(11)

i 1

In Eq. (11),  f    is an approximation of the original function f   . Later, we

CR IP T

will prove that the approximation is effective and accurate. f i are the constant coefficients of the truncated Legendre series and the Legendre series Li   are defined over    1,1 by the recurrence relation as

AN US

L0    1   L1       i +1 Li 1     2i  1 Li    iLi 1   , i  1, 2, For general multi-dimensional problems, the interval function f

(12)

 η can also be

M

expanded by multi-dimensional Legendre polynomials with interval variables: p

p

p

 f       fi ,i ,

ED

i1  0 i2  0

im  0

1 2

,im

Li1 ,i2 ,

,im

  ,im

  denote the

PT

where   =[-1,1]m denotes the m-dimensional interval variable; Li1 ,i2 ,

(13)

m-dimensional Legendre polynomials, which are defined as the tensor of each

CE

dimensional Legendre series

AC

Li1 ,i2 ,

,im

1 ,2 , ,m   Li 1  Li 2 

The constant coefficients fi1 ,i2 ,

1

,im

2

Lim m 

(14)

in Eq. (13) can be calculated by Gaussian-

Legendre interpolation integral formula. The interpolation points are the tensor product of each dimensional interpolation points, so the total number of interpolation points is

 p  1

m

. It is found that the total number of interpolation points is extremely large for

high dimensional problems, which will be computationally expensive. In order to reduce the computational efforts, Eq. (13) is changed into the form of

11

ACCEPTED MANUSCRIPT

 f    

k 1



0i1 i2   im  p

fi1 ,i2 ,

,im Li1 ,i2 ,

,im

        i 0

i

(15)

i

Eq. (15) is termed as Legendre metamodel, where  i corresponds coefficients fi1 ,i2 , one by one, and  i   corresponds to Li1 ,i2 ,

,im

,im

  one by one, i.e.

 0  f 0,0, ,0 ,  0    L0, ,0     1  f1,0, ,0 ,  1    L1, ,0       k 1    L0, , p    k 1  f 0,0, , p ,

CR IP T

(16)

In Eq. (13), the terms higher than p order are deleted, and the remaining terms equals to k   m  p !/  m! p ! . It is recommend that M  2  m  p !/  m! p !

 =  T T T

1

AN US

interpolation points are chosen to calculate the coefficients [37]:   0  1   T Y, T  0  M  T

 k 1  1  

   k 1  M  

M

where M denotes the number of interpolation points, Y   f  1  the model output vector at the interpolation points 1

(17)

f  M   denotes

M , and  =[ 0 , 1 ,

T

, k 1 ]T

ED

denotes the coefficients of Legendre polynomials. After obtaining the polynomial coefficients, the next step is to calculate the bounds k 1

PT

of Eq. (15) i.e.  f       i i   . Exactly, the key is to calculate the bounds of i 0

 i   with     1,1 . Using interval arithmetic, the bounds of the interval

CE

m

AC

function are calculated: k 1

 k 1



i 0

 i 1



 f       i i     0 +   i  i 

(18)

In Eq. (18),  i denotes the bounds of  i   with     1,1 . Note that  i is an m

interval number, so it can be expressed by  i   i , i  , where  i denotes the lower bound and i the upper bound. So the upper bound of Legendre metamodel is

12

ACCEPTED MANUSCRIPT  k 1  f L U   0 +   i  i   i 1 

(19)

and the lower bound of Legendre metamodel model is

 k 1  f L L   0 +   i  i   i 1 

(20)

In order to explain Eq. (18) clearly, firstly take the one-dimensional problem (m=1) for

CR IP T

example, the first five Legendre polynomials are:

AN US

  0   = L0    1     L     1  1  2  2    L2     3  1 / 2  3  3    L3     5  3  / 2  4 2  4    L4     35  30  3 / 8

As   -1,1 , then L1    -1,1 , L2    -1/ 2,1 , L3    -1,1 , L4     3 / 7,1 . In words

for

one-dimensional

problem  1 = -1,1  2 = -1/ 2,1  3  -1,1

M

other

ED

4 1 3  4   3 / 7,1 . So the upper bound is   i , and lower bound is  0   1   2   3   4 . 2 7 i 0

For the two-dimensional problem (m=2,     1,1 ), if the order of the Legendre 2

PT

polynomial is 4 (p=4), then the two-dimensional Legendre polynomials with

AC

CE

k   m  p !/  m! p !  15 terms are:

13

ACCEPTED MANUSCRIPT

 

AN US

 

CR IP T

 0   = L0 1  L0  2   1   1     L0 1  L1 2    2  2  2    L0 1  L2 2    32  1 / 2  3  3    L0 1  L3  2    5 2  3 2  / 2  4 2  4    L0 1  L4 2    35 2  30 2  3 / 8   5    L1 1  L0 2   1     L   L      1 1 1 2 1 2  6  7     L1 1  L2  2   1  3 2 2  1 / 2       L   L      5 3  3  / 2 1 1 3 2 1 2 2  8 2     L   L     3  1 / 2 2 1 0 2 1  9     L   L    3 2  1 / 2   1  2 2 1 1 2  10  2 2  11    L2 1  L2  2    31  1 / 2  3 2  1 / 2  3  12    L3 1  L0  2    51  31  / 2  3  13    L3 1  L1  2    51  31  / 2  2  4 2  14    L4 1  L0  2    351  301  3 / 8  



M





ED

The bounds of  1   ~  14   are products of corresponding bounds of the onedimensional Legendre polynomials.

 1 = -1,1 ;  2 = -1/ 2,1 ;  3  -1,1 ;  4   3 / 7,1  5  -1,1 1= -1,1 ;

PT

 6  -1,1  -1/ 2,1 = -1,1 ; 7  -1,1  -1,1 = -1,1 ; 8  -1,1  -3/7,1 = -1,1 ;

CE

 9  -1/ 2,1 1= -1/ 2,1 ; 10  -1/ 2,1  -1,1 = -1,1 ;

 11  -1/ 2,1  -1/ 2,1 = -1/ 2,1 ;  12  -1,1 1= -1,1 ;  13  -1,1  -1,1 = -1,1 ;

AC

 14   3 / 7,1 1=  3 / 7,1 8 13 1 3 1 1 3 So the lower bound is  0   1   2   3   4    i   9   10   11    i   14 , 2 7 2 2 7 i 5 i 12

14

and the upper bound is   i . i 0

From these two examples, we found that  i always equals 1, and  i equals -1 in most cases. For higher dimensional problems, the bounds of the multi-dimensional 14

ACCEPTED MANUSCRIPT Legendre polynomials can also be easily obtained in a similar way. It would be possible to evaluate the bounds information once and store it in the computer before interval analysis. During the implementation of Legendre metamodel, the major task is to calculate the coefficient vector  =[ 0 , 1 ,

, k 1 ]T by Eq. (17).

3.3. Polynomial-chaos-Legendre-metamodel method for hybrid uncertainty

CR IP T

The polynomial-chaos-Legendre-metamodel method (PCLM) is presented to solve Eq. (2) with random and interval parameters. Consider a function F  ,   in which =[1 ,2 ,

,n ]T is a random variable and     1,1 is an interval variable, the m

AN US

output of the function should have the characteristics of both random and interval variables. Firstly consider the random variable  only, use Eq. (6) to expand the function F  ,   .

s 1

F  ,      j j   

(21)

M

j 0

ED

where  j denotes the PC coefficients. The coefficients  j should be a function with respect to   , namely  j   . Using Legendre expansion Eq. (15), the coefficients

PT

 j   can be expanded by

k 1

CE

  j       j ,i i  

(22)

i 0

AC

where  j ,i denotes the element in the coefficients matrix  with k rows and s columns. Substituting Eq. (22) into Eq. (8), the mean value and variance can be obtained as follows: k 1

     0      0,i i   i 0



s 1

2

        j 1

2 j

15

2

 k 1       j ,i i     j2 j 1  i  0  s 1

2 j

(23)

ACCEPTED MANUSCRIPT As the expression of the mean and variance contains interval variables, the two statistical are also interval numbers: the interval mean IM    and interval variance IV

 2  . Based on Legendre polynomials, the IM and IV are expressed by  k 1

  k 1   i , i     0,0     0,i   1,1  i 1   i 1  2 2 s 1   s 1        k 1   k 1  2 2          j ,0     j ,i   i , i    j       j ,0     j ,i   1,1   j2   j 1    j 1    i 1   i 1       

CR IP T

      0,0    0,i

(24)

where  i   i , i   -1,1 is used to simplify computation. Since    and  2  are

AN US

functions of interval numbers  -1,1 , the above equations involve the overestimation according to the interval arithmetic. So Eq. (24) is not used in PCLM. Finally, the bounds of IM, IV and ISD   are calculated by 

k 1



i 0

k 1



M

      ,     min  0,i i   , max  0,i i   1 1 1 1 i 0

s 1  k 1       ,     min      j ,i i     j2  1 1 j 1   i 0   2

2

2



2 s 1  k 1     , max      j ,i i     j2  1 1 j 1   i 0   

   (25)  

PT

    ,  

ED

2

In interval analysis, the scanning method or global optimization algorithms are

CE

often applied to find the bounds [32,38]. In this case, the overestimation of the interval

AC

computation can be reasonably controlled and even avoided. If the number of interval variables is less than 3, the scanning method [39] can directly produce accurate bounds; otherwise, global optimization algorithm may be effective. The coefficients matrix  should be calculated before the numerical implementation of Eq. (25). Firstly  j   is considered as an interval function in terms of   . Using Eq. (17), its Legendre coefficients can be calculated as follows.

16

ACCEPTED MANUSCRIPT   j ,0  T      X1   X1     j ,k 1   



  0  1   where X1      0  M 



1

  j  1     X1       j  M     T

(26)

 k 1  1  

   k 1  M  

  0,0     0,k 1

 s 1,0 



T    X1   X1    s 1,k 1 

   s 1  M  

(27)

T



1

X2   

T

 F  1 , i        F   N , i  

(28)

s 1  1  

M

  s 1   N  

ED

 0  1   where X2      0   N 

 s 1  1  

,  s 1  i  can be calculated as

  0  i   T      X2    X2      s 1  i  



  0  1  T  X1     0  M 

AN US

Using Eq. (7), the vector  0  i  ,



1

CR IP T

Repeating above Eq. (26) from j  0 to j  s  1leads to

AC

CE

calculated by

PT

Applying Eq. (28) to all interpolation points from 1 to M , the matrix    can be

  0  1      0  M 

 s 1  1  



T    X 2   X 2    s 1  M  

where F is the function value matrix at the interpolation points.  F  1 , 1   F  F   N , 1 

F  1 , M     F   N , M  

Substituting equation (29) into Eq. (27) leads to

17



1

X2   F T

(29)

ACCEPTED MANUSCRIPT



  X1   X1   T



1



X1   FT X2    X2    X2    T

T



1

(30)

When the matrix  is obtained, the scanning method or the global optimization algorithm can be used to evaluate the bounds of the IM, IV and ISD.

Validation method

CR IP T

The Monte Carlo method and the scanning method are combined to validate against the hybrid uncertain analysis method [32]. Consider the function F  ,   with   N  0,1 ,

   [1,1]m ,

n

the number of Monte Carlo sampling points is N 0 , choosing the q

AN US

asymmetry scanning points in each dimension of interval variables, then the total number of experiment points is M 0 =N0  q m , the experiment points can be written in the form of

  ,   

  1 , 1   C  ,        N0 , 1

M

1





N0

, M 0



  

(31)

ED



M0

The output of the function at these experiment points can also be written in the

PT

form of a matrix.

CE

 F  1 , 1   F = F  C  ,        F  N0 , 1







(32)



AC





F 1 , M 0     F  N 0 , M 0  

The unbiased estimator of the mean value and variance of F  ,   at the scanning

points i will be

1 i  N0

N0



1 N0 F  j , i ,     F  j , i   i N 0  1 j 1 j 1 2 i



The lower bounds and upper bounds of i and  i2 are calculated by 18

(33)

ACCEPTED MANUSCRIPT i , max i      min i i 

(34)

    min  i2 , max  i2  i  i  2

4. Applications and discussions

CR IP T

4.1. 36-bar truss structure A 36-bar truss structure is shown in Fig.1. The structure is characterized by the following geometrical and mechanical properties: nominal cross-sectional areas of the bars A  0.001m2 ; nominal lengths of the bars are specified in Fig.1 where L  1m ;

AN US

nominal material Young’s modulus E  2.1e11Pa ; nominal material density

 =7830kg / m3 . The truss structure is subjected to a horizontal nominal harmonic sinusoidal external force P  t   a sin(2 ft ) at node 13, where a  20000 N and

M

f  20Hz , with initial conditions u0  x  0 and u0  x  0 . Rayleigh damping is

assumed so that the damping matrix is given by C   M   K , where  =20s1 and

AC

CE

PT

ED

 =0.00001s .

19

ACCEPTED MANUSCRIPT

P  t   a sin(2 ft )

L 1m

L 1m

(29)

(31)

13

14

(33)

15 (35)

(30)

(28) (34)

(36) (20)

(22)

11

(24)

12 (26)

(21)

(19) (25)

(27)

(13)

7

9

8 (17)

AN US

(15) (12)

(16)

(18)

(4)

5

6

(8)

M

(6)

L 1m

(14)

(2)

4

L 1m

(23)

(11)

(1)

CR IP T

10

(10)

L 1m

(32)

(3)

ED

(7)

1

(5)

L 1m

(9)

3

PT

2

CE

Figure 1. A schematic of truss structure

AC

4.1.1. Analysis of material property uncertainties modelled as random variables The Young’s modulus and density of the material are modelled as Gaussian random variables: E  E 1  1 N1  ,    1   2 N 2  , where E ,  are the Gaussian variables, N1 ~ N  0,1 and N 2 ~ N  0,1 are both standard Gaussian variables; 1 ,  2 are the

dimensionless constants and here 1   2  0.02 .

20

ACCEPTED MANUSCRIPT The 4th, 6th, 8th polynomial chaos are used in PC method, respectively. PC method requires 3.5s, 9.3s and 19.4s for 4th, 6th and 8th order polynomial, respectively. To validate the PC approach, a 1000-run Monte Carlo simulation is carried out. Monte Carlo requires 64s. The simulation results are presented in Fig.2-4. Figs.2 and 3 present the mean values of the horizontal displacement of node 13 and

CR IP T

of stress of bar 1, respectively. The results are indistinguishable. Figs.4 and 5 present the standard deviations of the horizontal displacement of node 13 and of stress of bar 1, respectively. The results are very close for standard deviations. Both Monte Carlo and PC simulation can produce the probability density function

AN US

(PDF) of the results. The PDFs of the horizontal displacement of node 13 and stress of bar 1 at time = 0.05s are presented in Figs.6 and 7. It is observed from Figs.4-7 that

AC

CE

PT

ED

M

increasing the order of polynomial chaos leads to more accurate results.

Figure 2. Mean values of horizontal displacement of node 13

21

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 3. Mean values of stress of bar 1

AC

CE

PT

Figure 4. Standard deviations of horizontal displacement of node 13

Figure 5. Standard deviations of stress of bar 1 22

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 6. Probability density of horizontal displacement of node 13

PT

Figure 7. Probability density of stress of bar 1

CE

4.1.2. Analysis of geometrical uncertainties modelled as random variables

AC

The bars of the structure are divided into two types: diagonal and non-diagonal bars. The nominal lengths of the non-diagonal bars equal L and the nominal lengths of the diagonal bars equal 2L . Considering the manufacturing errors, the lengths of bars are modelled as Gaussian random variables: Li  L 1+3 N3  for i  1-5,10-14,19-23,28-32;

Li  2 L 1+3 N3  for i  6-9,15-18,24-27,33-36, where i denotes the number of the bars; N3 ~ N  0,1 is the standard Gaussian variable;  3 is a dimensionless constant and 23

ACCEPTED MANUSCRIPT here 3  0.01 . The cross-sectional areas of the bars are also modelled as Gaussian variables: A  A 1+ 4 N 4  , where N 4 ~ N  0,1 is the standard Gaussian variable;  4 is a dimensionless constant and here 3  0.01 . The 4th, 6th, 8th polynomial chaos are used, respectively. PC method requires 3.5s,

CR IP T

8.4s and 19.5s for 4th, 6th and 8th order polynomial, respectively. A 1000-run Monte Carlo simulation requires 65s.

Figs.8 and 9 present the mean values of the horizontal displacement of node 13 and of stress of bar 1, respectively. The results are indistinguishable for the mean values.

AN US

Figs.10 and 11 present the standard deviations of the horizontal displacement of node 13 and of stress of bar 1, respectively. The results are very close for standard deviations.

Fig.12 and 13 presents the PDFs of the horizontal displacement of node 13 and

AC

CE

PT

ED

can get the more accurate results.

M

stress of bar 1 at time = 0.05s. It is observed from Figs.10-13 that the higher order PC

()

Figure 8. Mean values of horizontal displacement of node 13

24

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 9. Mean values of stress of bar 1

AC

CE

PT

Figure 10. Standard deviations of horizontal displacement of node 13

Figure 11. Standard deviations of stress of bar 1 25

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 12. Probability density of horizontal displacement of node 13

PT

Figure 13. Probability density of stress of bar 1

CE

4.1.3. Analysis of uncertain inputs modelled as interval variables The amplitude and frequency of the sinusoidal force are then modelled as interval

AC

variables:  a   a 1+1  I1  ,  f   f 1+ 2  I 2  , where  a  and  f  are the interval variables;

 I1    1,1

and

 I 2    1,1

are two independent interval variables;

1   2  0.02 . The 4th, 5th and 6th order Legendre polynomials are used to get the bounds of responses. The scanning method with 20 symmetrical points in each dimension of

26

ACCEPTED MANUSCRIPT interval parameters is used to produce the reference accurate bounds. The 4th, 5th and 6th order Legendre metamodel (LM) require 2s, 2.9s, 3.9s, respectively. The scanning method requires 26s. Fig.14 presents the bounds of horizontal displacement of node 13. Fig.15 presents the bounds of stress of bar1. It is observed that increasing the order of LM leads to more

M

AN US

CR IP T

accurate result and longer computation time.

AC

CE

PT

ED

Figure 14. Bounds of horizontal displacement of node 13

Figure 15. Bounds of stress of bar 1

27

ACCEPTED MANUSCRIPT 4.1.4. Analysis of hybrid uncertainties The uncertainties of material properties, geometrical properties and external force are considered simultaneously. Totally, there are four random variables N1 ~ N  0,1 , N 2 ~ N  0,1 , N3 ~ N  0,1 and N 4 ~ N  0,1 , and two interval variables  I1    1,1

CR IP T

and  I 2    1,1 . In addition, 1  2  3  4 =0.02 and 1   2  0.02 . Because 8th order PC and 6th order LM can provide accurate simulation results, in the PCLM, 8th order PC and 6th order LM are combined to compute IM and ISD of the structure. Monte-Carlo-Scanning test is carried out to provide the reference, in which

AN US

the number of Monte Carlo sampling points is 1000, and 20 scanning points for each interval parameters, so the total number is 1000x202=400,000. The CPU times are 2,377s for PCLM and 24,124s for Monte-Carlo-Scanning.

Figs.16 and 17 present the interval mean (IM) of horizontal displacement of node

M

13 and of stress of bar 1, respectively. The results show the IM of PCLM and Monte-

ED

Carlo-Scanning (reference) are very close.

Figs.18 and 19 present the interval standard deviation (ISD) of horizontal

PT

displacement of node 13 and of stress of bar 1, respectively. Small differences are observed between PCLM and Monte-Carlo-Scanning. It means that PCLM can also

AC

CE

provide good estimation for ISD.

28

CR IP T

ACCEPTED MANUSCRIPT

Figure 17. IM of stress of bar 1

AC

CE

PT

ED

M

AN US

Figure 16. IM of horizontal displacement of node 13

Figure 18. ISD of horizontal displacement of node 13 29

CR IP T

ACCEPTED MANUSCRIPT

Figure 19. ISD of stress of bar 1

AN US

4.2. Cantilever beam

The second example concerns a cantilever Euler-Bernoulli beam uniformly discretized into 10 elements, as shown in Fig.20. The beam is characterized by the following geometrical and mechanical nominal properties: total length L  2m ; Young’s modulus

P  t   5000sin(10 t )

ED

material density  =2500kg / m3 .

M

E  3.0e10Pa ; rectangular cross-section area with width b  0.1m and height h  0.15m ;

(1)

PT

(2)

CE

1

2

(9)

9

3

(10)

10

11

Fig.20 A schematic of cantilever beam

AC

The beam is subjected to harmonic sinusoidal excitations P  t   5000sin(10 t ) with initial conditions u0  x  0 and u0  x  0 . Rayleigh damping is assumed so that the damping matrix is given by C   M   K , where  =20s1 and  =0.00001s .

30

ACCEPTED MANUSCRIPT 4.2.1. Analysis of material property uncertainties modelled as random variables The Young’s modulus and density of the material are modelled as Gaussian random variables: E  E 1  1 N1  ，    1   2 N 2  ,where E ,  are the Gaussian variables, N1 ~ N  0,1 and N 2 ~ N  0,1 are both standard Gaussian variables; 1 ,  2 are the

CR IP T

dimensionless constants and here 1   2  0.01 . The 4th, 6th, 8th polynomial chaos are used in PC method, respectively. PC method requires 1.1s, 1.6s and 3.2s for 4th, 6th and 8th order polynomial, respectively. To validate the PC approach, a 1000-run Monte Carlo simulation is carried out. Monte

AN US

Carlo requires 21.3s. The simulation results are presented in Figs.21-24.

Figs.21 and 22 present the mean values of deflection of node 11 and of maximum normal stress of node 1, respectively. The results are indistinguishable for the means. Figs.23 and 24 present the standard deviations of deflection of node 11 and of

M

maximum normal stress of node 1, respectively. The results are very close for standard

AC

CE

PT

ED

deviations.

Figure 21. Mean values of deflection of node 11

31

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 22. Mean values of maximum normal stress of node 1

AC

CE

PT

Figure 23. Standard deviations of deflection of node 11

Figure 24. Standard deviations of maximum normal stress of node 1 32

ACCEPTED MANUSCRIPT 4.2.2. Analysis of geometrical uncertainties modelled as interval variables The total length, width and height of cross-section are modelled as interval variables:

 L  L 1+1  I1  , b  b 1+ 2  I 2  ,  h  h 1+3  I3  ,where  L  , b and  h are interval variables;  I1    1,1 ,  I 2    1,1 and  I3    1,1 are three independent

CR IP T

interval variables; 1   2  3  0.01. The 4th, 5th and 6th order Legendre polynomials are used to get the bounds of responses. The scanning method with 20 symmetrical points in each dimension of interval parameters is used to produce the reference accurate bounds. The 4th, 5th and 6th

AN US

order Legendre metamodel (LM) require 3.5s, 6.2s, 9.8s, respectively. The scanning method requires 212.2s.

Fig.25 presents the bounds of deflection of node 11. Fig.26 presents the bounds of maximum normal stress of node 1. It is observed that increasing the order of LM leads

AC

CE

PT

ED

M

to more accurate result and longer computation time.

Figure 25. Bounds of deflection of node 11

33

CR IP T

ACCEPTED MANUSCRIPT

4.2.3. Analysis of hybrid uncertainties

AN US

Figure 26. Bounds of maximum normal stress of node 1

The uncertainties of material and geometrical properties are considered simultaneously. N1 ~ N  0,1 and N 2 ~ N  0,1 are standard Gaussian variables.

,

are three independent interval variables. In addition,

M

 I 2    1,1 and  I3    1,1

 I1    1,1

ED

1  2  0.01 , 1   2  3  0.01.

The PCLM is used to solve this problem, in which 8th order PC and 6th order

PT

Legendre metamodel are used to compute IM and ISD of the structure. The Monte-

CE

Carlo-Scanning test is carried out to provide the reference, in which the number of Monte Carlo sampling points is 1000, and 20 scanning points for each interval

AC

parameters, so the total number is 1000x203=8,000,000. The CPU times are 21,117s for PCLM and 191,360s for Monte-Carlo-Scanning. Figs.27 and 28 present the interval mean (IM) of deflection of node 11 and of

maximum normal stress of node 1, respectively. The results show the IM of PCLM and Monte-Carlo-Scanning (reference) are very close. Figs.29 and 30 present the interval standard deviation (ISD) of deflection of node 11 and of maximum normal stress of node 1, respectively. Small differences are 34

ACCEPTED MANUSCRIPT observed between PCLM and Monte-Carlo-Scanning. It means that PCLM can also

AN US

CR IP T

provide good estimation for ISD.

CE

PT

ED

M

Figure 27. IM of deflection of node 11

AC

Figure 28. IM of maximum normal stress of node 1

35

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 29. ISD of deflection of node 11

CE

PT

Figure 30. ISD of maximum normal stress of node 1

4.3. Spring-mass-damper system

AC

As shown in Fig.31, a vehicle vibration model with five degrees of freedom is composed of one sprung mass that joints to three un-sprung masses. A similar model is also used in [40]. Symbols ms , kss , css , mb , Ib , mtf , mtr , l f , lr , r, ktf and ktr denote seat mass, seat stiffness coefficient, seat damping coefficient, sprung mass, pitch mass moment of inertia of the sprung mass, front tire mass, rear tire mass, front suspension position in relation to the centre of sprung mass, rear suspension position in relation to the centre of sprung mass, seat position in relation to the centre of sprung mass, front tire stiffness 36

ACCEPTED MANUSCRIPT coefficient and rear tire stiffness coefficient, respectively. Symbols ksf , csf , ksr and csr are defined as front suspension stiffness coefficient, front suspension damping coefficient, rear suspension stiffness coefficient and damping coefficient, respectively.

zs

css

k ss

lr

zbr



k sr

mtf

csr

m tr

ztr

AN US

zrf

zb

csf

k sf

ztf

r mb

lf

zbf

CR IP T

ms

ktf

k tr

zrr

Figure 31. Vehicle model

M

It is assumed that the vectors of generalized and excitation coordinates are given by T

ED

q1  [ zs , zb , , ztf , ztr ]T and q0   zrf , zrr  . The motion in vertical direction for sprung

mass and un-sprung masses and pitching motion for sprung mass, in terms of

PT

acceleration, velocity and displacement, are considered in the formulation of motion

CE

equations Mq1  Cq1  Kq1  K t q0 , where M , C , K are the ( 5  5 ) mass, damping,

AC

stiffness matrices, respectively, and K t is  5  2  forced stiffness matrix.  ms   M   

mb

Ib mtf

      mtr 

 css  c  ss C   rcss   0  0

css css  csf  csr lr csr  rcss  l f csf csf csr

37

css r css r  csf l f  csr lr l 2f csf  lr2csr  r 2css csf l f csr lr

0 csf l f csf csf 0

0  csr  lr csr   0  csr 

ACCEPTED MANUSCRIPT

 kss  k  ss K   rkss   0  0

kss kss  ksf  ksr l f k sf  lr k sr  rkss ksf ksr

kss r k ss r  k sf l f  k sr lr l 2f k sf  lr2 k sr  r 2 k ss l f ksf lr ksr

0 k sf l f k sf ksf  ktf 0

0 0  Kt   0   ktf  0

0  k sr  lr k sr   0  ksr  ktr 

0 0  0  0 ktr 

The nominal value of certain parameters and also, upper bound and lower bound of

CR IP T

uncertain parameters are given in Table 1. The main parameters are adopted from [40]. Meanwhile, the stiffness and damping coefficients of seat and suspensions are considered as uniformly distributed random variables. The mass of seat, the mass and pitch inertia of vehicle body are considered as interval variables. Note that the mass and

AN US

pitch inertia of vehicle body are not independent. It is assumed that the vehicle body pitch inertia ( ) is proportional to the vehicle body mass (

. Therefore, there are 2

independent interval variables and 6 independent random variables.

parameters name

Nominal value

Upper bound value

Type

_

_

Certain

1.803

_

_

Certain

40

_

_

Certain

and dimension

(m)

CE

(kg)

1.011

PT

(m)

Lower bound value

ED

Parameters

M

Table 1 Nominal value of certain parameters and upper bound and lower bound of uncertain

35.5

_

_

Certain

(N/m)

175,500

_

_

Certain

175,500

_

_

Certain

0.5

_

_

Certain

1,230

1,107

1,353

Interval*

(kg)

730

657

803

Interval*

(kg)

80

72

88

Interval

(N/m)

_

90,000

110,000

Uniform

AC

(kg)

(N/m)

r (m)

(kgm2)

38

ACCEPTED MANUSCRIPT (Ns/m)

_

2,250

2,750

Uniform

(N/m)

_

18,000

22,000

Uniform

(Ns/m)

_

900

1,100

Uniform

(N/m)

_

18,000

22,000

Uniform

(Ns/m)

_

900

1,100

Uniform

and

are dependent. It is assumed that

is proportional to

.

CR IP T

* Note:

This vehicle model is excited by a random terrain excitation derived from a roughness model of the form [41]

Gz       ;   0,   0

(35)

AN US

where Ω is the spatial frequency that is related to forward speed V and the temporal frequency f , such that f  V , and constants α and β are the roughness coefficient and waviness of the terrain, respectively. It is supposed that the vehicle moves at

M

constant speed V=20m/s over the terrain profile shown in Fig. 32, and it is further

AC

CE

PT

t   l f  lr  / V .

ED

assumed that the rear tire follows the same trajectory as the front tire with a delay of

Figure 32 Terrain profile for front tire Vehicle vibration response with hybrid uncertain parameters is analysed under random terrain excitations. The 4th order LM and 4th order PC are used in PCLM. The Monte39

ACCEPTED MANUSCRIPT Carlo-Scanning test is carried out to provide the reference, in which the number of Monte Carlo sampling points is 1000, and 20 scanning points for each interval parameters, so the total sampling number is 1000x202=400,000. The CPU times are 28,110s for PCLM and 110,992s for Monte-Carlo-Scanning. Some simulation results are presented in Figs.33-40. Fig.33 and 34 present IM and ISD

CR IP T

of the seat vertical displacement, respectively. Fig.35 and 36 present IM and ISD of the vehicle body pitch angle, respectively. Fig.37 and 38 present IM and ISD of seat suspension spring force, respectively. Fig.39 and 40 present IM and ISD of seat suspension damper force, respectively.

AN US

The results show that IM of PCLM and Monte-Carlo-Scanning (reference) are very close, and ISD differences between PCLM and Monte-Carlo-Scanning are very small. This example shows that the PCLM is also suitable for hybrid uncertain structures under

AC

CE

PT

ED

M

random excitations.

Figure 33. IM of vertical displacement of seat

40

CR IP T

ACCEPTED MANUSCRIPT

ED

M

AN US

Figure 34. ISD of vertical displacement of seat

AC

CE

PT

Figure 35. IM of pitch angle of vehicle body

Figure 36. ISD of pitch angle of vehicle body

41

ED

M

AN US

(a) Overall view

CR IP T

ACCEPTED MANUSCRIPT

PT

(b) Magnified view

AC

CE

Figure 37. IM of seat spring force

Figure 38. ISD of seat spring force 42

ED

M

AN US

(a) Overall view

CR IP T

ACCEPTED MANUSCRIPT

PT

(b) Magnified view

AC

CE

Figure 39. IM of seat damper force

(a) Overall view 43

(b) Magnified view

CR IP T

ACCEPTED MANUSCRIPT

AN US

Figure 40. ISD of seat damper force

5. Conclusions

The dynamic time response of structure with random or/and interval uncertainties is

M

analysed. Geometrical uncertainties, material properties uncertainties and external loads

ED

uncertainties are studied. Polynomial chaos (PC) theory is applied for structure with random uncertainty. Legendre metamodel (LM) method is presented for structure with

PT

interval uncertainty. PCLM method is presented for structure with hybrid random and interval uncertainty. To validate the PCLM method, the Monte-Carlo-Scanning test

CE

scheme is applied to calculate two types of evaluation indexes: the interval mean (IM)

AC

and interval variance/interval standard deviation (IV/ISD). Three engineering examples demonstrate the efficiency and accuracy of proposed methods. Numerical results show that the PCLM can provide accurate results for both IM and ISD in much shorter simulation time than Monte-Carlo-Scanning test.

44

ACCEPTED MANUSCRIPT Acknowledgements This research was supported by the National Natural Science Foundation of China (11772136, 11472112). This financial support is gratefully acknowledged.

References

AC

CE

PT

ED

M

AN US

CR IP T

[1] M. Beer, S. Ferson, V. Kreinovich, Imprecise probabilities in engineering analyses, Mechanical Systems and Signal Processing. 37 (2013) 4–29. doi:10.1016/j.ymssp.2013.01.024. [2] A. Sandu, C. Sandu, M. Ahmadian, Modeling multibody systems with uncertainties. Part I: Theoretical and computational aspects, Multibody System Dynamics. 15 (2006) 369–391. doi:10.1007/s11044-006-9007-5. [3] C. Sandu, A. Sandu, M. Ahmadian, Modeling multibody systems with uncertainties. Part II: Numerical applications, Multibody System Dynamics. 15 (2006) 241–262. doi:10.1007/s11044-006-9008-4. [4] A.S. Balu, B.N. Rao, High dimensional model representation based formulations for fuzzy finite element analysis of structures, Finite Elements in Analysis and Design. 50 (2012) 217–230. doi:10.1016/j.finel.2011.09.012. [5] D. Behera, S. Chakraverty, Fuzzy finite element analysis of imprecisely defined structures with fuzzy nodal force, Engineering Applications of Artificial Intelligence. 26 (2013) 2458–2466. doi:10.1016/j.engappai.2013.07.021. [6] J. Wu, Z. Luo, Y. Zhang, N. Zhang, L. Chen, Interval uncertain method for multibody mechanical systems using Chebyshev inclusion functions, International Journal for Numerical Methods in Engineering. 95 (2013) 608–630. doi:10.1002/nme.4525. [7] C. Jiang, R.G. Bi, G.Y. Lu, X. Han, Structural reliability analysis using nonprobabilistic convex model, Computer Methods in Applied Mechanics and Engineering. 254 (2013) 83–98. doi:10.1016/j.cma.2012.10.020. [8] C. Jiang, X. Han, G.Y. Lu, J. Liu, Z. Zhang, Y.C. Bai, Correlation analysis of nonprobabilistic convex model and corresponding structural reliability technique, Computer Methods in Applied Mechanics and Engineering. 200 (2011) 2528– 2546. doi:10.1016/j.cma.2011.04.007. [9] R.E. Moore, Methods and applications of interval analysis, Prentice Hall, London, 1979. [10] G. Alefeld, Introductions to interval computations, Academic Press, New York, 1983. [11] Z. Qiu, X. Wang, M.I. Friswell, Eigenvalue bounds of structures with uncertainbut-bounded parameters, Journal of Sound and Vibration. 282 (2005) 297–312. doi:10.1016/j.jsv.2004.02.051. [12] S.H. Chen, H.D. Lian, X.W. Yang, Interval eigenvalue analysis for structures with interval parameters, Finite Elements in Analysis and Design. 39 (2003) 419–431. doi:10.1016/S0168-874X(02)00082-3. [13] C. Wang, Z. Qiu, X. Wang, D. Wu, Interval finite element analysis and reliabilitybased optimization of coupled structural-acoustic system with uncertain parameters, Finite Elements in Analysis and Design. 91 (2014) 108–114. doi:10.1016/j.finel.2014.07.014.

45

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN US

CR IP T

[14] W. Gao, Interval natural frequency and mode shape analysis for truss structures with interval parameters, Finite Elements in Analysis and Design. 42 (2006) 471– 477. doi:10.1016/j.finel.2005.09.003. [15] Q. Li, Z. Qiu, X. Zhang, Eigenvalue analysis of structures with interval parameters using the second-order Taylor series expansion and the DCA for QB, Applied Mathematical Modelling. 49 (2017) 680–690. doi:10.1016/j.apm.2017.02.041. [16] A. Sofi, G. Muscolino, I. Elishakoff, Natural frequencies of structures with interval parameters, Journal of Sound and Vibration. 347 (2015) 79–95. doi:10.1016/j.jsv.2015.02.037. [17] A. Sofi, E. Romeo, A novel Interval Finite Element Method based on the improved interval analysis, Computer Methods in Applied Mechanics and Engineering. 311 (2016) 671–697. doi:10.1016/j.cma.2016.09.009. [18] G. Muscolino, A. Sofi, Bounds for the stationary stochastic response of truss structures with uncertain-but-bounded parameters, Mechanical Systems and Signal Processing. 37 (2013) 163–181. doi:10.1016/j.ymssp.2012.06.016. [19] B. Xia, D. Yu, X. Han, C. Jiang, Unified response probability distribution analysis of two hybrid uncertain acoustic fields, Computer Methods in Applied Mechanics and Engineering. 276 (2014) 20–34. doi:10.1016/j.cma.2014.03.014. [20] Z. Qiu, X. Wang, Comparison of dynamic response of structures with uncertainbut-bounded parameters using non-probabilistic interval analysis method and probabilistic approach, International Journal of Solids and Structures. 40 (2003) 5423–5439. doi:10.1016/S0020-7683(03)00282-8. [21] Z. Qiu, X. Wang, Parameter perturbation method for dynamic responses of structures with uncertain-but-bounded parameters based on interval analysis, International Journal of Solids and Structures. 42 (2005) 4958–4970. doi:10.1016/j.ijsolstr.2005.02.023. [22] Y. Xia, Z. Qiu, M.I. Friswell, The time response of structures with bounded parameters and interval initial conditions, Journal of Sound and Vibration. 329 (2010) 353–365. doi:10.1016/j.jsv.2009.09.019. [23] G. Muscolino, A. Sofi, Stochastic analysis of structures with uncertain-butbounded parameters via improved interval analysis, Probabilistic Engineering Mechanics. 28 (2012) 152–163. doi:10.1016/j.probengmech.2011.08.011. [24] B. Xia, Y. Qin, D. Yu, C. Jiang, Dynamic response analysis of structure under time-variant interval process model, Journal of Sound and Vibration. 381 (2016) 121–138. doi:10.1016/j.jsv.2016.06.030. [25] R.L. Muhanna, R.L. Mullen, Uncertainty in mechanics problems—interval–based approach, Journal of Engineering Mechanics. 127 (2001) 557–566. doi:10.1061/(ASCE)0733-9399(2001)127:6(557). [26] B. Xia, D. Yu, J. Liu, Hybrid uncertain analysis for structural–acoustic problem with random and interval parameters, Journal of Sound and Vibration. 332 (2013) 2701–2720. doi:10.1016/j.jsv.2012.12.028. [27] B. Xia, S. Yin, D. Yu, A new random interval method for response analysis of structural–acoustic system with interval random variables, Applied Acoustics. 99 (2015) 31–42. doi:10.1016/j.apacoust.2015.05.002. [28] N. Chen, D. Yu, B. Xia, Hybrid uncertain analysis for the prediction of exterior acoustic field with interval and random parameters, Computers & Structures. 141 (2014) 9–18. doi:10.1016/j.compstruc.2014.05.004. [29] J. Feng, D. Wu, W. Gao, G. Li, Uncertainty analysis for structures with hybrid random and interval parameters using mathematical programming approach,

46

ACCEPTED MANUSCRIPT

[35]

[36]

[37] [38]

[39]

[40]

CR IP T

AC

[41]

AN US

[34]

M

[33]

ED

[32]

PT

[31]

CE

[30]

Applied Mathematical Modelling. 48 (2017) 208–232. doi:10.1016/j.apm.2017.03.066. W. Gao, D. Wu, C. Song, F. Tin-Loi, X. Li, Hybrid probabilistic interval analysis of bar structures with uncertainty using a mixed perturbation Monte-Carlo method, Finite Elements in Analysis and Design. 47 (2011) 643–652. doi:10.1016/j.finel.2011.01.007. S. Yin, D. Yu, H. Yin, B. Xia, Interval and random analysis for structure–acoustic systems with large uncertain-but-bounded parameters, Computer Methods in Applied Mechanics and Engineering. 305 (2016) 910–935. doi:10.1016/j.cma.2016.03.034. J. Wu, Z. Luo, N. Zhang, Y. Zhang, A new uncertain analysis method and its application in vehicle dynamics, Mechanical Systems and Signal Processing. 50– 51 (2015) 659–675. doi:10.1016/j.ymssp.2014.05.036. D. Wu, W. Gao, Hybrid uncertain static analysis with random and interval fields, Computer Methods in Applied Mechanics and Engineering. 315 (2017) 222–246. doi:10.1016/j.cma.2016.10.047. J. Feng, D. Wu, W. Gao, G. Li, Hybrid uncertain natural frequency analysis for structures with random and interval fields, Computer Methods in Applied Mechanics and Engineering. 328 (2018) 365–389. doi:10.1016/j.cma.2017.09.004. D. Xiu, G.E. Karniadakis, Modeling uncertainty in flow simulations via generalized polynomial chaos, Journal of Computational Physics. 187 (2003) 137– 167. doi:10.1016/S0021-9991(03)00092-5. G. Kewlani, J. Crawford, K. Iagnemma, A polynomial chaos approach to the analysis of vehicle dynamics under uncertainty, Vehicle System Dynamics. 50 (2012) 749–774. doi:10.1080/00423114.2011.639897. S.S. Isukapalli, Uncertainty analysis of transport-transformation models, The State University of New Jersey, 1999. Z. Wang, Q. Tian, H. Hu, Dynamics of spatial rigid–flexible multibody systems with uncertain interval parameters, Nonlinear Dyn. 84 (2016) 527–548. doi:10.1007/s11071-015-2504-4. A.J. Buras, M. Jamin, M.E. Lautenbacher, A 1996 analysis of the CP violating ratio ε′/ε, Physics Letters B. 389 (1996) 749–756. doi:10.1016/S03702693(96)80019-0. A. Jamali, M. Salehpour, N. Nariman-zadeh, Robust Pareto active suspension design for vehicle vibration model with probabilistic uncertain parameters, Multibody Syst Dyn. 30 (2013) 265–285. doi:10.1007/s11044-012-9337-4. A. Pazooki, S. Rakheja, D. Cao, Modeling and validation of off-road vehicle ride dynamics, Mechanical Systems and Signal Processing. 28 (2012) 679–695. doi:10.1016/j.ymssp.2011.11.006.

47