Large-eddy simulations of flow normal to a circular disk at Re=1.5×105

Accepted Manuscript

Large-eddy simulations of flow normal to a circular disk at Re = 1.5 × 105 Xinliang Tian, Muk Chen Ong, Jianmin Yang, Dag Myrhaug PII: DOI: Reference:

S0045-7930(16)30323-1 10.1016/j.compfluid.2016.10.023 CAF 3311

To appear in:

Computers and Fluids

Received date: Revised date: Accepted date:

20 June 2016 10 October 2016 20 October 2016

Please cite this article as: Xinliang Tian, Muk Chen Ong, Jianmin Yang, Dag Myrhaug, Large-eddy simulations of flow normal to a circular disk at Re = 1.5 × 105 , Computers and Fluids (2016), doi: 10.1016/j.compfluid.2016.10.023

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

Highlights • Flow normal to a circular disk at Re = 1.5 × 105 is studied by LES.

CR IP T

• Two distinct vortex spreading modes in the far wake are observed.

• Three dominant frequencies characterizing the wake instability mechanisms are revealed.

AC

CE

PT

ED

M

AN US

• Time and azimuth averaged statistics of the wake flow are presented.

1

ACCEPTED MANUSCRIPT

Large-eddy simulations of flow normal to a circular disk at Re = 1.5 × 105

CR IP T

Xinliang Tiana,b,∗, Muk Chen Ongc , Jianmin Yanga,b , Dag Myrhaugd a

AN US

State Key Laboratory of Ocean Engineering, Shanghai Jiao Tong University Shanghai 200240, China b Collaborative Innovation Center for Advanced Ship and Deep-Sea Exploration, Shanghai 200240, China c Department of Mechanical and Structural Engineering and Materials Science, University of Stavanger, 4036 Stavanger, Norway d Department of Marine Technology, Norwegian University of Science and Technology NO-7491 Trondheim, Norway

Abstract

AC

CE

PT

ED

M

Numerical simulations of the flow normal to a circular disk have been carried out using the large-eddy simulation (LES) method with Smagorinsky subgrid scale model. The Reynolds number (Re) based on the free stream velocity and the diameter of the disk is 1.5 × 105 . The thickness to diameter ratio of the disk is 0.02. The instantaneous vortical structures in wake of the disk are revealed. Toroidal vortices are formed along the disk edges due to Kelvin-Helmholtz instability. The toroidal vortices break up as the flow goes downstream and then hairpin-like vortices are formed. Worm-like vortices parallel to the axis of the disk are observed in the center region of the medium wake. Two distinct vortex spreading modes in the far wake are observed: one has a main shedding plane; the other has not. Three dominant frequencies at the Strouhal number St1 = 0.01, Stn = 0.148 (corresponding to the natural frequency) and St3 ≈ 0.8–1.35, characterizing that the wake instability mechanisms are found through frequency analysis. The time and azimuth averaged statistics, e.g., streamlines, pressure, vorticity and profiles of velocity and kinetic energy, are also presented and discussed. ∗

Corresponding author Email address: [email protected] (Xinliang Tian)

Preprint submitted to Journal of Computers and Fluids

October 21, 2016

ACCEPTED MANUSCRIPT

Keywords: large-eddy simulation, circular disk, turbulent wake, axisymmetric flow

CR IP T

1. Introduction

AC

CE

PT

ED

M

AN US

The flow around around bluff bodies is of great interest for a wide range of engineering applications. As one of the typical cases, the flow normal to a circular disk may result in complex wake flow. Several experimental studies on this topic have been reported. One of the first published observations of the flow normal to a circular disk is the measurements of the gravity driven flow in a vertical tube by Marshall and Stanton (1931). The thickness ratio (χ = td /D, where td and D are the thickness and diameter of the disk, respectively) of the disk was 0.1 and the Reynolds number (Re = U∞ D/ν, where U∞ is the free stream velocity and ν is the kinematic viscosity of the fluid) was about 195 in their study. Kuo and Baldwin (1967) measured the mean velocity and turbulence intensity distributions in the far wake (downstream at 16D, 32D and 90D) behind a circular disk in a low-speed wind tunnel at Re = 7.1 × 104 . Roos and 2 Willmarth (1971) measured the drag coefficient [CD = 8Fx /(ρU∞ πD2 ), where Fx is the streamwise force component acting on the disk and ρ is the density of the fluid] for a normal circular disk by towing the model through a channel with glycerine-water mixtures. The Reynolds number ranges from 5.06 to 6.04 × 104 in their tests, and the mean drag coefficient shows relatively small variations when Re > 103 , i.e., hCD it ∈ [1.1, 1.3] (where hit represents average in time). Roberts (1973) measured the velocity variations in the wake of a stainless steel disk in a closed circuit wind tunnel at Re = 7.8 × 104 . The spectral and cross-spectral analyses were applied to study the correlation and phase relationships between the hot wire signals. Berger et al. (1990) investigated the near wake of a circular disk at 1.5 × 104 . Re . 3 × 105 in a wind tunnel. The non-dimensional form of the frequency is defined as the Strouhal number St = f D/U∞ , where f is the reference frequency. The wake was found to be dominated by three instability mechanisms: axisymmetric pulsation of the recirculation bubble at a low frequency St1 ≈ 0.05; antisymmetric fluctuations induced by helical vortex structures at a natural frequency Stn ≈ 0.135; and a high frequency instability St3 ≈ 1.62 of the separated shear layer. Experimental studies on the wake flow of a freely falling/rising disk in fluid have also 3

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN US

CR IP T

