Flow-Induced Vibration Analysis of Supported Pipes Conveying Pulsating Fluid Using Precise Integration Method

Dynamic analysis of supported pipes conveying pulsating fluid is investigated in Hamiltonian system using precise integration method PIM . First, symplectic canonical equations of supported pipes are deduced with state variable vectors composed of displacement and momentum. Then, PIM with linear interpolation formula is proposed to solve these equations. Finally, this approach’s precision is testified by several numerical examples of pinned-pinned pipes with different fluid velocities and frequencies. The results show that PIM is an efficient and rapid approach for flowinduced dynamic analysis o f supported pipes.


Introduction
As the pipes are widely used in many industrial fields, flow-induced vibration analysis of pipes conveying fluid has been one of the attractive subjects in structural dynamics.It is well known that pipeline systems may undergo divergence and flutter types of instabilities generated by fluid-structure interaction.Over the last sixty years, extensive studies have been carried out on dynamic analysis of pipeline systems subject to different boundary conditions and loadings.Notable contributions in this area include the works of Chen 1 and Paidoussis 2, 3 .At present, most of the research is concentrated on nonlinear dynamic analysis of pipes conveying pulsating fluid.A recent survey on bifurcations for supported pipes can be found in 4 .Folley and Bajaj 5 considered nonlinear spatial dynamic characteristics of cantilever pipes conveying fluid.In most cases, the corresponding ordinary differential motion equations of fluidconveyed pipes are deduced using Galerkin's method in Lagrange system.Then many Mathematical Problems in Engineering numerical methods, such as transfer matrix method, finite element method, perturbation method, Runge-Kutta method, and differential quadrature method, are applied to solve these equations.
For example, Jensen 6 analyzed dynamic behaviors of vibrating pipe containing fluid subject to lateral resonant base excitation using the perturbation method of multiple scales.Yang et al. 7 investigated the effect of fluid viscosity and mass ratio on instability regions of a Kelvin-type viscoelastic pipe conveying harmonically pulsating fluid using multiple scales method.Wang et al. 8 studied the nonlinear dynamics of curved fluid conveying pipe with differential quadrature method.
Jeong et al. 9 proposed a finite element model of pipes conveying periodically pulsating fluid and analyzed the influence of fluid velocities on pipe's stability.Stangl et al. 10 solved the extended version of Lagrange nonlinear equations for cantilevered pipes using implicit Runge-Kutta solver HOTINT.Wang 11 explored numerically the effect of the nonlinear motion constraints on dynamics of simply supported pipes conveying pulsating fluid via the fourth-order Runge-Kutta scheme.
Nikolić and Rajković 12 used Lyapunov-Schmidt reduction and singularity theory to investigate the behaviors of extensive fluid-conveying pipe supported at both ends around the neighborhood of the bifurcation points.Furthermore, Modarres-Sadeghi and Païdoussis 13 studied the possible postdivergence flutter instabilities of this complete nonlinear supported pipe's model with Houbolt's finite difference method 14 and AUTO Software package.Xu et al. 15 proposed the analytical expression of natural frequencies of fluid-conveying pipes with the help of homotopy perturbation method.Those calculated frequencies were in good agreement with experiment results.
Considering the effect of the internal and external fluids, the three-dimensional nonlinear differential equations of a fluid-conveying pipe undergoing overall motions were derived based on Kane's equation and the Ritz method 16 .Moreover, the time histories for the displacements were obtained using the incremental harmonic balance method.Based on Timoshenko beam model, Shen et al. 17 studied the band gap properties of the flexural vibration for periodic pipe system conveying fluid using the transfer matrix method.These methods have proved to be effective in analyzing flow-induced vibration of certain pipes.
It is well known that analysis of pipe dynamics could be conducted based on the energy-based approach according to Hamiltonian principle 18, 19 .However, these approximation methods mentioned above are not ideal for Hamiltonian systems 20 , because they are not structurally stable, which means that the Hamiltonian system will become dissipative.
Recently, many numerical algorithms, which can inherit the symplectic structure of Hamiltonian system, have been studied.Especially, Zhong and Williams 21 have proposed the precise integration method, which can give the highly precise numerical integration result and approach the full computer precision for these homogeneous equations.Moreover, this approach has been applied to solving complicated inhomogeneous problems with nonlinear time-variant item, for example, Floquet transition matrix, control problems, and so on 22-25 .In this paper, a Hamiltonian model of nonlinear flow-induced dynamics of supported pipes is analyzed numerically using precise integration method.Firstly, nonlinear equations of supported pipes conveying harmonically fluctuating fluid are deduced to two-order ordinary differential equations using the Galerkin's method.Then the equations are transformed into symplectic canonical equations composed of displacement and momentum.Moreover, PIM with linear interpolation formula is proposed.Finally, several numerical examples of pinned-pinned pipe conveying pulsating fluid are used to testify the precision of this approach.The results are compared with those using traditional Runge-Kutta method.The influence of different fluid parameters on nonlinear behaviors of supported pipes is also discussed.

