International Journal of Mechanical Sciences 127 (2017) 176–190

Contents lists available at ScienceDirect

International Journal of Mechanical Sciences

journal homepage: www.elsevier.com/locate/ijmecsci

2DOF CFD calibrated wake oscillator model to investigate vortex-induced MARK

vibrations
⁎
Andrey Postnikov, Ekaterina Pavlovskaia , Marian Wiercigroch
Centre for Applied Dynamics Research, School of Engineering, Aberdeen University, King's College, Aberdeen AB24 3UE, Scotland, UK

A BS T RAC T

In this study a new two degrees-of-freedom wake oscillator model is proposed to describe vortex-induced
vibrations of elastically supported cylinders capable of moving in cross-flow and in-line directions. The total
hydrodynamic force acting on the cylinder is obtained here as a sum of lift and drag forces, which are defined as
being proportional to the square of the magnitude of the relative flow velocity around the cylinder. The two van
der Pol type oscillators are then used to model fluctuating drag and lift coefficients. As the relative velocity
around the cylinder depends both on the fluid flow velocity and the velocity of the cylinder, the equations of
motions of the cylinder in cross-flow and in-line directions become coupled through the fluid forces. It is shown
that such approximation of the fluid forces allows to obtain the well known low dimensional models in the limit
case, and the model proposed by Facchinetti et al. [1] to describe the cross-flow vibrations is used as an
example. Existing experimental data and Computational Fluid Dynamics (CFD) results are used to calibrate the
proposed model and to verify the obtained predictions of complex fluid-structure interactions for different mass
ratios. The "super upper" branch phenomenon, exclusive for a two degrees-of-freedom motion at low mass
ratios, has been observed. The influence of the empirical parameters of the wake oscillators and fluid forces
coefficients on the dynamic responses is also discussed.

1. Introduction tical models. Analytical VIV models are represented by numerous

approaches to modelling both the structure and fluid, with some of
Slender marine structures such as risers, mooring cables, umbilicals them incorporating a van der Pol type equation as the governing
and tethers play crucial roles in global offshore exploration, installation equation for the fluid force acting on the structure. In present work, we
and production activities. As offshore oil and gas fields are moving into focus on this type of analytical models, known as wake oscillator
deeper waters, the nonlinearities in the system and the fluid-structure models.
interaction phenomena such as vortex induced vibrations (VIVs) The nature of the vortex shedding process behind cylindrical
become more and more important. Many of VIV aspects are far from structure suggests that the forces acting on the structure from the
being understood and advanced modelling is required to investigate the fluid can be modelled by a nonlinear, self-excited oscillator called the
impact of the phenomenon which significantly affects the service life of wake oscillator. This idea was first proposed by Bishop and Hassan [2]
marine structures. and then investigated by Skop and Griffin [3], and Blevins and Iwan
This work is motivated by the need of industry in effective toolkit [4]. In wake oscillator models the system is usually described by two
that would allow predicting loads and fatigue damage on riser systems, coupled ordinary differential equations. One of the equations is the
especially most common Top Tensioned Risers (TTRs) and Steel equation of transversal motion of the rigid cylinder. The second
Catenary Risers (SCRs), which represent a crucial part of offshore equation is a semi-empirical description of fluid: the nonlinear self-
facilities. Accurate prediction of VIVs can help to produce more robust excited fluid oscillator. A number of wake oscillator models were
structural design and lead to substantial savings in the offshore developed and applied to slender flexible structures undergoing VIV.
applications. The Balasubramanian and Skop model [5] proposed in 1997
The problem of vortex-induced vibrations is addressed by different included a van der Pol equation driven by the local transverse motion
approaches, which can be roughly categorized into one of the three of the cylinder as a governing equation for one component of the
major groups: experiments, computational fluid dynamics and analy- fluctuating lift force and a so-called stall term which is linearly

⁎
Corresponding author.
E-mail address: [email protected] (E. Pavlovskaia).

http://dx.doi.org/10.1016/j.ijmecsci.2016.05.019
Received 30 December 2015; Received in revised form 15 May 2016; Accepted 20 May 2016
Available online 24 May 2016
0020-7403/ © 2016 Elsevier Ltd. All rights reserved.
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 1. Wake oscillator model for the cylinder moving in transversal direction only (adopted from [1]).

proportional to the local transverse velocity of the cylinder and cylinder elastic support. These nonlinearities are included in the model
responsible for energy dissipation associated with motion of the by using so-called coupled Duffing oscillators, previously investigated
cylinder in the fluid. Because of the stall term they were able to capture by Raj and Rajasekar in 1997 [16]. Numerical results were tuned
an asymptotic, self-limiting structural response at zero structural against existing experimental data on two degrees-of-freedom cylinders
damping. The original form of van der Pol equation has also been obtained by several independent researchers. In particular, results used
reinterpreted by Krenk and Nielsen [6], Mureithi et al. [7], Plaschko for comparison were obtained by Jauvtis and Williamson in 2004 [17],
[8], and Skop and Luo [9] among others. who investigated the behaviour of low mass ratio cylinders and
A critical analysis on this class of low dimensional models in terms discovered the so-called “super-upper” branch of the response; by
of the fundamental behaviour has been done by Facchinetti et al. [1] by Stappenbelt et al. in 2007 [18], who examined behaviour of low mass
considering transverse VIV of a single degree-of-freedom structure in ratio cylinders at very low damping ratio ζ = 0.006 ; and by Blevins and
stationary uniform flow (see Fig. 1). Three different types of coupling Coughran in 2009 [19], who carried on investigations on the effect of
between structural and wake equations have been examined: via damping, varying ζ from 0.002 to 0.4 at a fixed mass ratio. The
displacement, velocity, and acceleration. It was shown that the accel- numerical results obtained by Srinil and Zanganeh demonstrated
eration coupling provides the best match to the experimental data. reasonable correspondence with experiments by relying on additional
The semi-empirical approach to vortex-induced vibrations has been structural nonlinearities instead of nonlinearities provided by the fluid.
further studied in the work by Keber and Wiercigroch [10]. The effect Another two degrees-of-freedom wake oscillator model called VIVTAS
of a weak structural nonlinearity on the dynamical behaviour of a was proposed by Furnes et al. [14] for a free span pipeline undergoing
vertical offshore riser undergoing VIVs was investigated, and the VIV. In their approach they used a complex, rather than a vector,
authors demonstrated that the structural nonlinearity has a stiffening representation of the total hydrodynamic force, with one equation of
effect on the oscillation of the riser, which becomes more pronounced motion for the structure and two nonlinear wake oscillators for real and
when the internal flow is incorporated into the model. It was shown imaginary parts of the fluid force. Preliminary model validation was
that additional nonlinearities in the structure affects the system performed using Marintek experimental data. In the work by Xu Bai
significantly, and thus, for the further insight, it is important to and Wei Qin in 2014 [20], another wake oscillator model was proposed
investigate new approaches to modelling the fluid. for two degrees-of-freedom VIV of elastically supported cylinder where
The new type of wake oscillator based on the Facchinetti model [1] a displacement variable related to the vortex strength was introduced.
with frequency dependent coupling was proposed by Ognik and Two-dimensional potential flow approach was used in the study, with
Metrikine [11], where an attempt has been made to introduce a wake fluctuating fluid forces acting on the cylinder simplified and quantified.
oscillator model that conforms to both the free and forced vibration Despite of the recent development of these semi-empirical models
experiments. Frequency dependent coupling in this case allows repro- focused on 2DOF motion, more work is still required to be done in the
duction of the measured frequency dependence of the fluid force on the area. Many of the empirical parameters used in the existing models are
cylinder. For the time domain representation of such coupling a evaluated using experimental studies on the stationary cylinders, where
convolutional integral is used. the fluid force acting on the cylinder is estimated using the velocity of
All models described so far focussed on the transversal vibrations of the fluid flow. When the cylinder is elastically supported and is allowed
the structure. However the experimental investigations show that the to move in the fluid flow, the force acting on the cylinder will depend on
in-line vibrations also play significant role (especially at the low mass the relative velocity of the flow around the cylinder as suggested by [21]
ratios), and some attempts to develop 2DOF wake oscillator model and later by [11]. However, none of the existing 2DOF wake oscillator
were made, where two wake equations are used for in-line and cross- models incorporates this dependence and therefore in this paper we
flow directions, see, for example, [12–14]. Ge et al. [12] proposed a two propose a new wake oscillator model which will fill this gap.
degrees-of-freedom wake oscillator model based on the approach by The rest of the paper is organised as follows. In the next section the
Wang et al. [15], which comes from considerations that the lift and new wake oscillator model is introduced and equations of motion are
drag forces do not coincide with X and Y axes (shown in Fig. 1), but act developed. The fluid forces are calculated using the relative velocity of
under a particular angle that indicates the direction of the cylinder's the flow around the cylinder and oscillating lift and drag coefficients
instantaneous velocity. Thus, a simple force decomposition can be described by two van der Pol equations. Here it is shown how the
performed and the forces acting in X and Y directions estimated. Srinil damping associated with the cylinder motion in the fluid is derived
and Zanganeh [13] adapted the same principle in evaluating the fluid from the suggested form of the fluid forces and how the proposed
forces but also introduced additional geometrical nonlinearities to the model could be reduced to the known model by Facchinetti [1] in the