been carried out by Willmarth et al. (1964) at 1.58 6 Re 6 1.9 × 104 and Fernandes et al. (2007) at 80 6 Re 6 330. Several numerical simulations of the flow normal to a circular disk have been performed. Michael (1966) calculated the drag force for the flow normal to an infinitely thin circular disk at 1.5 6 Re 6 50. The flow was assumed to be steady and axisymmetric, and the calculated drag results agreed well with the experimental data of Willmarth et al. (1964) in the considered Re range. Rimon (1969) calculated the unsteady flow past a thin oblate spheroid at Re = 10, 100, 300, 600, respectively. The flow was restricted to axial symmetry and any variations with azimuth angle were not allowed. At Re = 100, the calculated drag agreed well with the experimental results. As Re increased to higher values, the calculated drag was significantly lower than the experimental results; this was attributed to the strict axisymmetry assumption prohibiting the occurrence of periodic vortex shedding at Re > 100. Rivet et al. (1988) simulated the three-dimensional (3D) flow normal to a circular disk at Re = 190 using Lattice Gas Methods. The axial symmetry was observed to break down and the wake became fully 3D at this Re. More recently, Shenoy and Kleinstreuer (2008) studied the flow normal to a circular disk of χ = 0.1 at 10 6 Re 6 300 using direct numerical simulations (DNS). The wake flow was classified into five regimes: (I) steady axisymmetric flow (10 6 Re < 135); (II) steady asymmetric flow (135 6 Re < 155); (III) 3D periodic flow with regular rotation of the separation region (155 6 Re < 172); (IV) unsteady flow with a plane of symmetry (172 6 Re < 280) and (V) unsteady flow with loss of plane of symmetry (280 6 Re 6 300). The studies on the bifurcations in the wake of a circular disk at low Re were followed up by Fabre et al. (2008) at 115 . Re 6 150, Meliga et al. (2009) at 100 . Re . 150, Auguste et al. (2010) at 150 6 Re 6 218, Chrust et al. (2010) at 100 < Re < 500, Yang et al. (2014c) at 115 6 Re 6 300. The wake instabilities were found to be dependent on both χ and Re. Zhong et al. (2014) conducted LES studies for the extraction and recognition of large-scale structures in the near wake behind a circular disk at Re = 2.2 × 104 . The modes 1 and 2 were found to be the most energetic at all considered locations. Yang et al. (2014a) investigated the turbulent wake behind a circular disk of χ = 0.2 at Re = 104 . Yang et al. (2015) carried out numerical simulations for the flow normal to a circular disk. Due to the limited length of the data, the fourier spectra analysis of the velocity is conducted and the analysis revealed that the instabilities of the wake flow at very low frequencies, i.e., fp ≈ 0.03 and fr ≈ 0.02. Recently, the wake patterns and the paths of a freely falling/rising 4

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN US

CR IP T

disk in a fluid at rest have also been investigated extensively (see e.g., Chrust et al. (2014); Tchoufag et al. (2014); Auguste et al. (2013); Shenoy and Kleinstreuer (2010)). It was reported that the wake transition scenarios were similar to that of the flow patterns behind a fixed disk (Zhong and Lee, 2012). The characteristics of the flow around an oscillating thin circular disk have also been investigated (see e.g., Tao and Thiagarajan (2003a,b); Yang et al. (2014b)). A detailed review on the dynamics of the flow around a disk could be found in Ern et al. (2012). It is noted that ”thin” and ”axisymmetric” are the two key geometric features of a circular disk. The wake patterns behind another typical ”thin” body, a flat plate, and another typical ”axisymmetric” body, a sphere, have been investigated extensively, see e.g. Yang et al. (2012) for an inclined flat plate and Johnson and Patel (1999) for a sphere. The similarities in the geometry of the bluff bodies may result in some similarities in the underlying physics of their wake patterns. Based on the open literature, most of the numerical studies were confined to the cases at low and moderate Reynolds numbers. Characteristics of the turbulent wake behind a circular disk, especially in the near wake, have not been well documented in high Re region. Therefore, investigations on the turbulent wake flow at higher Reynolds numbers are required. In this study, the incompressible flow normal to a circular disk is investigated using the large-eddy simulation (LES) at a Reynolds number of Re = 1.5 × 105 . The reason for selecting of this particular Reynolds number is that the numerical simulations for the flow around a flat plate at the same Reynolds number has been reported recently (Tian et al., 2014). The comparisons between the present numerical results for a circular disk and those for a flat plate could bring some in-depth understanding on the fluid dynamics. The main objective of this work is to investigate the characteristics of both the instantaneous and the mean flow normal to the circular disk. The present study provides more physical insights of the flow which has not been discussed in the previous experimental and numerical studies. The paper is organized as follows. The mathematical formulation and numerical methods are given in Section 2. The computational overview, convergence and validation studies are presented in Section 3. The results and discussion are described in Section 4. Finally, the concluding remarks are given in Section 5.

5

ACCEPTED MANUSCRIPT

2. Mathematical formulation and numerical methods

AN US

∂ui =0 ∂xi

CR IP T

2.1. Mathematical formulation Large-eddy simulation is used in the present study, in which large-scale motions are explicitly computed and small-scale eddies are modeled with a subgrid scale (SGS) model. A filtering is introduced to decompose the variable φ into the sum of a filtered component, φ and a subgrid scale component, φ0 , i.e. φ = φ + φ0 . By applying the filtering operation to the Navier-Stokes (NS) equations for an incompressible flow, the filtered NS equations can be expressed as follows

∂ui ∂ 1 ∂p ∂ 2 ui ∂τij + (ui uj ) = − +ν − ∂t ∂xj ρ ∂xi ∂xj ∂xj ∂xj

(1)

(2)

M

where i, j = 1, 2, 3; x1 , x2 and x3 (or x, y and z) denote three directions in a Cartesian coordinate system; u1 , u2 and u3 (or ux , uy and uz ) are the corresponding resolved velocity components; p is the resolved pressure; t is the time; and the subgrid stress τij is given by (3)

ED

τij = ui uj − ui uj

CE

PT

The first term on the right-hand side of equation (3) can not be computed directly from equations (1) and (2), thus τij must be modelled in order to close the equations. The SGS model proposed by Smagorinsky (1963) is used in the present study, and τij is written as 1 τij − δij τkk = −2νt S ij 3

(4)

AC

where δij is the Kronecker delta; νt is the subgrid eddy viscosity; S ij = 0.5(∂ui /∂xj + ∂uj /∂xi ) is the strain rate in the resolved velocity field. Using dimensional analysis, the eddy viscosity νt can be calculated by the following formula

1

νt = (Cs ∆)2 |S|

(5)

where |S| = (2S ij S ij ) 2 ; the constant Cs = 0.2; ∆ is the filter-width. In this study, the Van Driest damping function (Van Driest, 1956) is applied; 6

ACCEPTED MANUSCRIPT

CR IP T

∆ = min{∆mesh , (κ/C∆ )n[1−exp(−n+ /A+ )]}, where ∆mesh is the cubic root of the mesh cell volume; κ = 0.41 is the von K´arm´an constant; C∆ = 0.158; A+ = 26; n is the normal distance of the first node to the wall; n+ is the nondimensionalized wall normal distance, taken as n+ = u∗ n/ν, where u∗ is the wall friction velocity. The reason for choosing this SGS model is that it has been validated and used to calculate the flow around a flat plate at high Reynolds numbers with success (Tian et al., 2014). It is expected that this SGS model could also be able to calculate the high Reynolds number flows around a circular disk.

ED

M

AN US

2.2. Numerical methods The open source Computational Fluid Dynamics (CFD) code OpenF OAM is used here to solve the equations described in Section 2.1. OpenF OAM is mainly applied for solving problems in continuum mechanics. It is built based on the tensorial approach and object oriented techniques (Weller et al., 1998). The PISO (Pressure Implicit with Splitting of Operators) scheme (pisoF oam) is used in the present study. The spatial schemes for interpolation, gradient, Laplacian and divergence are linear, Gauss linear, Gauss linear corrected and Gauss linear schemes, respectively. All these schemes are in second order. The second order Crank-Nicolson scheme is used for the time integration. Further details of these schemes are given in OpenFOAM (2009). 3. Computational overview, convergence and validation studies

