`Advances in Mechanical EngineeringVolume 2013 (2013), Article ID 890423, 14 pageshttp://dx.doi.org/10.1155/2013/890423`
Research Article

## Numerical Simulation of Vortex-Induced Vibration with Three-Step Finite Element Method and Arbitrary Lagrangian-Eulerian Formulation

1State Key Laboratory of Coastal and Offshore Engineering, Dalian University of Technology, Dalian 116024, China
2Center for Deepwater Engineering, Dalian University of Technology, Dalian 116024, China

Received 26 August 2013; Accepted 11 September 2013

Copyright © 2013 Guoqiang Tang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

Numerical simulations were performed in this paper to investigate an elastically mounted circular cylinder subjected to vortex-induced vibration (VIV). A three-step finite element method (FEM) is introduced for solving the incompressible fluid flow equations in two dimensions. The computational procedure is coupled with a mesh movement scheme by use of the arbitrary Lagrangian-Eulerian (ALE) formulation on account of the body motion in the flow field. On running the numerical simulations, the Reynolds number was kept constant of Re = 100 and the reduced velocity = U/(D) was varied from 3.0 to 10.2 by changing the natural frequency of the cylinder. The mass ratio m* = 4m/ρπD2 and damping ratio ξ are set to be 10.0 and 0.01, respectively, where U is free-stream velocity, D the diameter of the circular cylinder, the mass of the cylinder per unit length, and ρ the density of the fluid. Numerical results are examined for the response amplitude of transverse direction as well as the phase angle, φ, between the lift force and the transverse displacement of the cylinder. The numerical results reveal that the transverse amplitudes present only two branches, namely, initial branch and lower branch, rather than three branches as the results obtained from high-Re experiments with low . On the other hand, the phase angles present almost linear increase with the reduced velocity in the synchronization region. However, experiments concerned with high Re exhibit a sudden jump in phase angle of approximate . The difference between the present study and the high-Re experiment is attributed to no substantial vortex shedding mode transition at the present numerical results of laminar flow.

#### 1. Introduction

Vortex-induced vibrations of bluff bodies have been studied for a long time due to thier importance in both academic researches and engineering applications. The prediction of offshore structures subjected to VIV is one of the most challenging tasks since the complicated fluid-structure interaction phenomenon is involved. For the purpose of understanding this issue, the fluid flow over an elastically-mounted circular cylinder has received numerous attentions in the past decades. It is well known that the lock-in or synchronization will occur as the vortex shedding process is governed by the natural frequency of the cylinder and results in the attendant failure in fatigue. Comprehensive comments and reviews focused on the various aspects of VIV can be found in Bearman [1], Sarpkaya [2], Gabbai and Benaroya [3], and the recent review of Williamson and Govardhan [4].

Making an attempt at understanding the fundamental phenomena involved in VIV, a majority of experimental investigations have been conducted to concentrate on a self-excited cylinder exposed to the incoming flow. Feng [5] was among the first researchers to investigate the VIV problem with a single degree of freedom flexible cylinder at high Reynolds numbers (~). In the experimental results, the jump in the phase angle φ was observed when the lock-in occurred. In addition, the transverse amplitude with respect to reduced velocity showed the hysteresis effect by continuously increasing and decreasing the fluid speed. Brika and Laneville [6] investigated the VIV of a flexible cylinder through experimental study (Re = 3400~11800). Two different response amplitude branches were identified. The transition in the different branches was accompanied by the change in the vortex shedding mode. Khalak and Williamson [7] found that there are two distinct types of amplitude response for the elastically mounted cylinder limited to oscillate in transverse direction, depending on the values of the mass-damping parameter . For the low mass-damping parameter , three different amplitude branches were detected, defined as the initial branch, the upper branch, and the lower branch. However for the high mass-damping parameter , only initial branch and lower branch exist in the transverse displacement which is similar to the previous observations [5].

It is worth pointing out that all the above mentioned experimental works focused on high Reynolds numbers. The fluid flow possesses evidently three-dimensional characteristics. Apart from the experimental investigations, the computational fluid dynamics (CFD) method based on Navier-Stokes equation is an alternative choice for understanding the mechanism of VIV. As authors know, only few sets of numerical simulations have been conducted in three dimensions related to high Re due to the limitation of computer resources [8, 9]. Although the fluid flow is inherently three-dimensional, the two-dimensional computations still can provide useful insight on the various aspects associated with the VIV problem [10]. Several two-dimensional numerical simulations were conducted to study the VIV at high Re in conjunction with the necessary turbulence model [1113]. However, in the present study, our way of understanding VIV is performed at the low Re in consideration of few studies concentrating on the main concern of this paper. Mittal and Kumar [14] carried out numerical simulations of an elastically mounted circular cylinder subjected to VIV at . They reported a soft-lock-in phenomenon, defined as the vortex shedding frequency, that does not exactly match the structural frequency. Shiels et al. [15] suggested a new parameter, defined as effective elasticity , to demonstrate the dynamic response of a cylinder subjected to VIV at . The response displacements and hydrodynamic coefficients with respect to were examined. Following the above idea, Placzek et al. [16] also conducted numerical simulations of laminar fluid flow over a cylinder restricted to oscillate in transverse direction at . Similar to the work of Shiels et al. [15], they examined the response displacement and hydrodynamic coefficients. In addition, the vortex shedding mode was addressed in their work.

