Computational and Experimental Assessment of Turbulence Stimulation on Flow Induced Motion of a Circular Cylinder

In the MRELab of the University of Michigan, flow-induced motion (FIM) is studied as a means to convert marine hydrokinetic (MHK) energy to electricity using the vortex-induced vibrations for aquatic clean energy (VIVACE) energy harvester [1–5]. Turbulence stimulation by selectively distributed surface roughness, in the form of sand-strips, referred to as PTC, was added to oscillating cylinders in 2008 [5]. PTC enabled VIVACE to harness hydrokinetic energy from currents/tides over the entire range of FIM including VIV and galloping. In 2011, the MRELab produced experimentally the PTC-to-FIM Map defining the induced cylinder motion based on the location of PTC [6]. In 2013, the robustness of the map was tested and dominant zones were identified [7]. Even though the PTC-to-FIM Map has become a powerful tool in inducing specific motions of circular cylinders, several parameters remain unexplored. Experiments, albeit being the ultimate verification tool, are time consuming and hard to provide all needed information. A computational tool that could predict the FIM of a cylinder correctly would be invaluable to study the full parametric design space. A major side benefit of PTC was the fact that PTC enabled CFD simulations to generate results in good agreement with experiments by forcing the location of the separation point [8]. This valuable tool, along with experiments, is used in this paper to investigate PTC design parameters such as width and thickness and their impact on flow features with the intent of maximizing FIM and thus, hydrokinetic energy conversion.

Passive turbulence control is a means to altering flow kinematics around a cylinder by covering parts of the cylinder with strips that may be smooth or have surface roughness. Roughness was added to the whole cylinder surface before Bernitsas and Raghavan used selectively distributed surface roughness to enhance or minimize VIVs, depending on the location of the strips [5]. The strips change the boundary layer separation point and the flow separation characteristics may be controlled passively [9,10]. Chang et al. [11] and Park et al. [6] built the PTC-to-FIM Map, which defines zones of placement of PTC strips for enhancing or minimizing FIM. General information on VIV is available in review papers [12–14] and on galloping in Ref. [15]. VIVACE converts MHK energy to mechanical in the oscillations of the cylinder in FIM and subsequently converts it to electricity [3]. It uses VIV and galloping with a single cylinder. Those motions are further enhanced by properly implementing gap flow with multiple cylinders [16]. The higher and faster the FIM oscillations, the more power the cylinder will convert from hydrokinetic to mechanical. Therefore, from the point of view of MHK energy harnessing, vibration amplitude and frequency should be enhanced. When PTC is applied to the upstream part of the cylinder as established in the PTC-to-FIM Map, PTC will induce galloping. PTC induces galloping instability, which has been reported to result in oscillation amplitude up to 3 to 4 times the diameter of the cylinder [6,7].

This paper studies the effect of PTC particulars on the flow-induced motions of a circular cylinder, computationally and experimentally. The physical model is described in Sec. 2. The computational approach is explained in Sec. 3 and validated in Secs. 4 and 5. The PTC thickness is studied in Sec. 6, and the PTC width in Sec. 7. CFD complements experimental results due to the limited possibilities of visualizing the flow at these high Reynolds numbers (varying between 30,000 < Re < 120,000). Flow visualization is achieved by CFD as presented in Sec. 8. This facilitates discussion of the physics of the flow in VIV and galloping. Conclusions are provided in Sec. 9.

The cylinder properties used in this study both in the computational and experimental studies are provided in Table 1. Mass ratio, m*, is the oscillating mass over mass displaced mass.

The forward (upstream) stagnation point is referred to as the 0° point. Positioning of the PTC on the cylinder is measured from the forward stagnation point. An example of PTC positioning on the cylinder at $30−46deg$ is shown in Fig. 1.

The roughness strips are placed on the cylinder symmetrically with respect to the flow as shown in Fig. 1. The system has one degree-of-freedom. The in-flow movement is constrained and thus, the cylinder is allowed to move only in the direction perpendicular to the flow.