AC

CE

PT

3.1. Computational overview As shown in Figure 1, the computational domain is a cylindrical domain with circular cross-section. A Cartesian coordinate system (x, y, z) is used in the present study and the origin of the coordinates locates at the center of the disk. Here x, y and z represent the streamwise, vertical and horizontal crossstream directions, respectively. The inlet boundary locates at 5D upstream to the center of the disk. The distance from the outlet boundary to the center of the disk is 15D. The diameter of the cross-section of the computational domain is 10D, giving a blockage ratio of 1%, which is considered to give negligible blockage effects. The thickness ratio of the disk is χ = 0.02. In order to present the results of this axisymmetric body more clearly, (y, z) are also expressed by polar coordinates (r, θ), i.e., y = r cos θ and z = r sin θ, see Figure 1 (b). 7

ACCEPTED MANUSCRIPT

Free stream

Free stream

Disk

T

Disk 10 D

y

Outlet

r

x

CR IP T

D

Inlet

y

z

x

x

5D

15D

(b)

AN US

(a)

z

Figure 1: Computational domain (not to scale) and boundary conditions: (a) side view and (b) streamwise view.

AC

CE

PT

ED

M

A uniform flow, ux = U∞ , uy = uz = 0, is prescribed at the inlet boundary; the pressure p and the eddy viscosity νt are set as zero normal gradient at the inlet boundary. At the outlet boundary, the velocity ui and the eddy viscosity νt are specified as zero normal gradient boundary conditions; the pressure is set to zero. On the surface of the disk, no-slip boundary condition is prescribed, i.e. ui = 0; the pressure p and the eddy viscosity νt are set as zero normal gradients. On the free stream boundary, the free-slip and zero normal gradient pressure boundary conditions are applied. The computational domain is divided into blocks in order to control the grid distributions. Due to the difficulties in numerical modelling of the sharp edges (especially for the high Re flow), the plate corners are rounded with the circular arcs with a radius of 0.01D. A very fine grid resolution is employed near the disk surface and the averaged n+ of the first layer is 0.3 (with a local maximum of 3.7 along the disk corners). In the wake region, the grids are also refined in order to capture the wake flow structures more accurately. Figure 2 shows an example of the grid used in this study (for case 4 in Table 1). The computations are carried out on an IBM p575+ distributed sharedmemory system in a Notur project. The parallel computing technique is adopted and the computational domain is divided into 128 subdomains. The field mapping technique is used to map the fully developed flow field from a coarse mesh to a fine mesh, for the purpose of saving computational resources. 8

CR IP T

ACCEPTED MANUSCRIPT

AN US

(a)

y x

ED

M

z

(c)

PT

(b)

y

x

x

CE

z

y z

AC

Figure 2: Cutaway view of grid structures (a) whole computational domain, (b) grids around the disk and (c) grids on the disk surface.

9

ACCEPTED MANUSCRIPT

∆tU∞ /D 0.0025 0.0020 0.0015 0.0010 -

hCD it 1.135 1.139 1.126 1.124 1.15–1.28 -

hLw i/D 2.29 2.42 2.61 2.64 2.5

AN US

Case Elements 1 732,870 2 1,812,000 3 3,980,964 4 7,421,952 EXP-1 EXP-2 -

CR IP T

Table 1: Validation of the results with different grid resolutions and time steps, and comparison of the calculated results with experimental data. Here, EXP-1 and EXP-2 represent the experiments by Roos and Willmarth (1971) at Re=(1–7)×104 and Berger et al. (1990) at Re=2.1×105 , respectively.

ED

M

Then the new mapped case runs for at least 60D/U∞ to make sure that the flow is fully redeveloped on the fine mesh. Finally, the case runs for additional 600D/U∞ from which all the statistical quantities presented in this paper are obtained. The averaging methods based on the time-dependent resolved velocity ui and pressure p are defined as follows: the symbols hit , hiθ and hi represent the average in time, in azimuth direction, and in both time and azimuth direction, respectively. Then, the fluctuation of ui is obtained as(Garnier et al., 2002; Xu et al., 2010) u00i = ui − hui it .

AC

CE

PT

3.2. Grid resolution and validation studies The studies for evaluating the effects of the grid resolution and the time step on the calculated results are carried out. The main characteristics of four typical cases with different configurations are shown in Table 1, where ∆t and hLw i represent time step and the mean recirculation length, respectively. The mean recirculation length hLw i is defined as the streamwise distance from the center of the disk to the position where the time-averaged streamwise velocity changes its sign from negative to positive, see Figure 4(a). As shown in Table 1, the number of mesh elements has been validated from 732,870 to 7,421,952, while the time step has been considered from 0.0025 to 0.001 nondimensional time units. The results of hCD it and hLw i obtained from cases 3–4 with different grid resolutions and time steps show small differences. The present calculated mean drag coefficients agree well with the experimental results by Roos and Willmarth (1971) at Re=(1–7)×104 ; and the values of 10

ACCEPTED MANUSCRIPT

1.0

Cp

0.5

0.0

Case 1 Case 2 Case 3

AN US

back side

Case 4

-0.5 0.0

CR IP T

front side

0.1

0.2

rD

0.3

0.4

0.5

/

ED

M

Figure 3: Comparisons of the distributions of the time and azimuth averaged pressure coefficient (hCp i) over the surface of the disk for the simulations with different grid resolutions and time steps.

CE

PT

hLw i also agree well with the experimental results by Berger et al. (1990) at Re = 2.1 × 105 . Figure 3 shows the time and azimuth averaged distributions of the pressure coefficient hCp i over the disk surface. The definitions of Cp is given as Cp =

p − p∞ 1 2 ρU∞ 2

(6)

AC

where the reference pressure p∞ is taken as the mean pressure at the center of the inlet boundary and τw is the wall shear stress. Figure 3 indicates that the results of hCp i obtained from cases 3–4 almost coincide with each other. Further discussions of Figure 3 will be given in Section 4.3. Figure 4 shows the distributions of the time-averaged streamwise velocity hux it and the pressure coefficient hCp it along the centerline of the computational domain, and it appears that the results of hux it and hCp it obtained from cases 3–4 with agree well with each other. 11

ACCEPTED MANUSCRIPT

a

( )

1.5

CR IP T

/

ux t U

1.0

0.5

Lw

Case 1

0.0

AN US

Case 2

-0.5

-5

0

5

xD

Case 3 Case 4

10

15

/

b

( )

M

1.0

Case 1 Case 2 Case 3

ED

Case 4

PT

Cp t

0.5

CE

0.0

AC

-0.5

-5

0

5

xD

10

15

/

Figure 4: Comparisons of the results with different grid resolutions and time steps: distribution of the time-averaged (a) streamwise velocity hux it and (b) pressure coefficient hCp it along the centerline of the disk.