Formulation of Problem in Hamiltonian System
In this section, typical governing equations of supported pipe conveying fluid are deduced in Hamiltonian system.

Equation of Motion
We consider a straight supported pipe conveying the harmonically pulsating flow Figure 1 .It is assumed that the motion is planar, and the pipe is nominally horizontal.The crosssectional area of the flow is assumed constant.The effects of gravity and external tension are ignored.Moreover, the pipe behaves like an Euler-Bernoulli beam in transverse vibration and the fluid is assumed to be incompressible.
The transverse motion equation of the pipe is given by Païdoussis and Issid 26 , where x is the longitudinal coordinate, y the transverse deflection, M the moment of flexure of the pipe, L the pipe length, m 1 the mass of the fluid conveyed per unit length, m 2 pipe mass per unit length, u the fluid velocity, p the fluid pressure, and A the cross-sectional area of the flow.
Then the viscoelastic Kelvin-Voigt damping model is introduced, where EI is the flexural stiffness of the pipe material, and η is the coefficient of Kelvin-Voigt viscoelastic damping.Moreover, the velocity u t of pulsating fluid is assumed to be harmonically fluctuating, and has the following form: where u 0 is the mean flow velocity, μ the amplitude of the harmonic fluctuation assumed small , and ω the fluid pulsating frequency.This fluctuating flow velocity appears as parametric excitation term in the equation of motion and may lead to parametric instabilities.

Mathematical Problems in Engineering
Substituting 2.2 and 2.3 into 2.1 yields that

2.4
Incorporate the following dimensionless quantities: .

2.5
Then the equation of motion can be nondimensionalized as

2.6
The motion equation above is inhomogeneous, as the derivative coefficients of W are explicit functions of ξ and τ.
Then we discretize 2.6 using the Galerkin's method.Let where q i τ i 1, 2, . . ., n are generalized coordinates of the discretized pipe and φ i ξ are eigenfunctions of the beam with the same boundary conditions.It has been pointed out that instability boundaries for supported pipes could be determined with adequate precision using the two-mode expansion 2 .So the two-mode expansion of 2.7 is used in the analytical model for simplicity to investigate qualitative behaviors of supported pipes conveying fluid.Substitute 2.7 with n 2 into 2.6 .Then according to the orthogonal property of modal modes, the partial differential equation could be transformed into the second-order ordinary differential equation q HΛ 2β 0.5 vB q Λ v 2 T C q − 2β 0.5 vμ cos ωτB q − v 2 2μ cos ωτ C β 0.5 vμ sin ωτ D − C q, 2.8 q G q Kq f 1 τ q f 2 τ q, 2.9 where

2.10
In 2.9 , G and K denote the structural damping matrix and stiffness matrix, respectively.These two matrices are associated with systematic parameters, such as dimensionless flow velocity v and mass ratio β. λ i i 1, 2 are the ith eigenvalues of the supported pipe and Λ is the diagonal matrix with elements λ c rr s r.

2.12
Mathematical Problems in Engineering y x, t x u t

Figure 1:
The simply supported pipe conveying fluid.

Symplectic Canonical Equation
According to Hamiltonian principle, the nonlinear equation of supported pipe conveying pulsating fluid can be transformed into symplectic canonical equation with state variable vectors composed of displacement Q and corresponding momentum P , where