The PTC used in experiments and CFD in the MRELab have proven to result in agreement between CFD and experiments for the following reason [8,9]. Alike, two-dimenisional (2D)-Unsteady Reynolds-Averaged Navier-Stokes (URANS), three-dimensional (3D)-URANS, large eddy simulation, detached eddy simulation codes fail to simulate cylinder flows at Reynolds numbers greater than 10,000 [17]. The source of the error is the poor prediction of the motion of the flow separation point [8]. PTC in the form of straight sandpaper forces separation at the leading edge of the sandpaper strip. Actually, for Re > 10,000, the separation point in CFD simulations is stuck at 90 deg instead of oscillating around 81 deg for subcritical flow. This failure in 2D and 3D codes alike makes simulations useless even as complements to experimental results. PTC is the only way to achieve good agreement between experiments and CFD simulations. In the experiments, in addition to the PTC, models reach within 2 cm of the sidewalls of the recirculating channel, thus, limiting tip-flow effects and keeping the cylinder outside the boundary layer of the flow [18]. So, the 3D nature of flow kinematics in VIV is limited in the experiments. The effect of tip flow and its interference with the 3D nature of flow in VIV are investigated in detail in Ref. [19].

Flow analysis with CFD is performed using the 2D-URANS approach implementing the $k−ω$ shear stress transport (SST) turbulence model. $k−ωSST$ is a powerful model which can be used for flows around bluff bodies. $k−ωSST$ also takes advantage of the strengths of the $k−ε$ turbulence model in the outer flow region [20].

To avoid extended computational periods, the problem is simplified by using a 2D approach. Two-dimensional flow assumption to solve for VIV is a common practical approach implemented by many researchers in the field. A broad literature review on computational approach is provided by Sarpkaya [13]. In two dimensions, the oscillating mass is nondimensionalized by the total length of the cylinder. Due to the 2D structure of the computational approach, some flow characteristics like cross-flows, tip vortices, and cellular shedding are neglected. The grid structure used in this study is shown in Fig. 2. A close-up view of the grid structure in the vicinity of the cylinder is given in Fig. 3, where the roughness strips can also be seen located symmetrically with respect to the flow. Quadrilateral elements are used near the cylinder and its wake, while tetrahedral elements are used in the outer region. The inner domain, which has quadrilateral elements, moves together with the cylinder and is not deforming. The tetrahedral elements in the outer domain are deformed and remeshed in response to cylinder movement. The boundaries of the whole domain are constructed far away so that they do not have any effect on the cylinder's VIV and galloping response. The error limit is set to be on the order of $10−5$ and numerical results converge well.

where $Cmax$ is selected to be equal to one. Here, $uxi$⁠, $Δxi$⁠, and $Δt$ denote the velocity, height of the grid elements on the cylinder, and time-step size, respectively. $i$ denotes the degree-of-freedom of the system, and n is the number of degrees-of-freedom, which for the current case is n = 1.

The roughness strips used in experiments in this study are composed of backing paper with grits of different size on it. In CFD, however, thickness is modeled as a step on the cylinder surface with total thickness equal to the sum of the backing paper (h) and the average grit height (⁠$k$⁠) to make the complexity of grid generation manageable.

About 12,000 elements are used in the whole fluid domain at the start of the simulation. The number of elements, however, is subject to change as new elements may be formed during automatic remeshing. Depending on the response amplitude of the oscillating cylinder, the number of elements may reach up to 16,000. The computational method reaches steady-state at around 10 s of real-time simulation. So unlike experiments, the maximum point that the cylinder reaches at each oscillation does not change. Once URANS reaches a steady-state and finds its stable amplitude, the cylinder oscillates at that maximum amplitude for the duration of the simulation.

A real-time simulation of 10 s takes around 1.5 h of computation time with an Intel Xeon central processing unit E5-2630 at 2.30 GHz. The workstation that is used for these simulations has a 64-bit operating system with 64 GB of installed memory (RAM). The code uses 12 cores with two processors in each core.