177
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

limit case where only transversal vibrations of the cylinder are formulation, the rotation tensor can be written as
considered. In the following section, the model is calibrated using first ⎛ π →⎞ ⎛ π ⎞⎛ → →⎞ ⎛ π ⎞→ → →
the published experimental data [18,19,22] and then with our R⎜ , k ⎟ = cos⎜ ⎟⎜I − k ⊗ k ⎟ + sin⎜ ⎟ k × I + k ⊗ k =
⎝2 ⎠ ⎝ 2 ⎠⎝ ⎠ ⎝2⎠
Computational Fluid Dynamics (CFD) results. The CFD model used
→ → → → → →
for calculations is introduced and the results for transversal vibrations =j ⊗ i − i ⊗j + k ⊗ k, (5)
and combined transversal and in-line vibrations are used for compar- → → → → → →
ison with the proposed wake oscillator model predictions. Special where I = i ⊗ i + j ⊗ j + k ⊗ k is the unit tensor. Substituting
→ → →
attention is paid to the empirical coefficients of the wake oscillators UR = (U − x )̇ i − yj̇ into Eq. (4), we obtain
equations and their influence on the system response is investigated. 2 → →
Finally, some concluding remarks are given. ⎯→
⎯ 1 → y ̇ i + (U − x )̇ j
FL = ρf CLD UR → .
2 UR (6)
2. Two degrees-of-freedom wake oscillator model
Here the other parameters are ρf, the fluid density, CL, lift coefficient,
and CD, total drag, which can be represented as a sum of constant mean
In this work we consider an elastically supported cylinder experi-
sectional drag CD0 and fluctuating drag, CDfl
encing VIV, that is free to vibrate in cross-flow and in-line directions.
As mentioned in [23], for a cylinder capable of oscillating in both CD = CD0 + CDfl . (7)
directions, the equations of motion on an XY plane in terms of the
displacements in in-line and cross-flow directions, x and y, are To determine the values of total hydrodynamic force components
⎯→
⎯ ⎯→
⎯ ⎯→
⎯
FX and FY, the sum of lift and total drag forces F = FL + FD should be
m⋆x¨ + rsx ̇ + hx = FX , (1) projected on the appropriate axis (see Fig. 2) and therefore we have
m⋆y¨ + rsy ̇ + hy = FY , (2) ⎯→
⎯ ⎯→
⎯ →
FX = ( FL + FD )· i , (8)
where the total hydrodynamic force components in X and Y directions ⎯→
⎯ ⎯→
⎯ →
are FX and FY. Here m⋆ is mass per unit length including an added mass FY = ( FL + FD )· j , (9)
1 ⎯→⎯ ⎯→
⎯
m⋆ = ms + 4 πCM ρf D 2 , rs is structural damping, and h is stiffness of the where FL and FD are given by Eqs. (3) and (6), and thus the projections
support. ⎯→
⎯
⎯→
⎯ → → of total hydrodynamic force F are
This total hydrodynamic force, F = FX i + FY j , is the result of the
⎯→
⎯ ⎯→
⎯ ⎛1 → ⎛ → →⎞ 1 → ⎛ → → ⎞⎞ →
actions of the sectional vortex-induced drag FD and lift FL forces which FX = ⎜ ρf CLD UR ⎜y ̇ i + (U − x )̇ j +⎟ ρf CDD UR ⎜(U − x )̇ i − yj̇ ⎟⎟ · i
⎯→
⎯ ⎝2 ⎝ ⎠ 2 ⎝ ⎠⎠
are shown in Fig. 2. As can be seen from this figure, the drag force FD is
→ → → ⎛ ⎞
acting along the velocity, UR = U − V which is the fluid velocity relative 1 → 1 →
→ → → = ρf CLD UR y ̇ + ρf CDD UR ⎜U − x ⎟̇ ,
to the cylinder [21] (V is velocity of the cylinder and U = U i is the 2 2 ⎝ ⎠ (10)
⎯→
⎯
velocity of the flow). The lift force FL is then acting in the perpendicular
⎛1 → ⎛ → →⎞ 1 → ⎛ → → ⎞⎞ →
directions and the magnitudes of lift and drag forces depends on the FY = ⎜ ρf CLD UR ⎜y ̇ i + (U − x )̇ j ⎟ + ρf CDD UR ⎜(U − x )̇ i yj̇ ⎟⎟ · j
→ ⎝2 ⎝ ⎠ 2 ⎝ ⎠⎠
magnitude of relative velocity UR as [24]
1 → ⎛ ⎞ 1 →
2 → = ρf CLD UR ⎜U − x ⎟̇ − ρf CDD UR y .̇
⎯→
⎯ 1 → UR
FD = ρf CDD UR → , 2 ⎝ ⎠ 2 (11)
2 UR (3)
In general, Eqs. (10)–(11) provide the values of the fluid forces acting
2 → on the cylinder in X and Y directions and these forces should be
⎯→
⎯ 1 → U
FL = ρf CLD UR R· →R . substituted in the Eqs. (1)–(2) and solved together with wake equations
2 UR (4) Eqs. (14)–(15). By combining Eqs. (1)–(2) with Eqs. (10)–(11) we
Here to determine the direction of the lift force, the rotation tensor arrive with the system of equations describing motion of the cylinder
π →
R = R( 2 , k ) is used, which rotates the relative velocity vector in 1 → 1 → ⎛ ⎞
→ m⋆x¨ + rsx ̇ + hx = ρ CLD UR y ̇ + ρf CDD UR ⎜U − x ⎟̇ ,
counterclockwise direction for 90° around axis k . Using Rodrigues 2 f 2 ⎝ ⎠ (12)

→ → ⎯→
⎯ →
Fig. 2. Forces on a vibrating cylinder with instantaneous velocity V in the flow of velocity U . Drag force FD acts in line with relative stream velocity UR .

178
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 3. (a) Wake lift and (b) drag time histories; (c) 2D cylinder trajectory for m* = 1.275 at Ured = 6.68. CD 0 = 1.2, CDfl0 = 0.2, CL 0 = 0.3, Ax = Ay = 12, εx = εy = 0.3.

1 → ⎛ ⎞ 1 → 2.1. Full two degrees-of-freedom wake oscillator model in a non-