2.14
So, we can see that H is a 4 × 4 Hamiltonian symplectic matrix and F is a time-variant matrix related to state variable vectors.

Precise Integration Method with Linear Interpolation Approximation
In this section, the principle of precise integration method is briefly introduced.For a more detailed explanation, it is suggested that 21, 22 are consulted.
The precise integration method for homogeneous equations with initial value is fundamental, so it is described in the next subsection firstly.

Integration of Homogeneous Equation
The general solution of homogeneous equation V HV V 0 V 0 can be expressed as Suppose that the time step is τ t k 1 − t k , and then we have the following recursive steps: where T e Hτ .It is seen that how to compute the exponential matrix T is essential for the integration precision.
Then split the time interval τ into a smaller one.Define Δt τ/m and m 2 N .For example, m 1048576 when N 20.As τ is small, Δt is an extremely small time interval.
Assume that T a HΔt I HΔt /2, execute the cycle where I is the identity matrix.So T a is no longer a very small matrix.It could be computed by the following function: T I T a .

3.4
The algorithm given above is called precise integration method.It has no serious numerical round-off error and could approach full computer precision 20 .

Integration of Inhomogeneous Equation
In this subsection, PIM with linear interpolation formula would be proposed to solve inhomogeneous equations.
With the solution of homogeneous equation, 2.13 could been written as V t HV t F V, t .

3.5
Then its solution could be given by the Duhamel's integration as V t e Ht V 0 t 0 e H t−τ F V, τ dτ.

3.6
Similarly, the duration of structural dynamic response is also divided into small time intervals.The response between t k , t k 1 can be written as To solve this inhomogeneous equation, the analytical expression of the time-variant item F V, t is required.But it is not always available.
In this study, the linear interpolation formula is used to approximate this nonlinear item.

Mathematical Problems in Engineering
Assume that where 3.9 Substituting 3.8 into 3.7 gives the linear interpolation expression Thus, we have the numerical expression of symplectic canonical equation using PIM with linear interpolation formula.
In the next section, this method would be used to investigate the motion of supported pipes conveying pulsating fluid under different conditions.

Numerical Examples
In this section, several numerical examples of pinned-pinned pipes are used to testify the effectiveness of precise integration method.

Dynamic Response of Pipes Conveying Stable Fluid
In this subsection, this approach is used to analyze the dynamic response of stable fluidconveying pipes, especially for their computation stability after a long period.In this case, the pipe's dynamic function is a homogeneous equation.The results are compared with those using traditional forth-order Runge-Kutta method.
Consider that the dimensionless mean flow velocity v is 2.0, the mass ratio β is 0.32 and the fluid pressure T 1.The initial conditions are chosen to be q 1 q 2 q1 q2 T −0.1 0.2 0.1 0.4 T .Time increases from 0 to 1000 s and the time step Δt is selected as 0.2 s. Figure 2 illustrates time history of four state variables q 1 q 2 p 1 p 2 of pipe's middle point using the Runge-Kutta method, while Figure 3 shows the results calculated by precise integration method.
It can be found that there are evident differences for four state variables' amplitudes using two methods.The amplitudes in Figure 2 decrease gradually with time.When the simulation time is long enough, state variables may converge to zero.However, those in Figure 3 still keep constant with time, which are almost unaffected by the time step.
So we can conclude that there is the energy dissipation using traditional Runge-Kutta method, which cannot get the accurate numerical results.However, precise integration method is an energy conservative method and could maintain the stability of the numerical simulation in the long period of time.Table 1 lists the computation time needed for two methods, as two different time intervals are selected.It can be noted that precise integration method needs much less computing time than Runge-Kutta method.