The computational approach is first validated with the benchmark experiments of Khalak and Williamson for a smooth cylinder at low Reynolds numbers [21]. These measurements were taken at Re = 3800 and are valid till Re of about 10,000. These numbers are in the TrSL2 flow regime. The computational and experimental amplitude responses of the cylinder in the TrSL2 regime are given in Fig. 4. It may be said that the computational approach is successful to accurately predict pre-VIV range and initial branch while it can only partially capture the upper branch, lower branch, and desynchronization.

Computationally, capturing the upper branch is a challenging issue as many researchers in the field have tried different methods or approaches to accurately calculate the range of synchronization and the maximum amplitude achieved. One could find many computational results implementing Reynolds-Averaged Navier–Stokes Equations (RANSE) methods, large eddy simulation, or even direct numerical simulation, but none of these methods are completely accurate or computationally economical. Other computational results implementing RANSE methods are summarized in Refs. [8], [22], and [23]. A comparison of results implementing 2D-URANS is given in Ref. [24].

The effect of the PTC thickness is investigated for a fixed location of the two strips on the cylinder selected to cover the surface between $30 deg$ and $46 deg$ on both sides of the cylinder with respect to the forward point of symmetry. The properties of the roughness strips (commercial sandpaper) used are listed in Table 2.

Figure 5 shows the amplitude response for a smooth cylinder and for two single cylinders with PTC. One PTC-cylinder is appended with P60 commercial-grade sandpaper and the other with P180. Thickness and roughness are important parameters and are shown in Table 2 [6,7]. Figure 5 shows CFD simulations performed in this study for the two PTC-cylinders in comparison to experimental measurements. PTC forces the flow separation at leading edge of the sandpaper strip and, thus, alters the flow-induced motion of the cylinder. Figure 6 shows the corresponding data for the frequency response of the smooth and rough cylinders.

The form of the FIM response of the cylinders follows the typical response observed and measured consistently in the MRELab for $30,000<Re<120,000$⁠, which falls in the high-lift TrSL3 flow regime [2,4,11]. VIV response in TrSL3 differs dramatically from TrSL2 VIV response even though each branch is clearly identifiable in both flow regimes [14]. The major difference is that in TrSL3, the upper branch overtakes the lower branch and amplitude increases linearly rather than remaining constant [2,4]. On the basis of Figs. 5 and 6, we can make the following observations.

Pre-VIV range: For $U*<4,$ a smooth cylinder does not exhibit FIM. For PTC-cylinders, the pre-VIV range starts at $U*=4$⁠. Pre-VIV range occurs slightly earlier in CFD generated results and it starts at around $U*=3$ both with smooth and rough cylinders.

Initial branch: It is short characterized by a large jump in amplitude and frequency response. For a smooth cylinder, it is in the range $4<U*<5$⁠. For the PTC-cylinder, it is in the range $4<U*<4.5$ experimentally and $3<U*<4$ computationally.

Upper branch: In the TrSL3 flow regime, the upper branch is wide with range $5<U*<8$ for a smooth cylinder and characterized by a nearly steady increase in amplitude and frequency response leading to the lower branch. The upper branch is wider with experimental setups that can better approximate 2D flows. The gap between the tip of the cylinder and the channel wall is smaller in such setups. Thus, the flow induces lift more efficiently and cannot escape through the tips.

CFD and experiments show the same dependence on Reynolds number. For a PTC-cylinder, the upper branch starts at $U*=4.5−5$ and extends through the range of transition from VIV to galloping around $U*=11$⁠. Both computationally and experimentally, galloping occurs prior to the lower branch and desynchronization of VIV; thus, those two VIV branches are not observed. For cases of higher damping, the amplitude of the VIV branches reduces. Further, galloping initiates later, thus separating VIV from galloping [11]. PTC may initiate galloping but also results in reduced VIV response.