m⋆y¨ + rsy ̇ + hy = ρ CLD UR ⎜U − x ⎟̇ − ρf CDD UR y ,̇
2 f ⎝ ⎠ 2 (13) dimensional form
→
with UR = (U − x )̇ 2 + y 2̇ . The proposed full two degrees-of-freedom wake oscillator model in
Following the approach employing nonlinear oscillator equations of non-dimensional form is described by the following four second order
the van der Pol type [2,3,24], the fluctuating lift CL and drag CDfl coupled nonlinear differential equations
coefficients could be modelled by two wake oscillators using q and w
variables (q = 2CL / CL 0 and w = 2CDfl / CDfl0 ) ⎛U ⎞2 ⎛ ⎛ ⎞
1 1
x∼″ + 2ζx∼′ + x∼ = 8π 2St 2 ⎜ red − x∼′⎟ + y∼′ 2 ⎜⎜ MLqy∼′ +⎜MD + MDflw⎟
⎝ 2π ⎠ ⎝2 ⎝ 2 ⎠
w¨ + 2εxΩF (w 2 − 1)ẇ + 4ΩF2w = Sx , (14)
⎛ Ured ⎞⎞
q¨ + εyΩF (q 2 − 1)q ̇ + ΩF2q = Sy, (15) ⎜ − x∼′⎟⎟⎟ ,
⎝ 2π ⎠⎠
(18)
where εx and εy are van der Pol parameters, ΩF = 2πSt (U / D ) is the
frequency of vortex shedding, St is the Strouhal number, and Sx and Sy ⎛
⎯→
⎯ → → ⎛U ⎞2 1 ⎛U ⎞
are components of total structural force FS = Sx i + Sy j coupling wake y∼″ + 2ζy∼′ + y∼ = 8π 2St 2 ⎜ red − x∼′⎟ + y∼′ 2 ⎜⎜ MLq⎜ red − x∼′⎟
⎝ 2π ⎠ ⎝2 ⎝ 2π ⎠
equations with equations of cylinder motions, and CL0 and CDfl0 are lift
and fluctuating drag coefficients on a stationary cylinder. For in-line ⎛ ⎞ ⎞
1
vibrations, the frequency doubling is introduced to reflect an experi- −⎜MD + MDflw⎟y∼′⎟⎟ ,
⎝ 2 ⎠ ⎠
mentally observed phenomenon often mentioned in the literature, e.g. (19)
[25]. In next sections of this work this phenomenon will be discussed
further using the CFD approach. Here the acceleration coupling is w″ + 2εxΩ(w 2 − 1)w′ + 4Ω 2w = Ax x∼″, (20)
adopted as proposed by Facchinetti et al. [1], and therefore the Sx and
q″ + εyΩ(q 2 − 1)q′ + Ω 2q = Ay y∼″, (21)
Sy components are
Sx = (Ax / D )x¨, (16) where prime denotes differentiation with respect to non-dimensional
time τ and the following variables and system parameters are intro-
Sy = (Ay / D )y¨. (17) duced
As mentioned earlier, it was demonstrated by Facchinetti et al. [1] that τ = ωnt , x∼ = x / D, y∼ = y / D,
this type of coupling, in comparison to displacement coupling and
velocity coupling, provides best results when compared to experimental ωn = h / m⋆ , ζ = rs /(2ωnm⋆ ), Ω = ΩF / ωn , Ured = 2πU /(ωnD ),
data. Specifically, acceleration coupling captures lock-in domains at
low mass ratios with better accuracy than the other two types of ⎛ 1 ⎞
μ = ⎜ms + πCM ρf D 2⎟ / ρf D 2 , ML = CL 0 /16π 2St 2μ,
coupling. ⎝ 4 ⎠
The developed equations of motion describe the vibrations of the MD = CD0 /16π 2St 2μ, MDfl = CDfl0 /16π 2St 2μ,
cylinder in the fluid flow. However, a careful calibration of the model is
required and specifically empirical wake oscillators equations para- where μ is mass ratio in Facchinetti notation [1]. The first two
meters Ax, Ay, εx and εy need to be found. In case of a single degree-of- equations describe the dynamics of the structure whilst the remaining
freedom system, numerical results by Facchinetti [1] were fitted against two mimic the forces acting from the fluid. The relation between mass
experimental data, with Ay and εy estimated as 12 and 0.3 respectively. ratios in Williamson [26] and Facchinetti [1] notations, m⁎ and μ
However, further investigation and calibration for 2DOF models are correspondingly, is described by
essential and will be discussed in the Section 3. It should be noted that
m⁎ = 4μ / π − CM . (22)
in general these coefficients may be a function of various parameters of
the system such as mass ratio, damping ratio, reduced velocity, added Figs. 3 and 4 present two examples of system response calculated
mass coefficient, Reynolds number, etc. using this model for Ured = 6.68 and Ured = 8.01, respectively. Parts (a)
Before the model calibration is discussed, the model equations will be and (b) show the lift and drag coefficient time histories and parts (c)
transformed into the non-dimensional form. Then the approximations of demonstrate 2D cylinder trajectories computed for m⁎ = 1.275. In this
the fluid forces for small cylinder velocity and limit case of the case the wake oscillators coefficients are chosen to be the same as in [1]
transversal vibrations only will be discussed in the following subsections. for 1D case, i.e. Ax = Ay = 12 and εx = εy = 0.3

179
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 4. (a) Wake lift and (b) drag time histories; (c) 2D cylinder trajectory for m* = 1.275 at Ured = 8.01. CD 0 = 1.2, CDfl0 = 0.2, CL 0 = 0.3, Ax = Ay = 12, εx = εy = 0.3.