The objective of this work is to investigate the VIV of an elastically mounted circular cylinder freely oscillating in both cross-flow and in-line directions. A three-step finite element method coupled with ALE formulation was developed for solving the fluid-structure problem. In the present study, we concentrate our thoughts on laminar fluid flow with and mass-damping parameter . It is hoped that the present study can contribute to not only the various aspects involved in VIV at low Re but also to the ascertainment of the similarities and differences of VIV between the low Re and high Re. On the basis of present numerical solution, the dynamic responses of the oscillating cylinder are analyzed. At the same time, the timing of vortex shedding patterns is also investigated in this work.

#### 2. Numerical Model

##### 2.1. Governing Equations and Numerical Discretization

The governing equations of the laminar fluid flow, represented by the nondimensional continuity equation and Navier-Stokes equations, can be written in the following ALE formulation [17]: where is the fluid velocity components correlating with the Cartesian coordinate in the present two-dimensional numerical model (i.e., , corresponding to , ) and and denote pressure and time, respectively. The Reynolds number is defined as with being the dynamic viscosity of the fluid. Under the ALE description, the convection term in Navier-Stokes equations is revised so as to accommodate the moving boundary conditions; it reads as follows: where is the velocity components of the moving mesh occupied by the fluid.

In the present study, the standard finite element method is adopted to generate the solution of Navier-Stokes equations. As is well known, the standard finite element method shows weak performance for solving the fluid problem in comparison with the solution of solid mechanics due to the convection dominated characteristic in the fluid flow [18]. To overcome this problem, the upwind finite element methods are usually employed, such as the SUPG (Streamline Upwind/Petrov-Galerkin) method and Taylor-Galerkin method. In the present numerical investigation, a three-step finite element method affiliated to Taylor-Galerkin method is utilized to discretize the Navier-Stokes equations. The three-step finite element method is first proposed by Jiang and Kawahara [19] which exhibited good performance to produce stable numerical results for advection-diffusion problems because of its inherently higher-order upwind scheme [20]. However, the original three-step finite element method was developed for the fixed computational meshes from the Eulerian point of view. Therefore, it is extended to the ALE reference coordinates as described in (2) in order to consider the moving boundary conditions. Following the idea of three-step finite element method, the momentum equation can be discretized as follows: where the superscripts , , , and illustrate the time levels at , , , and () , respectively. It should be noted here that the calculations of the at different time levels are in terms of the same mesh moving velocity as the th time step.

Observation of (6) suggests that the pressure at time level has to be obtained beforehand in order to have the solutions of velocity . In the present numerical model, by means of applying the divergence operation on both sides of (6) and forcing the mass conservation , the Poisson-type pressure equation can be obtained as follows:

Under the description of three-step finite element method, the velocity and pressure in Navier-Stokes equations are decoupling based on the mentioned high-order Taylor-Galerkin scheme. Therefore, the solution of fluid velocity and pressure can be performed independently through the standard Galerkin finite element method. Consequently, the overall procedures of the three-step finite element algorithm for solving the present two-dimensional VIV can be summarized as the following steps.(1)Solving the momentum equation (4) to obtain the intermediate velocity .(2)Calculating the velocity based on (5) with the available , , and .(3)Solving the Poisson-type pressure equation with BI-CGSTAB method [21] to obtain the distribution of pressure at time level.(4)Correcting the fluid velocity with (6) and getting the final .(5)Calculating the structural response with the available external force.(6)Updating the computational meshes and getting the grid velocities for the next time step.

The time-dependent drag and lift coefficients, normalized by , are calculated by the integration operation of the instantaneous pressure and vorticity along the surface of the cylinder: where is the angle from the positive axis and is the local vorticity.

The Strouhal number St, denoting the nondimensional vortex shedding frequency, is defined as which can be evaluated through the Fast Fourier Transform (FFT) analysis of the time-dependent lift force with being the vortex shedding frequency.

##### 2.2. Computational Domain and Boundary Conditions.

The sketch definition for oscillating cylinder subjected to VIV in uniform flow is shown in Figure 1. The cylinder is placed in a rectangular computational domain. The boundaries of upstream and downstream are located at and relative to the center of the cylinder, respectively, while the half breadth in the lateral directions is set to be . For the purpose of reducing the computational efforts associated with the ALE mesh update, a rather smaller square subdomain enclosing the maximal possible oscillation of the cylinder is used, referring to the ALE domain in Figure 1, where the mesh deformation is carried out in order to accommodate the motion of the circular cylinder. Otherwise, outside the ALE domain, the fixed Eulerian meshes are adopted, allowing the elemental coefficient matrix to be formed beforehand and used repetitively for the computational advance.

Figure 1: Computational domain of incompressible viscous flow over elastically mounted circular cylinder.

The boundary conditions are imposed as follows: the no-slip boundary condition is specified on the surface of the cylinder, that is, , , which is the same as the velocities of the cylinder oscillation, and , where is the outward unit vector. At the inlet, the dimensionless velocity components , are prescribed together with . As for the two lateral boundaries, far away from the oscillating cylinder, the approximations of , , and are imposed. Along the outlet boundary, the free outflow boundary condition is applied, namely, , , and . The computations start from the still water state, implying the initial conditions with zero pressure and velocities.

##### 2.3. Grid Update Strategy

In the present study, the circular cylinder is assumed to be mounted on springs in both cross-flow and in-line directions with the same structural attribution. Hence, the motion of the body is governed by the following linear oscillator equations: where , , and represent the mass of the cylinder per unit length, structural damping, and the rigidity of the spring, respectively. , , and denote the acceleration, velocity, and displacement of the circular cylinder in the in-line direction, respectively, while , , are the cross-flow acceleration, velocity, and displacement, respectively. Employing the following relationships, , and , the oscillator equations can be rewritten as

The nondimensional forms of the variables in (11) can be defined by

Inserting (13) into (11), the final nondimensional governing equations can be obtained as

The reduced velocity is defined as . In the present study, the reduced velocity is varied from 3.0 to 10.2 by changing the natural frequency of the circular cylinder. Solving (14) can obtain the dynamic response of the cylinder due to VIV. Under the ALE description, the nodal coordinates inside the ALE domain should be updated to accommodate the motion of the oscillating cylinder. In the present model, the Laplace equation is employed to perform the grid update at each time step by assuming that the grids in the ALE domain are composed of linear elastic materials [22]. Therefore, the ALE Cartesian coordinates can be solved by where and are the horizontal and vertical coordinates after the ALE grid update, respectively, and is the two-dimensional Laplacian operator. It can be discretized by the standard Galerkin finite element method and solved by the same subroutine used for the pressure equation. During the numerical solution, the unchanged coordinates are prescribed along the outer boundaries of the ALE domain, while the new coordinates on the cylinder surface are imposed as the boundary conditions. The updated coordinates and are remapped to the Eulerian computational meshes, denoted by for the flow simulation at the time level. Therefore, the velocities of the moving mesh grids in the ALE domain can be evaluated as

#### 3. Numerical Validations

In this study, the numerical model is first verified by considering the grid independence. The main variables of the present study are the hydrodynamic coefficients as well as the displacements. Therefore, for the purpose of validating the accuracy of the present numerical model, two different cases are selected to compare the above mentioned parameters with the previous available work. The first case is the flow over a fixed cylinder at . After that, the case of an elastically mounted circular cylinder oscillating in transverse direction behind a stationary cylinder is considered at .

##### 3.1. Spatial and Time Convergence

For the purpose of establishing the grid independence of the present numerical model, numerical calculations were conducted by considering four different mesh sizes. The number of elements and nodes of different grids utilized is shown in Table 1. Here, the extreme situation is considered with , , , and which produces the largest response displacement in the transverse direction. The detailed numerical results with respect to the response displacements and hydrodynamic coefficients are also listed in Table 1, where , denote the number of grid points along the surface of the cylinder and in the radial direction. , , , , , and represent the root mean square (RMS) transverse displacement, mean in-line displacement, RMS in-line displacement, mean drag force, RMS drag, and lift force, respectively. In addition, during the numerical computations, the dynamic time step is utilized in this model, which is determined by where is the area of a computational cell, is the absolute velocity at the center of the mesh cell, and the function Min(*) denotes the minimal value of all the computational meshes. The numerical stability is further guaranteed by a safe coefficient considering that the present three-step finite element scheme proves to be consistent with the uniform CFL condition.

Table 1: Grid independence test for freely oscillating cylinder at Re = 100, = 5.2, m* = 10.0, and = 0.01.

It is evident in Table 1 that the numerical results of the concerned variations are not sensitive to the four different grid sizes. Therefore, Mesh 3 is selected throughout the numerical simulations from accurate and efficient points of view.

##### 3.2. Numerical Validations

The numerical code is first validated by considering laminar fluid flow over a fixed circular cylinder at . The corresponding numerical results are shown in Table 2. The present results exhibit good agreement with the numerical results obtained from the references as displayed in Table 2, suggesting that the present numerical model works well.

Table 2: Comparisons of hydrodynamic forces and St at Re = 100.

The present numerical model is further validated against an elastically mounted circular cylinder oscillating in transverse direction behind a stationary cylinder at . Two cylinders with the same diameter are arranged in tandem form and the center-to-center distance between them is = 3.0D. The comparisons of the nondimensional oscillating amplitude at various reduced velocities are shown in Figure 2. The mass ratio and damping ratio are set to be 2.0 and 0.007, respectively. It can be seen that the transverse displacement of the rear cylinder increases rapidly with the reduced velocity and attains the maximum around , which means that the lock-in happens. After that, the displacement response decreases slowly and approaches to constant value at the extremely large reduced velocity. Figure 2 shows that the numerical results of this work are in good agreement with those obtained from Carmo et al. [26], indicating that the present numerical model is adequate to produce enough accuracy to calculate the VIV problem.

Figure 2: Comparison of vortex-induced vibration of elastically mounted circular cylinder in steady flow.

#### 4. Numerical Results and Discussion

In the present study, numerical simulations were conducted for VIV of an elastically mounted circular cylinder freely oscillating in both cross-flow and in-line directions. All of the numerical results herein were obtained at , mass ratio , and structural damping ratio , leading to a low combined mass-damping parameter . The reduced velocity was varied from 3.0 to 10.2, and the change in the reduced velocity was achieved by altering the natural frequency of the cylinder.

##### 4.1. Dynamic Response of the Cylinder

The numerical results associated with the maximum transverse displacement , the response frequency , and the phase angle between lift force and displacement are shown in Figure 3.

Figure 3: Vortex excited dynamic response of the cylinder in transverse direction at different reduced velocities.

The variation of the maximum transverse displacement with respect to reduced velocity is shown in Figure 3(a). It shows that there is no evident upper branch as defined by Khalak and Williamson [7] in the response displacement. Instead, only initial branch and lower branch are clearly seen which coincides with experimental results reported by Anagnostopoulos and Bearman [27] for laminar fluid flow with very low . Khalak and Williamson [7] reported that the has a great influence on the response amplitude branch. For the low , the response amplitude should exhibit three different branches which is different from the present numerical investigations concerned with laminar flow. Therefore, by the collaboration of the present numerical simulation and those obtained from the experimental results of Anagnostopoulos and Bearman [27] and Khalak and Williamson [7], it is inferred that is not only the factors responsible for the number of branches in the amplitude response but also the Reynolds number. Further inspection reveals that the peak value of the transverse displacement is approximately 0.5D occurring at which illustrates the onset of the synchronization. In addition, the response displacement undergoes sudden transition from initial branch to lower branch at and without having an intermediate value of reduced velocity.

Figure 3(b) displays the frequency response against reduced velocity. In this picture, is the vortex shedding frequency obtained from the FFT analysis of the time-dependent lift force. One of the important phenomena associated with VIV is the synchronization. It has been pointed out that the synchronization can be defined as the vortex shedding frequency matches the natural frequency of the cylinder [18]. This is well demonstrated in Figure 3(b) where the two vertical dashed lines illustrate the boundaries of the synchronization region. Over the synchronization range, , it is observed that the vortex shedding frequency is modulated to the natural frequency of the elastically mounted circular cylinder and the oscillations are dominated by the oscillating frequency . Although only one frequency is plotted in Figure 3(b), it is not suggestive that a single frequency is involved in the lift force. In fact, the interaction between the oscillating frequency and the natural vortex shedding frequency usually introduces some additional frequencies for not only the synchronization cases but also for the nonsynchronization ones. Typical examples can be found over the range of reduced velocities . In this reduced velocity range, multiple frequencies are presented in the spectra of the response making the lift and displacement fluctuations behave large deviation from the pure-tone oscillations. A similar phenomenon was also observed by Blackburn and Henderson [10] which is divided into three different oscillating patterns according to the periodical characteristics in the oscillations. This will be demonstrated later.

The averaged phase angle between lift force and response displacement with reduced velocity is plotted in Figure 3(c). Before the onset of lock-in, , the phase angle is found to be around . After that, a linear increase in the phase angle is observed in the synchronization region rather than the sudden jump around 180° observed in the high-Re experiment. In other words, the phase angle shows great dependence on the reduced velocity in this region. As for the transition from the initial branch to lower branch, the phase angle jumps from 129.82° to 161.96°. For the lower branch, the averaged phase angle presents to be independent on the reduced velocity and keeps almost constant values of approximate 180°.

As the above mentioned numerical results state, the present study exhibits two distinct aspects in the dynamic response in comparison with the results obtained from the high-Re experiment. First, only two amplitude response branches are displayed in the transverse displacement with low . Secondly, the phase angle gradually ascends with the increasing reduced velocity in the synchronization rather than an abrupt jump was observed in the high-Re experiments. We now want to make a brief discussion about what is reason to result in the difference of dynamic response between the low-Re two-dimensional flow and the high-Re experiment. Inspections of the high-Re experiment [7] as well as the two-dimensional numerical simulations [12, 13] concerned with high Re reveal that the jump in the phase angle is usually accompanied by the transition of 2S vortex shedding mode to 2P. Here, 2S vortex shedding mode is defined as the vortices with opposite signs and shed alternately from each side of the cylinder, while in the 2P vortex shedding mode, two vortex pairs were shed from the cylinder as the interpretation of Williamson and Roshko [28]. However, in the present study, only 2S vortex shedding mode is detected. It is inferred that the strong dissipation effect for the laminar flow prohibits the formation of 2P vortex shedding mode, the consequent absence of the upper branch, and the jump in the phase angle.

In Figure 4, we present the numerical results of RMS displacement as well as the mean displacement of in-line direction. Again, the two vertical dashed lines demonstrate the synchronization region.

Figure 4: Variations of the RMS and mean displacement of in-line direction versus reduced velocity.

In Figure 4(a), compared with the transverse displacement, the RMS displacement of in-line direction behaves extremely small value. This is the reason why the in-line dynamic response is often neglected in the previous study. However, the oscillating frequency of in-line is usually two times of that of the transverse direction, so it is believed that the stream-wise oscillation may produce the similar effect on the fatigue damage as the transverse direction [29]. Further examination of Figure 4(a) shows that two peaks are clearly seen. They are adjacent to the lower and higher reduced velocity boundaries of the synchronization region. Figure 4(b) displays the mean displacement of in-line direction. From this figure, it is observed that, with the increasing of reduced velocity, the mean displacement increases continuously. However, a lower slope is found in the synchronization region than the other reduced velocity ranges.

##### 4.2. Hydrodynamic Forces

In Figure 5, the numerical results of hydrodynamic forces associated with RMS lift force, RMS drag force, and mean drag force are presented. It should be noted here that the hydrodynamic forces are normalized by the corresponding results obtained for a fixed cylinder at the same Reynolds number. The two vertical dashed lines represent the lock-in region as displayed in Figure 3.

Figure 5: Variations of hydrodynamic forces at different reduced velocities, where , , and represent the root mean square lift, drag force, and mean drag force obtained from the results of the fixed cylinder at , respectively.

Figure 5(a) displays associated with different reduced velocities. Compared with the fixed cylinder, at the beginning of initial branch, VIV makes the oscillating cylinder experiences slightly larger value of than that of the fixed cylinder until . After that, the undergoes rapidly increase with the increasing of reduced velocity and attains maximum value at which is adjacent to the onset of the lock-in region. At this point, we find that the of the freely oscillating cylinder is approximately four times that of the stationary cylinder. In the lock-in region, with the increasing of reduced velocity, the decreases continuously, but it still holds larger values than the fixed cylinder at the lower reduced velocity range of the synchronization. However, it is found that transition occurs at a critical reduced velocity for the . Beyond this critical reduced velocity, the is small in comparison with the result of the fixed cylinder for the higher reduced velocity range in the synchronization. For the lower branch with , it is evident that the switch from the initial branch to lower branch results in the change of trend in the . Over this reduced velocity range, increases slowly and approaches to almost constant value at the large reduced velocities. The for the fixed cylinder is found to be approximately 1.3 times that of the oscillating cylinder subjected to VIV.

Figures 5(b) and 5(c) show the variables of and , respectively. Contrasted to the variable , exhibits the similar trend with the reduced velocity. However, in the synchronization region, it is observed that almost linearly decreases with the increasing of the reduced velocity rather than the parabolic decrease in . Moreover, the maximum value of the is found to be 55 times of that for the fixed cylinder. As for the as presented in Figure 5(c), it is observed that, at the beginning of the initial branch, the is approximately equal to the . After a rapid increase, reaches its maximum value at the exact onset of the synchronization. For the lower branch, it is evident that the increasing in the reduced velocity almost has limited effect on the which is slightly smaller than 1.0.

##### 4.3. Time Traces of Displacement and Hydrodynamic Coefficient

In the above mentioned study, we report a multiple-frequency phenomenon over the range of reduced velocity . Here, three typical examples are selected to demonstrate this phenomenon; they are , 4.9, and 5.0, respectively. The corresponding numerical results with respect to hydrodynamic force and response displacement are presented in Figure 6. In this picture, the left column represents the time traces of lift force, while the response displacements for transverse direction are shown in the right column.

Figure 6: Time traces of hydrodynamic coefficients and displacement at three typical reduced velocities: (a) , (b) , and (c) .

Time traces of lift force and transverse displacement with respect to are displayed in Figure 6(a). It is evident that both the lift force and transverse displacement clearly behave the beating oscillations and share the same beating period. Further, increasing to , the absence of periodical repetition discloses the characteristic of random oscillation. As for Figure 6(c) associated with , long-period oscillations are evident for both the lift force and transverse displacement giving rise to the classification of this case as weakly chaotic defined by Blackburn and Henderson [10] but concentrating on a cylinder forced oscillating in transverse direction.

One of the most important aspects involved in VIV is the synchronization. Hence the numerical results associated with the time traces of transverse displacement and lift force at various reduced velocities in the synchronization region are presented in Figure 7.

Figure 7: Time history of lift force and transverse displacement in the synchronization.

As shown in Figure 7, the amplitude of displacement is observed to decrease with the increase of reduced velocity which is coincident with the observations of Figure 3(a) in the synchronization. All the displacements are in purely harmonic fluctuations. As for the lift forces, it is evident that the oscillating cylinder introduces an additional frequency resulting in that the lift forces are slight deviation from the pure-tone oscillations. Figure 8 shows the spectra analysis of the lift force time series corresponding to Figure 7. In Figure 8, the oscillating frequency of the lift force is normalized by the natural frequency of the cylinder. It is found that the primary frequency, which exhibits the highest amplitude in the spectra, peaks at the natural frequency of the cylinder. However, the oscillating cylinder gives rise to not only the primary frequency but also the frequency component of three times of the primary frequency . Both the two frequencies result in the lift force is slight deviation from the simple harmonic oscillation.

Figure 8: Spectra of lift forces at various reduced velocities.
##### 4.4. Timing of Vortex Shedding

From the above discussion, we find that the vortex shedding mode has a great influence on the dynamic response of cylinder subjected to VIV. To support this claim, in Figure 9, the instantaneous vorticity contours are shown at various reduced velocities. The vorticity contours cover three different phase angle ranges. The first range is at and 4.6 with the phase angle near 0°. The second range corresponds to the gradual increase of phase angle in the synchronization with , 5.8, 6.2, 6.6, 7.1, respectively. At last, and 9.0 are selected for the lower branch where the phase angle is approximately 180°. It should be mentioned here that the vorticity contours are taken at the instant where the cylinder is at its maximum transverse displacement. The negative vorticity is represented by blue color, while the red color illustrates the positive vorticity.

Figure 9: Vorticity contours at different reduced velocities.

It is seen in Figures 9(a) and 9(b) that, corresponding to the beginning of the initial branch as presented in Figure 3(a) with and 4.6, the classical Von Kármán vortex sheet is displayed for these two cases, that is, 2S vortex shedding mode. Further increasing to which is responsible for the onset of synchronization and the maximum transverse displacement, compared with the vorticity contour at and 4.6, the vortex shedding develops into two rows and coalesces in the far wake instead of a single vortex suggesting the occurrence of the (2S) vortex shedding mode as reported by [30]. In the synchronization, the oscillations grow gradually as presented in Figure 3(a); the evolution of vortex shedding with the increase of reduced velocity is displayed by Figures 9(b)9(g). The vortex shedding transition from (2S) to 2S shedding mode is clearly seen with the increasing of reduced velocity. As the above discussion, the vortex coalesces in the far wake at . Although two-row configuration of vortex shedding is shown at , there is no evident coalescence effect in the far wake suggesting again the occurrence of the 2S vortex shedding mode. For the values of between 6.2 and 7.1, with the decrease of response displacement shown in Figure 3(a), the trend of vortex shed in a single row aligned in the medial plane is clearly exhibited. Finally, the classical Von Kármán vortex sheet is presented again at which is associated with the low amplitude response of VIV. Further investigations reveal that the length of the vortex attached to the upper side of the oscillating cylinder continuously increased with the increasing . The increasing negative vorticity interacts with the positive vorticity generated on the base of the oscillating cylinder making the vortex affiliated to the lower side of the cylinder roll up with a long distance relative to the cylinder. The vortex shedding mode which is accompanied by the transition from initial branch to lower branch is displayed in Figure 9(h) by increasing the reduced velocity from to . From the above discussion, it is found that there is no substantial vortex shedding mode transition between such two different response branches which is alien from the experimental results of high Re [7]. In the high-Re experiment, it is observed that the transition in the response branch results in the change in the vortex shedding pattern from 2S to 2P as well as the jump in phase angle. Although the same 2S vortex shedding mode is observed for and , there is still difference between the vorticity contours. At , the amount of both the positive and negative vorticity generated in the shear layer is less than that of resulting in the reduction in the length of the opposite-sign vortex pair attached to the cylinder. Moreover, compared with vortex at , the lower vortex has rolled up tightly making the upper vortex split into two vortices. As for the vortices far away from the cylinder, an interesting phenomenon is found that the vorticity for presents almost mirror image of that for . For the lower branch cases with and , the same 2S vortex shedding mode is shown for these two reduced velocities. From the above discussion, we find that there are two different vortex shedding patterns, namely, (2S) and 2S, respectively. Observations suggest that the (2S) vortex shedding pattern is found at the high-amplitude oscillation, while 2S is responsible for the low-amplitude oscillations.

In Figure 6, the numerical results of both the hydrodynamic forces and response displacements display instantaneous characteristics which should be closely related to the timing of vortex shedding. Hence, in Figure 10, we display the vorticity contours corresponding to the different instants of transverse displacements at .

Figure 10: Vorticity contours at different instantaneous transverse displacements with .

In Figure 10(b) corresponding to rather small response displacement, the single vortex is shed from the upper side of the cylinder, while the vortex coalesces between each other at the lower side of the cylinder suggesting that this is a competitive shedding mode between 2S and (2S). The vorticity contours presented in Figures 10(c) and 10(d) are obtained at the moderate response displacements. Both the two cases show the 2S vortex shedding mode and slight difference between each other. At the large response displacement, the (2S) vortex shedding patterns are observed as shown in Figures 10(e) and 10(f). However, Figure 10(f) presents longer length in the vortex behind the oscillating cylinder than that in Figure 10(e). As time progresses to , the 2S shedding mode is observed in Figure 10(g) which is associated with small response displacement. Because of the long-period characteristic between the instantaneous displacement of Figures 10(b) and 10(h), the vorticity contour presented in Figure 10(h) shows the approximately same shedding pattern as Figure 10(b). From the above discussions, we find that the periodic vortex shedding process now is governed by the long period as presented in the displacement of Figure 10(a) instead of the oscillating frequency of the cylinder. On the other hand, it is inferred that the competition between the 2S and vortex shedding mode resulted in the random and weakly chaotic oscillations as displayed in Figures 6(b) and 6(c).

#### 5. Conclusions

Numerical results have been presented for a freely oscillating cylinder subjected to VIV with low combined mass-damping parameter in the regime of laminar flow . Computations were conducted by the three-step finite element method in conjunction with ALE formulation on account of the body motion in the fluid field. Although the present study focuses on the laminar fluid flow, two amplitude branches, namely, initial branch and lower branch, are identified which have been found by the previous research at high-Re experiment with high combined mass-damping parameter. However, with low combined mass-damping parameter , three different amplitude branches are exhibited in the high-Re experiment which is different from the present study. We attribute the absence of the 2P vortex shedding mode to the lack of the upper branch. The phase angle between the lift force and transverse displacement is also addressed in the present work. It is found that the phase angle gradually climbs with the increasing of the reduced velocities in the synchronization rather that an abrupt phenomenon at high-Re experiment. High-Re experiments reveal that the jump in the phase angle is usually accompanied by the substantialtransition of vortex shedding mode. However, very similar vortex shedding modes of (2S) and 2S were found in the present study. The corresponding numerical results show that the 2S vortex shedding mode is observed as the cylinder oscillates in low amplitude, while the high-amplitude oscillations are responsible for the occurrence of (2S) vortex shedding mode. It is found that the VIV can lead to large response hydrodynamic forces on the cylinder especially for the RMS drag force. The maximum value of the is found to be 55 times of that for the fixed cylinder at the same Reynolds number.

#### Acknowledgments

The financial supports from the Natural Science Foundation of Liaoning Province, China, with Grant no. 2013020075 and the Fundamental Research Funds for the Central Universities with Grant no. DUT13JS05 are acknowledged. This work is also partially supported by NSFC of China with Grant no. 51279029.

#### References

1. P. W. Bearman, “Vortex shedding from oscillating bluff bodies,” Annual Review of Fluid Mechanics, vol. 16, pp. 195–222, 1984.
2. T. Sarpkaya, “A critical review of the intrinsic nature of vortex-induced vibrations,” Journal of Fluids and Structures, vol. 19, no. 4, pp. 389–447, 2004.
3. R. D. Gabbai and H. Benaroya, “An overview of modeling and experiments of vortex-induced vibration of circular cylinders,” Journal of Sound and Vibration, vol. 282, no. 3–5, pp. 575–616, 2005.
4. C. H. K. Williamson and R. Govardhan, “A brief review of recent results in vortex-induced vibrations,” Journal of Wind Engineering and Industrial Aerodynamics, vol. 96, no. 6-7, pp. 713–735, 2008.
5. C. C. Feng, The measurement of vortex-induced effects in flow past a stationary and oscillating circular cylinder and D-section cylinders [M.S. thesis], University of British Columbia, 1968.
6. D. Brika and A. Laneville, “Vortex-induced vibrations of a long flexible circular cylinder,” Journal of Fluid Mechanics, vol. 250, pp. 481–508, 1993.
7. A. Khalak and C. H. K. Williamson, “Motions, forces and mode transitions in vortex-induced vibrations at low mass-damping,” Journal of Fluids and Structures, vol. 13, no. 7-8, pp. 813–851, 1999.
8. D. Lucor, H. Mukundan, and M. S. Triantafyllou, “Riser modal identification in CFD and full-scale experiments,” Journal of Fluids and Structures, vol. 22, no. 6-7, pp. 905–917, 2006.
9. N. Kondo, “Three-dimensional computation for flow-induced vibrations in in-line and cross-flow directions of a circular cylinder,” International Journal for Numerical Methods in Fluids, vol. 70, no. 2, pp. 158–185, 2012.
10. H. M. Blackburn and R. D. Henderson, “A study of two-dimensional flow past an oscillating cylinder,” Journal of Fluid Mechanics, vol. 385, pp. 255–286, 1999.
11. H. Al-Jamal and C. Dalton, “Vortex induced vibrations using Large Eddy Simulation at a moderate Reynolds number,” Journal of Fluids and Structures, vol. 19, no. 1, pp. 73–92, 2004.
12. E. Guilmineau and P. Queutey, “Numerical simulation of vortex-induced vibration of a circular cylinder with low mass-damping in a turbulent flow,” Journal of Fluids and Structures, vol. 19, no. 4, pp. 449–466, 2004.
13. J. B. V. Wanderley, G. H. B. Souza, S. H. Sphaier, and C. Levi, “Vortex-induced vibration of an elastically mounted circular cylinder using an upwind TVD two-dimensional numerical scheme,” Ocean Engineering, vol. 35, no. 14-15, pp. 1533–1544, 2008.
14. S. Mittal and V. Kumar, “Finite element study of vortex-induced cross-flow and in-line oscillations of a circular cylinder at low reynolds numbers,” International Journal for Numerical Methods in Fluids, vol. 31, no. 7, pp. 1087–1120, 1999.
15. D. Shiels, A. Leonard, and A. Roshko, “Flow-induced vibration of a circular cylinder at limiting structural parameters,” Journal of Fluids and Structures, vol. 15, no. 1, pp. 3–21, 2001.
16. A. Placzek, J.-F. Sigrist, and A. Hamdouni, “Numerical simulation of an oscillating cylinder in a cross-flow at low Reynolds number: forced and free oscillations,” Computers & Fluids, vol. 38, no. 1, pp. 80–100, 2009.
17. W. Dettmer and D. Perić, “A computational framework for fluid-structure interaction: finite element formulation and applications,” Computer Methods in Applied Mechanics and Engineering, vol. 195, no. 41-43, pp. 5754–5779, 2006.
18. M. R. H. Nobari and J. Ghazanfarian, “A numerical investigation of fluid flow over a rotating cylinder with cross flow oscillation,” Computers & Fluids, vol. 38, no. 10, pp. 2026–2036, 2009.
19. C. B. Jiang and M. Kawahara, “A three-step finite element method for unsteady incompressible flows,” Computational Mechanics, vol. 11, no. 5-6, pp. 355–370, 1993.
20. L. Lu, J.-M. Qin, B. Teng, and Y.-C. Li, “Numerical investigations of lift suppression by feedback rotary oscillation of circular cylinder at low Reynolds number,” Physics of Fluids, vol. 23, no. 3, Article ID 033601, 2011.
21. H. A. van der Vorst, “BI-CGSTAB: a fast and smoothly converging variant of BI-CG for the solution of nonsymmetric linear system,” SIAM Journal on Scientific and Statistical Computing, vol. 13, no. 2, pp. 631–644, 1992.
22. M. Sahin and K. Mohseni, “An arbitrary Lagrangian-Eulerian formulation for the numerical simulation of flow patterns generated by the hydromedusa Aequorea victoria,” Journal of Computational Physics, vol. 228, no. 12, pp. 4588–4605, 2009.
23. P. K. Standsby and A. Slaouti, “Simulation of vortex shedding including blockage by the random-vortex and other methods,” International Journal for Numerical Methods in Fluids, vol. 17, no. 11, pp. 1003–1013, 1993.
24. C. Norberg, “Fluctuating lift on a circular cylinder: review and new measurements,” Journal of Fluids and Structures, vol. 17, no. 1, pp. 57–96, 2003.
25. L. Baranyi and R. I. Lewis, “Comparison of a grid-based CFD method and vortex dynamics predictions of low Reynolds number cylinder flows,” Aeronautical Journal, vol. 110, no. 1103, pp. 63–70, 2006.
26. B. S. Carmo, S. J. Sherwin, P. W. Bearman, and R. H. J. Willden, “Flow-induced vibration of a circular cylinder subjected to wake interference at low Reynolds number,” Journal of Fluids and Structures, vol. 27, no. 4, pp. 503–522, 2011.
27. P. Anagnostopoulos and P. W. Bearman, “Response characteristics of a vortex-excited cylinder at low reynolds numbers,” Journal of Fluids and Structures, vol. 6, no. 1, pp. 39–50, 1992.
28. C. H. K. Williamson and A. Roshko, “Vortex formation in the wake of an oscillating cylinder,” Journal of Fluids and Structures, vol. 2, no. 4, pp. 355–381, 1988.
29. A. D. Trim, H. Braaten, H. Lie, and M. A. Tognarelli, “Experimental investigation of vortex-induced vibration of long marine risers,” Journal of Fluids and Structures, vol. 21, no. 3, pp. 335–361, 2005.
30. T. K. Prasanth and S. Mittal, “Vortex-induced vibrations of a circular cylinder at low Reynolds numbers,” Journal of Fluid Mechanics, vol. 594, pp. 463–491, 2008.