12

a

0.05

b

0.05

( )

CR IP T

ACCEPTED MANUSCRIPT

( )

Case 2 Case 3

0.04

0.04

AN US

0.04

Case 4

0.03

0.02

0.02

0.03

0.02

M

ux"ux" /U

2

0.03

0.01

ED

0.01

0.01

0.00

0.00 0

1

2

c

0.05

( )

Case 1

0.00 0

1

r/ D

rD

PT

/

2

0

1

2

r/ D

AC

CE

Figure 5: Comparisons of the results with different grid resolutions and time steps: distribution of the time and azimuth averaged streamwise Reynolds stress component hu00x u00x i at (a) x/D = 1, (b) x/D = 5 and (c) x/D = 8.

13

ED

M

AN US

CR IP T

ACCEPTED MANUSCRIPT

PT

Figure 6: Instantaneous flow visualizations behind the disk (a) present result (case 4, identified by Q = 2) and (b) reprint from the smoke visualization in a wind tunnel by Berger et al. (1990) at Re = 1.5 × 104 .

AC

CE

Figure 5 shows the distribution of the time and azimuth averaged streamwise Reynolds stress component hu00x u00x i at x/D = 1, 5 and 8, respectively, for the cases with different grid resolutions and time steps. The results indicate an favorable convergence for the spatial and temporal resolutions. One example of the instantaneous wake flow observation obtained from case 4 is presented in Figure 6(a) and the smoke visualization in a wind tunnel by Berger et al. (1990) at Re = 1.5 × 104 is presented in Figure 6(b) for comparison. Here the vortical structures in Figure 6(a) are identified by the iso-surface of the Q criterion (Hunt et al., 1988)

14

ACCEPTED MANUSCRIPT

ED

M

AN US

CR IP T

1 Q = − (kSk2 − kΩk2 ), (7) 2 where S and Ω denote the strain and the rotation tensor, respectively. As shown in Figure 6, the wake flow structures in these two pictures show similar features for the larger-scale structures. It is indicated that the dominate flow characteristics have been obtained from the present LES study. However, it is noted from Figure 6 that the smaller-scale structures are not well-captured by the present simulations. This is because that the smallest eddies in the wake are not resolved completely in the present LES study. More detailed discussions will be given in Section 4.1. In brief, the results obtained from cases 3–4 agree well with each other, indicating a good convergence for the grid resolution and time step. In order to make an accurate prediction, the results presented in the following are obtained from case 4, i.e., the case with 7 421 952 grid elements and a time step of 0.001D/U∞ . The results calculated in the present study compare favorably with the experimental results (Berger et al., 1990; Roos and Willmarth, 1971). Moreover, the present numerical approach has already been applied successfully to calculate the incompressible flow normal to an infinitely long plate at the same Reynolds number(Tian et al., 2014). 4. Results and discussion

AC

CE

PT

4.1. Instantaneous flow visualizations 4.1.1. Vortical structures in the near and medium wake Figure 6(a) shows the instantaneous flow structures in the near and medium wake (x/D 6 4). Toroidal vortices (VT in Figures 6a and 7) resulting from the Kelvin-Helmholtz instability are formed along the disk edges. The Kelvin-Helmholtz instability here is due to the presence of the velocity shear in the fluid. However, these toroidal vortices were not seen in the low-Re cases (Re 6 300) by Shenoy and Kleinstreuer (2008). The toroidal vortices break up as the flow goes downstream and hairpin-like vortices (VH in Figure 6a) are formed instead. The legs of these hairpin vortices are immersed into the inner part of the wake. Similar hairpin vortices are also observed in the case of flow normal to a circular disk case by Shenoy and Kleinstreuer (2008) at Re = 180, but are more ordered compared with the vortical structures in the present case. Figure 7 shows 15

AN US

CR IP T

ACCEPTED MANUSCRIPT

Figure 7: Grid-overlaid view of the vortical structures (Q = 20) in the (x, y)-plane.

AC

CE

PT

ED

M

the vortical structures identified with the iso-surface of Q = 20 in the near wake. Figure 7 is plotted with the same view and spatial region as in Figure 6(a). Worm-like vortical structures (VW in Figure 7) are clearly seen in the inner region of the wake at x/D > 1.5. In general, these worm-like vortices are in parallel with the streamwise direction. To further explain the vortex evolution mechanism in the wake of the disk, Figures 8(a) and 8(b) show the instantaneous contour plots of the streamwise and vertical vorticity components in the plane z = 0, respectively. The streamwise vorticity component is defined as ωx = ∂uz /∂y − ∂uy /∂z, while the vertical vorticity component is defined as ωz = ∂uy /∂x − ∂ux /∂y. As shown in Figure 8(a), the streamwise vorticity ωx with large magnitudes are mainly distributed in the medium wake region, i.e., x/D > 1.5. The instantaneous contour of the vertical vorticity ωz in Figure 8(b) shows that the large magnitude vorticities mainly exist in the separated shear layer along the edge of the disk. Following the procedure by Gallardo et al. (2014), the inviscid approximation of the vorticity equation is introduced, Dω = (ω · ∇)u Dt

16

(8)

PT

ED

M

AN US

CR IP T

ACCEPTED MANUSCRIPT

AC

CE

Figure 8: Instantaneous contours of the streamwise (a) and the vertical (b) components of vorticity in the plane z = 0.

17

ACCEPTED MANUSCRIPT

Then the generation of the streamwise vorticity ωx could be written as Dωx ∂ux ∂ux ∂ux = ωx + ωy + ωz Dt ∂x ∂y ∂z

(9)

AN US

CR IP T

It is noted that the generation of streamwise vorticity ωx consists of two parts, i.e., the stretching term ωx ∂ux /∂x and the tilting term ωy ∂ux /∂y + ωz ∂ux /∂z. The tilting term represents the interaction between the crossstream gradient of streamwise velocity ∂uz /∂z (or ∂ux /∂y) and the crossstream vorticity ωz (or ωy ). Therefore, as shown in Figure 8(a), the increase of the streamwise vorticity ωx in the medium wake of the disk is partially attributed to the tilting of the cross-stream vorticities by the cross-stream gradients of streamwise velocity.

AC

CE

PT

ED

M

