# [

###### Abstract

Direct numerical simulations of Taylor-Couette flow (TC), i.e. the flow between two coaxial and independently rotating cylinders were performed. Shear Reynolds numbers of up to , corresponding to Taylor numbers of , were reached. Effective scaling laws for the torque are presented. The transition to the ultimate regime, in which asymptotic scaling laws (with logarithmic corrections) for the torque are expected to hold up to arbitrarily high driving, is analysed for different radius ratios, different aspect ratios and different rotation ratios. It is shown that the transition is approximately independent of the aspect- and rotation- ratios, but depends significantly on the radius-ratio. We furthermore calculate the local angular velocity profiles and visualize different flow regimes that depend both on the shearing of the flow, and the Coriolis force originating from the outer cylinder rotation. Two main regimes are distinguished, based on the magnitude of the Coriolis force, namely the co-rotating and weakly counter-rotating regime dominated by Rayleigh-unstable regions, and the strongly counter-rotating regime where a mixture of Rayleigh-stable and Rayleigh-unstable regions exist. Furthermore, an analogy between radius-ratio and outer-cylinder rotation is revealed, namely that smaller gaps behave like a wider gap with co-rotating cylinders, and that wider gaps behave like smaller gaps with weakly counter-rotating cylinders. Finally, the effect of the aspect ratio on the effective torque versus Taylor number scaling is analysed and it is shown that different branches of the torque-versus-Taylor relationships associated to different aspect ratios are found to cross within of the Reynolds number associated to the transition to the ultimate regime. The paper culminates in phase diagram in the inner vs outer Reynolds number parameter space and in the Taylor vs inverse Rossby number parameter space, which can be seen as the extension of the Andereck et al. (J. Fluid Mech. 164, 155-183, 1986) phase diagram towards the ultimate regime.

Exploring the phase diagram of fully turbulent Taylor-Couette flow]Exploring the phase diagram of fully turbulent Taylor-Couette flow R. Ostilla Mónico and others]RODOLFOOSTILLAMÓNICO,ERWINP.VANDERPOEL,ROBERTOVERZICCO,SIEGFRIEDGROSSMANN,ANDDETLEFLOHSE

Key words:

## 1 Introduction

Taylor-Couette flow (TC), i.e. the flow between two independently rotating concentric cylinders, has for long been used as a model system in fluid dynamics. Couette (1890) was the first to investigate it, and he pioneered its usage as a viscometer. But it was Mallock (1896) who, by rotating the inner cylinder, and not the outer as Couette had done, found the first indications of turbulence in the system. Taylor (1923, 1936) further studied the system, finding that it was linearly unstable, unlike pipe-flow and other studied systems to the date. Wendt (1933) expanded the study of the turbulent regime, measuring torques and velocities in the system. Since then, and due to its simplicity, TC has been used as a model system for studying shear flows. For a broader historical context, we refer the reader to Donnelly (1991).

Recently, a mathematically exact analogy between TC and Rayleigh-Bénard flow (RB), i.e. the convective flow between two parallel plates heated from below and cooled from above was found by Eckhardt, Grossmann and Lohse (2007), (here referred to as EGL07) . Within this context, TC can be viewed as a convective flow, driven by the shear between both cylinders where angular velocity is transported from the inner to the outer cylinder. As explained by Grossmann et al. (2014), as long as the driving of the system is small, the transport is limited by the laminar boundary layers. But if the driving becomes strong enough the boundary layers become turbulent and the system enters the so-called “ultimate” regime. The study of the transition to this regime, expected to be also present in RB, has attracted recent interest, as most applications of TC and RB in geo- and astro-physics are expected to be in this ultimate regime.

For RB flow, the transition to an ultimate regime was first qualitatively predicted by Kraichnan (1962), and later quantitatively by Grossmann & Lohse (2000, 2001, 2011) and then experimentally found by He et al. (2012b, a); Ahlers et al. (2012); Roche et al. (2010). It lies outside the present reach of DNS. The analogous boundary layer transition to an ultimate regime in TC flow was first found in the experiments by Lathrop et al. (1992b, a), and analysed more precisely in Lewis & Swinney (1999), even though earlier work by Wendt (1933) already showed some transition in the torque scaling around the same Reynolds number. The transition was not related to the transition to the ultimate regime until later (van Gils et al. 2011; Paoletti & Lathrop 2011; Huisman et al. 2012; Grossmann et al. 2014). In DNS, it was observed for the first time in (Ostilla-Monico et al. 2014b).

In TC flow this transition is easier to achieve as the mechanical driving is more efficient than the thermal one, and thus the frictional Reynolds numbers in the boundary layer are much larger. By using the analogy between both systems, better understanding of the transition in TC can thus also lead to new insight in RB, where it is more elusive.

Ostilla-Monico et al. (2014b) numerically studied the transitions in TC for pure inner cylinder rotation for a radius ratio of , where and are the outer and inner radii respectively, and an aspect ratio , where is the axial period in the DNS. In that study, the flow transitions and boundary layer dynamics were revealed in the range of Taylor numbers between and , where the Taylor number is defined as:

(1.0) |

with and the angular velocities of the outer and inner cylinder, respectively, the gap width, and the kinematic viscosity of the fluid. can be considered as a geometric quasi-Prandtl number (EGL07).

We now describe the series of events when increasing . For small enough , the flow is in the purely azimuthal, laminar, state. When the system is driven beyond a critical driving, one passes the onset of instability and the purely azimuthal, laminar, flow disappears and large-scale Taylor rolls form. Further increasing of the driving breaks up these rolls, causing the onset of time-dependence as the system transitions from the stationary Taylor vortex regime to the modulated Taylor vortex regime and finally the breakup of these into chaotic turbulent Taylor vortices. These changes of the flow are reflected in transitions of the local scaling laws for the torque versus driving, i.e. versus Taylor number . All this has been studied extensively and summarized e.g. in Andereck et al. (1983); Lathrop et al. (1992b, a); Lewis & Swinney (1999). The mentioned breakup of the rolls leads to the existence of a transitional regime, where the large-scale coherent structures still can be identified when looking at the time-averaged quantities. Looking at the details of the flow, a mixture of turbulent and laminar boundary layers is present.

In this transitional regime, hairpin vortices, which, in the context of RB, can also be viewed as plumes, are ejected from both inner and outer cylinders, and these contribute to large-scale bulk structures. These structures in turn cause an axial pressure gradient, which couples back to the boundary layers, causing plumes to be ejected there. But this only happens from preferential spots in the boundary layers. Once the driving is strong enough, the large-scale structures slowly vanish, and the plumes no longer feel an axial pressure gradient. The boundary layers now become fully turbulent and the flow transitions to the “ultimate” regime. As the flow enters the ultimate regime, and the boundary layer become turbulent, a logarithmic signature in the angular velocity boundary layers is expected, which indeed has been found experimentally (Huisman et al. 2013) and numerically (Ostilla-Monico et al. 2014b).

In the ultimate regime, an effective scaling relation between the Nusselt number , i.e. the non-dimensional torque where is the torque, and the torque in the purely azimuthal state, and the Taylor number is expected, with an effective scaling exponent which exceeds that for the laminar-type boundary layer case (Malkus 1954), for which . I.e. in the ultimate regime, we expect an effective scaling law with . In fact, for that regime, the relation law with logarithmic corrections was suggested, (Kraichnan 1962; Spiegel 1971; Grossmann & Lohse 2011). The logarithmic corrections are quite large, and lead to an effective scaling law with for (Grossmann & Lohse 2011; van Gils et al. 2012). We note that this scaling law is analog to the scaling of the friction factor with Reynolds number in fully turbulent pipes Prandtl (1933).

For the largest drivings, remnants of the larger rolls, which can be seen as a large scale wind, are still observed at even the largest Reynolds numbers studied numerically (Ostilla-Monico et al. 2014b), and experimentally, even up to (Huisman et al. 2014). In Ostilla-Monico et al. (2014b), the remnants of the large scale structures played a crucial role in the transition to the ultimate regime. However, large scale structures are not present in the whole parameter space of TC. Andereck et al. (1986) showed how rich a variety of different states exists at low Reynolds number when the outer cylinder is also rotated. Brauckmann & Eckhardt (2013b) reported that the strength of the large scale wind was most pronounced at the position of optimal transport. However, if the outer cylinder is counter-rotated past the position of optimal transport, bursts arise from the outer cylinder. The flow is very different outside and inside the neutral surface, which separates Rayleigh-stable from Rayleigh-unstable regions of the gap, changing completely the dynamics of the system. The Taylor vortices no longer penetrate the whole gap, extending thus the unstable region effectively somewhat outside the neutral surface of laminar type flow (Ostilla et al. 2013).

The geometry of the system can be expected to play an important role in determining the strength of the large scale wind, and how the transition takes place. In the context of understanding the radius-ratio dependence of the transition to the ultimate regime, Merbold et al. (2013) reported a higher transitional Reynolds numbers for than what was seen for by Ostilla-Monico et al. (2014b) and for by Ravelet et al. (2010). Also the aspect-ratio plays a role. Although different vortical states were known to coexist at low Reynolds number (Benjamin 1978), it was previously thought that if the driving was sufficiently large, only one branch of the torque versus Taylor number curve would survive (Lewis & Swinney 1999). Brauckmann & Eckhardt (2013a) found that the difference in the global response between different vortical states becomes smaller with increasing Reynolds number. Recently, Martinez-Arias et al. (2014) reported on the existence of different vortical states associated to different global torques at a given Taylor number for , and that there is a crossing between those torque-versus- curves around the transition to the ultimate regime. Furthermore, Huisman et al. (2014) showed that different vortical states survive up to Reynolds number of , corresponding to Taylor numbers of order . Furthermore, by combining measurements of global torque and local velocity, Huisman et al. (2014) found that the optimal transport is connected to the existence of the large-scale coherent structures at high Taylor numbers.

Therefore, some questions arise which we want to address in the present paper: How does the transition in the boundary layers take place across the full parameter space of TC? Is the vanishing of the large-scale wind a necessary and/or a sufficient condition for the boundary layer transition? Why does the transition occur later for than for larger values of ? And finally, what is the effect of the vortical wavelength and why do different branches of the torque versus Taylor number scaling curves cross near the transition to the ultimate regime?

## 2 Explored parameter space

### 2.1 Control parameters

To answer these questions, direct numerical simulations (DNS) of TC have been performed across all dimensions of the parameter space, not only adding outer cylinder rotation, but also varying both geometrical parameters and . To do this, the rotating-frame formulation of Ostilla et al. (2013) was used. In that paper, TC was formulated in a frame rotating with the outer cylinder, such that it looks like a system in which only the inner cylinder is rotating, but with a Coriolis force term, which represents the original presence of the outer cylinder rotation. The shear driving of the system is non-dimensionally expressed as a Taylor number, introduced previously:

(2.0) |

is the analog to the Rayleigh number in RB, as elaborated in EGL07. The outer cylinder rotation reflects in a Coriolis force, characterized by a Rossby number . The Rossby number or rather is the parameter which appears in the equations of motion for the fluid:

(2.0) |

where , a geometrical parameter. The Rossby number is related to the frequency ratio via

(2.0) |

Thus fixed means fixed and vice versa. describes co-rotation or , while means counter-rotation. The radius ratio is presented by the geometrical amplitude factor , being small for small gap () and large for large gap (). A resting outer cylinder is described by .

There are also other ways of choosing the control parameters. Classically, they have been expressed as two non-dimensional Reynolds numbers corresponding to the inner and outer cylinders: , where are the azimuthal velocities of the inner and outer cylinders. The classical flow control parameters can be transformed to the parameter space by:

(2.0) |

and

(2.0) |

Viceversa, we have

(2.0) |

and

(2.0) |

The driving can also be expressed as a shear Reynolds number .

### 2.2 Numerical scheme

A second–order finite–difference code was used with fractional time integration. The code was parallelized using hybrid OpenMP and MPI-slab decomposition. Simulations were run on local clusters and on the supercomputer CURIE (Thin nodes) using a maximum of 8192 cores. Details about the code can be found in Verzicco & Orlandi (1996) and in Ostilla et al. (2013). The explored parameter space from previous work (Ostilla et al. 2013; Ostilla-Monico et al. 2014b) was extended through further simulations. Figure 1 shows the parameter space explored in this manuscript. Circles show simulations of a “full” geometry, i.e. a complete cylinder and with . Following the work of Brauckmann & Eckhardt (2013a), the simulations with the largest were performed on “reduced” geometries to reduce computational costs, and these are indicated as squares in the plots. The idea is that instead of simulating the whole cylinder, a cylinder wedge with rotational symmetry of order is considered. The aspect ratio was also reduced to , accommodating a single vortex pair with the wavelength . The vortical wavelength remains the same, although there is a single vortex instead of the three vortex pairs having also the wavelength . Other vortical wavelengths were also simulated using reduced geometries for . We note that the aspect ratio is a geometrical control parameter, but is a response of the system, which depends both on and on the amount of vortex pairs which fit in the system. They are related by , where is the amount of vortex pairs which fit in the system. For all simulations axially periodic boundary conditions were used. Its consequences on the vortex wavelength are analyzed in section 4. Further details on the numerical resolution can be found in Table 1 in the appendix.

### 2.3 Explored parameter space

The top two panels of figure 1 show the parameter space explored for in both and to study the effects of outer cylinder rotation. For , reduced geometries simulate one sixth of the cylinder, i.e. as used in Ostilla-Monico et al. (2014b). The chosen values of include a co-rotating outer cylinder (), a weakly counter-rotating outer cylinder (), counter-rotation near the asymptotic position of optimal transport, (), and two values of in the strongly counter-rotating regime ( and ). No simulations were run in the Rayleigh-stable regime (i.e. when ) as no evidence of turbulence was found in that regime up to in Ostilla-Monico et al. (2014c).

In addition, to study the effects of geometry, i.e. both the radius ratio , and the vortical wavelength (controlled through the aspect ratio ), additional simulations were performed. The bottom left panel shows that two additional radius ratios were simulated up to , one with a larger gap () and one with a smaller gap (). For , one third of the cylinder () was simulated for larger than . This value of for was shown not to affect the values of the torque obtained in the simulations in Brauckmann & Eckhardt (2013b). For , one twentieth of the geometry () was used.

The bottom right panel shows the simulations with varying vortical wavelength done for and pure inner cylinder rotation. was chosen as we expect the effects of the coherent structures, and thus of , to be stronger for larger (see later sections 4 and 5 for an explanation). The values of simulated are around the range where Martinez-Arias et al. (2014) have experimentally observed the crossing of different branches in and also coincides with the onset of the “ultimate” regime.

### 2.4 Non–dimensionalization

The following non-dimensionalizations will be used: as the flow is simulated in a rotating frame, the outer cylinder is stationary, and the system has an unique velocity scale, equal to in the laboratory frame. All velocities are non-dimensionalized using , i.e. . The gap width is the characteristic length scale, and thus used for normalizing distances.

We define the normalized (non-dimensional) distance from the inner cylinder and the normalized height . We furthermore define the time- and azimuthally-averaged velocity fields as:

(2.0) |

where indicates averaging of the field with respect to . As mentioned previously, the torque is non-dimensionalized as an angular velocity “Nusselt” number (EGL07), defined as , where is the torque in the purely azimuthal flow. The torque is calculated from the radial derivative of at the inner and outer cylinders. The simulations are run in time until the respective values are equal within . The torque is then taken as the average value of the inner and outer cylinder torques. Therefore, the error due to finite time statistics is smaller than .

From here on, for convenience we will drop the overhead tilde on all non-dimensionalized variables.

### 2.5 Structure of paper

The organization of the paper is as follows. In section 3, we analyze the effect of rotating the outer cylinder. This is followed by section 4, where we study the influence of , and notice an analogy between the effects of smaller and larger . In section 5, we consider the effects of the last parameter, the vortical wavelength . We finish in section 6 with a summary of the results and an outlook for future work.

## 3 The effect of outer cylinder rotation or the inverse Rossby number dependence

In this section we will study the effect of the Coriolis force (), originating from the rotation of the outer cylinder, on the scaling of with and, more specifically, the effect of on the transition to the ultimate regime. Depending on the value of , two distinct regimes will be identified: First a co- and weakly counter-rotating range, denoted from here on as CWCR regime, and second the strongly counter-rotating range, denoted from here on as SCR regime. The CWCR regime is found when the outer cylinder either is at rest, co-rotates with the inner cylinder, or slowly counter-rotates. The counter-rotation must be slow enough such that no Rayleigh-stable zones are generated in the bulk of the flow. In this CWCR regime the Coriolis force is balanced through the bulk gradient of . This can be derived from a large scale balance in the -component of the velocity in equation (2.1). In summary, the non–linear term and the Coriolis force term balance each other out on average (cf. Ostilla et al. (2013) for the full derivation). This results in a linear relationship between and (Ostilla-Monico et al. 2014a).

Taylor-Couette flow can be considered as being in the SCR regime, if the outer cylinder strongly counter-rotates and generates a Coriolis force which exceeds what the -gradient can balance. The threshold value of corresponds to the flattest profile. This also is the value of , for which is found to be largest (van Gils et al. 2012; Ostilla et al. 2013), denoted henceforth as . In this regime the turbulent plumes originating from the inner cylinder are not strong enough to overcome the stabilizing effect of the outer cylinder, and the flow is divided into two regions, a Rayleigh-stable region in the outer gap region, which plumes do not reach, and a Rayleigh-unstable region in the inner parts of the gap. For given Coriolis force, the relative sizes of these spatial regions depend on , as for a stronger driving (i.e. larger ), the turbulence originating from the inner cylinder “pushes” these zones more towards the outer cylinder. This may lead to switching between vortical states and jumps in global quantities as seen in Ostilla et al. (2013). The boundary between both regimes is at . Of course, depends on too, due to effect of viscosity in the Coriolis force balance (Ostilla et al. 2013), and only saturates to for sufficiently high drivings of and more (cf. both panels of figure 2 and Brauckmann & Eckhardt (2013a)).

Figure 2 shows both and the compensated Nusselt number versus for and the six values of studied. For the largest drivings (i.e. ) all values of reach the effective scaling law (with a different amplitude), similar to what was reported in the experiments by van Gils et al. (2011). However, very different behavior can be seen for , i.e. before the onset of the ultimate regime.

In the CWCR regime (), the Coriolis force is reflected in the flow structure through the bulk gradient of , making it either flatter as in the case of weak-counter rotation, or steeper, as in the case of co-rotation (if the driving is sufficiently large). A consequence of the angular velocity gradient in the bulk is that large scale structures can be weakened or even completely dissapear in the CWCR regime. These changes in -gradient strongly affect the capability of plumes to “coordinate” and form a large-scale wind, which in turn leads to an earlier (or later) onset of the sharp decrease in the local exponent in the scaling law associated with the breakdown of coherence, and the onset of time dependence in (Ostilla-Monico et al. 2014b).

For the case of co-rotating cylinders (), this happens when the system enters the so-called “wavelet” regime, characterized by moving waves in the boundary regions between a pair of Taylor vortices (Andereck et al. 1983, 1986). These waves move with different speeds, and as a consequence this regime is not stationary in any reference frame. This regime only persists for a small range of , and eventually all remnants of Taylor vortices vanish. Axial dependence of the flow structure is almost completely lost, even at as low as . Unlike the case of studied in Ostilla-Monico et al. (2014b), however, in this transitional regime, the large–scale rolls already completely vanished, but for this does not immediately lead to the transition to the ultimate regime. After its sharp decrease, does not exceed . Instead, at a driving strength around (coinciding with the disappearance of the structures), the local effective scaling exponent has increased to , and then stops growing. Only if increases further and the shear in the boundary layers grows past a threshold, a shear-instability takes place, and the system transitions to the ultimate regime.

For the case of counter-rotating cylinders, (i.e. ), can locally grow beyond in the classical regime. This is unexpected, as values of larger than one third have been associated to the transition to turbulence of the boundary layers in the context of both Rayleigh-Bénard convection (He et al. 2012b), and TC with a stationary outer cylinder (Ostilla-Monico et al. 2014b). However, in this case, the shear in the boundary layers is too low so the boundary layers still stay laminar.

For counter-rotating cylindres, a wide range of flow configurations is available in the low- regime (Andereck et al. 1986). We can relate local steps in to the switching between such flow configurations. The interplay between Rayleigh-stable and -unstable regions can also play a role. Larger drivings cause the Rayleigh-unstable region to grow, and thus to increase the transport. These two effects lead to larger increases in the non-dimensional torque than what is expected for pure inner cylinder rotation, and explain the large values of seen.

To further illustrate the effect of the Coriolis force on the large-scale structures, figure 3 presents a contour plot of in the CWCR regime , around the optimum and in the SCR regime . Figure 4 shows the axially-averaged angular velocity profiles for and the six values of simulated here. The large–scale structures cannot be seen in the left panel of figure 3, which corresponds to (co-rotating cylinders), but they are pronounced for the other two panels ( and ). As shown in figure 4, in the CWCR regime, the bulk sustains a large gradient, and to accomodate for this, there is smaller jump across the boundary layers. Plumes ejected from both cylinders can now mix easier when entering the bulk. As a consequence, the large-scale structures, which essentially consist of unmixed plumes, break up easier and thus do that for lower values of . For this reason they have completely vanished in the left panel of figure 3.

If we now decrease , the profile becomes flatter. The effect of this is visible in the middle panel of Figure 3 showing for . It can be seen from figure 4 that this value of corresponds to the flattest -profile available, and it is also the closest to the experimental optimum transport . A very marked signature of the large-scale structure on can be seen. This is because a very flat profile will sustain a large jump across the boundary layer, and thus plumes detach less violently into the bulk, thus stabilizing the large-scale structures. Therefore, we can relate the flatness of the -profile to the strength of the large-scale circulation, and this in turn can be related to the optimum in . As mentioned in Brauckmann & Eckhardt (2013b), optimum transport coincides with the strongest mean circulation. Plumes travel faster from one cylinder to the other when the large-scale circulation is strongest, and thus more angular momentum is transferred. We also highlight that the signature of the large-scale structures on the mean azimuthal flow remains even in the ultimate regime, and is also seen in experiment at (Huisman et al. 2014). Thus in general the vanishing of the rolls appears to be independent from the transition to the ultimate regime. Only in the special case of pure inner cylinder rotation these two effects coincidentally occur at the same .

In the right panel of figure 3, we can see that once the Coriolis force is sufficiently large, the vortices cannot fully penetrate the domain. Near the outer cylinder, the flow is predominantly Rayleigh-stable. Rayleigh-stable zones are well mixed, as transport here happens through intermittent turbulent bursts, instead of convective transport by plumes and vortices (Brauckmann & Eckhardt 2013b). Thus, in Rayleigh-stable regions, no rolls can be seen in the averaged fields. The effect of the neutral surface can also be observed in the averaged profiles (cf. figure 4). The two simulated cases in the SCR regime, ( and ) show an outer cylinder boundary layer which with more and more negative extends deeper into the flow, and the distinction from the bulk is blurred away.

To further disentangle the effect of axial dependence and the transition to the ultimate regime we show the loss of axial dependence characterized by a special spread measure as a function of the driving in figure 5. is defined as , with , the mid-gap, defined as , the arithmetic mean of the inner and outer cylinder radii. When measuring the axial spread, the velocity is averaged in time, and azimuthally, as the flow is homogeneous in the azimuthal direction. As stated previously, for co-rotating cylinders, the axial dependence disappears for low drivings corresponding to those in the transitional regime, and associated to the appearance of the “wavelet” states. For counter-rotating cylinders, a sharp jump in can be noticed. This is due to being measured at the mid-cylinder . For low drivings, is located in the Rayleigh-stable zones, and the flow is mixed better. As the driving increases, turbulence from the inner cylinder pushes the neutral surface, which divides the stable and unstable zones further towards the outer cylinder. As a consequence of this pushing, is no longer in the Rayleigh-stable zone, but instead in the Rayleigh-unstable zone. This zone is dominated by large-scale structures. This makes the axial dependence increase and provides more evidence that the vanishing of the Taylor-rolls is only coincidental with the transition to the ultimate regime for pure inner cylinder rotation.

As mentioned previously, the value of , and thus of the border between the CWCR and the SCR regimes depends on . This is summarized in figure 6, which shows the approximate division between the different flow regimes explored in this paper in both the and the parameter spaces, both for . , and thus the division between the regims can be seen to saturate for , when driving is large enough, and the mean profile at is completely flat.

Finally, to further justify the division of the flow into the CWCR and the SCR regimes with decreasing inverse Rossby number , we can quantify the distribution of Rayleigh-stable and unstable zones as a function of . This is done by looking at the PDF of , i.e. the collection of points outlining the neutral surface . This is, the border between Rayleigh-stable outer gap range and the Rayleigh-unstable inner gap parts, and given as the points for which in the laboratory (non-rotating) frame. For counter-rotating cylinders, the neutral surface defines the instantaneous border between Rayleigh-stable and Rayleigh-unstable zones. For co-rotating cylinders, the neutral line does not exist, and the whole flow is either Rayleigh-stable or Rayleigh-unstable. In principle, the neutral surface might be fragmented, and thus the position of multivalued. However, this is usually not the case. When taking the ensemble, all values are considered, as this does not change the PDFs significantly.

Figure 7 shows the PDFs of calculated for the four negative values of at the largest driving simulated here. The difference between the two regimes can clearly be noticed. In the CWCR regime and near the optimum, the border between the zones is located very closely to the outer cylinder, which means that almost all the domain is Rayleigh-unstable and dominated by plumes or rolls. In the SCR regime, the border between the zones is pushed closer towards the inner cylinder, and Rayleigh-stable zones appear all over the gap. For the most negative simulated value of , i.e. , the areas near the outer cylinder are permanently Rayleigh-stable, and transport occurs in intermittent bursts which mix this zone well. This causes the partial dissappearance of axial dependence seen in the right panel of Figure 3.

## 4 The effect of radius ratio or the -dependence

In the previous section we showed that for the transition to the ultimate regime and the vanishing of the rolls only (incidentally) co-occur at the same for pure inner cylinder rotation. Flatter bulk -profiles result in stronger large-scale structures, and steeper bulk -profiles result in weaker large-scale structures which vanish at . Now we will show that we can modify the profile in the bulk not only by varying the Coriolis force, but also by changing the radius ratio (or the gap width). In this section, we will thus analyze the influence of , to understand whether the co-ocurrence of the vanishing large scales and the boundary layer transition observed for pure inner cylinder rotation is just a coincidence seen in the case .

Figure 8 shows both the Nusselt number and the compensated Nusselt number plotted as a function of for the three values of simulated. As seen in Ostilla-Monico et al. (2014a) for (and now also for ), the flow undergoes a structural transition at around , where the local exponent of the effective scaling law rapidly decreases. This is associated with the breakdown of coherence in the flow and the onset of time-dependence in the Nusselt number. For and , the effective exponent begins to increase again after this breakdown. We can say that the flow transitions to the ultimate regime once , and this happens at about . This value coincides with the experimentally observed value for the transition to the ultimate regime for , cf. Ravelet et al. (2010).

For a different behavior can be seen. After the breakdown of coherence, the transitional regime with goes on for three decades in , up to (last three data points of the panel). An increase in only happens for the last three data points, with . This might be the beginning of the transition to the ultimate regime, observed at about that value of in the experiments by Merbold et al. (2013). We emphasize that the behavior of the curve for is very similar to the one seen for and (cf. figure 2), while the curve for is similar to the one for and .

We thus can draw an analogy between the effects of varying and those of changing . The larger the gap or the smaller is, the more the flow feels the curvature. This is reflected in an asymmetry between inner and outer cylinder, since the inner cylinder curvature becomes increasingly stronger relative to the outer cylinder curvature. Also the exact relationship (cf. van Gils et al. (2012)) must hold in both boundary layers due to the -independence of the angular velocity current (EGL07). For we have and the -slope at the inner cylinder is eight-fold steeper than the outer cylinder -slope. Thus the inner-outer asymmetry is expected to become much more dominant for in comparison to () as well as (), for which it is hardly visible anymore.

While the inner and outer cylinder boundary layers extend into the bulk equally for pure inner cylinder rotation (cf. (Ostilla-Monico et al. 2014a)), the jump of in the boundary layers is much larger in the inner cylinder as compared to the outer cylinder due to the different slopes and equal extents. Therefore, the plumes are highly asymmetric, and smaller drivings break up the “plume conveyor belts”, which form the large-scale structures seen in the time-averaged azimuthal velocity. On top of this plume asymmetry, originating from the boundary layers, a larger curvature has an effect on the bulk. The underlying profile is less flat, and thus the drop in angular velocity inside the bulk is the larger the smaller the value of is.

Both effects can be appreciated in figure 9, which shows contour plots of the azimuthally- and time-averaged angular velocity at for the three simulated values of . This also explains the left panel of figure 10, where the now also axially averaged angular velocity is shown for the same three values of . For comparison, the right panel of 10 shows three profiles of in the CWCR regime for .

The analogy between the effect of and the effect of on is also demonstrated in figure 10. The rolls are weak for , as they are weak for co-rotating cylinders, and the rolls are strongest for and for . This also explains why, for large enough , is highest at a given for the largest . However, the analogy is not perfect. For pure inner cylinder rotation, i.e., for the wide variety of flow states seen in Andereck et al. (1983) and Andereck et al. (1986) is greatly reduced. The system essentially goes from Taylor vortex flow to modulated Taylor vortex flow to finally turbulent Taylor vortex flow. It does not undergo transitions to different states (such as e.g. the “wavelet” state), and thus the rolls do not vanish for the lower drivings at which this happens in co-rotating cylinders. This can be seen in figure 11, which shows the measure for the axial velocity spread as function of . With increased driving, the rolls progressively lose importance until reaches a value of . However, the effect of , and thus of the cylinder wall curvature on the profiles can be clearly noticed in the residual axial dependence and behaves as expected from the analogy. The behaviour of the transition to the ultimate regime and associated sub–regimes is summarized in figure 12, which is analogous to figure 6, but now for the parameter space explored.

Finally, one may ask the question of why the onset of the ultimate regime happens at a much higher for than for the two other values of studied. For , the transition seems to set in for the same value of independently of . A factor ten increase in shear in the boundary layers is required for the boundary layer instability to occur and the ultimate regime to set in. Convex curvature is known to produce a stabilizing effect on boundary layers (Görtler 1940a; Muck et al. 1985), and this will have a more significant effect on the inner cylinder for than for the larger . On the other hand we might expect that the destabilizing effect of concave curvature (Görtler 1940b; Hoffmann et al. 1985) would also play a role in accelerating the transition. Due to the boundary layer asymmetry however, the outer boundary layer is much more “quiet”, and has less fluctuations. This also delays the transition, and can be seen in figure 13, which shows the rms-fluctuations of the angular velocity , for and the three values of simulated. The levels of fluctuations at the outer cylinder are significantly reduced for when compared to the other values of . Finally, the large gradient of angular velocity sustained in the bulk will also reduce the shear in the outer cylinder, as the bulk angular velocity is smaller for . Thus, a combination of reduced fluctuations, stabilizing effect due to curvature at the inner cylinder, and reduced shear due to bulk angular velocity gradients is causing the delayed transition.

## 5 Dependence on number and size of rolls

Finally, we will quantify how the torque depends on the number and the size of the rolls, i.e. the vortical wavelength. The wavelength of a roll is restricted to the values , where is a strictly positive integer. For all simulations in this paper, , and thus . This is not necessarily always the case, is a response of the system, and if is large enough, i.e. the system can accomodate more than one vortex pair, can take several values depending on how the final state of the system is reached. Brauckmann & Eckhardt (2013a) showed that for , the “optimal” vortex wavelength, i.e. the vortex wavelength which corresponds to a maximum , increased when comparing for two Taylor numbers, one in the Taylor vortex regime and another in the turbulent Taylor vortex regime. For the higher , the dependence of on was quite weak. Martinez-Arias et al. (2014) showed that for , different branches in the relationship, associated to distinct vortical states cross around . This corresponds to a driving of , around the value at which the transition to the ultimate regime occurs for . The large-scale circulation could still be seen to play a role in determining the system response after the transition to the ultimate regime. Furthermore, large scale patterns were observed in Ostilla-Monico et al. (2014b) when looking at the correlation at , even though they are absent when looking only at .

Figure 14 shows the compensated torque as function of for the four values of the vortical wavelength studied. Experimental data by Martinez-Arias et al. (2014) and DNS data by Ostilla-Monico et al. (2014a) is also plotted. It is worth noting that experimental data will have some end-plate effects, even if the aspect ratio of the experiments is larger than , while the DNSs have periodic axial boundary conditions. Even so, very similar behaviour can be seen. The transition to the asymptotic scaling laws of the ultimate regime seem to occur around the same value of , but are less pronounced the smaller the vortical wavelength is.

The change in behaviour of the curves can be associated to the change of behaviour of the wind-sheared regions in the ultimate regime. As seen in Ostilla-Monico et al. (2014b), plume ejection is supressed outside the ultimate regime in regions of the flow, the so called “wind-sheared” regions due to the sweeping by the large scale rolls. This reduction in plume ejection results in a reduced transport of angular velocity (torque). A similar reduction in the torque caused by a mean flow was also seen when forcing the flow with an axial pressure gradient by Manna & Vacca (2009). Vortices with a smaller wavelength have smaller wind-sheared regions and thus result in a larger , if this suppression is taking place. After the transition to the ultimate regime, the suppression ceases, and these regions become active ejectors of plumes, leading to increased transport.

The difference between , and the is very small for , of the order of , but for the branch the difference is almost . Only at , when the distinction between wind-sheared and ejection regions is completely blurred away, and the whole inner cylinder can emit plumes (or hairpin vortices), loses its dependence, within the error bars of the numerics. This sudden transition of wind-sheared regions to ejection regions causes the jump we see in the curve at around for .

Note that for the largest drivings axially periodic boundary conditions have been used, with only one vortex pair. This does not prevent the creation of two pairs of vortices with wavelength by a breakup of one pair of vortices of in a domain, which has . And indeed this is seen to happen for the lower drivings both in DNS and experiment. On the other hand, this axial periodicity affects the stability of one pair of vortices of wavelength in a domain of . Therefore, vortices with might be an artifact due to the numerical constraintment, and not be stable if a system with large at large is considered. States with are not reported in Martinez-Arias et al. (2014).

Even if we do not expect a quantitative agreement of the present DNS results with those of Brauckmann & Eckhardt (2013a) and experimental data by Huisman et al. (2014), as we simulate a different , the results reported in this section even do not agree qualitatively. Brauckmann & Eckhardt (2013a) see a maximum in torque for in the turbulent Taylor vortex regime (), while in the present simulations for at the same , this maximum is clearly at , and not near . In the experiments of Martinez-Arias et al. (2014), states with smaller than one are not reported, and a direct comparison cannot be made.

We also note that the relationship between larger vortices and larger torque in the ultimate regime is the inverse of what was recently reported by Huisman et al. (2014). Huisman et al. (2014) found multiple states, with different in highly turbulent TC flow. For different states they found that the torque differs less than , although they note that this might be due to the fact the torque is only measured on part of the inner cylinder, not on the entire inner cylinder. Furthermore their results are for , for higher , and for different , as compared to the current research.

## 6 Summary and conclusions

Numerical simulations of turbulent Taylor-Couette flow in the range were performed to explore the transition of TC flow to the (fully turbulent) ultimate regime. The four dimensions of the parameter space were explored, including the dependence of the transition on the radius ratio , the vortex wavelength and Coriolis force or rotation ratio .

First, the effect of the outer cylinder rotation, in the equations of motion in the frame co-rotating with the outer cylinder, present as a Coriolis force, was analyzed for . Depending on the value of two regimes were identified, (i) the co-rotating and weakly counter-rotating cylinder regime (CWCR) and (ii) the strongly counter-rotating cylinder regime (SCR), both with their respective sub–regime. Our findings of that chapter culminate in the phase diagram fig. 6, in the regime (fig. 6a) and in the regime (fig. 6b & c). The transition to the ultimate regime could be observed for all values of around . However, for these two regimes a rather different behavior in the scaling laws was found before the transition. We also found very different flow structures in the respective ultimate regimes in accordance with the description by Brauckmann & Eckhardt (2013b). An explanation why the Coriolis force, proportional to stabilizes the large-scale structures was illustrated; the large-scale structures were found to not vanish at the transition to the ultimate regime for , unlike what was seen in Ostilla-Monico et al. (2014b) for resting outer cylinder.

After this, the transition was analyzed for various gap widths, namely for , , and without Coriolis forces, i.e., for . The transition was found to occur at about the same for and . However, the transition was considerably delayed to for , due to the combined effects of stabilizing curvature of the inner cylinder, and the reduced shear as well as smaller fluctuations in the vicinity of the outer cylinder. An analogy between the effect of in the CWCR regime and the effect of on the large scale rolls was described: Decreasing was found to have the same effect as adding a positive –corresponding to co-rotating cylinders– , while increasing behaved like (weakly) counter-rotating the outer cylinder.

Finally, as the large-scale structures were found to be strongest for , the effect of varying the vortical wavelength was analysed for this value of . As in Martinez-Arias et al. (2014), different branches of the curve were found to cross around the transition to the ultimate regime. Before this transition, the influence of the vortical wavelength (and thus of the aspect ratio) on was quite noticeable. After the ultimate range transition, this effect decreased drastically. The results of our DNS agree qualitatively with those in the experiments by Martinez-Arias et al. (2014) for even though the axial boundary conditions are different. However, they are qualitatively different from those reported for by Brauckmann & Eckhardt (2013a) and by Huisman et al. (2014)

In this work, the vortical wavelength by using periodic boundary conditions was fixed. Some of these states might not be accessible in experiment or might be a product of the periodic boundary conditions. Studying the coexistence of different states for large , like done in Martinez-Arias et al. (2014) or Huisman et al. (2014) with DNS requires a large amount of computational resources for high . Switches between two and three vortex pairs were seen at lower for (Ostilla-Monico et al. 2014a). Switching between states might also occur at high , although they are not captured in the DNS presented in this work. In the future, additional DNS for with large at high should be run to improve the understanding of the switching between different states.

Our ambition also is to further understand why the transition is delayed at , but also the curvature effects on the -profiles in the boundary layers along the ideas of Grossmann et al. (2014). Curvature effects at and are too small to be appreciated, and the flow for is still in the transition to the “ultimate” regime. Thus, higher simulations for will provide further understanding on how curvature makes the boundary layers of TC flow different from those of channel and pipe flow.

Acknowledgements: We would like to thank H. Brauckmann, J. Peixinho, M. Salewski, and C. Sun for various stimulating discussions during the years. We would also like to thank the Dutch Supercomputing Consortium SurfSARA for technical support, FOM, COST from the EU and ERC for financial support through an Advanced Grant. We acknowledge that these results come from computational resources at the PRACE resource Curie, based in France at GENCI/CEA.

- Ahlers et al. (2012) Ahlers, G., He, X., Funfschilling, D. & Bodenschatz, E. 2012 Heat transport by turbulent Rayleigh-Bénard convection for and : aspect ratio . New J. Phys.. 14, 103012.
- Andereck et al. (1983) Andereck, C. D., Dickman, R. & Swinney, H. L. 1983 New flows in a circular Couette system with corotating cylinders. Phys. Fluids 26 (1395).
- Andereck et al. (1986) Andereck, C. D., Liu, S. S. & Swinney, H. L. 1986 Flow regimes in a circular Couette system with independently rotating cylinders. J. Fluid Mech. 164, 155–183.
- Benjamin (1978) Benjamin, T. B. 1978 Bifurcation phenomena in steady flows of a viscous liquid. Proc. R. Soc. London A 359, 1–43.
- Brauckmann & Eckhardt (2013a) Brauckmann, H. & Eckhardt, B. 2013a Direct Numerical Simulations of Local and Global Torque in Taylor-Couette Flow up to Re=30.000. J. Fluid Mech. 718, 398–427.
- Brauckmann & Eckhardt (2013b) Brauckmann, H. J. & Eckhardt, B. 2013b Intermittent boundary layers and torque maxima in Taylor-Couette flow. Phys. Rev. E 87 (3), 033004.
- Couette (1890) Couette, M. 1890 Études sur le frottement des liquides. Gauthier-Villars et fils.
- Donnelly (1991) Donnelly, R. 1991 Taylor-Couette flow: the early days. Physics Today pp. 32–39.
- Eckhardt et al. (2007) Eckhardt, B., Grossmann, S. & Lohse, D. 2007 Torque scaling in turbulent Taylor-Couette flow between independently rotating cylinders. J. Fluid Mech. 581, 221–250.
- van Gils et al. (2011) van Gils, D. P. M., Huisman, S. G., Bruggert, G. W., Sun, C. & Lohse, D. 2011 Torque scaling in turbulent Taylor-Couette flow with co- and counter-rotating cylinders. Phys. Rev. Lett. 106, 024502.
- van Gils et al. (2012) van Gils, D. P. M., Huisman, S. G., Grossmann, S., Sun, C. & Lohse, D. 2012 Optimal Taylor-Couette turbulence. J. Fluid Mech. 706, 118–149.
- Görtler (1940a) Görtler, H. 1940a Über den Einflusss der Wandkrümmung auf die Entstehung der Turbulenz. Z. Angew. Math. Mech. 20, 138–147.
- Görtler (1940b) Görtler, H. 1940b Über eine dreidimensionale Instabilität laminarer Grenzschichten an konkaven Wänden. Z. Angew. Math. Mech. 21, 250–252.
- Grossmann & Lohse (2000) Grossmann, S. & Lohse, D. 2000 Scaling in thermal convection: A unifying view. J. Fluid. Mech. 407, 27–56.
- Grossmann & Lohse (2001) Grossmann, S. & Lohse, D. 2001 Thermal convection for large Prandtl number. Phys. Rev. Lett. 86, 3316–3319.
- Grossmann & Lohse (2011) Grossmann, S. & Lohse, D. 2011 Multiple scaling in the ultimate regime of thermal convection. Phys. Fluids 23, 045108.
- Grossmann et al. (2014) Grossmann, S., Lohse, D. & Sun, C. 2014 Velocity profiles in strongly turbulent Taylor-Couette flow. Phys. Fluids 26, 025114.
- He et al. (2012a) He, X., Funfschilling, D., Bodenschatz, E. & Ahlers, G. 2012a Heat transport by turbulent Rayleigh-Bénard convection for and : ultimate-state transition for aspect ratio . New J. Phys.. 14 (6), 063030.
- He et al. (2012b) He, X., Funfschilling, D., Nobach, H., Bodenschatz, E. & Ahlers, G. 2012b Transition to the ultimate state of turbulent Rayleigh-Bénard convection. Phys. Rev. Lett. 108, 024502.
- Hoffmann et al. (1985) Hoffmann, P. H., Muck, K. C. & Bradshaw, P. 1985 The effect of concave surface curvature on turbulent boundary layers. J. Fluid Mech. 161, 371–403.
- Huisman et al. (2012) Huisman, S. G., van Gils, D. P. M., Grossmann, S., Sun, C. & Lohse, D. 2012 Ultimate turbulent Taylor-Couette flow. Phys. Rev. Lett. 108, 024501.
- Huisman et al. (2013) Huisman, S. G., Scharnowski, S., Cierpka, C., Kähler, C., Lohse, D. & Sun, C. 2013 Logarithmic boundary layers in strong Taylor-Couette turbulence. Phys. Rev. Lett. 110, 264501.
- Huisman et al. (2014) Huisman, S. G., van der Veen, R. C. A., Sun, C. & Lohse, D. 2014 Multiple states in ultimate Taylor-Couette turbulence. Nature Comm. 5 (3820).
- Kraichnan (1962) Kraichnan, R. H. 1962 Turbulent thermal convection at arbritrary Prandtl number. Phys. Fluids 5, 1374–1389.
- Lathrop et al. (1992a) Lathrop, D. P., Fineberg, Jay & Swinney, H. S. 1992a Transition to shear-driven turbulence in Couette-Taylor flow. Phys. Rev. A 46, 6390–6405.
- Lathrop et al. (1992b) Lathrop, D. P., Fineberg, Jay & Swinney, H. S. 1992b Turbulent flow between concentric rotating cylinders at large Reynolds numbers. Phys. Rev. Lett. 68, 1515–1518.
- Lewis & Swinney (1999) Lewis, G. S. & Swinney, H. L. 1999 Velocity structure functions, scaling, and transitions in high-Reynolds-number Couette-Taylor flow. Phys. Rev. E 59, 5457–5467.
- Malkus (1954) Malkus, M. V. R. 1954 The heat transport and spectrum of thermal turbulence. Proc. R. Soc. London A 225, 196–212.
- Mallock (1896) Mallock, A. 1896 Experiments on fluid viscosity. Phil. Trans. R. Soc. Lond. A 187, 41–56.
- Manna & Vacca (2009) Manna, M. & Vacca, A. 2009 Torque reduction in Taylor–Couette flows subject to an axial pressure gradient. J. Fluid Mech. 639, 373–401.
- Martinez-Arias et al. (2014) Martinez-Arias, B., Peixinho, J., Crumeyrolle, O. & Mutabazi, I. 2014 Effect of the number of vortices on the torque scaling in Taylor-Couette flow. J. Fluid Mech. 748, 756–767.
- Merbold et al. (2013) Merbold, S., Brauckmann, H. & Egbers, C. 2013 Torque measurements and numerical determination in differentially rotating wide gap Taylor-Couette flow. Phys. Rev. E. 87 (023014).
- Muck et al. (1985) Muck, K. C., Hoffmann, P. H. & Bradshaw, P. 1985 The effect of convex surface curvature on turbulent boundary layers. J. Fluid Mech. 161, 347–369.
- Ostilla et al. (2013) Ostilla, R., Stevens, R. J. A. M., Grossmann, S., Verzicco, R. & Lohse, D. 2013 Optimal Taylor-Couette flow: direct numerical simulations. J. Fluid Mech. 719, 14–46.
- Ostilla-Monico et al. (2014a) Ostilla-Monico, R., Huisman, S. G., Jannink, T.J.G., van Gils, D. P. M., Verzicco, R., Grossmann, S., Sun, C. & Lohse, D. 2014a Optimal Taylor-Couette flow: Radius ratio dependence. J. Fluid Mech. 747, 1–29.
- Ostilla-Monico et al. (2014b) Ostilla-Monico, R., van der Poel, E. P., Verzicco, R., Grossmann, S. & Lohse, D. 2014b Boundary layer dynamics at the transition between the classical and the ultimate regime of Taylor-Couette flow. Phys. Fluids 26 (015114).
- Ostilla-Monico et al. (2014c) Ostilla-Monico, R., Verzicco, R., Grossmann, S. & Lohse, D. 2014c Turbulence decay towards the linearly-stable regime of Taylor-Couette flow. J. Fluid Mech. 747, 1–29.
- Paoletti & Lathrop (2011) Paoletti, M. S. & Lathrop, D. P. 2011 Angular momentum transport in turbulent flow between independently rotating cylinders. Phys. Rev. Lett. 106, 024501.
- Prandtl (1933) Prandtl, L. 1933 Neuere ergebnisse der turbulenzforschung. Zeit. Vereines Deutscher Ingenieure 77 (5), 105–114.
- Ravelet et al. (2010) Ravelet, F., Delfos, R. & Westerweel, J. 2010 Influence of global rotation and Reynolds number on the large-scale features of a turbulent Taylor–Couette flow. Phys. Fluids 22 (5), 055103.
- Roche et al. (2010) Roche, P. E., Gauthier, G., Kaiser, R. & Salort, J. 2010 On the triggering of the ultimate regime of convection. New J. Phys. 12, 085014.
- Spiegel (1971) Spiegel, E. A. 1971 Convection in stars. Ann. Rev. Astron. Astrophys. 9, 323–352.
- Taylor (1923) Taylor, G. I. 1923 Experiments on the motion of solid bodies in rotating fluids. Proc. R. Soc. Lond. A 104, 213–218.
- Taylor (1936) Taylor, G. I. 1936 Proc. R. Soc. London A 157, 565–578.
- Verzicco & Orlandi (1996) Verzicco, R. & Orlandi, P. 1996 A finite-difference scheme for three-dimensional incompressible flow in cylindrical coordinates. J. Comput. Phys. 123, 402–413.
- Wendt (1933) Wendt, F. 1933 Turbulente Strömungen zwischen zwei rotierenden Zylindern. Ingenieurs-Archiv 4, 577–595.

APPENDIX: NUMERICAL DETAILS

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||

- | - | |||