2nd WSEAS Int. Conf. On NON-LINEAR ANALYSIS,
NON-LINEAR SYSTEMS and CHAOS (WSEAS NOLASC 2003)
Vouliagmeni, Athens, Greece, December 29-31, 2003
Hybrid Finite Element in Non-Linear Structural Dynamics of Anisotropic Tubes Conveying Axial Flow
M.H. Toorani1 and A.A. Lakis[2]
1. Nuclear Engineering Department, Babcock & Wilcox Canada
P.O. Box 310, Cambridge, Ontario, Canada N1R 5V3
2. Mechanical Engineering Department, Ecole Polytechnique of Montreal
C.P. 6079, Station Centre-ville, Montreal (QC), Canada H3C 3A7
Tel.:(514) 340-4711, x4906, Fax: (514) 340-4176
Abstract: A semi-analytical approach has been developed in the present theory to determine the geometrical non-linearity effects on the natural frequencies of anisotropic cylindrical tubes conveying axial flow. Particular important in this study is to obtain the natural frequencies of the coupled system of the fluid-structure, taking the geometrical non-linearity of the structure into consideration and estimating the critical flow velocity at which the structure loses its stability. The displacement functions, mass and stiffness matrices, linear and non-linear ones, of the structure are obtained by exact analytical integration over a hybrid element developed in this work. Linear potential flow theory is applied to describe the fluid effect that leads to the inertial, centrifugal and Coriolis forces. Stability of tubes subjected to the flowing fluid is also discussed. Numerical results are given and compared with those of experiment and other theories to demonstrate the practical application of the present method.
Key-Words: Non-Linear, Dynamics, Flowing Fluid, Anisotropic, Cylindrical Tubes
1
1Introduction
Component failures due to excessive flow-induced vibrations continue to affect the performance and reliability of nuclear components, piping system and tube heat exchangers. Fluid-elastic vibrations have been recognized as a major cause of failure in shell-and-tube-type heat exchangers. Fluid-elastic vibrations result from coupling between fluid-induced dynamic forces and the motion of structure. Depending on the boundary conditions, static (buckling) and dynamic (flutter) instabilities are possible in the structures at sufficiently high flow velocities. The nature of fluid-elastic instability can be illustrated as a feedback mechanism between structural motion and the resulting fluid forces. A small structural displacement due to fluid forces or whatever alerts the flow pattern, inducing a change in the fluid forces; this in turn leads to further displacement, and so on. When the flow velocity becomes higher, the vibration amplitude becomes larger and the impact phenomenon occurs that can lead to unacceptable tube damage due to fatigue and /or fretting wear in critical process equipment. Therefore, the evaluation of complex vibrational behavior of these structures is highly desirable in nuclear industry to avoid such problems.
While several studies have been conducted on the dynamic stability of circular cylindrical shells conveying fluid, but in contrast, non-linear studies of shells subjected to either internal or external axial flow are few. Particularly interesting, in case of internal axial flow, are the studies of Païdoussis &Dennis [1], Lakis & Païdoussis [2, 3], Weaver & Unny [4], Pettigrew [5], Selmane &Lakis [6] and Amabili et al. [7]. These works are performed based on the classical shell theory by neglecting the shear deformation effect while this later plays a very important role in reducing the effective flexural stiffness of composite shells and also the moderately thick structures. The present work addresses the question of stability of anisotropic cylindrical shells, based on a shearable shell theory, subjected to internal and external axial flow. The non-linearities due to large amplitude shell motion are taken into account, by using the modal coefficient approach, while the amplitude of shell displacement remains within the linear range from the fluid point of view.
2Mathematical Model
The analytical solution involves the following steps:
a) The strain-displacement relations expressed in an arbitrary orthogonal curvilinear coordinate are inserted into the equations of motion, obtained based on shearable shell theory, of anisotropic cylindrical shells. The mass and linear stiffness matrices are determined for an empty finite element, Fig. 1, and assemble the matrices for the complete shell.
b) The coefficients of the modal equations are derived using the non-linear part of the kinematics relations.
c) A finite fluid element bounded by two nodal lines, Fig. 1-B, is considered to account the effect of the fluid on the structure.
d) The linear and non-linear natural vibration frequencies are then obtained and compared with the available results.
2.1Structure Model
The equations of motion of anisotropic cylindrical shells in terms of U, V and W (axial, tangential and radial displacements of the mean surface of the shell), and (the rotations of the normal about the coordinates of the reference surface), see Fig. 1, and in terms of Pij’s elements are written as follows:
(1)
where anisotropic elasticity matrix (Pij’s elements), and five linear differential operators (Lm), are fully given in [8].
The finite element developed is shown in Fig. 1. It is a cylindrical panel segment defined by two nodal lines i and j. Each node has five degrees of freedom, three displacements and two rotations.
Figure 1: (A) Finite element discretization.
(B) Nodal displacements at node i.
The displacement functions associated with the axial wave number are assumed to be:
(2)
where m is the axial mode, and is a complex number. A system of five homogeneous linear functions is obtained by substituting (2) into equations of motion (1). For the solution to be non-trivial, the determinant of this system must be equal to zero. This brings us to the following characteristic equation (see References [8, 9] for more detail):
(3)
Each roots of this equation yields a solution to the linear equations of motion (1). The complete solution is obtained by adding the ten solutions independently. After carrying out the some intermediate manipulations, that are not displayed here (see Reference [9]), the following equations are obtained:
(4)
where is the displacement vector at the boundaries and [N] represents the displacement function matrix. The constitutive relation between the stress and deformation vector of cylindrical shells is given as [9]:
The matrices [P], as a function of geometrical and mechanical parameters of anisotropic cylindrical shells, and [B] are given in [9].
2.1.1 Mass and Linear Stiffness Matrices
Using the procedure of the classical finite element, the mass and stiffness matrices are then calculated. For one finite element, they may be written as follows:
(6)
where is the density of the shell, h its thickness, dA a surface element, [P] the elasticity matrix and the [N] and [B] are derived from equations (4, 5). The matrices [ms] and [ks(L)] are obtained analytically by carrying out the necessary matrix operations over in equation (6). The global matrices [Ms] and [Ks(L)] may be obtained, respectively, by superimposing the mass and stiffness matrices for each individual panel finite element [9].
2.1.2Non-Linear Stiffness Matrices of Structure
Theexact Green strain-displacement relations are used in order to describe the non-linear behavior, including large displacements and large rotations, of anisotropic cylindrical shells. In common with linear theory, it is based on refined shell theory in which the shear deformations and rotary inertia effects are taken into account. The approach developed by Radwan and Genin [10] is used with particular attention to geometric non-linearities.
The coefficients of the modal equations are obtained through the Lagrange method. Thus, the non-linear stiffness matrices of second and third order are then calculated by precise analytical integration and superimposed on the linear part of equations to establish the non-linear modal equations. The main steps of this method are as follow:
2.1.2.a Shell displacements are expressed as generalized product of coordinate sums and spatial functions:
(7)
where the qi (t)’s functions are the generalized coordinates and the spatial functions U, V, W,x and are given by equation (2).
2.1.2.bThe deformation vector is written as a function of the generalised coordinates by separating the linear part from non-linear one:
(8)
This vector is given in [11]. The subscripts “L” and “NL” mean “linear” and “non-linear”, respectively. In general, these terms can be expressed in the following form:
(9-a)
(9-b)
Note: AAij=AAji, BBij=BBji and etc.
2.1.2.cLagrange’s equations of motion in the generalized coordinates qi (t) is defined as:
(10)
Where T is the total kinetic energy, V the total elastic strain energy of deformation and the Qi’s are the generalized forces. Assuming , the strain energy V can be defined as follow:
(11)
Where:
a =1if i = j or k = l(i, j =1,2,...,10),(k,l=3,6) and i , j 3,6
a =2if ij or kl(i ,j=1,2…,10)(k,l=3,6)
2.1.2.dAfter developing the total kinetic energy and strain energy, using definitions (9), and then substituting into the Lagrange equation (10) and carrying out a large number of the intermediate mathematical operations, while are not given here due to the complexity of the manipulations, the following non-linear modal equations are obtained. These non-linear modal equations are used to study the dynamic behavior of an empty anisotropic cylindrical shell.
(12)
Where mij, kij(L) are the terms of mass and linear stiffness matrices given by equation (6). The terms of kijk(NL2) and kijks(NL3) , which represent the second and third-order non-linear stiffness matrices, are given by the following integrals in the case of anisotropic laminated cylindrical shell based on the refined shell theory in which the shear deformation and rotary inertia effects are considered:
(13-a)
and
(13-b)
Where dA=R dx d and:
m =1,2,...,9I=1,3,5,...,55andm, n3,6
n=m+1 to 10J=I+1
The Pij’s are the terms of the elasticity matrix [P] and the terms AAijk, BBijk, ..., AUXijk58 and AAijks, BBijks,..., AUXijks58 represent the coefficients of the modal equations in step (2.1.2.d). Follow are the expressions for the coefficients ai, AAij, AAijk and AAjkrs, the others coefficients are obtained in the same way, details are given in [11].
(14)
where U,V and W are spatial functions determined by equations (2). In equation (14), the subscript “i,j”, “i,j,k” and “i,j,k,s” represent the coupling between two; three and four mode, respectively. Substituting equation (2) into equations (14), we obtain
(15)
where ( i =1,...,10) are the roots of characteristic equation (3) and m is the axial mode number. The same definitions, as relation (15), are obtained for other parameters and given in [11]. The constants Ci(i=1,...10) can be determined using ten boundary conditions for each element. The axial, tangential and radial displacements as well as the rotations have to be specified for each node.
Substituting these definitions into equation (13) and then integrating over , the two expressions for the second- and third-order non-linear matrices are obtained, as given in equation (12).
2.2Fluid Model
Linear potential flow theory is applied to describe the fluid effects that lead to the inertial, centrifugal and Coriolis forces. The mathematical model is based on the following hypothesis: i) the fluid flow is potential; ii) the fluid is irrotational, incompressible and non-viscous.
2.2.1Dynamic Pressure
Based on the previous hypothesis, the potential function must satisfy the Laplace equation. This relation is expressed in the cylindrical coordinate system by:
(16)
is the potential function that represents the velocity potential. The components of the flow velocity are given by:
(17)
where are respectively the axial, tangential and radial components of the fluid velocity; Uxu is the velocity of the liquid through the shell section. The Bernoulli equation is given by:
(18)
A full definition of the flow requires that a condition be applied to the structure-fluid interface. The impermeability condition ensures contact between the shell and the fluid. This should be:
(19)
From the theory of shells, we have:
(20)
Assuming then,
(21)
The function is explicitly determined by applying the impermeability condition (19) and using the radial displacement (20). Substituting the assumed function into equation (16) leads to the following differential Bessel equation:
(22)
where “i” is the complex number, i2=-1 and is the complex solution of the characteristic equation for the empty shell. The general solution of equation (22) is given by:
(23)
where are, respectively, the Bessel functions of the first and second kind of complex order “”. For inside flow, the solution (23) must be finite on the axis of shell (r=0); this means we have to set the constant “B” equal to zero. For outside flow (); this means that the constant “A” is equal to zero. When the shell is simultaneously subjected to internal and external flow, we have to take the complete solution (23). We carry the Bessel equation solution back into (21) to obtain the final expression of velocity potential evaluated at the shell wall:
(24)
where
(25)
and
(26)
where are the roots of the characteristic equation of the empty shell; are, respectively, the Bessel functions of the first and second kind of order “”; “m” is the axial mode number; “R” is the mean radius of the shell; “L” its length; the subscript “u” is equal to “i” for internal flow and is equal to “e” for external flow.
Substituting relation (24) into (18), we obtain the equation for the pressure on the shell wall.
(27)
By introducing the displacement function (20) into pressure expression (27), performing the matrix operation and thereafter integration over the fluid element required by the finite element method, the linear matrices (mass [mf], damping [cf] and stiffness [kf]) of moving fluid are obtained. Finally, the global linear matrices [Mf], [Cf] and [Kf] may be obtained, respectively, by superimposing the different matrices for each individual fluid finite element.
3Non-Linear Differential Relations
The structural and fluid mass and stiffness matrices, either linear or non-linear, as well as the fluid damping matrix, obtained in the previous sections, are only determined for one element. The global mass, stiffness and damping matrices are obtained by assembling the matrices for each element. Assembling is done in such way that all the equations of motion and the continuity of displacements at each node are satisfied. These matrices are designated as [M], [K] and [C], respectively.
(28)
where {} is the displacement vector and [Ms], [KsL], [KsNL2] and [KsNL3] are, respectively the mass, linear and second- and third-order non-linear stiffness matrices of the structure, respectively, and [Mf], [Kf] and [Cf] are the inertial, centrifugal and Coriolis forces, respectively, due to the fluid effect.
Setting:
(29)
where represents the square matrix for eigenvectors of the linear system and is a time related vector. Numerical solution of the coupled system (28) is difficult and costly. Here, we limit ourselves to solving the uncoupled system. In this case, equation (28) is reduced to the following equation:
(30)
where
(31)
where “h” represents the shell thickness. The square root of coefficient represents the ith linear vibration frequency of system. The solution of the non-linear differential equations (30), which satisfies the conditions in (29), is calculated by a fourth order Runge-Kuta numerical method. The linear and non-linear natural frequencies are evaluated by a systematic search for the roots as a function of time. The ratio of linear and non-linear frequency is expressed as a function of non-dimensional ratio where is the vibration amplitude.
4Numerical Results and Discussions
This research work is focused on the shear deformation and geometrically non-linear effects on the dynamic behavior of anisotropic cylindrical shells conveying fluid. Non-linearity effects produce either hardening or softening behavior in circular cylindrical shells. Considering the shear deformation effects leads to reducing the flexural stiffness of the structures. The developed method in this work is a hybrid finite element based on a combination of refined shell theory, modal expansion approach and potential flow theory. This method is capable of obtaining the high as well as low frequencies with high accuracy. The values of shear correction factors used in the calculation have been taken .
Figure 2: Variation of non-dimensional natural frequencies in conjunction with variation of m.
a) Linear Vibrations of empty and liquid filled isotropic and anisotropic cylindrical shells- It should be noted that in the two first examples, the natural frequencies of the structures are also obtained using Sanders’ theory (non-shearable shell theory), by authors.
In the first example, the different longitudinal vibration modes () as a function of the circumferential wave number are drawn in Fig. 2. This figure shows the results for four symmetric layers cross-ply (0o/90o/90o/0o) laminated shell whose mechanical properties are given as:
E1=25E2; G23=0.2E2; G13=G12=0.5E2, ν12=0.25; ρ=1
Figure 3: Variation of non-dimensional natural frequencies in terms of mn.
In the next example (Fig. 3), the effect of axial mode number on the non-dimensional natural frequencies () of an isotropic cylindrical shell is studied and the results are compared with the obtained corresponding values based on the Sanders’ theory.
Fig. 4 shows the natural frequencies computed for closed simply supported, circular cylindrical shell for m=2 and compared with the experimental results, given in [13]. To determine natural frequencies with the developed program, based on the present theory, only 10 elements are required to provide acceptable accuracy. As can be seen, there is good agreement between the present theoretical results and those of experimental. Dimensions and material properties are given as follow: