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

Abstract

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 flow-induced dynamic analysis o f supported pipes.

1. 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 fluid-conveyed pipes are deduced using Galerkin’s method in Lagrange system. Then many 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.

2. Formulation of Problem in Hamiltonian System

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

2.1. 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 cross-sectional 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.

Figure 1 

The simply supported pipe conveying fluid.

The transverse motion equation of the pipe is given by Païdoussis and Issid [26], 𝜕²𝑀𝜕𝑥²+𝑚₁(𝐿−𝑥)𝜕𝑢𝜕𝜕𝑡²𝑦𝜕𝑥²+(𝑚₁𝑢²𝜕+𝑝𝐴)²𝑦𝜕𝑥²+2𝑚₁𝑢𝜕²𝑦𝜕𝑥𝜕𝑡+(𝑚₁+𝑚₂)𝜕²𝑦𝜕𝑡²=0, (2.1) where 𝑥 is the longitudinal coordinate, 𝑦 the transverse deflection, 𝑀 the moment of flexure of the pipe, 𝐿 the pipe length, 𝑚₁ the mass of the fluid conveyed per unit length, 𝑚₂ pipe mass per unit length, 𝑢 the fluid velocity, 𝑝 the fluid pressure, and 𝐴 the cross-sectional area of the flow.

Then the viscoelastic Kelvin-Voigt damping model is introduced,𝜕𝑀=(𝐸+𝜂𝜕𝑡)𝐼𝑦^′′, (2.2) where 𝐸𝐼 is the flexural stiffness of the pipe material, and 𝜂 is the coefficient of Kelvin-Voigt viscoelastic damping.

Moreover, the velocity 𝑢(𝑡) of pulsating fluid is assumed to be harmonically fluctuating, and has the following form:𝑢(𝑡)=𝑢₀(1+𝜇cos𝑢𝜔𝑡),₂(𝑡)≈𝑢²₀(1+2𝜇cos𝜔𝑡), (2.3) where 𝑢₀ 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.

Substituting (2.2) and (2.3) into (2.1) yields that𝜕(𝐸+𝜂𝜕𝜕𝑡)𝐼⁴𝑦𝜕𝑥⁴+[−𝑚₁𝑢₀(𝐿−𝑥)𝜇𝜔sin𝜔𝑡+𝑚₁𝑢²₀(1+2𝜇cos𝜕𝜔𝑡)+𝑝𝐴]²𝑦𝜕𝑥²+2𝑚₁𝑢₀(1+𝜇cos𝜕𝜔𝑡)²𝑦𝜕𝑥𝜕𝑡+(𝑚₁+𝑚₂)𝜕²𝑦𝜕𝑡²=0. (2.4)

Incorporate the following dimensionless quantities:𝑥𝜉=𝐿𝑦,𝑊=𝐿𝑡,𝜏=𝐿²(𝐸𝐼𝑚₁+𝑚₂)^0.5,𝑣=𝑢₀𝑚𝐿(¹)𝐸𝐼^0.5𝐿,𝑇=𝑝𝐴²,𝑚𝐸𝐼𝛽=¹𝑚₁+𝑚₂𝜂,𝐻=𝐿²(𝐼𝐸(𝑚₁+𝑚₂))^0.5,𝜔=𝜔𝐿²(𝑚₁+𝑚₂)𝐸𝐼^0.5. (2.5)

Then the equation of motion can be nondimensionalized as𝐻𝜕⁵𝑊𝜕𝜉⁴+𝜕𝜕𝜏⁴𝑊𝜕𝜉⁴+[𝑣²(1+𝜇cos𝜔𝜏)²−𝑣𝛽^0.5𝜕(1−𝜉)𝜇𝜔sin𝜔𝑡+𝑇]²𝑊𝜕𝜉²+2𝑣𝛽^0.5𝜕(1+𝜇cos𝜔𝜏)²𝑊+𝜕𝜕𝜉𝜕𝜏²𝑊𝜕𝜏²=0. (2.6)

The motion equation above is inhomogeneous, as the derivative coefficients of 𝑊 are explicit functions of 𝜉 and 𝜏.

Then we discretize (2.6) using the Galerkin’s method. Let𝑊(𝜉,𝜏)=^𝑛_𝑟=1𝜙_𝑟(𝜉)𝑞_𝑟(𝜏), (2.7) where 𝑞_𝑖(𝜏) (𝑖=1,2,…,𝑛) are generalized coordinates of the discretized pipe and 𝜙_𝑖(𝜉) 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 𝑛=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 𝑞̈+(𝐻Λ+2𝛽^0.5𝑣𝐵)𝑞̇+[Λ+(𝑣²+𝑇)𝐶]𝑞=−(2𝛽^0.5𝑣𝜇cos𝜔𝜏𝐵)𝑞̇−[𝑣²2𝜇cos𝜔𝜏)𝐶+𝛽^0.5̈̇𝑣𝜇sin𝜔𝜏(𝐷−𝐶)]𝑞,𝑞+𝐺𝑞+𝐾𝑞=𝑓₁(𝜏)𝑞̇+𝑓₂(𝜏)𝑞, (2.8)(2.9) where 𝐺=𝐻Λ+2𝛽^0.5𝑣𝐵,𝐾=Λ+(𝑣²𝑓+𝑇)𝐶,₁(𝜏)=−2𝛽^0.5𝑓𝑣𝜇cos𝜔𝜏𝐵,₂(𝜏)=−2𝑣²𝜇cos𝜔𝜏𝐶−𝛽^0.5𝑣𝜇𝜔sin𝜔𝜏(𝐷−𝐶). (2.10)

In (2.9), 𝐺 and 𝐾 denote the structural damping matrix and stiffness matrix, respectively. These two matrices are associated with systematic parameters, such as dimensionless flow velocity 𝑣 and mass ratio 𝛽. 𝜆_𝑖 (𝑖=1,2) are the 𝑖th eigenvalues of the supported pipe and Λ is the diagonal matrix with elements 𝜆⁴_𝑖.

Moreover, 𝐵, 𝐶, and 𝐷 are matrices with elements 𝑏_𝑠𝑟, 𝑐_𝑠𝑟, and 𝑑_𝑠𝑟 (𝑠,𝑟=1,2), respectively. They are defined as𝑏_𝑠𝑟=¹₀𝜙_𝑠𝜙_𝑟d𝜉,𝑐_𝑠𝑟=¹₀𝜙_𝑠𝜙_𝑟^′′d𝜉,𝑑_𝑠𝑟=¹₀𝜙_𝑠𝜉𝜙𝜂_𝑟^′′d𝜉. (2.11)

Different value should be taken for those three parameters depending on different boundary conditions of the pipe. For the pinned-pinned pipe, we have 𝑏_𝑠𝑟=2𝜆_𝑟𝜆_𝑠𝜆²_𝑟−𝜆²_𝑠{(−1)^𝑟+𝑠𝑐−1}𝑠≠𝑟,0𝑠=𝑟,_𝑠𝑟=0𝑠≠𝑟,−𝜆²_𝑟𝑑𝑠=𝑟,_𝑠𝑟=4𝜆³_𝑟𝜆_𝑠(𝜆²_𝑟−𝜆²_𝑠)²{1−(−1)^𝑟+𝑠1}𝑠≠𝑟,2𝑐_𝑟𝑟𝑠=𝑟. (2.12)

2.2. 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 𝑄 and corresponding momentum 𝑃, 𝑉̇=𝐻𝑉+𝐹, (2.13) where 𝑉=𝑄𝑃^𝑇,𝑄=𝑞₁𝑞₂^𝑇,𝑃=𝑝₁𝑝₂^𝑇=𝑄̇+𝐺𝑄2,−𝐺𝐻=2𝐼𝐺−𝐾−^𝑇𝐺4(𝐺2)^𝑇,𝐹=00𝑓₁𝑃+(𝑓₂−𝑓₁𝐺)𝑄^𝑇. (2.14)

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

3. 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.

3.1. Integration of Homogeneous Equation

The general solution of homogeneous equation 𝑉̇=𝐻𝑉(𝑉₀=𝑉(0)) can be expressed as𝑉(𝑡)=𝑒^𝐻𝑡⋅𝑉₀. (3.1)

Suppose that the time step is 𝜏=𝑡_𝑘+1−𝑡_𝑘, and then we have the following recursive steps:𝑉₁=𝑇𝑉₀,𝑉₂=𝑇𝑉₁,…,𝑉_𝑘+1=𝑇𝑉_𝑘, (3.2) where 𝑇=𝑒^𝐻𝜏. It is seen that how to compute the exponential matrix 𝑇 is essential for the integration precision.

Then split the time interval 𝜏 into a smaller one. Define Δ𝑡=𝜏/𝑚 and 𝑚=2^𝑁. For example, 𝑚=1048576 when 𝑁=20. As 𝜏 is small, Δ𝑡 is an extremely small time interval.

Assume that 𝑇_𝑎=(𝐻Δ𝑡)(𝐼+𝐻Δ𝑡)/2, execute the cycle For(𝑖=0;𝑖<𝑁;𝑖++){𝑇_𝑎=2𝑇_𝑎+𝑇_𝑎×𝑇_𝑎}, (3.3) where 𝐼 is the identity matrix.

So 𝑇_𝑎 is no longer a very small matrix. It could be computed by the following function:𝑇=𝐼+𝑇_𝑎. (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].

3.2. 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 𝑉̇(𝑡)=𝐻𝑉(𝑡)+𝐹(𝑉,𝑡). (3.5)

Then its solution could be given by the Duhamel’s integration as𝑉(𝑡)=𝑒^𝐻𝑡𝑉₀+^𝑡₀𝑒^{𝐻(𝑡−𝜏)}𝐹(𝑉,𝜏)𝑑𝜏. (3.6)

Similarly, the duration of structural dynamic response is also divided into small time intervals. The response between (𝑡_𝑘,𝑡_𝑘+1) can be written as𝑉(𝑡_𝑘+1)=𝑇𝑉_𝑘+^𝑡𝑘+1_𝑡𝑘𝑒^{𝐻(𝑡𝑘−𝜏)}𝐹(𝑉,𝜏)𝑑𝜏. (3.7)

To solve this inhomogeneous equation, the analytical expression of the time-variant item 𝐹(𝑉,𝑡) is required. But it is not always available.

In this study, the linear interpolation formula is used to approximate this nonlinear item.

Assume that𝐹_𝑘+1=𝐹_𝑘+(𝑡−𝑡_𝑘)𝐹̇_𝑘, (3.8) where𝐹_𝑘=𝐹(𝑉,𝑡_𝑘),𝐹̇_𝑘=(𝜕𝐹𝜕𝑡+_𝑛𝜕𝐹𝜕𝑉_𝑖𝜕𝑉_𝑖)𝜕𝑡_{𝑡=𝑡𝑘}. (3.9)

Substituting (3.8) into (3.7) gives the linear interpolation expression 𝑉(𝑡_𝑘+1)=𝑇[𝑉_𝑘+𝐻⁻¹(𝐹_𝑘+𝐻⁻¹𝐹̇_𝑘)]−𝐻⁻¹[𝐹_𝑘+𝐻⁻¹𝐹̇_𝑘+𝜏𝐹̇_𝑘]. (3.10)

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.

4. Numerical Examples

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

4.1. Dynamic Response of Pipes Conveying Stable Fluid

In this subsection, this approach is used to analyze the dynamic response of stable fluid-conveying 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 𝑣 is 2.0, the mass ratio 𝜷 is 0.32 and the fluid pressure 𝑇=1. The initial conditions are chosen to be [𝑞₁𝑞₂𝑞̇₁𝑞̇₂]^𝑇=[−0.10.20.10.4]^𝑇. Time increases from 0 to 1000 s and the time step Δ𝑡 is selected as 0.2 s.

Figure 2 illustrates time history of four state variables (𝑞₁𝑞₂𝑝₁𝑝₂) of pipe’s middle point using the Runge-Kutta method, while Figure 3 shows the results calculated by precise integration method.

(a)

(b)

(c)

(d)

Figure 2 

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

(a)

(b)

(c)

(d)

Figure 3 

Time history of four state variables using 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.

Table 1 

Computation time needed using two methods.

4.2. 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 𝑣 is 1, the amplitude 𝜇=0.4, the frequency 𝜔=2.5, the mass ratio 𝜷=0.32, and the fluid pressure 𝑇=1. The initial conditions are chosen to be [𝑞₁𝑞₂𝑞̇₁𝑞̇₂]^𝑇=[−0.0010.00100]^𝑇. Time increases from 0 to 100 s, and two time steps involved in this example are selected as Δ𝑡=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.

(a)

(b)

Figure 4 

Displacement response and phase diagram of pipe’s middle point using Runge-Kutta $(\mathrm{\Delta }t=0.01)$.

(a)

(b)

Figure 5 

Displacement response and phase diagram of pipe’s middle point using PIM $(\mathrm{\Delta }t=0.01)$.

(a)

(b)

Figure 6 

Displacement response and phase diagram of pipe’s middle point using Runge-Kutta ($\mathrm{\Delta }t=0.05$).

(a)

(b)

Figure 7 

Displacement response and phase diagram of pipe’s middle point using PIM $(\mathrm{\Delta }t=0.05)$.

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 long-term dynamic analysis of fluid-conveyed pipes.

Table 2 

Computation time needed using two methods.

4.3. 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 𝑣=1.0, 1.5, and 2.0. The fluid frequency step and the time step are selected as Δ𝜔=0.4 and Δ𝑡=0.01s, respectively. Others parameters are the same with the preceding Section 4.2.

The displacement response 𝑞₀ 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.

Figure 8 

Displacement response of pipe’s middle point as the function of the dimensionless pulsating fluid frequency $(v=1.0)$.

(a)

(b)

Figure 9 

Phase diagram of pipe’s middle point as $v=1.0$. (a) $\omega =5$; (b) $\omega =17$.

Figures 10 and 11 show the displacement response variations of the middle point versus the pulsating fluid frequency as 𝑣=1.5 and 2.0, respectively. It can be shown that the instable zone is changing with fluid frequencies. As 𝑣=1.5, the instable zones are at (14,18), (44,51), and (68,70). As 𝑣=2.0, the instable zones lie at (6,7), (11,16), (42,50), and (65,70).

Figure 10 

Displacement response of pipe’s middle point as the function of the dimensionless pulsating fluid frequency $(v=1.5)$.

Figure 11 

Displacement response of pipe’s middle point as the function of the dimensionless pulsating fluid frequency $(v=2.0)$.

It can be seen that with the increasing of fluid velocity, the critical fluid frequency gets smaller and the pipe shows complicated nonlinear vibration.

5. 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.