4.1.2. Vortical structures in the far wake The wake flow behind the circular disk at the high Reynolds number is quite different from those observed in the low Reynolds number cases. Figure 9 shows two typical snapshots of the vortical structures in the far wake identified by the iso-surfaces of Q = 0.1. As shown in Figure 9, the aforementioned hairpin vortices are observed in a wide range of the wake. The vortical structures shown on the left-hand side of Figure 9 appear to spread uniformly in every azimuthal direction. By contrast, those shown on the right-hand side of Figure 9 exhibit a pronounced major vortex shedding plane. As shown in Figure 9(b), the width of the wake perpendicular to the major plane does not vary very much as x increases and is approximately confined to the region of |y| < 1.5D, similar to the features of the wake observed in Figures 9(a) and 9(c). However, in the major vortex shedding plane [(x, z)-plane], the spread of the wake increases as the flow goes downstream and reaches z ≈ ±2D till x = 10D, see Figure 9(d). Despite the modulation due to the small-scale turbulence, the large-scale flow topology in the present LES results resembles that of the M Mπ solution (Mixed Mode), analytically reconstructed from the superposition of one stationary mode and two counter-rotating spiral modes for the low Re number case, see Meliga et al. (2009) and Fabre et al. (2008) for more details. 4.2. Frequency analysis To demonstrate the existence of the so-called axisymmetric pulsation (Berger et al., 1990) in the present simulation, the instantaneous 18

(b)

AN US

(a)

CR IP T

ACCEPTED MANUSCRIPT

(c)

M

(d)

ED

(f )

CE

PT

(e)

AC

Figure 9: Two typical snapshots of the vortical structures: left-hand side for one snapshot and right-hand side for the other snapshot. (a, b) (x, y)-plane, (c, d) (x, z)-plane and (e, f ) (y, z)-plane.

19

ACCEPTED MANUSCRIPT

1.0

Lw

0.5

2

Lw

/

ux U

1

0.0

Lw Lw Lw = (

-0.5

0

1

2

3

xD

4

/

3.5

3.0

/

Lw D

+

1

)/2

2

5

6

AN US

b

( )

CR IP T

a

( )

2.5

M

2.0

1.5

100

ED

0

(c)

St

200

300

tU /D

400

500

600

=0.01

AC

CE

PSD

PT

1

0.01

0.1

1

f D/U

Figure 10: (a) Definition, (b) the time-dependent variation and (c) the corresponding power spectrum density (PSD) of the instantaneous recirculation length Lw

20

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN US

CR IP T

recirculation length Lw is calculated. Figure 10 (a) shows one example definition of Lw , where Lw1 and Lw2 denote the streamwise distances from the center of the disk to the positions where the resolved instantaneous streamwise velocity changes its sign from negative to positive, respectively. Lw is calculated as the averaged value of Lwi , where i = 1, 2, 3, · · ·. The time-dependent variation of Lw is shown in Figure 10 (b). The instantaneous recirculation length Lw directly reflects the variations of the recirculation bubble size. Therefore, the frequency of the axisymmetric pulsation could be obtained based on the frequency analysis for Lw . The power spectrum density of Lw is shown in Figure 10 (c). It is clearly observed that the spectrum of Lw has a dominant frequency at St1 = 0.01. Therefore, the low-frequency axisymmetric pulsation of the wake is mainly dominated by the frequency at St1 = 0.01. The value of St1 is lower than the results of St1 ≈ 0.05 by Berger et al. (1990) and St1 = 0.041 by Shenoy and Kleinstreuer (2008). A similar peak at St = 0.01 was also found by Roberts (1973), but the author did not explain the reason for this peak frequency. The lack of capture of the dominant low-frequency around St = 0.01 in the experiment by Berger et al. (1990) is probably due to that the disk was forced to vibrate in a nutating mode. The duration of the calculation by Shenoy and Kleinstreuer (2008) was not long enough to capture such a low-frequency mode. The time traces of the pressure coefficient (Cp ) sampled along the line of x/D = 0.5, θ = 0◦ from r/D = 0 to 2 with 20 equal steps are plotted in Figure 11. Appropriate vertical offsets (without rescaling) are made for clarity. The sampled line of x/D = 0.5, θ = 0◦ origins from the axis of the wake region, goes through the shear layer and eventually reaches the far field region. The pressure fluctuations at the locations of r/D = 0.7 and 0.8 are more energetic than those at other locations. Unlike the case of flow around a cylindrical structure, the von K´arm´an vortex shedding is not obvious for this axisymmetric body. Therefore, the pressure fluctuation in the wake is dominated by the high-frequency instability of the separated shear layer which is associated with the Kelvin-Helmholtz instability. Furthermore, a low-frequency unsteadiness appears to exist in the time traces of Cp . Figure 12 shows the spectra of the time-dependent variations of the pressure coefficient shown in Figure 11. Appropriate vertical offsets and rescaling are made in Figure 12 for clarity. As shown, the previously observed dominant frequency at St1 = 0.01 is identified and indicated in the region of r/D 6 1.5. The pressure spectra in the region of r/D > 1 21

ACCEPTED MANUSCRIPT

rD /

2.0

1.5

CR IP T

1.0

Cp

1.0

0

0

100

200

AN US

0.5

300

Non-dimensional time,

400

tU D

500

600

/

M

Figure 11: Time-dependent variations of the pressure coefficient Cp obtained at 21 locations along the line of x/D = 0.5, θ = 0◦ with equal steps from r/D = 0 to 2.

CE

PT

ED

have a peak at approximately Stn = 0.148 due to the vortex shedding, which is slightly higher than the previous DNS result for Stn = 0.122 at Re = 300 by Shenoy and Kleinstreuer (2008) and the measurement for Stn = 0.135 at Re = 104 -105 by Berger et al. (1990). For the spectra obtained at r/D = 0.7–0.9, the dominant frequency St3 associated with the separated shear layer instability are observed at approximately 0.8–1.35, see Figure 12. Berger et al. (1990) reported that St3 = 1.62, which is higher than the value of the present calculated St3 .

AC

4.3. Reynolds-averaged statistics 4.3.1. Streamlines Figure 13 shows the streamlines of the time and azimuth averaged flow around the disk. As shown in Figure 13(a), the steady incoming flow is normal to the disk, separates at the center of the front side of the disk and then approaches the edge of the disk along the radial direction; the flow separates again at the edge of the disk and eventually a recirculation bubble is formed behind the disk. The length of the recirculation bubble in the streamwise direction corresponds to the mean recirculation length of hLw i = 22

CR IP T

ACCEPTED MANUSCRIPT

r/D 2.0

St

1.0

St

PT

0

= 0.01

St

3

0.8-1.35

ED

1

M

1.5

0.5

= 0.148

AN US

n

0.1

1

10

f D/U

CE

0.01

AC

Figure 12: Spectra of the pressure coefficient (Cp ) along the line of x/D = 0.5, θ = 0◦ from r/D = 0 to 2 with 20 equal steps.

23

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN US

CR IP T

2.64D, see Table 1. The maximum diameter of the recirculation bubble in the cross-stream direction (denoted as Dbub in Figure 13a) is 1.6D in this study. For the case of the flow normal to a circular disk by Shenoy and Kleinstreuer (2008) at Re = 100, in which the flow is steady and axisymmetric, the size of the recirculation region is smaller than in the present high Re case, and Lw ≈ 1.97D and Dbub ≈ 1.3D. One can imagine that, in a 3D view, that the recirculation region is ringlike in shape and the vortex core C1 (indicated in Figure 13a) is a circle at the position of approximately (x/D, r/D) = (1.55, 0.65). It appears that the vortex core C1 is squeezed to the margin of the recirculation region, unlike the case of the flow around a cylindrical structure, in which the vortex core is approximately located at the center of the half recirculation bubble, see Figure 5 in Breuer (1998) for a circular cylinder and Figure 5(a) in Narasimhamurthy and Andersson (2009) for an infinitely long plate, respectively. On the cross-stream section where C1 locates, the streamwise velocities in and out of C1 are negative and positive, respectively. If it is simply assumed that the streamwise velocities in and out of C1 are close in magnitudes, the area of the cross-stream section in and out of C1 should also be similar. In this case the vortex core C1 is squeezed to be closer to the outer margin of the recirculation bubble due to the present axisymmetric flow configuration. However, the absolute value of the streamwise velocity inside the vortex core C1 is smaller than that out of C1 (as will be shown in Figure 16), because the flow out of C1 is accelerated by the main flow. The position of C1 is squeezed again to be much closer to the outer margin of the recirculation bubble. It is interesting to note the existence of a flat secondary vortex adhering to the back side of the disk, see Figure 13(b). The rotating direction of this secondary vortex is opposite to the main recirculation vortex. Therefore, on the average, the direction of the velocity in close vicinity the back side of the disk is from the disk edge to the center. This secondary vortex is also observed in some other separating and reattaching flows, see the examples of the forward-facing step case in Figure 4 of Lamballais et al. (2010) and the circular cylinder case in Figure 4 of Breuer (2000), respectively. However, this secondary vortex is not observed near the back side of a flat plate with infinite span, see Narasimhamurthy and Andersson (2009). 4.3.2. Pressure distributions Figure 3 shows the distributions of the time and azimuth averaged pressure coefficient (hCp i). The pressure coefficient on the back side of the 24

ACCEPTED MANUSCRIPT

1 0.4

P1

Dbub

0.2

C1 -1

-1

0

1

2

3

4

xD

CR IP T

/

0

rD

/

rD

P2

0.0 -0.2

0.0

0.2

xD

/

/

(b)

AN US

(a)

Figure 13: Streamlines of the time and azimuth averaged flow field in (x, r)-plane: (a) around the disk and (b) zoomed view of (a).

AC

CE

PT

ED

M

disk is found to be constant, i.e., hCp i = −0.43. The flat distribution of pressure has also been reported widely for the flow around bluff bodies, see, e.g., Xu et al. (2010) for a circular cylinder and Najjar and Balachandar (1998) for a plate with infinite span. Note that the base pressure for the circular disk in this study is much higher than in the case of an infinitely long plate(Fage and Johansen, 1927) in which the base pressure coefficient is about -1.36. The reason for this difference in the base pressure is that the cross flow are free in both y- and z-directions for the circular disk, while the cross flow is restricted in the spanwise direction for the infinitely long plate. Figure 14 shows the contours of the time and azimuth averaged pressure coefficient hCp i in the (x, r)-plane. The pressure difference between the front and back sides of the disk dominates the drag force, see Figure 3. It is observed that the pressure gradient is relatively small near the back side of the disk with values in the range −0.45 to −0.4 (which has also been shown in Figure 3). The position of the minimum pressure coefficient is located at (x/D, r/D) ≈ (1.5, 0.55), which is close to the position of the core of the main recirculation bubble (C1 in Figure 13a), and the minimum value of the pressure coefficient is lower than −0.6.

4.3.3. Azimuthal vorticity contours Figure 15 shows the time and azimuth averaged azimuthal vorticity contours at the four levels of ωθ = −0.5, −1, −2 and −5. The azimuthal vorticity is defined as ωθ = ∂ur /∂x − ∂ux /∂r. As shown, the mean 25

ACCEPTED MANUSCRIPT

1

0.2

-0.45

/

rD

0.4

-0.5 -0.55

-0.6

0.8

0

-0.6

-1

-1

0

1

2

xD

3

4

AN US

/

CR IP T

-0.2 -0.4

AC

CE

PT

ED

M

Figure 14: Contours of the time and azimuth averaged pressure coefficient at 9 levels of hCp i = −0.6, −0.55, −0.5, −0.45, −0.4, −0.2, 0.2, 0.4 and 0.8, respectively.

Figure 15: Contours of the time and azimuth averaged azimuthal vorticity at 4 levels of hωθ i = −0.5, −1, −2 and −5, respectively.

26

ACCEPTED MANUSCRIPT

CR IP T

azimuthal vorticity remains low in the near wake, and the shape of the vorticity contours is qualitatively similar to the result by Shenoy and Kleinstreuer (2008) at Re = 100. This similarity indicates that, on the average, the present high Re = 1.5 × 105 case has some of the basic flow features as observed in the low Re = 100 case; although the wake flow at this high Re is actually turbulent.

PT

ED

M

AN US

4.3.4. Velocity profiles The time and azimuth averaged velocity profiles at 8 cross-stream sections from x/D = 0.5 to 4 with equal steps of 0.5 are plotted in Figures 16(a) and 16(b) for ux and ur , respectively. As shown in Figure 16(a), at all the presented cross-sections, the minimum hux i values are obtained at the center of the wake, i.e., r = 0. The minimum hux i in the entire recirculation bubble is obtained at the position (x/D, r/D) = (1.49, 0); see also Figure 4(a). The U-shaped hux i profile is observed in the near wake, e.g., at x/D = 0.5. The hux i profiles transform gradually from U-shaped profiles in the near wake to the V-shaped profiles in the far wake. Moreover, the hux i profiles become more smooth in the far wake than in the near wake. For the radial velocity (ur ), positive values represent that the flow direction is from the wake center (r = 0) to the outer region, and vice versa. As shown in Figure 16(b), the hur i profiles at x/D = 0.5 and 1 are positive indicating that the averaged flow goes from the inner wake toward the outer region; while the hur i profiles at x/D > 2 are negative, indicating that the averaged flow there goes from the outer region toward the inner wake. The hur i profile at x/D = 1.5 crosses the zero level near r/D = 0.7; this is because x/D = 1.5 is near the recirculation core C1, see Figure 13(a).

AC

CE

4.3.5. Kinetic energy profiles Figure 17(a) shows the time and azimuth averaged profiles of the fluctuating kinetic energy, hki = (hu00x u00x i + hu00y u00y i + hu00z u00z i)/2, at 8 cross-stream sections from x/D = 0.5 to 4 with equal steps of 0.5. The distributions of hki as well as its contributors in the streamwise direction and in the cross-stream plane, denoted as hkx i and hky + kz i, respectively, along the wake centerline are shown in Figure 17(b), where hkx i = hu00x u00x i/2 and hky + kz i = (hu00y u00y i + hu00z u00z i)/2. It appears that the contributions of hky + kz i to hki is about four times to that of hkx i except in the near wake region, i.e., x/D < 1. In Figure 17(b), it is observed that there are two peaks in hki. Careful examination shows that the locations of these two 27