There is an obvious bias of computational and experimental results at this branch both with smooth and rough cylinders. Computationally, the amplitude response increases linearly in the upper branch for a PTC-cylinder while it is nearly constant experimentally. The oscillation frequency, on the other hand, is higher as compared to computational results. CFD cannot capture the flow physics of the late upper branch and therefore it may be said that computational results are not in good accordance with experiments at this flow range of VIV for the PTC-cylinder. For the smooth cylinder case, the large amplitudes occurring in this phase of flow were not captured computationally. This is a typical error encountered in many computational approaches. A discussion was held relevant to this in Ref. [24] and in the same article, a comparison of computational results for the benchmark experiments conducted in Ref. [21] for flows in the TrSL2 regime are given. For flows at higher Reynolds numbers like in the TrSL3 regime, similar problems were observed in this study and also in Ref. [24].

Lower branch: For smooth cylinders at this flow regime (TrSL3), the lower branch is not present. This significant VIV response as compared to TrSL2 flow regime is first discussed in Ref. [4]. For PTC-cylinders, the lower branch is not present. The beginning of the lower branch is overtaken by the upper branch while the end of it is overtaken by the onset of galloping. Although the form of the curve generated by CFD does not capture the amplitude and frequency responses of the cylinder quite well, it may be said that the lower branch is not present in CFD calculations as well.

VIV to galloping transition: Smooth cylinders do not go into galloping unless there is an asymmetry introduced by geometry or upstream nonuniform flow. PTC-cylinders go into galloping as represented by the high amplitudes after $U*=12$⁠. In the range $11<U*<12$⁠, both mechanisms of VIV and galloping work together to induce sudden amplitude increase and frequency change. After certain amplitude is reached, the VIV mechanism desynchronizes and galloping takes over. Transition from VIV to galloping is more gradual computationally than experimentally. Nevertheless, the characteristic amplitude increase in transition from VIV to galloping is shown by CFD as well. Using CFD flow visualization in Sec. 5, the underlying hydrodynamic mechanism is also demonstrated. At this point, however, it should still be highlighted that computationally, the transition to galloping is not as sharp as it is in the experiments. Onset of galloping cannot be pinpointed from computational results.

Fully developed galloping: For $12.5<U*$⁠, galloping is fully developed and the amplitude of oscillation reaches values that result in more complex vortex structures that do not synchronize with the galloping motion. The galloping amplitude continues to increase albeit slower than shown experimentally. The underlying hydrodynamic mechanism is captured by CFD as discussed in Sec. 5 and is consistent with experimental measurements by Park et al. [6].

Based on the computational and experimental results in Sec. 5, we can make the following observations regarding the effect of PTC thickness on FIM.

VIV amplitude: As shown in Fig. 5, the amplitude response of the cylinder remains almost steady and PTC thickness has only a slight effect in the upper branch in VIV.

VIV frequency: As shown in Fig. 6, frequency response increases with PTC thickness throughout the entire VIV range considered in this work.

Galloping onset: PTC thickness affects the onset of galloping. Thicker PTC induces galloping slightly earlier than the thinner one. The thickness affects the point of flow separation and possible reattachment. It has been reported by Chang et al. that the thicker PTC promotes earlier galloping [11]. The same applies here where the thicker PTC (P60) induces galloping earlier than the thinner one (P180).The PTC thickness has negligible effect in the galloping region, where $A*$ has the same form for different values of the PTC thickness. Galloping is triggered by geometric asymmetry, resulting in negative damping and thus instability. The instability is slightly stronger for thicker PTC, which is consistent with the conclusions by Blevins [15].

Galloping amplitude: The PTC thickness does not affect the amplitude response significantly in the galloping region.

Galloping frequency: The thicker PTC (P60) causes f* = f_osc/f_n,water to have smaller values throughout the whole VIV range. The shape of the $f*$ versus U* curve does not change because galloping is an instability. The PTC thickness affects the magnitude of f*. When compared to the smooth cylinder, PTC increases the frequency of the cylinder.

