F Y Houdroge

T Leweke

K Hourigan

M C Thompson

Two-and three-dimensional wake transitions of an impulsively started uniformly rolling circular cylinder

. The main focus here is for Reynolds numbers beyond the transition to unsteady flow at Re c,2D = 88. From impulsive start up, the wake almost immediately undergoes transition to a periodic two-dimensional wake state, which, in turn, is three-dimensionally unstable. Thus, the previous threedimensional stability analysis based on the two-dimensional steady flow provides only an element of the full story. Floquet analysis based on the periodic two-dimensional flow was undertaken and new three-dimensional instability modes were revealed. The results suggest that an impulsively started cylinder rolling along a surface at constant velocity for Re 90 will result in the rapid development of a periodic two-dimensional wake that will be maintained for a considerable time prior to the wake undergoing threedimensional transition. Of interest, the mean lift and drag coefficients obtained from full three-dimensional simulations match predictions from two-dimensional simulations to within a few percent.

Introduction

Many previous studies have focused on the flow around a circular cylinder in an unbounded flow. In the Stokes range (Re = U d/ν 1), viscous effects dominate the flow. The flow around a stationary cylinder remains attached and symmetrical about the spanwise and streamwise axes through the centre point of the cylinder. As the Reynolds number is increased, the flow loses its upstream/downstream symmetry as the fluid separates at the rear of the cylinder. This results in the formation of two closed recirculation zones, first occurring at 5 Re 7 [START_REF] Taneda | Experimental investigation of the wakes behind cylinders and plates at low Reynolds numbers[END_REF][START_REF] Dennis | Numerical solutions for steady flow past a circular cylinder at Reynolds numbers up to 100[END_REF]. The length of these recirculation regions was found to increase linearly with Re, until at Re 46, the wake becomes absolutely unstable and undergoes a transition to a periodic flow state [START_REF] Taneda | Experimental investigation of the wakes behind cylinders and plates at low Reynolds numbers[END_REF][START_REF] Roshko | On the development of turbulent wakes from vortex streets[END_REF][START_REF] Henderson | Nonlinear dynamics and pattern formation in turbulent wake transition[END_REF][START_REF] Provansal | Bénard-von Kármán instability: transient and forced regimes[END_REF]. This transition is the result of a Hopf bifurcation (i.e., a steady to unsteady transition) of the steady flow that occurs as the flow becomes globally absolutely unstable [START_REF] Provansal | Bénard-von Kármán instability: transient and forced regimes[END_REF][START_REF] Henderson | Nonlinear dynamics and pattern formation in turbulent wake transition[END_REF]. The saturated state of vortex shedding in the wake of the cylinder takes the form of a Bénard-von Kármán vortex street [START_REF] Bénard | Formation de centres de giration à l'arrière d'un obstacle en mouvement[END_REF][START_REF] Von Kármán | Über den Mechanismus des Widerstandes, den ein bewegter Körper in einer Flüssigkeit erfährt[END_REF] that is characterised by a periodic, repeating pattern of swirling vortices of opposite sign that are shed from the rolled-up shear layers.

As the Reynolds number is increased further, the now-periodic wake undergoes a further transition to three-dimensional flow. [START_REF] Williamson | The existence of two stages in the transition to three-dimensionality of a cylinder wake[END_REF] found that two clearly identifiable transitions take place sequentially that are distinguished by the development of distinct spatio-temporal three-dimensional wake states designated mode A and mode B. The first of these transitions is also accompanied by a discontinuity in the Strouhal-Reynolds number curve. The mode A instability appears beyond Re ≈ 190 (Williamson 1996b;[START_REF] Henderson | Nonlinear dynamics and pattern formation in turbulent wake transition[END_REF], resulting in pairs of counter-rotating streamwise vortices forming along the span of the cylinder. This three-dimensional mode has a spanwise wavelength of λ ≈ 4d, where d is the diameter of the cylinder. The second transition to Mode B becomes fully developed at Re = 260 and has a preferred spanwise wavelength of λ 0.8d [START_REF] Henderson | Nonlinear dynamics and pattern formation in turbulent wake transition[END_REF]Williamson 1996b). The remnants of the streamwise Mode B vortical structures can be seen at much higher Reynolds numbers, well beyond the development of fully turbulent flow (e.g. [START_REF] Wu | Three-dimensional vortex structures in a cylinder wake[END_REF].

Imposing a rotation on a cylinder in an unbounded flow has a strong influence on the wake structure and transitions. The degree of rotation is often quantified by the nondimensional rotation rate, α = ωd/(2U ), defined as the ratio of the tangential surface speed (ωd/2, with ω the angular velocity) and the free-stream speed U . Many authors, including [START_REF] Tang | On steady flow past a rotating circular cylinder at Reynolds numbers 60 and 100[END_REF], have shown that imposing a rotation on the body renders the wake asymmetrical and, at Re 60, depending on the rotation rate, the elimination of one or both of the recirculation regions in the wake can be observed. For larger Re, the imposed rotation may also suppress or delay the transition to unsteady flow in comparison to the case of a non-rotating body.

For the non-rotating cylinder, as the Reynolds number is increased, the wake becomes unsteady. At low rotation rates, a Bénard-von Kármán vortex street is observed [START_REF] Jaminet | Experiments on vortex shedding from rotating circular cylinders[END_REF], also known as Mode I shedding. For higher values of α, the unsteady wake narrows and is displaced in the direction of motion of the cylinder surface [START_REF] Díaz | Vortex shedding from a spinning cylinder[END_REF][START_REF] Mittal | Flow past a rotating cylinder[END_REF]. As the rotation rate increases beyond a critical value of α c 2 [START_REF] Mittal | Flow past a rotating cylinder[END_REF][START_REF] Díaz | Vortex shedding from a spinning cylinder[END_REF], the unsteady flow is completely suppressed. Instead of vortex shedding, the surrounding fluid is entrained by the rotation of the body and creates a layer around the cylinder that thickens as α increases [START_REF] Díaz | Vortex shedding from a spinning cylinder[END_REF][START_REF] Mittal | Flow past rotating cylinders: effect of eccentricity[END_REF]. Perhaps surprisingly, a second shedding regime is observed over a specific range of α at much higher α [START_REF] Mittal | Flow past a rotating cylinder[END_REF][START_REF] Kumar | Flow past a rotating cylinder at low and high rotation rates[END_REF], where single-sided vortex shedding occurs with a period much longer than that of Mode I shedding. This wake state is referred to as Mode II shedding.

Limited investigations have been carried out on the development of three-dimensional wakes for a rotating cylinder. At low rotation rates, α < 1, [START_REF] Akoury | The threedimensional transition in the flow around a rotating cylinder[END_REF] found that Mode A becomes unstable at higher Reynolds numbers as α is increased. At higher rotation rates, the flow becomes increasingly unstable to perturbations at Re = 200 in the range 3 α 5 [START_REF] Meena | Three dimensional instabilities in flow past a spinning and translating cylinder[END_REF]. Recent numerical and experimental investigations [START_REF] Mittal | Flow past a rotating cylinder[END_REF]Rao et al. 2013b,a;[START_REF] Radi | Experimental evidence of new three-dimensional modes in the wake of a rotating cylinder[END_REF]) have identified several new three-dimensional transitions for Re < 400. In the Mode I shedding regime, five three-dimensional modes were found to become unstable and, in the steady regime of flow α 2, four three-dimensional modes were observed.