a

CR IP T

ACCEPTED MANUSCRIPT

1.5

( )

xD xD xD xD xD xD xD xD

0.5

0.0

-0.5 0

1

/

= 0.5

/

= 1.0

/

= 1.5

AN US

/

ux U

1.0

2

rD

3

/

= 2.0

/

= 2.5

/

= 3.0

/

= 3.5

/

= 4.0

4

5

/

b

0.3

xD xD xD xD xD xD xD xD

M

( )

0.2

/

ur U

ED

0.1

/

= 0.5

/

= 1.0

/

= 1.5

/

= 2.0

/

= 2.5

/

= 3.0

/

= 3.5

/

= 4.0

0.0

PT

-0.1

-0.2

1

2

rD

CE

0

3

4

5

/

AC

Figure 16: Time and azimuth averaged profiles of the velocity components: (a) hux i and (b) hur i.

28

ACCEPTED MANUSCRIPT

a

0.15

CR IP T

( )

0.10

/

k U

2

xD xD xD xD xD xD xD xD

0.00 0

=0.5

/

=1.0

/

=1.5

/

=2.0

/

=2.5

AN US

0.05

/

1

2

3

/

=3.0

/

=3.5

/

=4.0

4

5

rD /

b

( )

M

0.02

/

k U

2

ED

/

k U

2

0.10

CE

xD

0.00 0.0

PT

0.05

0.2

/

0.4

k U kx U ky+kz U /

/

2

2

/

2

0.00

3

xD

6

/

AC

0

9

12

15

Figure 17: Time and azimuth averaged profiles of the kinetic energy hki: (a) at 8 crossstream sections from x/D = 0.5–4 with equal steps of 0.5 and (b) along wake centerline (r = 0).

29

ACCEPTED MANUSCRIPT

CR IP T

peaks correspond to the two stagnation points in the mean velocity field behind the disk (marked with P1 and P2 in Figure 13a, respectively). By denoting these two peaks as the major and the minor peaks, respectively, Figure 17(b) shows that the major peak has contributions from hkx i and hky + kz i. However, the minor peak near the back side of the disk is mainly due to hky + kz i, and hkx i is suppressed near the disk surface. Furthermore, 2 the peak value of hki/U∞ for the circular disk case (≈ 0.12) is lower than in the flat plate case with an infinite span(Narasimhamurthy and Andersson, 2009), in which the peak values varies up to 0.4.

AN US

5. Concluding remarks

ED

M

The flow normal to a circular disk at Re = 1.5 × 105 has been investigated by using the large-eddy simulations (LES). The Smagorinsky subgrid scale model is used to account for the turbulence scales filtered by the grid. The thickness to diameter ratio of the disk is 0.02. The calculated mean drag force agrees well with the previous experimental results. The instantaneous flow visualizations obtained in the present study show qualitative agreement with the experimental observations. Frequency analysis is applied to investigate the frequency characteristics in the turbulent wake. Moreover, the instantaneous vortical structures and the time and azimuth averaged quantities are presented in the paper. The main conclusions are:

AC

CE

PT

• Toroidal vortices along the disk edges induced by the Kelvin-Helmholtz instability is observed, which is contrary to that in the lower Reynolds number cases by Shenoy and Kleinstreuer (2008) (Re 6 300), in which the toroidal vortices are not observed. The toroidal vortices break up, and the hairpin-like vortices are formed as the flow goes downstream. Worm-like vortex structures exist in the inner region of the medium wake. Two distinct vortex modes in the far wake are identified: the spread of the vortex in the low-drag phase is approximately uniform in azimuth direction while the spread of the vortex in the high-drag phase has a major shedding plane. • Three dominant frequencies at St1 = 0.01, Stn = 0.148 and St3 ≈ 0.8–1.35 characterizing the wake instability mechanisms are revealed through frequency analysis. The latter two frequencies agree well with 30

ACCEPTED MANUSCRIPT

CR IP T

the values reported in the experimental studies by Berger et al. (1990). Frequency analysis of the instantaneous recirculation length proves that the low-frequency pulsation of the recirculation bubble occurs at St0 = 0.01.

AN US

• A flat secondary vortex adhering to the back side of the disk is observed in the time and azimuth averaged flow. The constant pressure distribution on the back side of the disk is observed. In the present case with Re = 1.5 × 105 , the time and azimuth averaged azimuthal vorticity contours show similar features as the low Re case by Shenoy and Kleinstreuer (2008). The streamwise velocity profiles along the radial direction changes from U-shaped in the near wake region to V-Shaped in the far wake region. Moreover, the double peaks in the mean fluctuating kinetic energy profile along the wake centerline are also identified. Acknowledgment

PT

ED

M

The authors thank Dr. Lars Erik Holmedal for the technical discussions. The authors also thank the support of the National Natural Science Foundation of China (Grant nos. 51509152, 51239007 and 11632011), Shanghai Yang Fan Program (Grant no. 15YF1406100) and Professor Bjørnar Pettersen’s support to use the IBM p575+ distribute shared-memory system at the Norwegian University of Science and Technology. The first author also acknowledges the support of the China Scholarship Council (CSC).

CE

References

Auguste, F., Fabre, D., Magnaudet, J., 2010. Bifurcations in the wake of a thick circular disk. Theor. Compu. Fluid Dyn. 24, 305–313.

AC

Auguste, F., Magnaudet, J., Fabre, D., 2013. Falling styles of disks. J. Fluid Mech. 719, 388–405. Berger, E., Scholz, D., Schumm, M., 1990. Coherent vortex structures in the wake of a sphere and a circular disk at rest and under forced vibrations. J. Fluids Struct. 4, 231–257.

31

ACCEPTED MANUSCRIPT

Breuer, M., 1998. Large eddy simulation of the subcritical flow past a circular cylinder: numerical and modeling aspects. Int. J. Numer. Methods Fluids 28, 1281–1302.

CR IP T

Breuer, M., 2000. A challenging test case for large eddy simulation: high Reynolds number circular cylinder flow. Int. J. Heat Fluid Flow 21, 648– 654.

Chrust, M., Bouchet, G., Duˇsek, J., 2010. Parametric study of the transition in the wake of oblate spheroids and flat cylinders. J. Fluid Mech. 665, 199–208.

AN US

Chrust, M., Bouchet, G., Duˇsek, J., 2014. Effect of solid body degrees of freedom on the path instabilities of freely falling or rising flat cylinders. J. Fluids Struct. 47, 55–70. Ern, P., Risso, F., Fabre, D., Magnaudet, J., 2012. Wake-induced oscillatory paths of bodies freely rising or falling in fluids. Ann. Rev. Fluid Mech. 44, 97–121.

M