2.2. Approximate two degrees-of-freedom wake oscillator model As can be seen from Eqs. (10) and (11) and also Eqs. (26) and (27),
components of the total hydrodynamic force on the cylinder, contain-
In order to compare the proposed model with the models available ing contributions from both the lift and drag, form the right hand side
in the literature we will use some commonly accepted assumptions. of the cylinder equations representing the motion in XY plane. By
Assuming that the horizontal and vertical velocities of the cylinder, ẋ omitting all terms with the order of ϵ3 in Eqs. (26) and (27), equations
→
and ẏ, are smaller than the magnitude of the flow velocity U = U [15], for this model can be written as
→
and introducing notation UR = UR , we can approximate the value of the 1 1
m⋆x¨ + rsx ̇ + hx = ρ DU2CD0 + ρf DU2CDfl0w − ρf DUCD0x ̇
relative velocity using the following expansion 2 f 4
⎛ ⎞ 1 1
⎛ ẋ ⎞
2 ⎛ y ̇ ⎞2 ẋ 1 y 2̇ ⎟ + ρf DUCL 0qy ̇ + ρf DCD0y 2̇
UR = (U − x )̇ 2 + y 2̇ = U ⎜1 − ⎟ + ⎜ ⎟ ≈ U ⎜⎜1 − + . 4 4
⎝ U⎠ ⎝U ⎠ ⎝ U 2 U ⎟⎠
2
1 1
+ ρf DCD0x 2̇ − ρf DUCDfl0wx ,̇
2 2
(23)
(28)
Substituting UR into the equations for FX and FY we obtain
1 1 1
⎛ m⋆y¨ + rsy ̇ + hy = ρf DU2CL 0q − ρf DUCD0y ̇ + ρf DCD0x ̇
1 xẏ ̇ ⎞ 1 ⎛ y 2̇ x 2̇ ⎞ 4 2 2
FX = ρf CLDU ⎜y ̇ − ⎟ + ρf CDDU ⎜U − 2x ̇ + + ⎟,
2 ⎝ U⎠ 2 ⎝ 2U U⎠ (24) 1 1
ẏ − ρf DUCDfl0wy ̇ − ρf DUCL 0qx ,̇
4 2
1 ⎛ y 2̇ x 2̇ ⎞ 1 ⎛ xẏ ̇ ⎞
FY = ρf CLDU ⎜U − 2x ̇ + + ⎟ − ρf CDDU ⎜y ̇ − ⎟. (29)
2 ⎝ 2U U⎠ 2 ⎝ U⎠ (25)
with wake equations as follows
fl
Now we can substitute x ̇ = ϵx ̇ , y ̇ = ϵy ̇, CL = = ϵCL , CDfl
where x ̇ ,
ϵCD ,
fl
w¨ + 2ΩF εx (w 2 − 1)ẇ + 4ΩF2w = (Ax / D )x¨, (30)
y ̇, CL , CD are no longer small, and rewrite expressions for forcing terms
as follows q¨ + εyΩF (q 2 − 1)q ̇ + ΩF2q = (Ay / D )y¨. (31)
1 ⎛ ⎛ fl ⎞ ⎛ 1 If the motion of the cylinder in the in-line direction is omitted
FX = ρf D⎜CD0U2 + ϵ⎜CD U2 − 2CD0Ux+̇⎟ ϵ2⎜CLUy ̇ + CD0y ̇2 + CD0x ̇2
2 ⎝ ⎝ ⎠ ⎝ 2 (ẋ = 0 ) and the fluctuating drag is ignored, keeping only terms of
⎞ ⎛ ⎞ ⎞ order of ϵ we can obtain one degree-of-freedom wake oscillator
fl 1 fl fl
− 2CD Ux+̇⎟ ϵ3⎜ CD y ̇2 + CD x ̇2 − CLx ẏ ̇⎟⎟ , equation which is widely used (see for example [1]) with FY including
⎠ ⎝2 ⎠⎠ (26) the wake force term and the stall term (damping associated with the
⎛ fluid motion)
1 fl
FY = ρ D⎜ϵ(CLU2 − CD0Uy ̇) + ϵ2(CD0x ẏ ̇ − CD Uy ̇ − 2CLUx ̇) 1 1
2 f ⎝ FY = ρf CLDU2 − ρf CD0DUy ,̇
2 2 (32)
⎛1 fl ⎞⎞
+ϵ3⎜ CLy ̇2 + CLx ̇2 + CD x ẏ ̇⎟⎟ . where CL = CL 0q /2 and CL0 is usually taken as 0.3. Then one degree-of-
⎝2 ⎠⎠ (27)
freedom wake oscillator model will look as follows
Analysing the obtained equations, we can see that the largest force
1 ⎛ 1 ⎞ 1 1
acting on the cylinder is the constant horizontal force 2 ρf DCD0U2 ⎜ms + πCM ρf D 2⎟y¨ + rsy ̇ + hy = ρf DU2CL 0q − ρf DUCD0y ,̇
⎝ 4 ⎠ 4 2 (33)
associated with the constant mean sectional component of the drag.
It is also clear from these equations that there are forces of the same
q¨ + εyΩF (q 2 − 1)q ̇ + ΩF2q = (Ay / D )y¨. (34)
magnitude of order of ϵ acting in both horizontal and vertical direction.
It can be noted that equations of motion contain a number of The analysis presented in this sub-section demonstrates that by using
1
nonlinear terms, and terms ρf DUCD0ẋ and 2 ρf DUCD0ẏ representing widely accepted assumption of the small velocity of the cylinder, the
damping from the fluid in in-line and cross-flow equations of motion proposed model could be reduced to the existing models in the limit
respectively. For the in-line equation of motion this damping coeffi- case of transversal vibrations only, naturally obtaining both well-
cient is of a double magnitude in comparison to the one in cross-flow known fluid force damping term (e.g. stall term) and the appropriate
equation, which is an interesting feature of the proposed model, and lift force term.
this effect has not been previously captured by existing two degrees-of- In the next section, the calibration of the model will be considered
freedom models [15,13]. and the responses obtained using the new model will be compared with

180
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 5. Resonance curves (dependence of response amplitude on the reduced velocity). New 2DOF wake oscillator model compared to experimental data by Stappenbelt et al. [18].
m* = 2.36 , ζ = 0.006 . Wake parameters are Ay=5, Ax=12, εx = 0.3. Reference fluid force parameters are CL 0 = 0.3, CDfl0 = 0.2 . (a) Influence of mean sectional drag CD 0 , εy = 0.008 ; (b)
Influence of wake parameter εy at CD 0 = 2.0 .

published experimental data and also with the results from the CFD Overall, tuning εy has given best results in fitting the response curve to
simulations. experimental data, although generally a combined Ay / εy parameter
must be considered.
The best fit at CD0 = 1.75 (Fig. 6a) was chosen to illustrate the in-
3. Model calibration using CFD and experimental results line amplitude predictions, which are presented in Fig. 6b. It can be
seen that in-line amplitudes are higher than those observed experi-
3.1. Calibration based on experimental data mentally. Our analysis indicates that adjusting parameters of the wake
oscillator equation in in-line direction Ax and εx does not have the same
First, the published experimental results are utilised to calibrate the effect on the in-line vibration amplitudes as parameters Ay and εy on
proposed wake oscillator model. Three sets of experimental data were the transversal amplitudes. It is clear that more in-depth investigation
considered, and comparisons are made for different mass-damping is required to achieve a better match between the model predictions
ratio parameters. As observed by Jauvtis and Williamson in experi- and the experimental data. However, even for this non-optimised
ments carried out in 2004 [17], an additional branch of response choice of parameters, the overestimation of the in-line amplitude could
appears for m⁎ < 6 . As low mass ratios are of a particular interest in two be acceptable if the model is to be used for the design calculations as
degrees-of-freedom case, two sets of data are chosen to specifically satisfying the design criteria in this case will improve overall safety
capture the “super-upper” branch, at m⁎ = 2.36 by Stappenbelt et al. factor.
[18] and m⁎ = 2.6 by Jauvtis and Williamson [17,22]. Fig. 7 presents comparisons with experiments by Jauvtis and
Fig. 5 presents the amplitude of the transversal vibrations as Williamson [22] at a slightly higher mass ratio m⁎ = 2.6 and lower
function of the flow velocity (Ured = 2πU /(Dωn )). The amplitude of damping ratio ζ = 0.0015. It can be noted here that the same
transversal vibration is calculated as a maximum value of cylinder recommendations apply for both amplified sectional drag CD0 and
displacement at a given value of Ured. Here the system responses are wake parameter εy since the presented results demonstrate similar
calculated using the proposed model for different values of the drag trends for both sets of experimental data. The main issue, common for
coefficient CD0 in part (a) and different values of wake oscillator both cases, is that the upper branch of the calculated response always
coefficients εy in part (b). In general, choosing the wake oscillator happens to be shifted to the left giving overestimation of the amplitude
parameters is a challenging task and in the future this could be done values at lower reduced velocities in the lock-in region.
using a comprehensive optimisation procedure. However, some pre- Finally, a comparison with experimental data by Blevins and
liminary tuning of the parameters can be done without it, and the Goughran [19] for mass ratio m⁎ = 5.4 is presented in Fig. 8. A notable
results are shown in Fig. 5 where amplitudes evaluated using the tuned difference from other presented experimental data is the absence of a
wake oscillator model are compared to experimental data by lower branch after an amplitude drop around Ured ≈ 8. The response
Stappenbelt et al. [18] (m⁎ = 2.36 , ζ = 0.006 ). In order to achieve a predicted by the wake oscillator model in this case shows a better
reasonable match, both the coupling coefficient Ay and the wake correspondence with experiments, although higher values of amplitude
oscillator coefficient εy were reduced to Ay=5 and εy = 0.008. CD0 is starting at Ured ≈ 3.5 were obtained indicating, similar to previously
normally taken as a constant value, estimated the same way as other presented results, that according to the model the entrance to the lock-
reference parameters, however this assumption neglects the fact that in region occurs at lower values of reduced velocity than in experi-
value of CD0 depends on the transverse amplitude of vibration. As can ments.
be seen in Fig. 5a, mean sectional drag affects the branch of response It has to be noted that experimental data presented here were
significantly. Recommended value for this parameter is CD0 = 2.0 if it is obtained using test rigs that inevitably differ from each other in their
not chosen to be modelled as a function of transversal amplitude for the technical characteristics and may include structural nonlinearities
sake of simplicity. It can be noted here, that alternation of this which are not incorporated in the proposed wake oscillator model.
parameter can vastly change the shape of response branch, specifically Further experimental studies would be useful for more refined model
to nullify the jump to the lower branch of response as illustrated in calibration to achieve a better match and to formulate clear recom-
Fig. 5a for the value CD0 = 2.5. mendations on the selection of the empirical wake oscillator coeffi-
Fig. 5b demonstrates the best fit at tuned van der Pol parameter εy. cients. However, from the conducted analysis, it could be concluded
A fairly good fit to experimental data is observed at 0.007 < εy < 0.009.

181
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 6. Resonance curves (dependence of response amplitude on the reduced velocity). New 2DOF wake oscillator model compared to experimental data by Stappenbelt et al. [18].
m* = 2.36 , ζ = 0.006 . Wake parameters are Ay=5, Ax=12, εy = 0.0008 , εx = 0.3. Reference fluid force parameters are CD 0 = 1.75, CL 0 = 0.3, CDfl0 = 0.2 . (a) Cross-flow response; (b) In-line
response.

that much lower values of the parameter εy should be selected 3.2.1. CFD model
(0.007 < εy < 0.009) than 0.3 value identified in the case of transversal In this section, we consider the behaviour of elastically supported
vibrations only [1], whereas for the lower mass ration coupling cylinder capable of moving in in-line and transversal directions (often
coefficient Ay should be reduced, but for the higher mass ratio, it referred in literature as XY motion).
should be increased in comparison with 12 value from [1]. CFD model [27] has been implemented in ANSYS Fluent 12.0.16
utilizing User Defined Functions (UDFs) in order to compute the
displacement of the cylinder on each time step based on the forces
obtained from the dynamic pressure. Relatively low Reynolds numbers
3.2. Calibration using CFD modelling (600–2000) were considered for the sake of simplicity. However, even
for these Reynolds numbers (Re > 300) the vortex street is turbulent,
As was mentioned earlier, experimental facilities are different from and a high quality grid is required for solution to converge.
each other and having access to these facilities is not always possible. To couple the motions of the cylinder and the fluid, the forces acting
In order to use the experimental data for the model calibration, it on the cylinder have been calculated by integrating the total wall
would be useful to refine the proposed generic model in order to pressure on the cylinder surface obtained from the CFD solver, and the
accommodate the specific features of the chosen experimental rig. In drag and lift coefficients have been obtained as non-dimensional
such cases, it is inevitable that the additional effects would complicate components of these forces.
the main phenomenon and might make it challenging to separate the The computational domain used in this study is shown in Fig. 9a,
core vortex induced vibrations from structural nonlinear vibrations of where a cylinder of a unit diameter was considered [27]. The domain
the rig, for example. Therefore, we will explore the CFD approach to consists of an upstream of 11.5 times the diameter to downstream of 20
calibrate the proposed 2DOF wake oscillator model. The appropriate times the diameter of the cylinder and 12.5 times on each cross stream
CFD model could be set up without extra structural complications and direction. The data grid contains 15380 quadrilateral cells and 15659
if it is properly verified, such a CFD model could be a very valuable tool nodes (Fig. 9a). A finer grid is created near the boundary layer around
for the wake oscillator model calibration.

Fig. 7. Resonance curves (dependence of response amplitude on the reduced velocity). New 2DOF wake oscillator model compared to experimental data by Jauvtis and Williamson [22].
m* = 2.6 , ζ = 0.0036 . Wake parameters are Ay=5, Ax=12, εx = 0.3. Reference fluid force parameters are CL 0 = 0.3, CDfl0 = 0.2 . (a) Influence of mean sectional drag CD 0 , εy = 0.008 ; (b)
Influence of wake parameter εy at CD 0 = 2.0 .

182
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 8. Resonance curves (dependence of response amplitude on the reduced velocity). New 2DOF wake oscillator model compared to experimental data by Blevins and Coughran [19].
m* = 5.4 , ζ = 0.002 . Wake parameters are Ay=24, Ax=12, εx = 0.3. Reference fluid force parameters are CL 0 = 0.3, CDfl0 = 0.2 . (a) Influence of mean sectional drag CD 0 , εy = 0.01; (b)
Influence of wake parameter εy at CD 0 = 2.0 .

Fig. 9. (a) Computational mesh showing a zoom-up of the area near a cylinder; (b) a single degree-of-freedom system [17] and (c) a two degrees-of-freedom system [17].

the cylinder and it gets coarser at far flow field, particularly outside of of moving in transverse direction only, and amplitudes of the oscilla-
the wake region. The inlet region is on the left side and the outlet is on tions were recorded under varying the flow velocity. In this case the
the right side of the grid and periodic boundary conditions were chosen equation of motion of the cylinder is coupled with the fluid simulations
for the upper and lower boundaries of the computational domain. The via lift coefficient which is calculated by CFD on each time step of the
PISO (Pressure Implicit solution by Split Operator method) pressure- simulation process. The equation of motion for the transversal
velocity coupling scheme [28] was used as a solution method to allow a displacement y is
larger time step size without compromising the stability of the solution.
Dynamic mesh operated by the UDF allows motion of the cylinder 1 ⎛⎞
msy¨ + rsy ̇ + hy = ρ DU2CY ⎜t ⎟ ,
in a two-dimensional plane. The spring-based smoothing method [28] 2 f ⎝⎠ (35)
is applied to all cells of the dynamic mesh. The UDF receives full
feedback from the fluid and controls the displacement of the cylinder where ms is mass per unit length, CY(t) is the non-dimensional lift
by incorporating the ‘Compute Force and Moment’ function into the coefficient obtained from CFD solver using the transversal (lift) fluid
equation of motion, providing lift and drag forces directly from the force as CY (t ) = 2FY (t )/(ρf DU2 ). A dot denotes the differentiation with
solver on each time step. The time step size was chosen with regards to respect to dimensional time t.
the Strouhal number appropriate for shedding frequency (defined as The second set of simulation was run for the cylinder moving in
the frequency of a complete vortex shedding cycle), stream velocity and both transversal and in-line directions. The equation of motion for the
diameter of the cylinder. To capture the vortex shedding correctly, at in-line displacement x is
least 50 time steps were performed in one shedding cycle which was
calculated as Tcycle = D / St U . At the start of the simulations, the cylinder 1 ⎛⎞
msx¨ + rsx ̇ + hx = ρ DU2CX ⎜t ⎟ ,
is at rest in its initial equilibrium position, and the initial conditions for 2 f ⎝⎠ (36)
cylinder's transverse and in-line displacements and velocities are
y(0) = 0 and ẏ(0) = 0 , and x(0) = 0 and ẋ(0) = 0 . where CX(t) is drag coefficient obtained from the CFD solver in a
The first set of simulations was carried out for the cylinder capable similar manner.

183
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 10. (a) Amplitude response at two different mass ratios m* = 1.275 and m* = 6.375 calculated using CFD, with a “super-upper” branch appearing at a lower mass ratio m* = 1.275,
[29]; (b) “Super-upper” branch obtained experimentally [26].

3.2.2. CFD results directly related to the vortex shedding frequency and that when the
In this sub-section, we present the comparison of the results for flow velocity generates the vortex shedding at the frequency close to the
vibrations of a single degree-of-freedom system (transversal vibration cylinder's natural frequency ωn, the shedding frequency locks onto it,
of the rigid cylinder) and a two degrees-of-freedom system (transversal and the lock-in resonant oscillations occur at or near the natural
and in-line vibrations of the cylinder). While the natural frequency of frequency of the structure and tend to have significantly greater
the structure is kept constant in order to fix structural properties, a amplitude. In the current study the entry to and exit from the lock-in
variation of reduced velocity, Ured = 2πU /(ωnD ), is undertaken by condition were considered in the range of reduced velocities from
altering the flow velocity, U. In this case Reynolds number is not a approximately 3 to 10 (or Reynolds numbers from 600 to 2000).
constant value and the solution quality is heavily dependent on and The system behaviour at various mass ratios was investigated, and
constrained by the quality of the grid. The calculations were performed particular attention was paid to low-mass ratio cylinders. Fig. 10a
for the damping ratio ζ = rs /(ms 2ωn )=0.01 and various mass ratios presents the amplitude of the transversal vibration as a function of
m⁎ = 4ms /(ρf πD 2 ). The results are presented for the displacement which reduced velocity for two different values of cylinder's mass ratio
is normalised with respect to the diameter of the cylinder, D. m⁎ = 1.275 and 6.375 where cylinder is allowed to move in both
It is well known that the frequency of the cylinder's vibrations is directions. As can be seen in the figure, lower mass ratio cylinder

Fig. 11. (a) Comparison of 1DOF and 2DOF amplitudes of non-dimensional cylinder displacement at various mass ratios [29]; (b) Trajectories of the cylinder motion on the phase plane
for m* = 1.275 for four different values of Ured = 5.34, 6.68, 8.01 and 10.69.

184
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 12. (a) Lift and drag history for m* = 1.275 at Ured = 6.68 (lock-in condition); the numbers indicate the time moments where the pressure distribution presented in Fig. 12b are
calculated; dash lines show the beginning and the end of one cycle; (b) Vorticity contours calculated at different times during one cycle. Grey arrows indicate the direction of vortex-
induced lift and drag forces.

demonstrates significantly larger peak value of dimensionless ampli- the in-line component of total pressure is switching its direction every
tude. In this case, the CFD simulations reveal that the drag force quarter of a vortex-shedding cycle as the flow around the cylinder
influences cross-flow amplitude, and it is observed that drag's impact at generates vortex-induced drag force along the stream (before passing
resonance becomes more significant as mass ratio of a cylinder the middle of the cylinder) and against the stream (after passing the
decreases. middle of the cylinder). The cross-flow component of total pressure
Experiments by Williamson and Jauvtis [17,26] demonstrated a always acts in one direction per half of a vortex-shedding cycle.
large increase of amplitude when mass ratio is decreased below m⁎ = 6 Components of the total hydrodynamic force acting on the cylinder,
as shown in Fig. 10b. A new branch of response appeared in case of two containing contributions from both the lift and drag, cause the cylinder
degrees-of-freedom motion of the structure, and it was defined as a motion in XY plane. This XY motion is schematically illustrated in
“super-upper” branch. As can be seen from Fig. 10, our computational Fig. 13. For each position of the cylinder in this figure, a corresponding
results agree well with the existence of this branch, although the peak snapshot of vorticity contours obtained with the CFD simulations of the
amplitude exhibits lower value than in physical experiments. cylinder vibrating under lock-in condition is shown. It is important to
Comparison between peak amplitudes of the single degrees-of- note here that ANSYS Fluent decomposes the total hydrodynamic force
⎯→
⎯ ⎯→
⎯
freedom cylinder and the two degrees-of-freedom cylinder at various into two projections ( FX and FY ) on X and Y axes, which do not
mass ratios is made and represented in Fig. 11a. It can be seen that, for ⎯→
⎯ ⎯→
⎯
normally coincide with the directions of lift and drag forces ( FL and FD )
the ratios below m⁎ < 3, the two degrees-of-freedom cylinder demon- of an oscillating cylinder (although they do for a fixed one). These
strates the presence of “supper-upper” branch. In this case in-line components CX(t) and CY(t) of dimensionless total force computed with
vibrations have significant impact on cross-flow amplitude, and this ⎯→
⎯
the CFD were used to determine the direction of FTotal , the total
effects becomes more pronounced as mass ratio decreases. To demon- hydrodynamic force vector, and to identify for each particular drag
strate this, the trajectories of the cylinder on XY plane are shown in ⎯→
⎯ ⎯→
⎯
force FD (see positions 1–4 in Fig. 13) the direction of lift force FL .
Fig. 11b for various values of the reduced velocity at m⁎ = 1.25. Here ⎯→
⎯
While the direction of FD coincides with the direction of the relative
two values Ured = 5.34 and 6.68 belong to “super-upper” branch and the → ⎯→⎯
other two values Ured = 8.01 and 10.69 are chosen from outside region. fluid flow UR , FL acts perpendicularly in one of two possible directions.
As can be seen in Fig. 11b, the trajectory of motion becomes It can be seen that the total hydrodynamic force vector always falls in I
asymmetric to the Y axis and the original “eight” shape is bent. The and IV quadrants. In this particular illustration, CY(t) is in phase with
increase of amplitude at illustrated reduced velocity values is signifi- dimensionless cylinder displacement y / D .
cant and estimates about 50–60%. At the same time, at higher values of
reduced velocities, outside of the lock-in region, the influence of in-line
vibrations becomes minimal as can be seen in Fig. 11b where the cross-
flow dimensionless amplitude is decreased to 0.6 in the presented case. 3.2.3. Model calibration using CFD results
In case of XY motion a frequency doubling phenomenon in in-line The CFD results are used in this sub-section to calibrate the
response is well-known and observed in experiments [25], and Fig. 12 proposed wake oscillator model. Specifically, empirical wake para-
demonstrates that is also captured properly by the CFD simulations. meters Ax, Ay, εx, εy are considered and determined using parametric
The time histories of non-dimensional drag and lift forces are shown in analysis. Also, fluid force (reference) parameters CL 0 , CD0 , CDfl0 , taken
Fig. 12a where a steady state response is presented for m⁎ = 1.275 and from a fixed cylinder lift and drag measurements, are discussed here.
Ured = 6.68. The vorticity contours shown in Fig. 12b are calculated at It is also important to note, that all presented comparisons are
the times marked by numbers 1 to 4 in Fig. 12a, and here the directions made using mass ratio notation as it appears in the work by Jauvtis and
of the fluctuating drag and the lift forces are given by the grey arrows. It Williamson [17], where m⁎ = 4ms / ρf πD 2 , i.e. the mass ratio does not
is observed that this frequency doubling phenomenon can be explained include added mass.
exclusively through the fluid motion and it is not a property of the For comparisons with data obtained from the CFD simulations, full
structure. As can be seen from Fig. 12b, during one period of motion wake oscillator model, introduced in Section 2 and described by Eqs.
(37)–(40), is used in its dimensional form

185
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 13. (Color online) Full vortex shedding cycle schematic. The snapshots of the vorticity contours of the cylinder vibrating in the lock-in conditions are presented at four different
moments of time during one cycle. In those moments of time the total force acting from the fluid on the cylinder calculated using CFD results for m* = 1.275 is shown by thick blue arrow.
⎯→
⎯ ⎯→
⎯ ⎯→
⎯
The components of this force FTotal acting along the relative fluid flow ( FD ) and in perpendicular directions ( FL ) are shown by red thick arrows. Velocities of the cylinder are shown by
→
black solid arrows and schematics to determine the relative flow velocity UR are also included.

186
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 14. Two degrees-of-freedom: CFD vs Wake Oscillator Model cross-flow response branches for m* = 1.275, ζ = 0.01. Wake parameters are Ay=3, Ax=12, εx = 0.3. Reference
parameters are CL 0 = 0.3, CDfl0 = 0.2 . (a) Effect of mean sectional drag CD 0 , εy = 0.0014 ; (b) Influence of wake parameter εy , CD 0 = 1.3.

⎛ 1 ⎞ ⎛1 Therefore, they should be adjusted first when the model calibration is