When the presence of a wall is considered, [START_REF] Taneda | Experimental investigation of vortex streets[END_REF] showed that a stationary wall near a cylinder stabilises the flow. [START_REF] Arnal | Vortex shedding from a bluff body adjacent to a plane sliding wall[END_REF] found that the onset of unsteady flow in the presence of a wall is shifted to Re ≈ 100, from Re ≈ 46 for an isolated cylinder in a free stream. This finding is correct as long as the gap ratio between the cylinder and the wall does not exceed a certain critical value, which depends on the Reynolds number. The steady flow around a stationary cylinder on a wall is similar to that of a backward-facing step [START_REF] Armaly | Experimental and theoretical investigation of backward-facing step flow[END_REF]; it is characterised by a single recirculation region, surrounded by fluid which separates from the body and reattaches on the wall downstream. However, the behaviour of the flow depends on many parameters: the distance between the cylinder and the wall, the motion of the wall relative to the cylinder, and the imposed rotation of the body. [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF] investigated the case of a rotating cylinder translating next to a moving wall at rotation rates in the range -1 < α < 1. In the steady flow regime, two recirculation zones are observed in the wake, and the drag and lift forces decrease as Re increases. Prograde rolling (α > 0) was found to destabilise the flow whereas retrograde rotation (α < 0) delayed the onset of unsteady flow. This was later confirmed/extended in a study by [START_REF] Rao | Flows past rotating cylinders next to a wall[END_REF] for which higher values of α were considered.

For unsteady flows in the wake of a cylinder placed near a wall, the strength of the vortex shedding decreases as the cylinder is placed closer to the wall [START_REF] Lei | Re-examination of the effect of a plane boundary on force and vortex shedding of a circular cylinder[END_REF]. At Re = 170, [START_REF] Taneda | Experimental investigation of vortex streets[END_REF] observed a single row of vortices for a cylinder moving near a wall. In the unsteady regime, vortex pairs with a net rotation appear in the wake as a result of the interaction between the shear layer shed from the top of the cylinder and induced secondary shedding from the wall boundary layer vorticity downstream [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF][START_REF] Rao | Flows past rotating cylinders next to a wall[END_REF]. Unlike the case of a cylinder placed in a free stream but similar to the flow over a backward facing step, the wake undergoes a transition to threedimensionality before the onset to two-dimensional vortex shedding [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF][START_REF] Rao | Flows past rotating cylinders next to a wall[END_REF].

In the current study, the results of [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF] for a cylinder rolling at α = 1 without slipping on a wall are extended to identify and characterise the different saturated flow states and modes of the three-dimensional instability, and their influence on the underlying two-dimensional structure of the flow. A more detailed description of the problem, the corresponding equations and the numerical methods are given in the first section of this paper. Then follow the results from the linear stability analysis and Floquet analysis, highlighting the appearance of new modes of the three-dimensional instability. Thorough comparisons of the similarities and differences in the flow structures and fluid forces show how the two-dimensional simulations can be used confidently to approximate this problem during the initial stage of flow development. The essential findings of this work are then summarised in the last section, along with some discussion of future research extending from this paper.

Problem description and methodology

Figure 1 illustrates the problem setup and parameter definitions: a cylinder of diameter d is rolling along a wall at a rotation rate ω. For computational simplicity, the frame of reference is placed at the centre of the cylinder, this being equivalent to the fluid and wall moving past the fixed, rotating cylinder at a speed U .

Background to new stability analysis studies

Consider a stationary cylinder placed in a freestream. In the steady, laminar regime below the critical Reynolds number Re c,2D = 46 at which transition to periodic shedding begins, it is well known that two mirror-symmetrical recirculation regions form in the wake whose length increases with the Reynolds number. For Reynolds numbers above the critical value, the steady wake is unstable, causing a decrease in the mean formation length of the recirculation region of the time averaged-flow (Williamson 1996a,b). Experiments and numerical simulations show that, for a given Reynolds number after the background flow is impulsively started, the evolution of the recirculation region is characterised by a linear, steady increase of its length followed by slowly growing waviness, and eventually, a fully developed Bénard-von Kármán vortex street. Thus in this case, the growth of perturbations leading to a fully developed periodic wake occurs after the steady symmetric wake state has fully (or almost fully) formed (e.g., see [START_REF] Thompson | The Stuart-Landau model applied to wake transition revisited[END_REF]. Moreover, the two-dimensional steady wake is not unstable to three-dimensional perturbations; the first three-dimensional transition occurs on the two-dimensional periodic wake [START_REF] Williamson | The existence of two stages in the transition to three-dimensionality of a cylinder wake[END_REF].

Conversely, when the cylinder is uniformly rolling along a wall, two-dimensional numerical simulations indicate that well below the critical Reynolds number for sustained two-dimensional vortex shedding (Re c,2D = 88), vortex shedding still occurs almost from startup. This can be seen from figure 2, which shows the time evolution of the lift coefficient for Re = 90 and 160, above the critical Reynolds number. The oscillations in the curves, observed from the initial starting time t 0 = 0, are the result of the immediate shedding of vortices into the wake. For Re = 90, it takes approximately six shedding cycles to reach the asymptotic periodic state, whilst for Re = 160, well beyond the critical Reynolds number, the transition to the fully developed 2D periodic state is essentially complete after just two shedding periods. Interestingly, even when the flow is laminar below Re c,2D , it undergoes immediate vortex shedding in the first instances of the flow development before settling down to its steady state. This effect was also made visible by [START_REF] Gal | Visualization of the space-time impulse response of the subcritical wake of a cylinder[END_REF] who studied experimentally the impulse response of the subcritical wake of a cylinder. These preliminary oscillations are visible on the plots of the lift force on figure 2 for Re = 60 and 80.

Thus, while [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF] documented that the onset of three-dimensional flow occurs at a Reynolds number of approximately 40 on a steady recirculating base flow, the actual flow transitions and dynamics observed in practice, for a cylinder rolling along a surface at constant velocity from an impulsive start, may be different. In particular, above the critical Reynolds number for 2D vortex shedding, a two-dimensional periodic wake seems likely to essentially fully develop prior to any observable development of threedimensionality. This hypothesis will be tested in the following sections. It has several consequences. It suggests that two-dimensional simulations have validity well beyond the critical Reynolds number at which three-dimensionality first occurs for the steady wake, at least for a non-negligible period after impulsive startup. It also suggests that the previous-documented three-dimensional stability analysis [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF] based on a steady two-dimensional base flow may not have a strong relevance for the overall wake dynamics for Reynolds numbers above Re c,2D . Finally, it is not clear what this means for the fully saturated wake state at different Reynolds numbers. Thus, in this paper, the stability analysis is extended to examine the three-dimensional stability of the two-dimensional periodic wake state. This is supplemented by direct three-dimensional simulations to examine the longer-term wake evolution.

Governing equations

The governing equations are the continuity and incompressible Navier-Stokes equations for the motion of the fluid. Let u(x, y, z, t) = (u, v, w) represent the velocity of the fluid in Cartesian coordinates. In the case of an incompressible flow, the continuity equation is:

∇ • u = 0, (2.1)
and the general form of the incompressible Navier-Stokes equation is:

∂u ∂t + u • ∇u = - 1 ρ ∇P + ν∇ 2 u, (2.2)
where ρ is the density of the fluid and P is the static pressure. The drag and lift coefficients per unit length are defined in the usual way:

C D = D 1 2 ρU 2 d and C L = L 1 2 ρU 2 d , (2.3)
where D and L represent the drag and lift forces, respectively.

Numerical scheme

The solver is based on a code that has been used extensively for similar problems, so it will only be described briefly here. Overall, the time-dependent incompressible Navier-Stokes equations for the fluid are solved in Cartesian coordinates using a spectral-element approach, in a cross-sectional plane for two-dimensional simulations or a combined spectral-element/Fourier approach in three dimensions. The spectral-element method is a formulation of a high-order finite-element method that uses high-order Lagrangian interpolants to approximate the solutions of partial-differential equations. It has the advantage of converging much faster than a typical finite-element method, considering that the error decreases exponentially with the order of the approximating polynomial, all the while retaining some of the flexibility for modelling complex geometries that finite-element methods provide. The (nodal) approach adopted is described in detail in [START_REF] Karniadakis | Spectral/HP element methods for CFD[END_REF]. The spatially discretised equations are then integrated forward in time using a three-step time-splitting approach, where the advection, pressure and diffusion terms are treated separately and sequentially [START_REF] Chorin | Numerical solution of the Navier-Stokes equations[END_REF][START_REF] Karniadakis | High-order splitting methods for the incompressible Navier-Stokes equations[END_REF][START_REF] Thompson | Hydrodynamics of a particle impact on a wall[END_REF]. The advection step is carried out using the third-order Adams-Bashforth approach. The pressure and viscous substeps are both solved directly using LU decomposition, which factors the matrices into a lower triangular matrix L and an upper triangular one U [START_REF] Turing | Rounding-off errors in matrix processes[END_REF], invoking the second-order Adams-Moulton (Crank-Nicholson) approximation for the linear viscous step. Whilst higher-order timestepping methods could be employed, because typically several hundred time-steps are required per shedding cycle to guarantee stability of the iterative approach, there is no discernible improvement in accuracy of the overall solution. The solver is explained in more detail by [START_REF] Thompson | Hydrodynamics of a particle impact on a wall[END_REF], and has widely been tested, validated and used for studies of flows around bluff bodies such as cylinders (Thompson et al. 2001b;[START_REF] Ryan | Three-dimensional transition in the wake of bluff elongated cylinders[END_REF][START_REF] Rao | Flows past rotating cylinders next to a wall[END_REF] and spheres (Thompson et al. 2001a;[START_REF] Rao | Transition to chaos in the wake of a rolling sphere[END_REF][START_REF] Thompson | Hydrodynamics of a particle impact on a wall[END_REF]. This code has also been modified to determine the linear stability of base flows, as explained in section 2.4 below.

In addition to the time-dependent solver, a steady solver is required to produce the steady base flows for the linear stability analysis. This is a modified version of the spectralelement code based on the penalty formulation (see [START_REF] Zienkiewicz | The Finite Element Method[END_REF]. This has been validated for a number of similar flow problems (e.g., [START_REF] Jones | A study of the geometry and parameter dependence of vortex breakdown[END_REF][START_REF] Thompson | The sensitivity of steady vortex breakdown bubbles in confined cylinder flows to rotating lid misalignment[END_REF].

Linear stability analysis

Linear stability analysis was undertaken in order to quantify flow transitions leading to vortex shedding, and to three-dimensional flow. For a two-dimensional steady or periodic base flow U, the velocity and pressure perturbation fields (u , P ) satisfy the continuity and linearised Navier-Stokes equations:

∇ • u = 0, (2.4) ∂u ∂t + U • ∇u + u • ∇U = - 1 ρ ∇P + ν∇ 2 u . (2.5)
Because the coefficients are independent of z, the perturbation fields can be further decomposed, representing z-dependence of variables as the sum of complex terms of a Fourier expansion:

u (x, y, z, t) = k ûk (x, y, t) sin(2πz/λ k ), (2.6) v (x, y, z, t) = k vk (x, y, t) sin(2πz/λ k ), (2.7) w (x, y, z, t) = k ŵk (x, y, t) cos(2πz/λ k ), (2.8) P (x, y, z, t) = k Pk (x, y, t) sin(2πz/λ k ).
(2.9)

Using these expansions, equations (2.4) and (2.5) give for each of the Fourier modes:

∂ ûk ∂t = -ûk ∂U ∂x + vk ∂U ∂y + U ∂ ûk ∂x + V ∂ ûk ∂y - 1 ρ ∂ Pk ∂x +ν ∂ 2 ûk ∂x 2 + ∂ 2 ûk ∂y 2 -(2π/λ k ) 2 ûk , (2.10) ∂v k ∂t = -ûk ∂V ∂x + vk ∂V ∂y + U ∂v k ∂x + V ∂v k ∂y - 1 ρ ∂ Pk ∂y +ν ∂ 2 vk ∂x 2 + ∂ 2 vk ∂y 2 -(2π/λ k ) 2 vk , (2.11) ∂ ŵk ∂t = -U ∂ ŵk ∂x + V ∂ ŵk ∂y -(2π/λ k ) 1 ρ Pk +ν ∂ 2 ŵk ∂x 2 + ∂ 2 ŵk ∂y 2 -(2π/λ k ) 2 ŵk , (2.12) ∂ ûk ∂x + ∂v k ∂y -(2π/λ k ) ŵk = 0. (2.13)
These perturbation field modes (û k , vk , ŵk , Pk ) can be further expressed as a sum of eigenmodes, each with its own growth rate. After choosing a wavelength and integrating these equations from initially random fields for sufficient time, the velocity perturbation fields will be dominated by the eigenmodes with the largest growth rates. Using a Krylov subspace together with Arnoldi decomposition allows a sequence of evolved fields to be decomposed into the dominant eigenmodes together with their corresponding growth rates (e.g., [START_REF] Mamun | Asymmetry and Hopf-bifurcation in spherical Couette flow[END_REF][START_REF] Barkley | Three-dimensional Floquet stability analysis of the wake of a circular cylinder[END_REF]. If  represents any of the perturbation fields ûk , vk , ŵk and Pk , then this method gives Â(x, y, t + T ) = Â(x, y, t) exp(σT ), if the eigenmode spatial distribution is not a function of time. Here σ is the growth rate and T is the time interval over which the growth of the mode is recorded.

It is also possible to get pairs of eigenmodes that have complex conjugate growth rates, providing the possibility of solutions with a periodic component on top of the exponential time variation. These pairs can also be extracted directly from the Arnoldi decomposition together with the complex growth rates σ = σ r + iσ i . For a three-dimensional transition on a two-dimensional periodic base flow, the procedure is the same, with the sequence of perturbation fields forming the Krylov subspace taken at full baseflow period intervals T .

In that case, the approach is called Floquet analysis, and the growth of each eigenmode is often expressed as a Floquet multiplier µ = exp(σT ), i.e., the amplitude of the mode after it has evolved for one period relative to the initial state. Again, it is possible to have pairs of eigenmodes with complex conjugate Floquet multipliers that are resolved through Arnoldi decomposition. For the three-dimensional analysis, the eigenmodes depend on the selected spanwise wavelength λ k . If |µ| > 1 (or σ r > 0), then the perturbation field will be amplified over time, while |µ| < 1 (or equivalently σ r < 0) for all λ implies that any perturbation will decay and, hence, the flow is linearly stable. Transition occurs when |µ| = 1 or σ r = 0. For the three-dimensional case, this condition has to be tested for every possible spanwise wavelength. The (eigen)modes that are obtained with Floquet analysis reported in this paper comprise periods equal to that of the periodic two-dimensional base flow, twice the period (subharmonic modes) and modes with different periods (quasiperiodic modes). The first two cases are characterised by a growth rate σ that is real. The quasi-periodic modes have a period that isn't commensurate with the base flow period

Domain size and resolution studies

For this study, the Reynolds number range is restricted to be Re 300, covering flow transitions to vortex shedding and to three-dimensionality.

Two domain sizes were investigated to quantify blockage effects. All positions are nondimensional and are scaled by the diameter of the cylinder. The first mesh M 1 is shown on figure 3 and consists of 1472 macro-elements. The upstream, downstream and upper boundaries are positioned at (x 1,u , x 1,d , y 1 ) = (-25, 25, 50), respectively. The second mesh M 2, also shown on figure 3, consists of 1906 elements and has the dimensions (x 2,u , x 2,d , y 2 ) = (-50, 50, 100). Both meshes have increased resolution in the vicinity and downstream of the cylinder that is located at x = 0 and y = 0. The differences in the drag and lift forces do not exceed 1% at the highest Reynolds number considered (Re = 300).

To ensure that the solution is converged with the chosen timestep ∆t = 0.0030, the latter was halved to 0.0015. This produced a variation in the body forces of less than 0.2% at Re = 150.

A last resolution study was carried out by increasing the number of the internal node points within each macro-element from N = 4 (×4) to N = 5 (×5), N = 6 (×6) and N = 7 (×7), which is taken as the reference value. We found that the drag force, the lift force and the period of oscillation differ by less than 1% respectively for N 5 at Re = 150. At the highest Reynolds number Re = 300, the drag and Strouhal number are well within the 0.5% at N = 5 (×5) whereas the error in the lift force reached 2.5%. From these results, and, considering that this study focuses mostly on a range of Reynolds number around the transition values Re c,3D and Re c,2D and up to 160, we can safely conclude that the mesh with 5 nodes per macro-element is converged and use it throughout our simulations.

The point of contact between the cylinder and the wall leads to a mesh singularity, and hence a small gap G is imposed between the cylinder and the wall. It has been shown previously that the flow structures visualised in the experiments and those observed numerically are in good agreement with G/d = 0.005 [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF][START_REF] Stewart | Flow dynamics and forces associated with a cylinder rolling along a wall[END_REF][START_REF] Rao | Flows past rotating cylinders next to a wall[END_REF], while reducing the gap ratio has little effect on flow quantities of interest [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF]) Thus a gap ratio G/d of 0.005 was used throughout this investigation.

For the full three-dimensional non-linear simulations, a Fourier expansion is used to represent the spanwise dependence of the flow variables (e.g., see [START_REF] Karniadakis | Three-dimensional dynamics and transition to turbulence in the wake of bluff objects[END_REF][START_REF] Thompson | Three-dimensional instabilities in the wake of a circular cylinder[END_REF][START_REF] Karniadakis | Spectral/HP element methods for CFD[END_REF]. The spanwise dimension is set to 54d for simulations for Re Re c,2D and 25.5d for Re < Re c,2D . These domains were chosen to fit approximately three wavelengths of the longest wavelength instability mode predicted by Floquet analysis (λ ∼ 8.5d and 18d for the steady and periodic base states). The choice of such large spanwise domains was to allow non-linear interactions between all relevant modes as the wake evolves to a fully saturated state. Another option would be to restrict the domain to wavelengths corresponding to each dominant mode to look at the super-/sub-critical nature of the each transition from a linear to a non-linear state (but not the fully saturated state). This could have been done, but it is not clear it would contribute much to the physical picture. For instance, for the Mode A transition of a circular cylinder in freestream, beyond saturation the wake evolves to allow dislocations (e.g., Williamson 1996a,b) that cannot be represented on a single wavelength spanwise domain. Also relevant, [START_REF] Karniadakis | Three-dimensional dynamics and transition to turbulence in the wake of bluff objects[END_REF] used a spanwise width that only allowed Mode B to grow. Because of this unphysical restriction, the authors observed period-doubling as the route to a fully-turbulent flow. However, this does not seem to be the situation for the real wake, or for computations using a sufficiently wide spanwise domain [START_REF] Henderson | Nonlinear dynamics and pattern formation in turbulent wake transition[END_REF]). Thus, it was decided to use a wider domain that would not put unphysical restrictions on mode development. Typically 48 and 96 Fourier modes are used for these simulations for the steady and periodic regimes, respectively. Since the shortest wavelength mode corresponds to λ 2.5d, this corresponds to approximately 10 Fourier planes to resolve the smallest important scales that develop in the wake. Tests with 144 Fourier modes confirm that this resolution accurately captures the wake evolution for the Reynolds number range considered.

Table 1 reports the results from the different resolution studies mentioned above.

Results

Linear stability analysis and flow transitions

Experimental studies by [START_REF] Taneda | Experimental investigation of vortex streets[END_REF] showed that the presence of the wall has a stabilising effect on the flow as long as the gap ratio does not exceed a certain critical value [START_REF] Lei | Re-examination of the effect of a plane boundary on force and vortex shedding of a circular cylinder[END_REF]. [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF] investigated the wake behind rolling cylinders at various rotation rates α, and found that as α varies from prograde (α > 0) to retrograde (α < 0) rolling, the critical Reynolds numbers for three-dimensional (Re c,3D ) and unsteady (Re c,2D ) transitions both increase. For comparison with that previous study, these critical Reynolds numbers were again predicted for the reference case of α = 1 (pure rolling without slipping). These transitions and the resulting flow states are depicted in figure 4. As indicated above, in contrast to the situation of a cylinder placed in an unbounded flow, the transition to three-dimensional flow for a cylinder placed near a wall occurs directly from a steady two-dimensional flow, similar to the situation for a backward- facing step [START_REF] Barkley | Three-dimensional instability in flow over a backward-facing step[END_REF]. As with that flow, its perturbation mode takes the form of periodic cells evenly distributed along the span of the cylinder. Figure 5 shows the perturbation velocity field projected into the plane just touching the top of the cylinder for Re = 60 and λ/d = 11 taken when the three-dimensional instability begins. This clearly shows the rotating cells associated with this instability mode. Stability analysis on the steady base flow shows that the onset of three-dimensional flow first occurs at Re c,3D = 36.8 for a spanwise wavelength λ c of 8.6d. Figure 6 shows that the maximum observed growth rate saturates by Re ∼ 60, only increasing further beyond Re ∼ 150, well beyond the onset of shedding. This is further confirmed in figure 7, where the growth rate is reported over the entire range of Reynolds numbers and spanwise wavelengths. The white regions in this figure represent negative growth rates and therefore stable wakes. Growth rates that are not real and are instead composed of a complex conjugate pair (i.e. the period of the mode is different from that of the base flow) were only detected in the blue region comprised between 90 Re 150 and 5 λ/d 8. The base flows for this analysis were generated using a steady version of the spectral-element code, allowing the stability to be investigated well beyond the transition to unsteady flow. The preferred wavelength can be seen to increase from λ = 8.6d at onset to reach values in excess of 20d at Re = 150, before suddenly dropping back on further increasing the Reynolds number. This is associated with the growth rate versus wavelength curves for Re 150 developing two peaks, with the lower wavelength peak dominating for Re 180.

Transition to two-dimensional unsteady flow

The transition from steady to unsteady two-dimensional flow occurs when the recirculation bubble at the rear of the cylinder becomes unstable and starts to shed vortices (figure 4b). These vortices interact with the wall through the no-slip condition, generating secondary vorticity as they advect downstream. In turn, this secondary vorticity is pulled away from the wall to combine with the primary generating vortex to form a vortex pair. The self-induced velocity of the pair causes it to move upwards and away from the wall, but because the primary vortex is stronger, the movement is along a curved path. The essential features of this process can be seen in figure 4b.

Figure 8 shows the form of the instability mode visualised by the perturbation spanwise vorticity, indicating that the mode has large amplitude where the base vorticity field is strong as well as close to the ground plane. Linear stability analysis of the steady base flow shows that this transition, which is characterised by a Hopf bifurcation, occurs at Re c,2D

88. This is shown in figure 9, which gives the growth rate and the preferred oscillation frequency as a function of the Reynolds number. Although not shown in the paper, it was verified that the fluctuating lift oscillation amplitude varied as (Re -Re c,2D ) close to the transition, as expected for a Hopf bifurcation. The oscillation frequency decreases with Reynolds number from St 0.066 at onset to 0.054 at Re = 200. Interestingly, the frequency of the fully saturated 2D wake stays relatively close to the perturbation mode frequency over this entire range. This is perhaps surprising given that the saturated periodic state deviates considerably from the steady wake base state, but it is probably an indication that the frequency selection is based on the separating shearlayer properties rather than the near-wake field. Discussion on frequency selection for the related case of a cylinder in freestream can be found in [START_REF] Pier | On the frequency selection of finite-amplitude vortex shedding in the cylinder wake[END_REF]; [START_REF] Barkley | Linear analysis of the cylinder wake mean flow[END_REF]; [START_REF] Sipp | Global stability of base and mean flows: a general approach and its applications to cylinder and open cavity flows[END_REF]; [START_REF] Leontini | A numerical study of global frequency selection in the time-mean wake of a circular cylinder[END_REF].

Stability of the fully developed 2D periodic wake

Given the discussion in section 2.1, it seems likely that the initial three-dimensional development of the wake at a particular Reynolds number above Re c,U = 88 will be determined by the stability of the 2D periodic flow. That will be tested in later sections through direct simulations. In this section, the 3D linear stability of the 2D periodic base flow is characterised first. Beyond the transition to unsteady flow, a number of different modes contribute to the wake becoming three-dimensional. The occurrence and growth rates of these modes are also strongly dependent on the Reynolds number, presumably because the structure of the time-dependent two-dimensional wake is also a strong function of the Reynolds number. Figure 10 summarises the situation by showing the growth rate corresponding to the dominant mode as a contour map over a wide range of spanwise wavelengths and

F. Y. Houdroge et al

Reynolds numbers beyond Re c,2D

. There are a few regions of substantial growth, notably corresponding to λ/d ∼ 4 and 8-9 covering different Reynolds number ranges. The picture is a little more complicated than indicated by this map, with local peaks corresponding to different mode types: synchronous modes (i.e., with the same period as the base flow); subharmonic modes (with twice the base period); and quasi-periodic modes (with periods different from the base period). To show this in more detail, Floquet multiplier variations as a function of wavelength at three different Reynolds numbers Re = 100, 130 and 160 are given in figure 11. Just above the transition at Re = 100, the fastest growing mode, marked (1) in the figure, reaches a maximum growth at λ 3.2d. This corresponds to a real or synchronous Floquet mode, i.e., the period is the same as the base flow period. At this Reynolds number, there are several other contributing modes with positive growth rates covering the wavelength range studied: another real mode (2) with a wavelength of λ ∼ 6-7d, a subharmonic mode (3) starting at λ 8d, and a quasi-periodic mode (4) for λ 15d. These all have strongly positive growth rates, although less than the short wavelength mode (1). At Re = 130, the short wavelength mode (1) is still present, but now the most dominant mode is the subharmonic mode (2) at a wavelength of λ 9d. The longer wavelength quasi-periodic mode is still present, although it gives way to another real mode at still higher wavelengths (λ 25d). At Re = 160, there are again several changes to the picture. The short wavelength mode (1) and subharmonic mode (3) are still present, with the subharmonic mode relatively much more dominant. At higher wavelengths (λ 15), a real mode (5) becomes more dominant than the quasi-periodic mode (4) over that wavelength range. The fact that all these modes are amplified and they cover a wide wavelength range suggests that the wake is likely to become chaotic quickly after the initial growth of the most dominant mode begins to saturate. (The evolution of the modal amplitude until saturation is shown in figure 20, as discussed later). This is investigated further using direct numerical simulations in the following section, but prior to this, the vorticity structure of the modes is examined.

Figure 12 shows the evolution of the perturbation spanwise vorticity field for mode (1) at Re = 100, where it is the fastest growing mode, and at Re = 160, where it is less dominant. Especially in the lower Reynolds number case, the structure of the perturbation field inside the newly formed and shed vortex cores clearly shows the characteristics of elliptical instability [START_REF] Bayly | Three-dimensional instability of elliptical flow[END_REF][START_REF] Pierrehumbert | Universal short-wavelength instability of two-dimensional eddies in an inviscid fluid[END_REF][START_REF] Landman | The three-dimensional instability of strained vortices in a viscous fluid[END_REF][START_REF] Waleffe | On the three-dimensional instability of strained vortices[END_REF][START_REF] Leweke | Three-dimensional instabilities in wake transition[END_REF]Thompson et al. 2001b;[START_REF] Kerswell | Elliptical instability[END_REF]. In particular, the perturbation vorticity shows two lobes of positive and negative vorticity, whose extrema align at ∼ 45 • to the main axes of the elliptically shaped vortex cores (marking regions with elliptic streamlines in the reference frame moving with the advection velocity at the centres of these cores. Also importantly, the orientation of the lobes is approximately maintained as the vortices advect downstream, allowing the perturbation to grow and allowing feedback from one shedding cycle to the next. Although somewhat far from the idealised cases for which the theory was developed [START_REF] Waleffe | On the three-dimensional instability of strained vortices[END_REF], for finite-sized vortices, the preferred wavelength is dependent on the core size. Le [START_REF] Dizès | Theoretical predictions for the elliptic instability in a twovortex flow[END_REF] showed that for Gaussian vortices under strain, the spanwise wavelength is given by λ = 2.78a, with a the Gaussian length scale. For highly strained vortices, such as is the case here, the appropriate length scale (a) is given by Le [START_REF] Dizès | Viscous interaction of two co-rotating vortices before merging[END_REF] as a 2 = (a 2 M + a 2 m )/2, with a M and a m corresponding to the semi-major and minor axis lengths, respectively. At Re = 100, figure 12 shows that the approximately invariant vorticity tube grows in size as the vortex cores advect downstream. For the first image, where the elliptical instability pattern is first recognisable, the length scales from the just formed and downstream cores obtained by fitting Gaussian profiles to the major and minor axes are 0.95d and 1.27d, giving preferred spanwise wavelengths of 2.7 and 3.5d, respectively. Figure 11 shows that Floquet analysis indicates that the maximum growth for this mode corresponds to a wavelength of λ = 3.2d, near the centre of the range of the theoretical prediction. At Re = 160, Floquet analysis shows the preferred wavelength of mode (1) drops to 2.5d. This is in line with the prediction of the more compact shed cores due to lower viscous diffusion at the higher Reynolds number. In this case, an asymmetrical counter-rotating vortex pair is formed before the newly formed vortex advects very far downstream. This composite structure also shows evidence of a perturbation pattern consistent with elliptic instability, as it advects away from the wall in an approximately circular arc. Various studies have identified elliptic instability in an isolated counter-rotating vortex pair (e.g., [START_REF] Leweke | Cooperative elliptic instability of a vortex pair[END_REF].

Figure 13 shows the perturbation streamwise vorticity structure of modes ( 3) and ( 5) at Re = 160. Mode (3), the shorter wavelength mode, is subharmonic, repeating every two base flow periods. Mode (5) only becomes dominant for Re 160. Below this Reynolds number, a quasi-periodic mode occupying this wavelength range has a higher growth rate. Interestingly, the amplitude distributions of these two modes appear similar. Inside the newly forming vortex, the distributions broadly match between the two modes, in terms of both the sign and distribution of perturbation vorticity. The vorticity distributions within the downstream vortex pairs are also similar, but of opposite sign, as is the case with the layer of vorticity near the ground between the two main vortical structures.

The physical nature of the instability in this case is more difficult to discern. The developing wake does not form a series of relatively long-lived elliptical-shaped vortices in this case, but rather a newly formed vortex generates secondary vorticity beneath it, which is subsequently drawn away from the boundary to form an unequal strength counter-rotating vortex pair. The process is shown in figure 14. This indicates that there are a number of different identifiable vortex structures and groupings that combine to lead to the observed flow instability modes.

The evolution of the circulation in the primary and secondary vortices as they advect downstream is shown in figure 15. The primary (clockwise) vortices that directly form from the separating shear layer, grow quickly in strength prior to the shear layer pinching off and releasing the vortices to move downstream. The circulation then slowly decays. During its initial growth and soon after its release, this primary vortex combines with the secondary vorticity generated at the boundary, to form the vortex pair. At a point in time when this pair becomes unambiguously defined, i.e., between images 2 and 3 of figure 14, the ratio of circulations between the component vortices of the pair is approximately -2:1. [START_REF] So | Short-wave instabilities on a vortex pair of unequal strength circulation ratio[END_REF] analysed the stability of unequal strength Lamb-Oseen vortices, which have a Gaussian vorticity distribution, examining how the growth rate and preferred wavelength varied with circulation ratio. Such a system is subject to both the shortwavelength elliptic instability and the longer-wavelength Crow instability. The results of [START_REF] So | Short-wave instabilities on a vortex pair of unequal strength circulation ratio[END_REF] can be used to obtain an estimate for the most unstable wavelength. A recent review of the Crow instability is given in [START_REF] Leweke | Dynamics and instabilities of vortex pairs[END_REF]. Examining the spanwise perturbation vorticity field of mode 5 (not shown) shows the characteristic lobe structure of Crow instability in the vortex pair as it moves away from the cylinder and the wall. By approximating the vorticity distributions within the pair at formation in terms of a sum of Gaussian distributions to extract the length scale for each vortex, together with the overall circulation ratio, a preferred wavelength of approximately λ 20d can be predicted from the work of [START_REF] So | Short-wave instabilities on a vortex pair of unequal strength circulation ratio[END_REF]. This is close to the observed preferred wavelength of mode 5 at Re = 160 of λ 18d. However, in this case, the perturbation field does not grow as the vortex pair advects, so the Crow instability of the pair alone cannot be responsible for the maintenance of the instability mode from one cycle to the During the formation and evolution of the wake vortices, it is also possible to take into account the image vortices, linked to the presence of the wall and symmetrically located with respect to it. The near wake vortex pair and its image form a symmetric four-vortex system, a configuration analyzed previously in the context of aircraft trailing wakes (e.g., [START_REF] Crouch | Airplane trailing vortices and their control[END_REF][START_REF] Winckelmans | Vortex methods and their application to trailing wake vortex simulations[END_REF][START_REF] Jacquin | Unsteadiness, instability and turbulence in trailing vortices[END_REF]. The existence of short-wave (elliptic) and Crow-type long-wave instabilities was also found in these systems. Although the identification of such systems is transitory for this wake, is seems plausible that the Crow instability would play a role in the global three-dimensional instability and wavelength selection for the wake.

Saturated three-dimensional state

Computed forces

In this section the force coefficients after the flow has reached its fully saturated threedimensional state are compared with predicted force coefficients from 2D steady and periodic simulations. The time-mean forces computed obtained from 3D direct numerical simulations are in good agreement with the two-dimensional ones. Figure 16 shows plots of the mean drag and lift coefficients, as defined in equations (2.3), versus the Reynolds number. For Reynolds numbers up to the Hopf bifurcation leading to a periodic two-dimensional state, the 2D and 3D curves are effectively indistinguishable. This is consistent with the relatively weak effect on the wake of the steady three-dimensional instability even as it saturates, as shown in figure 17. However, even beyond the transition to periodic flow, the difference between the 2D periodic and 3D predictions remains small, and is limited to be less than 5% for the mean drag coefficient and 4% for the lift coefficient at the highest Reynolds number considered here of Re = 180. In this case, the saturated three-dimensional wake is distinctly different from the two-dimensional periodic wake, as is explored further below. The figure also shows the lift and drag coefficients based on the steady flow for Re > Re c,2D obtained from the steady solver. These curves deviate considerably from the other two sets as the Reynolds number increases. This is consistent with the increasingly elongated recirculation region of the steady flow deviating further from the near-wake vortex shedding of the 2D-periodic and 3D flows as the Reynolds number is increased.

Figure 18 shows the temporal evolution of the lift coefficient obtained from full threedimensional simulations at Re = 60, 80, 90 and 160. The initial evolution follows the one observed from the two-dimensional simulations (figure 2) until the three-dimensional instability grows sufficiently to change the two-dimensional structure of the flow. This effect can be seen in the temporal development of the lift coefficient: below the Hopf bifurcation at Re c,2D = 88, the three-dimensional transition disturbs the otherwise constant lift coefficient at approximately t = 400 -600U/d (upper two plots), whilst beyond the 2D transition, the periodic oscillations in the curves die out at approximately t = 200 -400U/d (bottom two plots) as a result of the three-dimensional instability reaching a sufficient amplitude to substantially alter the otherwise two-dimensional periodic flow.

It is of interest why the oscillations in the lift signal are substantially suppressed once the wake reaches its saturated three-dimensional state. To investigate this, sectional lift coefficient signals, i.e., the lift coefficient per unit span at a particular spanwise position, were examined at different points across the span. Figure 19 shows the evolution of the lift signals at two points separated by half the span width (dashed lines) together with the mean lift signal (solid line), for Re = 160 in the saturated state. Clearly, there is a significant variation in the local lift coefficient across the span, indicating that the underlying two-dimensional vortex shedding is uncorrelated. In addition to this, even the sectional lift coefficients aren't very periodic. This is consistent with a change from strong, regular 2D shedding of vortices initially to a much more disordered 3D wake without a strong underlying 2D periodic vortical wake structure. Note that for these simulations, a low-level white-noise perturbation of amplitude 10 -4 U was added to each velocity component at startup to accelerate the development of the three-dimensionality.

Development and saturation of three-dimensional flow

Within the steady regime (Re < Re c,2D ), the evolution to a fully evolved 3D wake for Re > Re c,3D leads to the two-dimensional spanwise vorticity isosurfaces becoming wavy in the spanwise direction, with little alteration to the main underlying two-dimensional structure of the flow. As identified above, this effect can be seen on figure 17(a), and the extent of this deformation on figure 17(b). These predictions are consistent with those made by [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF][START_REF] Stewart | The dynamics and stability of flows around rolling bluff bodies[END_REF], who in addition conducted experiments in a water channel. Their experiments showed that the results of the experimental streaklines and of the predicted two-dimensional flow are in good qualitative agreement, at least while the three-dimensionality is developing. Note that the final saturated state in this case is weakly unsteady. For instance, at Re = 45 there is a low-level oscillation leading to the weak shedding of vorticity into the wake, while the global 3D structure shown in figures 5 and 17 is maintained. For Reynolds numbers beyond the transition to vortex shedding, which is the main focus here, direct numerical simulations of the flow indicate that, prior to its settling into its final state, the flow initially undergoes a rapid transition to two-dimensional vortex shedding, as indicated by the lift trace curves of figure 18. The numerical method involves representing the spanwise variation through a Fourier decomposition, hence the evolution of the spanwise modes can be easily extracted. A convenient measure of the amplitude of each mode is provided by the RMS amplitude of the spanwise velocity component of each complex Fourier mode, since that velocity component is zero prior to three-dimensional flow development. Specifically, the evolution of the modal amplitudes (A n ), computed as a RMS spatial average over all 2D node points (N xy ) of the moduli of the z-velocity complex Fourier coefficients (a n ) corresponding to mode index n, i.e.,

A n (t) =   1 N xy Nxy i=1 |a i n (t)| 2   ,
are shown in figure 20 for Re = 100 and Re = 160. The two figures in the left column show the development of modes corresponding to key wavelengths identified by the global stability analysis. Indeed, after an initial period over which the dominant mode for each wavelength emerges, the growth rates as measured by the slopes of the curves over many oscillation periods have values consistent with the linear stability analysis predictions. At some point in time, the modes grow sufficiently to begin to saturate non-linearly, leading the flow to reach its asymptotic state. After saturation, the figures in the right column show that the final state is influenced by many modes of different wavelengths, suggesting a rapid transition to fully chaotic flow. Figure 21 shows time-mean RMS amplitudes of each Fourier mode taken over the last 100 time units (in the fully saturated state) at (a) Re = 100 and (b) Re = 160. Here, the horizontal axis is the non-dimensional wavenumber kd. These spectra can be compared with figure 21(c) at Re = 45, where the saturated state shows a single spectral peak corresponding to kd = 0.74 (or λ/d = 8.5) together with harmonics accounting for the distortion of the saturated final state from the pure sinusoidal linear instability mode. For the two higher Reynolds number cases, the spectra are continuous as a result of the nonlinear interactions between modes, and this is indicative of a chaotic final wake state. Indeed, the modes corresponding to the dominant linear instability mode numbers do not dominate the spectra at saturation. Isosurface plots taken at key points in the flow development help elucidate the wake development. Figure 22(a) shows an isosurface of Q = 0.01 at t = 350d/U . The Qcriterion is a vortex identification method defined initially by [START_REF] Hunt | Eddies, streams, and convergence zones in turbulent flows[END_REF]. This isosurface is merged with isosurfaces of positive/negative streamwise vorticity to highlight the dominant spanwise mode at a time when the three-dimensionality is beginning to modify the otherwise two-dimensional wake structure. In this case, Re = 100. The wavelength of the streamwise vorticity pattern extracted from this image is consistent with the short wavelength mode (1) instability prediction (figure 22a) from stability analysis. Soon after, the wake develops non-linearly, with a typical snapshot shown in figure 22(b).

At a higher Reynolds number of Re = 160, the development is somewhat different. Figure 23 shows a sequence of wake states from the time that three dimensionality is beginning to develop. The first three images show the evolution at three consecutive shedding cycles. The three-dimensionality develops quickly, with the spanwise wavelength of the perturbation corresponding to that of mode (3) of figure 11(c). The second image shows substantial distortion of the previously shed two-dimensional vortex pair, whilst the third image, one cycle later, shows that the subsequently shed vortex pair is virtually destroyed. The final image is a plot taken a few cycles later, after the final asymptotic state is reached. This is similar to the final state at Re = 100 shown in figure 22, except it shows an even more complex finer-scale structure. The previous two-dimensional periodic wake state is no longer visible at all.

Conclusions

The stability analysis of the steady two-dimensional flow past a cylinder rolling at constant speed along a rigid surface has shown that the key transitions to steady threedimensional flow and to periodic vortex shedding occur at Re c,3D = 36.8 and Re c,2D = 88, respectively. These results are mainly confirmation of findings from a previous study by [START_REF] Stewart | The wake behind a cylinder rolling on a wall at varying rotation rates[END_REF]. However, the main emphasis of this paper is concerned with three-dimensional wake development in a more realistic configuration, i.e. after a cylinder starts rolling impulsively at a constant velocity, and how this evolution is related to stability theory.

Two main cases can be distinguished. The first is when the Reynolds number is lower than the critical Reynolds number leading to two-dimensional transition (Re c,2D ) and above the critical Reynolds number for three-dimensional steady transition (Re c,3D ). In this case, the asymptotic flow state is a three-dimensional flow that is not too far from the prediction based on assuming two-dimensional flow. Indeed, the drag and lift coefficients are practically unaffected by the three-dimensionality. When the Reynolds number is close to Re c,2D , the initially stationary flow develops a few cycles of shedding prior to settling towards a steady state. The second case is observed for Reynolds numbers above Re c,2D . After an impulsive start, the flow undergoes a rapid transition to two-dimensional periodic shedding. Within a few cycles, e.g., about three at Re = 90 and two at Re = 160, the wake evolves to be close to the periodic state predicted by two-dimensional simulations. It then continues in this near-periodic state for several cycles, depending on the background noise level. For the cases considered here, this period of evolution was about 15 and 10 cycles at Re = 90 and 160, respectively. Stability analysis of the two-dimensional periodic state, which has not been undertaken previously, then determines the subsequent development of the wake three-dimensionality. At Re = 100, the wake appears to initially undergo a short wavelength instability (mode 1 shown in figure 8a), consistent with an elliptic instability of the shed vortex cores. At longer times, many more spanwise modes come into play and interact non-linearly, leading to a chaotic final flow state. At a higher Reynolds number of Re = 160, the initial development of three-dimensionality is different. Here, mode 3 of figure 8(c) is the mode to break two-dimensionality. Again, the wake undergoes a rapid transition to a chaotic final state soon afterwards. With the noise levels used to initiate the three-dimensional flow development in the three-dimensional simulations, the fully saturated wake states take approximately 400 and 200 non-dimensional time units to develop, for Re = 90 and 160, respectively. These values are equivalent to the number of diameters the cylinder rolls whilst maintaining a two-dimensional state. Although experimental noise levels are likely to be higher, it is still an indication that, after an impulsive start of the cylinder, a two-dimensional periodic flow state will be maintained for a considerable rolling distance prior to the evolution to a complex three-dimensional wake. The two-and three-dimensional simulations also show that the mean lift and drag coefficients of the fully saturated three-dimensional flow are very close to predictions based on two-dimensional simulations.

It is interesting to speculate whether a similar scenario would apply to a sphere rolling at constant velocity along a wall. In that case, even in freestream, a non-axisymmetric steady transition occurs prior to the periodic transition. The presence of the wall seems likely to cause the premature generation of shedding on impulsive startup, bypassing the slow transition associated with a Hopf bifurcation of a steady flow. However, we will leave this as an open question at this stage.

The numerical model here is essentially an infinite two-dimensional cylinder forced to roll at constant speed. End effects may play a strong part in the wake evolution, just as it can with the flow past a cylinder away from a boundary. Additionally, if the cylinder is free to roll without any constraints on its movement and velocity, vortex-induced vibrations are likely to occur with an unsteady wake. The simulations and results presented in this paper aim to provide a reference study for the idealised case, and constitute an essential element of an ongoing study concerning the fluid-structure interaction of uniformly and freely rolling bodies translating along a wall.

Figure 1 :

 1 Figure1: A two-dimensional schematic of the problem under consideration: a cylinder of diameter d is rolling along a wall. The translational and angular velocities are represented by U and ω, respectively. In the simulations, the inertial frame of reference (x, y) is attached to the centre of the body.

Figure 2 :

 2 Figure 2: Time evolution (scaled by d/U ) of the lift coefficient from impulsive startup. (a) & (b) Evolution for Reynolds number below the transition to 2D shedding. (c) & (d) Evolution for Reynolds numbers beyond the transition to vortex shedding.

Figure 3 :

 3 Figure 3: View of the cylinder mesh M 1: (x 1,u , x 1,d , y 1 ) = (-25, 25, 50) (left image) and M 2: (x 2,u , x 2,d , y 2 ) = (-50, 50, 100) (right image). The cylinder is placed in the middle of the x-axis, and near the wall at a gap ratio of G/d = 0.005 in order to avoid numerical singularities from arising. The flow is from left to right, and the resolution in the vicinity and downstream of the cylinder is increased in order to accurately capture the flow structures in the wake.

  (a) 2D steady (left, Re = 30) to 3D steady at (right, Re = 40). (b) 2D steady (left, Re = 30) to 2D unsteady (right, Re = 200).

Figure 4 :

 4 Figure 4: Visualisation of the initial flow transitions under consideration. The images showing 2D flow depict the spanwise vorticity field, whereas the 3D flow pattern is depicted using isosurfaces of streamwise vorticity. Red and blue represent positive and negative vorticity, respectively.

Figure 5 :

 5 Figure 5: Top down view of the three-dimensional steady flow at Re = 60 and λ/d = 11, depicted through the projected velocity field in the plane grazing the top of the cylinder.

Figure 6 :

 6 Figure 6: Left: (a) Maximum growth rate of the most unstable mode for transition from 2D steady to 3D steady flow. Right: (b) Variation of the wavelength of the fastest growing mode with the Reynolds number. Beyond Re ∼ 150 there are two peaks in the growth rate curve with the shorter wavelength peak developing the higher amplitude for Re 180.

Figure 7 :

 7 Figure 7: Contour map of the growth rate σ as a function of the Reynolds number Re and the wavelength λ/d for the 2D steady to 3D transition. For most of the domain, the dominant 3D mode is steady, except for a small (blue) region centred around Re = 120 and λ/d = 5 and extending to higher Reynolds numbers.

Figure 8 :Figure 9 :

 89 Figure 8: The structure of the perturbation field at Re = Re c,2D depicted using perturbation spanwise vorticity with overlaid base flow vorticity contours at ±0.1U/d.

Figure 10 :Figure 11 :

 1011 Figure 10: Contour map of the dominant growth rate σ as a function of the Reynolds number Re and the wavelength λ/d for the 2D periodic base flow.

Figure 12 :

 12 Figure 12: Evolution of the spanwise perturbation vorticity for the mode (1) of figure 11 at Re = 100 (left column) and Re = 160 (right column), overlaid with the base flow vorticity contours at ±0.1U/d. The images are 1/5 of a period apart. The spatial distribution of the instability fields show strong signs of elliptic instability of the vortex cores, as discussed in the text.

Figure 13 :

 13 Figure 13: Visualisations of the streamwise perturbation vorticity for subharmonic mode (3) (left) and real mode (5) (right) at Re = 160. These images depict streamwise perturbation vorticity coloured contours with base flow vorticity contour line at ±0.1d/U overlaid to highlight the locations of the vortices.

Figure 14 :

 14 Figure 14: Evolution of the spanwise wake vorticity at Re = 160 showing the formation of new vortices from the separating shear layer, generation of secondary vorticity at the boundary and release into the wake, and the formation of counter-rotating vortex pairs that self-propel away from the wake as the pair moves downstream. Images are separated in time by 1/5 of a period.

Γ

  Figure 15: Evolution of the circulation Γ of the primary and secondary vortices over a shedding cycle as they form and advect downstream. Here, Re = 160.

Figure 16 :

 16 Figure 16: Comparison of time-averaged body force coefficients, drag C D and lift C L , between the two-and three-dimensional simulations for 30 Re 180. Note that beyond Re c,2D the 2D-steady predictions shown by the filled circles are based on flows calculated with the steady solver.

Figure 17 :Figure 18 :

 1718 Figure17: The fully developed wake state at Re = 50 visualised over a spanwise distance of three characteristic wavelengths. The wake state is depicted using an isosurface of the spanwise vorticity (ω z = ±0.5).

Figure 19 :

 19 Figure19: Evolution of the sectional lift coefficient at two spanwise locations (dashed lines) separated by half the span width at Re = 160 in the fully saturated state. The spanwise averaged lift coefficient is also shown by the solid line. The period corresponding to two-dimensional vortex shedding is approximately 20 units.

Figure 20 :Figure 21 :

 2021 Figure 20: Evolution of the amplitude of spanwise Fourier modes at different Reynolds numbers. Top row: Re = 100; bottom row: Re = 160. The subfigures on the left show the evolution of the three dominant wavelengths as predicted by Floquet analysis, and those on the right show the evolution of the amplitudes corresponding to the first 48 modes.The simulations are started impulsively, with a low level white noise to accelerate the development of the three-dimensionality. Measured slopes of the evolution curves in the linear regime give estimated growth rates (σ) in agreement with growth rate predictions from Floquet analysis given in figure11.

Figure 22 :

 22 Figure 22: Evolution of the wake at Re = 100 from direct simulations. (a): Isosurface of Q = 0.01 (blue) highlighting the predominantly 2D vortices, with isosurfaces of streamwise vorticity at ω x = ±0.1 (red/yellow) superimposed. At this time (t = 350d/U ), the three-dimensional instability modes have grown sufficiently to begin to affect the wake. (b) & (c): Later, at t = 500d/U , after the flow has fully saturated. The wake is complex and chaotic with many three-dimensional wavelength components contributing. The iso-surface corresponding to Q = 0.01 is plotted on figure (b), and the isosurfaces of streamwise vorticity at ω x = ±0.1U/d on figure (c).

Figure 23 :

 23 Figure 23: Stages in the evolution of the wake at Re = 160 as depicted by isosurfaces of Q = 0.01. The first three images show the wake structure for three consecutive cycles just after the three-dimensionality is starting to modify the flow. The final image is a typical image after the flow has reached its asymptotic state.

Table 1 :

 1

	Mesh Re N × N	∆t	CD	CL	Period
	M1 150 4 × 4 0.0030 3.2827(-2.45) 1.48186(2.81) 17.525(-0.51)
	M1 150 5 × 5 0.0030 3.3631(-0.06) 1.45020(0.62) 17.612(-0.02)
	M1 150 6 × 6 0.0030	3.3650	1.44133	17.615
	M1 300 5 × 5 0.0020 3.4175(-0.12) 0.8532(2.47) 19.646(-0.13)
	M1 300 6 × 6 0.0020 3.4200(-0.05) 0.8364(0.46) 19.670(-0.01)
	M1 300 7 × 7 0.0020 3.4205(-0.03) 0.8325(-0.01) 19.671(-0.005)
	M1 300 8 × 8 0.0020	3.4216	0.8326	19.672
	M1 150 5 × 5 0.0015 3.3637(0.02) 1.45179(0.11) 17.611(-0.01)
	M2 150 5 × 5 0.0030 3.3977(1.03) 1.45282(0.18) 17.400(-1.20)

Domain, and temporal and spatial resolution study. The numbers in parentheses show the error relative to the highest resolution (or number of nodes N ) used for the comparison. The mesh M 1 has a blockage ratio of 1% and the mesh M 2 a blockage ratio of 2%. Most of the simulations were undertaken with mesh M 2 , noting typical blockage-induced errors in the Strouhal and force coefficients of ∼ 2%. The comparisons indicate that the timestepping error is negligible, whilst the error induced using 5 × 5 nodes/element is typically ∼ 1% for a Reynolds number of 150.

† Email address for correspondence: farah.houdroge@monash.edu

Acknowledgements

This research was supported under Australian Research Council, Discovery Projects funding scheme DP130100822 and DP150102879. We also acknowledge computing time support through National Computational Infrastructure projects D71 and N67.