The effect of the PTC width is experimentally and computationally investigated in this section by comparing the results of the $16 deg$ PTC P60 in Sec. 6 with $8 deg$ PTC P60 strips. The starting points of the roughness strips are set at $30 deg$⁠. The narrower strip (⁠$8 deg$⁠) ends at $38 deg$ while the broader one (⁠$16 deg$⁠) ends at $46 deg$⁠. The amplitude and frequency responses for 8 deg of coverage and their comparison with 16 deg of coverage are given in Figs. 7 and 8, respectively. The following observations can be made:

VIV region: As shown in Fig. 7, the PTC width does not have a major impact on the oscillation amplitude and frequency on VIV. PTC helps alter the point of flow separation from its smooth cylinder location. Its width (coverage) plays a secondary role in VIV as long as the entire sand strip remains inside the same zone of the PTC-to-FIM Map [6,7].

Galloping region: Narrower PTC generates lower amplitudes in the galloping range but has negligible effect in the VIV range. Similarly, Chang et al. [11] observed that the bigger PTC coverage leads to increased amplitudes in the galloping range. They also observed that the bigger PTC coverage leads to decreased amplitudes in the VIV range, which is consistent with the PTC-to-FIM Map, as their PTC was closer to the dominant suppression zone [6]. Their PTC coverage started at $40 deg$⁠. When the robustness of the PTC-to-FIM map was studied [6], it was concluded that the strong suppression zone is the dominant zone when PTC covers multiple zones. That is, the conclusion in this paper is consistent with Ref. [11].

Similar results were observed regarding the effect of PTC thickness using a computational approach. For some part of the galloping range, CFD returns higher amplitudes for the narrower PTC although the difference is really small and hard to distinguish in Fig. 7. However, it is hard to generalize this result as this is only valid for two data points in the range covered in this study. More data points are needed in the galloping phase to establish a solid conclusion.

Frequency response: It is shown in Fig. 8. It is observed that the PTC thickness has a negligible effect on the frequency response of the cylinder. This is also captured computationally, although it is found out that the narrower PTC induces slightly higher frequency in the VIV range.

Additional CFD results with relevant flow visualization are provided in Sec. 8.

Upon verification and experimental validation of the CFD code applied to a specific problem, the code can become a very useful tool in investigating the FIM under study. To gain a greater insight in the flow details, the flow has to be visualized. It is difficult to visualize the flow experimentally at high Reynolds numbers using particle image velocimetry. Therefore, to visualize the flow characteristics in VIV and galloping, CFD can provide a unique tool.

In this section, pressure distributions as in Fig. 9 and shear stress around the cylinder are used to integrate pressure and calculate force. Velocity vector fields in the body vicinity in Figs. 10–13 are used to visualize separation of the flow. The following observations can be made based on Figs. 9–13.

Amplitude response: The difference in amplitudes can be confirmed by the pressure differences between the upper and lower parts of the cylinder as well. The pressure coefficients at the lowest position of the cylinder for reduced velocity $U*=14$ for the cylinders with PTC P60 and P180 are given in Fig. 9. The thicker PTC generates higher amplitudes experimentally and computationally (see Fig. 5) in the transition from VIV to galloping region. It also generates higher amplitudes in the galloping region but with smaller differences as the vortex shedding mechanism acts out of synchronization with the galloping instability driving mechanism.

Vortex shedding: The cylinder flow is studied in four quadrants as shown in Fig. 9. It may be noticed directly that there is a higher pressure gradient at the rear part of the cylinder (between quadrants II and III), which is where vortex shedding occurs. Higher pressure gradients lead to stronger vortex shedding which, in turn, results in higher amplitudes. The protrusions in the pressure coefficient graph are the places of the roughness strips.

Separation: Figures 10–13 present the velocity vectors around the cylinder for $16 deg$ and $8 deg$ coverage at $U*=7$ and $U*=14$⁠. These figures capture the flow speed and vector field at the lowest position of the cylinder oscillation. The boundary layer separation on the lower side of the cylinder occurs earlier at $U*=7$ when compared to the separation at $U*=14$⁠. There is a reattachment of the flow at the upper part of the cylinder, which is not observed at the lower part independent of the flow speed. The formation of the flow characteristics shown in Figs. 10–13 is periodic and at the upper end of the oscillation, they are approximately reversed. The boundary layer separation and flow reattachment at the upper part of the cylinder occur later in galloping compared to the VIV range.