⎜ms + πCM ρf D 2⎟x¨ + rsx ̇ + hx = (U − x )̇ 2 + y 2̇ ⎜ ρf CL 0qDy ̇
⎝ 4 ⎠ ⎝4 performed.
Comparisons of the results of the simulations obtained using the
1 ⎛ 1 ⎞ ⎛ ⎞⎞
+ ρf ⎜CD0 + CDfl0w⎟D⎜U − x ⎟̇ ⎟ 2DOF wake oscillator model described by Eqs. (37)–(40) with the
2 ⎝ 2 ⎠ ⎝ ⎠⎠ (37) 2DOF CFD results are presented in Figs. 14–17. Here the reference
parameters for lift and drag coefficients, CL 0 and CDfl0 are kept at their
⎛ 1 ⎞ ⎛1 ⎛ ⎞
⎜ms + πCM ρf D 2⎟y¨ + rsy ̇ + hy = (U − x )̇ 2 + y 2̇ ⎜ ρf CL 0qD⎜U − x ⎟̇ initial values [1], with CL 0 = 0.3 and CDfl0 = 0.2 .
⎝ 4 ⎠ ⎝4 ⎝ ⎠ From Fig. 14a it can be seen that the variation of drag coefficient
1 ⎛ 1 ⎞ ⎞ CD0 affects the amplitude of response during lock-in and has almost no
− ρf ⎜CD0 + CDfl0w⎟Dy⎟̇ effect after desynchronization. At the same time, εy primarily affects
2 ⎝ 2 ⎠ ⎠ (38)
exit from lock-in when desynchronization takes place as can be seen in
with wake equations Fig. 14b. Increasing CD0 leads to lower amplitudes in the lock-in region
and higher jump, while reducing εy leads to wider lock-in region, but
w¨ + 2ΩF εx (w 2 − 1)ẇ + 4ΩF2w = (Ax / D )x¨ (39)
has little to no influence on the height of the jump. Ax and εx were kept
at Ax=12 and εx = 0.3 as appear in [1].
q¨ + εyΩF (q 2 − 1)q ̇ + ΩF2q = (Ay / D )y¨ (40)
In Fig. 15a the case where parameters Ay, εy are tuned to
In general, tuning reference fluid parameters may be hard to justify, demonstrate the best fit with the CFD results for m⁎ = 1.275 is
and would require a clear explanation since CDst0 , CDfl0 and CL 0 represent presented. Although a fairly good correspondence with the CFD results
experimentally observed force components acting on a fixed cylinder: is shown in Fig. 15a for this mass ratio m⁎ = 1.275, with the same values
non-amplified mean sectional drag, fluctuating drag coefficient and of reference and wake parameters used for mass ratio m⁎ = 6.375 in
fluctuating lift coefficient respectively. All these parameters are usually Fig. 15b the model fails to predict the amplitudes of the response.
taken from literature and modifying their values could be questionable. Better results can be achieved when these parameters are treated as a
On the other hand, the empirical wake parameters do not have clear function of mass ratio, or, more generally, as a function of mass-
physical meaning and they do not necessarily have to be fixed and/or damping ratio, instead of being kept constant as shown in Fig. 16.
remain the same for different values of mass or damping of the system. These conclusions have also been drawn by other researches (see [13]).

Fig. 15. Two degrees-of-freedom: CFD and wake oscillator model cross-flow response branches for (a) m* = 1.275 and (b) m* = 6.375. Parameters
CD 0 = 1.3, CDfl0 = 0.2, CL 0 = 0.3, Ay = 3, Ax = 12, εy = 0.0014, εx = 0.3 are fixed for both mass ratios. Damping ratio ζ = 0.01.

187
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 16. Two degrees-of-freedom: CFD and wake oscillator model cross-flow response branches at initial reference parameters [1] CDfl0 = 0.2, CL 0 = 0.3. (a)
m* = 1.275, CD 0 = 1.3, Ay = 3, Ax = 12, εy = 0.0014, εx = 0.3; (b) m* = 6.375, CD 0 = 2.0, Ay = 26, Ax = 12, εy = 0.05, εx = 0.3. Damping ratio ζ = 0.01.

In Fig. 16 two sets of tuned wake parameters at different values of Finally, the in-line response predictions are shown in Fig. 18. All
amplified drag CD0 demonstrate a good match with the CFD data for parameters are the same as for the previous cross-flow responses, and
two different mass ratios. It is clear that a more detailed study should are also presented for m⁎ = 1.275 (Fig. 18a) and m⁎ = 6.375 (Fig. 18b).
primarily focus on the combined Ay / εy and Ax / εx parameters as a As was mentioned before, it has been found that wake parameters Ax
function of the mass-damping parameter (m⁎ + CM )ζ . As can be seen in and εx have much smaller effect on in-line amplitude predictions than
Fig. 16a, tuned wake parameters for m⁎ = 1.275 are Ay=3 and Ay and εy have on cross-flow amplitudes. As can be seen from both
εy = 0.0014 , which give us Ay / εy = 2143. For a higher value of mass figures, predictions by wake oscillator model do not capture the in-line
ratio m⁎ = 6.375, as illustrated in Fig. 16b, a good match with CFD response amplitudes as good as they do the cross-flow ones. The values
results is achieved at Ay=26 and εy = 0.05, with Ay / εy = 520 . However, of the in-line amplitudes, are, however, significantly lower than those
at the higher mass ratio m⁎ = 6.375 the wake oscillator model does not for cross-flow motion, and often negligible, especially for higher mass
capture a jump from the upper to the lower branch of response. ratios. It should be noted here that overestimation of the in-line
Different tuning is shown in Fig. 17 with the same value of amplitudes would be acceptable from the design point of view as the
amplified drag CD0 = 1.3 for both masses. As can be seen in Fig. 17a, safe solution will be developed based on the prediction of this model.
wake parameters are the same as they were in the previous figure for
m⁎ = 1.275 with Ay=3 and εy = 0.0014 . For a higher value of mass ratio
m⁎ = 6.375 shown in Fig. 17b, wake parameters are changed to Ay=24 4. Concluding remarks
and εy = 0.2 , with Ay / εy = 120 . Although amplified drag would generally
be different for different masses as it depends on the value of In this work, vortex-induced vibrations of elastically supported
amplitude, for simplicity, when it is taken as constant, it is taken so cylinders capable of moving in cross-flow and in-line directions were
for all values of mass ratio. In can also be noted that for this set of investigated. A new wake oscillator model for two degrees-of-freedom
parameters the jump from the upper to the lower branch is now vortex-induced vibrations has been proposed, where vortex-induced lift
captured for m⁎ = 6.375 (Fig. 17b), although amplitude values predicted and drag forces were modelled with two nonlinear self-excited oscilla-
by wake oscillator model after the jump are less than those predicted by tors of van der Pol type.
CFD. Total hydrodynamic force was represented as a sum of lift and drag
forces acting on a vibrating cylinder in a uniform flow. The lift and drag

Fig. 17. Two degrees-of-freedom: CFD and wake oscillator model cross-flow response branches at initial reference parameters [1] CDfl0 = 0.2, CL 0 = 0.3. (a)
m* = 1.275, CD 0 = 1.3, Ay = 3, Ax = 12, εy = 0.0014, εx = 0.3; (b) m* = 6.375, CD 0 = 1.3, Ay = 24, Ax = 12, εy = 0.2, εx = 0.3. Damping ratio ζ = 0.01.

188
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

Fig. 18. Two degrees-of-freedom: CFD and wake oscillator model in-line response branches at initial reference parameters [1] CDfl0 = 0.2, CL 0 = 0.3. (a)
m* = 1.275, CD 0 = 1.3, Ay = 3, Ax = 12, εy = 0.0014, εx = 0.3; (b) m* = 6.375, CD 0 = 1.3, Ay = 24, Ax = 12, εy = 0.2, εx = 0.3. Damping ratio ζ = 0.01.

forces were defined as being proportional to the square of the to be done and is also planned for the future work. Instead of taking it
magnitude of the relative flow velocity around the cylinder with the as constant, amplified drag can be represented in the form
drag force acting in the direction of relative velocity and lift force acting (1 + 2Ay / D )CDst0 [30], where Ay is the cross-flow amplitude of vibration
perpendicularly. Equations of motion of the cylinder in cross-flow and and CDst0 is the drag coefficient for a stationary cylinder in the subcritical
in-line directions are coupled through the fluid forces calculated from range of Reynolds numbers. This will introduce an additional non-
the instantaneous relative flow velocity around the cylinder that linearity into structural equations and contribute to refinement of the
depends both on the fluid flow velocity and the instantaneous velocity proposed model.
of the cylinder; the type of coupling for a two degrees-of-freedom wake
oscillator model that has not been done before. An acceleration Acknowledgements
coupling was used based on the recommendations by Facchinetti
et al. [1]. This work is supported by the National Subsea Research Institute
It was demonstrated that the approximation of the fluid forces (NSRI) UK.
allows to obtain the well-known low-dimensional models in the limit
case, and the model proposed by Facchinetti et al. [1] was used an References
example of the model reduction to a single degree-of-freedom system.
A two-dimensional CFD model developed as a part of this work was [1] Facchinetti ML, de Langre E, Biolley F. Coupling of structure and wake oscillators
used together with the published experimental data to calibrate the in vortex-induced vibrations. J Fluids Struct 2004;19:123–40.
[2] Bishop RED, Hassan AY. The lift and drag fores on a circular cylinder in a flowing
wake oscillator model parameters. field. Proc Roy Soc Ser A 1964;277:51–75.
It was understood that to achieve a reasonable match with the [3] Skop RA, Griffin OM. A model for the vortex-excited resonant vibrations of bluff
experiments much lower values of the parameter εy should be selected bodies. J Sound Vib 1973;27:225–33.
[4] Iwan WD, Blevins RD. A model for vortex-induced oscillation of structures. J Appl
(for the proposed model, as presented in this work, 0.007 < εy < 0.009) Mech 1974;41:581–6.
than 0.3 value identified in the case of transversal vibrations only [1], [5] Balasubramanian S, Skop RA. A new twist on an old model for vortex-excited
whereas for the lower mass ratios coupling coefficient Ay should be reduced, vibrations. J Fluids Struct 1997;11:395–412.
[6] Krenk S, Nielsen SRK. Energy balanced double oscillator model for vortex-induced
and for the higher mass ratios, it should be increased in comparison with 12 vibrations. ASCE J Eng Mech 1999;125:263–71.
value from [1]. It was also found that the variation of drag coefficient CD0 [7] Mureithi NW, Kanki H, Nakamura T. Bifurcation and perturbation analysis of some
affects the amplitude of response during lock-in and have almost no effect vortex shedding models. In: Proceedings of the seventh international conference on
flow-induced vibrations. Luzern, Switzerland. Balkema, Rotterdam; 2000. p. 61–8.
on the predicted length of the lock-in region, whereas εy primarily affects
[8] Plaschko P. Global chaos in flow-induced oscillations of cylinders. J Fluids Struct
this length. Increasing CD0 leads to lower amplitudes in the lock-in region 2000;14:883–93.
and higher jump, while reducing εy leads to wider lock-in region, but has [9] Skop RA, Luo G. An inverse-direct method for predicting the vortex-induced
vibrations of cylinders in uniform and nonuniform flows. J Fluids Struct
little to no influence on the hight of the jump.
2001;15:867–84.
Our analysis has shown that the wake oscillators parameters [10] Keber M, Wiercigroch M. Dynamics of a vertical riser with weak structural
depend on the mass ratio as the best agreement with the experiments nonlinearity excited by wakes. J Sound Vib 2008;315:685–99.
is observed for different values of these parameters at low and high [11] Ogink RHM, Metrikine AV. A wake oscillator with frequency dependent coupling
for the modeling of vortex-induced vibration. J Sound Vib 2010;329:5452–73.
mass ratios. In the considered cases for the chosen system parameters, [12] Ge F, Long X, Wang L, Hong YSh. Flow-induced vibrations of long circular
the model slightly overestimated the in-line amplitudes of the vibra- cylinders modeled by coupled nonlinear oscillators. Sci China Ser G: Phys Mech
tions and also amplitudes of the transversal vibrations in the beginning Astron 2009;52(7):1086–93.
[13] Srinil N, Zanganeh H. Modelling of coupled cross-flow/in-line vortex-induced
of the lock-in region. It was observed that adjusting parameters of the vibrations using double duffing and Van Der pol oscillators. Ocean Eng
wake oscillator equation in the in-line direction Ax and εx does not have 2012;53:83–97.
the same effect on the in-line vibration amplitudes as parameters Ay [14] Furnes GK, Sorensen K. Flow induced vibrations modeled by coupled non-linear
oscillators. In: Proceedings of the 17th international offshore and polar engineering
and εy on the transversal amplitudes. conference. Lisbon, Portugal; 2007. p. 2781–7.
It should be noted here that although general recommendations on [15] Wang XQ, So RMC, Chan KT. A nonlinear fluid force model for vortex-induced
the choice of empirical wake parameters can be given based on this vibration of an elastic cylinder. J Sound Vib 2003;260:287–305.
[16] Raj SP, Rajasekar S. Migration control in two coupled Duffing oscillators. Phys Rev
study, in the future a development of some type of optimisation
E 1997;55(5):6237–40.
procedure for in-depth model calibration would be beneficial. [17] Williamson CHK, Jauvtis N. A hight-amplitude 2T mode of vortex-induced
A more accurate estimation of the amplified drag coefficient CD0 has vibration for a light body in XY motion. Eur J Mech B/Fluids 2004;23:107–14.

189
A. Postnikov et al. International Journal of Mechanical Sciences 127 (2017) 176–190

[18] Stappenbelt B, Lalji F, Tan G. Low mass ratio vortex-induced motion. The 16th [24] Hartlen RT, Currie IG. Lift-oscillator model of vortex-induced vibration. J Eng
Australian fluid mechanics conference. Gold Coast, Australia; 2007. p. 1491–7. Mech Division 1970;EM5:577–91.
[19] Blevins RD, Coughran CS. Experimental investigation of vortex-induced vibration [25] Vandiver JK, Jong J-Y. The relationship between in-line and cross-flow vortex-
in one and two dimensions with variable mass, damping, and Reynolds number. J induced vibration of cylinders. J Fluids Struct 1987:381–99.
Fluids Eng 2009;131(10):101202–7. [26] Williamson CHK, Govardhan R. Vortex-induced vibrations. Ann Rev Fluid Mech
[20] Bai Xu, Qin Wei. Using vortex strength wake oscillator in modelling of vortex 2004.
induced vibrations in two degrees of freedom. Eur J Mech B/Fluids [27] Postnikov A. Wake oscillator and CFD in modelling of VIVs. [Ph.D. thesis].
2014;48:165–73. Aberdeen: School of Engineering, University of Aberdeen; 2016.
[21] Vandiver JK. Dimensionless parameters important to the prediction of vortex- [28] ANSYS FLUENT 12.0 “Theory Guide”, ANSYS Inc. April 2009.
induced vibration of long, flexible cylinders in ocean currents. J Fluids Struct [29] Postnikov A, Pavlovskaia EE, Wiercigroch M. Modelling of vortex-induced vibra-
1993;7:423–55. tions of slender marine structures. In: Meskell C, Bennett G, editors. Proceedings of
[22] Williamson CHK, Jauvtis N. The effect of two degrees of freedom on vortex-induced the 10th international conference on flow-induced vibration & flow-induced noise
vibration at low mass and damping. J Fluid Mech 2004;509:23–62. (FIV2012). Flow-induced vibration; 2012. p. 569–76.
[23] Pavlovskaia E, Keber M, Postnikov A, Reddington K, Wiercigroch M. Multi-mode [30] Pantazopoulos MS. Vortex-induced vibration parameters: a critical review. OMAE -
modelling of vortex-induced vibration. Int J Nonlinear Mech 2016;80:40–51. Volume I, Offshore Technology, ASME; 1994. p. 199–255.

190