Fabre, D., Auguste, F., Magnaudet, J., 2008. Bifurcations and symmetry breaking in the wake of axisymmetric bodies. Phys. Fluids 20, 051702.

PT

ED

Fage, A., Johansen, F.C., 1927. On the flow of air behind an inclined flat plate of infinite span. Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character 116, 170– 197.

CE

Fernandes, P.C., Risso, F., Ern, P., Magnaudet, J., 2007. Oscillatory motion and wake instability of freely rising axisymmetric bodies. J. Fluid Mech. 573, 479–502.

AC

Gallardo, J.P., Andersson, H.I., Pettersen, P., 2014. Turbulent wake behind a curved circular cylinder. J. Fluid Mech. 742, 192–229. Garnier, E., Sagaut, P., Deville, M., 2002. Large eddy simulation of shock/homogeneous turbulence interaction. Comput. Fluids 31, 245–268. Hunt, J.C.R., Wray, A.A., Moin, P., 1988. Eddies, streams, and convergence zones in turbulent flows, in: Center for Turbulence Research Report CTRS88, NASA, USA. pp. 193–208. 32

ACCEPTED MANUSCRIPT

Johnson, T., Patel, V., 1999. Flow past a sphere up to a Reynolds number of 300. J. Fluid Mech. 378, 19–70.

CR IP T

Kuo, Y.H., Baldwin, L.V., 1967. The formation of elliptical wakes. J. Fluid Mech. 27, 353–360.

Lamballais, E., Silvestrini, J., Laizet, S., 2010. Direct numerical simulation of flow separation behind a rounded leading edge: Study of curvature effects. Int. J. Heat Fluid Flow 31, 295–306.

AN US

Marshall, D., Stanton, T.E., 1931. On the eddy system in the wake of flat circular plates in three dimensional flow. Proceedings of the Royal Society of London. Series A, Containing Papers of a Mathematical and Physical Character 130, 295–301. Meliga, P., Chomaz, J.M., Sipp, D., 2009. Global mode interaction and pattern selection in the wake of a disk: a weakly nonlinear expansion. J. Fluid Mech. 633, 159–189.

M

Michael, P., 1966. Steady motion of a disk in a viscous fluid. Phys. Fluids 9, 466.

ED

Najjar, F.M., Balachandar, S., 1998. Low-frequency unsteadiness in the wake of a normal flat plate. J. Fluid Mech. 370, 101–147.

PT

Narasimhamurthy, V.D., Andersson, H.I., 2009. Numerical simulation of the turbulent wake behind a normal flat plate. Int. J. Heat Fluid Flow 30, 1037–1043.

CE

OpenFOAM, 2009. The Open Source CFD Toolbox, Programmer’s Guide, Version 1.6. OpenCFD Limited, Boston, MA, USA.

AC

Rimon, Y., 1969. Numerical solution of the incompressible time-dependent viscous flow past a thin oblate spheroid. Phys. Fluids 12, II–65–75. Rivet, J.P., Henon, M., Frisch, U., d’Humieres, D., 1988. Simulating fully three-dimensional external flow by lattice gas methods. Europhys. Letters 7, 231. Roberts, J.B., 1973. Coherence measurements in an axisymmetric wake. AIAA J. 11, 1569–1571. 33

ACCEPTED MANUSCRIPT

Roos, F.W., Willmarth, W.W., 1971. Some experimental results on sphere and disk drag. AIAA J. 9, 285–291.

CR IP T

Shenoy, A.R., Kleinstreuer, C., 2008. Flow over a thin circular disk at low to moderate Reynolds numbers. J. Fluid Mech. 605, 253–262.

Shenoy, A.R., Kleinstreuer, C., 2010. Influence of aspect ratio on the dynamics of a freely moving circular disk. J. Fluid Mech. 653, 463–487. Smagorinsky, J., 1963. General circulation experiments with the primitive equations. Mon. Weather Rev. 91, 99–164.

AN US

Tao, L., Thiagarajan, K., 2003a. Low KC flow regimes of oscillating sharp edges I. Vortex shedding observation. Appl. Ocean Res. 25, 21–35.

Tao, L., Thiagarajan, K., 2003b. Low KC flow regimes of oscillating sharp edges II. Hydrodynamic forces. Appl. Ocean Res. 25, 53–62.

M

Tchoufag, J., Fabre, D., Magnaudet, J., 2014. Global linear stability analysis of the wake and path of buoyancy-driven disks and thin cylinders. J. Fluid Mech. 740, 278–311.

ED

Tian, X., Ong, M.C., Yang, J., Myrhaug, D., 2014. Large-eddy simulation of the flow normal to a flat plate including corner effects at a high Reynolds number. J. Fluids Struct. 49, 149–169.

PT

Van Driest, E.R., 1956. On turbulent flow near a wall. J. Aero. Sci. 23, 1007–1011.

CE

Weller, H.G., Tabor, G., Jasak, H., Fureby, C., 1998. A tensorial approach to computational continuum mechanics using object-oriented techniques. Comput. Phys. 12, 620–631.

AC

Willmarth, W.W., Hawk, N.E., Harvey, R.L., 1964. Steady and unsteady motions and wakes of freely falling disks. Phys. Fluids . Xu, C.Y., Chen, L., Lu, X., 2010. Large-eddy simulation of the compressible flow past a wavy cylinder. J. Fluid Mech. 665, 238–273. Yang, D., Pettersen, B., Andersson, H.I., Narasimhamurthy, V.D., 2012. Vortex shedding in flow past an inclined flat plate at high incidence. Phys. Fluids 24, 084103. 34

ACCEPTED MANUSCRIPT

Yang, J., Liu, M., Wu, G., Liu, Q., Zhang, X., 2015. Low-frequency characteristics in the wake of a circular disk. Phys. Fluids 27, 064101.

CR IP T

Yang, J., Liu, M., Wu, G., Zhong, W., Zhang, X., 2014a. Numerical study on coherent structure behind a circular disk. J. Fluids Struct. 51, 172–188.

Yang, J., Tian, X., Li, X., 2014b. Hydrodynamic characteristics of an oscillating circular disk under steady in-plane current conditions. Ocean Eng. 75, 53 – 63. URL: http://www.sciencedirect.com/science/article/pii/S0029801813003995, doi:http://dx.doi.org/10.1016/j.oceaneng.2013.11.011.

AN US

Yang, J., Wu, G., Zhong, W., Liu, M., 2014c. Numerical study on bifurcations in the wake of a circular disk. Int. J. Comput. Fluid Dyn. 28, 187–203. Zhong, H.J., Lee, C.B., 2012. The wake of falling disks at low Reynolds numbers. Acta Mech. Sin. 28, 1–5.

AC

CE

PT

ED

M

Zhong, W., Liu, M., Wu, G., Yang, J., Zhang, X., 2014. Extraction and recognization of large-scale structures in the turbulent near wake of a circular disc. Fluid Dyn. Res. 46, 1–21.

35