Dynamic Response of Pipes Conveying Pulsating Fluid
In this subsection, PIM with linear approximation is used to analyze the dynamic response of supported pipes conveying harmonically pulsating fluid.Similarly, the results are compared with those using forth-order Runge-Kutta method.
Consider the dimensionless mean flow velocity v is 1, the amplitude μ 0.4, the frequency ω 2.5, the mass ratio β 0.32, and the fluid pressure T 1.The initial conditions are chosen to be q 1 q 2 q1 q2 T −0.001 0.001 0 0 T .Time increases from 0 to 100 s, and two time steps involved in this example are selected as Δt 0.01 s and 0.05 s, respectively.By the way, the relative stable dynamic response after 20 s is considered.
Figure 4 illustrates time history of displacement response and phase planes of pipe's middle point using the forth-order Runge-Kutta method, while Figure 5 shows the results calculated by PIM, when the time step is 0.01 s.Figures 6 and 7 show the results with time step 0.05 s using two methods, respectively.
These figures show that phrase planes of two modes would shrink gradually with time and converge to a point if time is long enough.The precise integration method shows nearly the same precision during calculating the dynamic response of supported pipes conveying pulsating fluid.
Furthermore, Table 2 lists the computation time needed for two methods.Similarly, PIM with linear interpolation formula need much less computing time than Runge-Kutta method.This approach is very quick to obtain dynamic response because of running a number of cycles during the computation, which is shown in 3.3 .So it is suitable for longterm dynamic analysis of fluid-conveyed pipes.

Stability Analysis of Pipes under Different Fluid Velocities and Frequencies
In this subsection, the influence of different fluid parameters on nonlinear behaviors of pinned-pinned pipes is discussed using PIM with linear interpolation formula.The dimensionless fluid frequency ω increases from 0 to 70, and three fluid velocities in this example are selected as v 1.0, 1.5, and 2.0.The fluid frequency step and the time step The displacement response q 0 of the pipe's middle point under different fluid parameters is calculated.Figure 8 shows the displacement responses of the middle point versus the pulsating fluid frequency.It can be seen that the pipe keeps stable at most frequencies domain.For example, as ω 5, the pipe behaves stable on the limit loop condition, which is shown as Figure 9 a .However, when the frequencies lie between 16, 18 and 46, 49 , the pipes are unstable.For example, the pipe will be divergent as Figure 9 b shows when ω 17.So, it is very dangerous for pipes' operating safety.
Figures 10 and 11 show the displacement response variations of the middle point versus the pulsating fluid frequency as v 1.5 and 2.0, respectively.It can be shown that the instable zone is changing with fluid frequencies.As v 1.5, the instable zones are at 14, 18 , 44, 51 , and 68, 70 .As v 2.0, the instable zones lie at 6, 7 , 11, 16 , 42, 50 , and 65, 70 .
It can be seen that with the increasing of fluid velocity, the critical fluid frequency gets smaller and the pipe shows complicated nonlinear vibration.

Conclusion
In this study, PIM with linear interpolation formula is presented to analyze nonlinear dynamics of Hamiltonian model of supported straight pipe conveying pulsating fluid.
Several numerical examples are used to testify the effectiveness of this approach.The results show this approach could keep stable even for long period of time, and is much more rapid than traditional Runge-Kutta method.Moreover, the pipe's nonlinear behaviors under the condition of different fluid parameters are discussed.The work presented here provides an alternative approach for investigating the nonlinear dynamic response of the pipes conveying fluid.However, it should be pointed out that linear interpolation formula is a rough approximation method, and more accurate methods should be studied to analyze nonlinear flow-induced dynamics in Hamiltonian system.
4  i .Moreover, B, C, and D are matrices with elements b sr , c sr , and d sr s, r 1, 2 , respectively.They are defined as value should be taken for those three parameters depending on different boundary conditions of the pipe.For the pinned-pinned pipe, we have b sr

Figure 2 :
Figure 2: Time history of four state variables using Runge-Kutta method.

Figure 3 :
Figure 3: Time history of four state variables using precise integration method.

Figure 8 :Figure 9 :
Figure 8: Displacement response of pipe's middle point as the function of the dimensionless pulsating fluid frequency v 1.0 .

Figure 10 :Figure 11 :
Figure 10: Displacement response of pipe's middle point as the function of the dimensionless pulsating fluid frequency v 1.5 .

Table 1 :
Computation time needed using two methods.

Table 2 :
Computation time needed using two methods.