Trapped fluid: The trapped dead fluid in the VIV region (⁠$U*=7$⁠) is greater when compared to the galloping region (⁠$U*=14$⁠). The higher flow momentum in galloping, compared to VIV, pushes the trapped dead fluid zone to a downstream part of the cylinder. In addition, recall that in galloping, vortex shedding is not the driving mechanism. On the contrary, vortex shedding is desynchronized and may act against the driving force in galloping. In the VIV region, the trapped dead fluid restricts the amplitude of the cylinder when compared to the cylinder in the galloping region. In the galloping region, the dead fluid zone is not carried with the cylinder in oscillation. At $U*=7$⁠, which is in the VIV region, the vortex that will be shed once the cylinder reaches its highest position is built at the upper part of the cylinder and waiting to be shed. At $U*=14$⁠, which is in the galloping region, the vortex formed at the upper part of the cylinder is about to be shed. In the galloping region, the number of shed vortices per half cycle is higher when compared to the VIV region. The trapped dead fluid in Figs. 10 and 12 increases the added mass at $U*=7$⁠. This CFD result is consistent with the added mass curve in Ref. [25], showing that the added mass decreases with respect to increasing fluid velocity [25].

Selectively distributed local roughness was applied to the surface of an otherwise smooth circular cylinder in the form of sand strips to control passively the flow kinematics around the cylinder in FIMs in a steady uniform flow. This method of altering the flow is called PTC. Application of PTC has been shown in previous work to result in good agreement between CFD and experiments even in the TrSL3 flow regime. On the contrary, CFD results for smooth cylinders agree with experiments for lower Reynolds numbers only. Nevertheless, CFD results in the TrSL3 flow regime can be used for comparative study of cases. Some of the effects of PTC are studied as function of PTC parameters computationally and are compared to experiments. Computational results were validated quantitatively in comparison to experimental results. Then, CFD visualization of the flow was used to study flow properties that would be difficult to study experimentally. The major conclusions of this paper are as follows:

• CFD visualization enables studying flow kinematic features at the boundary layer scale, such as flow separation and reattachment, which would be difficult to visualize experimentally in the TrSL3 flow regime.

• With PTC, the amplitude response of the cylinder remains steady in the VIV range.

• Galloping is initiated earlier with a thicker PTC.

• The onset of galloping using CFD with the $k−ω SST$ turbulence model is at the same location as shown in experiments, albeit not as steep.

• CFD results reveal that the width of the PTC does not have a major impact on the amplitude and the frequency response of the cylinder. This is in agreement with the experimental conclusion in Ref. [6] as long as the roughness strip remains in one zone of the PTC-to-FIM Map [6,7].

• CFD results reveal the difference between near-cylinder flow characteristics in the VIV and galloping phases. The dead fluid zone attached to the cylinder in the VIV region restricts the maximum amplitude achieved. The dead fluid zone is at a later part of the cylinder in the galloping phase and it breaks loose quicker when compared to the VIV phase. The dead fluid zone is not carried with the cylinder during oscillation, which increases the maximum amplitude in the galloping phase.

Future studies using CFD will cover broader ranges of PTC width and location to compare to experimental results that defined the PTC-to-FIM map [6,7]. The numerical methodology adopted in this paper was insufficient at some specific phases of fluid induced motions like the VIV phase of the cylinder. In future studies, denser grids in the near vicinity of the cylinder and particularly the PTC location will be implemented. The complementary information that CFD visualization contributes is most valuable in interpreting numerical and experimental information based on the location of separation and reattachment of the flow, which is difficult to visualize in the high Reynolds number range of this study.

The financial support of TUBITAK for the first author is gratefully acknowledged. This paper was prepared under Cooperative Agreement No. DE-EE0006780 between Vortex Hydro Energy, Inc., and the U.S. Department of Energy. The Marine Renewable Energy Laboratory is a subcontractor through the University of Michigan.