Issue 
A&A
Volume 619, November 2018



Article Number  A44  
Number of page(s)  17  
Section  Astrophysical processes  
DOI  https://doi.org/10.1051/00046361/201833024  
Published online  01 November 2018 
Subcritical transition to turbulence in accretion disc boundary layer
Sternberg Astronomical Institute, Moscow M.V. Lomonosov State University, Universitetskij pr., 13, Moscow, 119234 Russia
email: zhuravlev@sai.msu.ru
Received:
14
March
2018
Accepted:
1
August
2018
Context. Enhanced angular momentum transfer through the boundary layer near the surface of weakly magnetised accreting star is required in order to explain the observed accretion timescales in lowmass Xray binaries, cataclysmic variables, or young stars with massive protoplanetary discs. The accretion disc boundary layer is locally represented by incompressible homogeneous and boundless flow of the cyclonic type, which is linearly stable. Its nonlinear instability at the shear rates of the order of the rotational frequency remains an issue.
Aims. We put forward a conjecture that hydrodynamical subcritical turbulence in such a flow is sustained by the nonlinear feedback from essentially threedimensional vortices, which are generated by quasitwodimensional trailing shearing spirals grown to high amplitude via the swing amplification. We refer to those threedimensional vortices as crossrolls, since they are aligned in the shearwise direction in contrast to streamwise rolls generated by the antiliftup mechanism in rotating shear flow on the Rayleigh line.
Methods. Transient growth of crossrolls is studied analytically and further confronted with direct numerical simulations (DNS) of the dynamics of nonlinear perturbations in the shearing box approximation.
Results. A substantial decrease of transition Reynolds number R_{T} is revealed as one changes a cubic box to a tall box. DNS performed in a tall box show that R_{T} as a function of shear rate accords with the line of constant maximum transient growth of crossrolls. The transition in the tall box has been observed until the shear rate is three times higher than the rotational frequency, when R_{T} ∼ 50 000.
Conclusions. Assuming that the crossrolls are also responsible for turbulence in the Keplerian flow, we estimate R _{T} ≲ 10^{8} in this case. Our results imply that nonlinear stability of Keplerian flow should be verified by extending turbulent solutions found in the cyclonic regime across the solidbody line rather than entering a quasiKeplerian regime from the side of the Rayleigh line. The most favourable shear rate to test the existence of turbulence in the quasiKeplerian regime may be subKeplerian and equal approximately to 1/2.
Key words: hydrodynamics / accretion, accretion disks / instabilities / turbulence / protoplanetary disks
© ESO 2018
1. Introduction
Rotating boundary layers are believed to exist in the vicinity of weakly magnetised stars accumulating the material from an accretion disc. A slowly spinning star surface causes transformation of the rotational energy of the disc into thermal energy inside the boundary layer. Thus, the boundary layer is expected to have a brightness comparable to that of the accretion disc and be the main source of radiation in a relatively higher frequency band than the disc itself. The structure of the boundary layer and the adjacent disc inner part must be determined together (see BisnovatyiKogan 1994). A number of boundary layer models have been developed for various classes of objects. Usually they employ slim disc equations (Abramowicz et al. 1988) in order to match consistently the disc inner part to the boundary layer. Among them the accreting premain sequence stars were studied by Popham et al. (1993), whereas the accreting white dwarfs and neutron stars were studied by Narayan & Popham (1993) and Popham & Sunyaev (2001). In the latter case, the formation of a spreading layer is also possible in systems with high enough accretion rate, as was suggested by Inogamov & Sunyaev (1999).
Just as for the formation of an accretion disc, the cornerstone of boundary layer formation is the origin of effective shear viscosity providing enhanced angular momentum and mass transfer onto the star surface. Conventionally, it is parametrised by α, which is the kinematic viscosity coefficient scaled by the speed of sound and pressure radial scaleheight (see Shakura & Sunyaev 1988). This variant of viscosity prescription is similar to what is done in the disc models, where α is scaled by the disc thickness (see Shakura & Sunyaev 1973). However, the physical process veiled by an effective viscosity as well as magnitude of α remains a matter of debate, probably, to a greater extent, in the case of boundary layer rather than in the case of accretion disc. Turbulence is known to be a natural solution to this issue, but its simple (supercritical) variant is discarded by the centrifugal stability of the rotating shear flows, which represent both accretion discs and boundary layers. A “magic wand” working in hot magnetised discs, where angular velocity Ω decreases with distance to the rotation axis r, is a magnetorotational instability. However, it does not work if dΩ^{2}/dr > 0, which is the case for the boundary layers. Such flows are locally linearly stable even in the presence of the magnetic field. In this context, an alternative mechanism of angular momentum transfer has been proposed by Belyaev et al. (2013) and more recently by Philippov et al. (2016). It was shown that global sonic instability of the flow associated with the presence of the star surface excites global nonaxisymmetric acoustic modes, which are responsible for effective viscosity even in the absence of turbulence.
At the same time, the possibility that astrophysical boundary layers acquire the effective viscosity through the turbulence is not ruled out, since flows with shear rates including that of the order of Ω may be locally nonlinearly unstable. Along with the hydrodynamical, nonlinear stability of homogeneous Keplerian shear flow, which locally represents the accretion and protoplanetary discs, this is one of the major unresolved problems in astrophysical fluid dynamics. Below we consider both problems jointly. In order to be nonlinearly unstable, such flows first require a candidate for transiently growing perturbations, which could be amplified by background shear by orders of magnitude for the Reynolds number R typical in astrophysical situations; second, in phase of high amplitude, those perturbations have to give a positive nonlinear feedback to the basin of small perturbations subject to transient growth. Two of these conditions together introduce what is usually called a bypass scenario for transition to turbulence.
The bypass scenario is best developed in application to plane shear flows and, in particular, to the Couette flow (see Hamilton et al. (1995) and Waleffe (1997), who proposed a nonlinear mechanism sustaining the turbulence). This mechanism is called the selfsustaining process (SSP) and is sketched, for example in Fig. 1 by Waleffe (1997). The transient growth of initially small streamwise rolls produces high amplitude streamwise streaks, the process usually called liftup effect after the work by Ellingsen & Palm (1975). As streamwise streaks modify the background velocity yielding an inflexion point in a spanwise direction, the flow becomes linearly unstable with respect to nonaxisymmetric modes, which in turn regenerate streamwise rolls via nonlinear interactions with each other.
However, the situation becomes less clear in centrifugally stable rotating shear flows. Let us define rotation number, which characterises a dynamical contribution of rotation to shear, R_{Ω} = −2/q, where q = −(r/Ω)dΩ/dr is the local dimensionless shear rate. Steady, nonlinear selfsustaining solutions obtained in the nonrotating case, for which R_{Ω} = 0, can be nonlinearly continued into both the centrifugally unstable regime, −1 < R_{Ω} < 0 (q > 2), and the cyclonic regime with high shear rate, R_{Ω} ≪ 1 (−q ≫ 1) (see Rincon et al. 2007). Nevertheless, this operation cannot be extended either until the centrifugally stable anticyclonic (or quasiKeplerian) regime with R_{Ω} ≤ −1 (0 < q < 2), or far into the cyclonic regime with R_{Ω} ≲ 1 (q ≲ −1). Presumably, this takes place because the liftup mechanism ceases to work in both cases. On the other hand, direct numerical simulations (DNS) allow us to obtain turbulent solutions until, at least, R_{Ω} ≈ 0.3 (see Lesur & Longaretti 2005; LL05 hereafter). It is interesting to note that such a value of rotation number is far beyond the interval of small positive R_{Ω} corresponding to a substantial transient growth of streamwise rolls (see Sect. 5.1.1 of this paper). Another notable fact is that the quasiKeplerian regime demonstrates a strong asymmetry in comparison with the cyclonic regime. Namely, DNS show dramatic stabilisation of the quasiKeplerian shear, as one goes away from the Rayleigh line, R_{Ω} = −1, to smaller rotation numbers (cf. Figs. 4 and 7 of LL05; see also the earlier DNS of the quasiKeplerian regime made by Hawley et al. 1999). We have already noted that when R_{Ω} = −1 exactly, the liftup mechanism does not work, however, there is an inverse process of small streamwise streaks growing into high amplitude streamwise rolls, which is called the antiliftup effect (see Sect. 2.3 in Rincon et al. 2008 for more detail). The regime R_{Ω} = −1 is neutrally stable with respect to local modal perturbations, that is, the growth of streaks can be regarded as epicyclic oscillations with an infinite period. Because finite amplitude streamwise rolls of the antiliftup are quite different from finite amplitude streamwise streaks of the liftup, Rincon et al. (2007) put forward several arguments that the turbulent solution should be represented by a bypass scenario other than SSP. In any case, since a bypass transition in the regime R_{Ω} = −1 should involve the antiliftup mechanism, it looks natural that LL05 observed a sharp increase of the transition Reynolds number R_{T} as one goes to R_{Ω} < −1, in contrast to the situation with the cyclonic regime. Indeed, the streak’s growth is strongly suppressed at R_{Ω} < −1, as it turns into epicyclic oscillations with nearly rotational frequency (see also Sect. 5.1.2 of this paper). LL05 observed a decay of initial turbulence in the ranges −∞ < R_{Ω} < −1.03 and 0.3 < R_{Ω} < ∞ for R ≲ 10^{5} within a local model of the shearing box with shearing boundary conditions. A bypass scenario proposed in this work naturally resolves the issues mentioned above.
There are a number of studies that examine properties of centrifugally stable flows either experimentally or employing DNS in the framework of the TaylorCouette problem (the flow between the coaxial corotating cylinders, see Sect. 6 of the review by Grossmann et al. 2016). Among the recent ones are Burin & Czarnocki (2012), who observed the turbulisation of the cyclonic flow at R_{T} ∼ 10^{5} in the case of R_{Ω} ≈ 0.8, and OstillaMónico et al. (2016), who considered the nonstationary dynamics in the same regime numerically excluding the influence of axial boundaries, which is unavoidable in the laboratory setup. At the same time, the quasiKeplerian regime was investigated in the experiment by Schartman et al. (2012) subsequently continued by Edlund & Ji (2014), and numerically by OstillaMónico et al. (2014) subsequently continued by Shi et al. (2017). All those studies demonstrated an outstanding stability of quasiKeplerian flow with respect to finite amplitude perturbations from R ∼ 10^{5} up to higher than R ∼ 10^{6}. Yet, this result should be accepted with a caution, because the presence of walls may prevent the onset of selfsustaining solutions inherent to the anticyclonic subcritical regime as was argued by Rincon et al. (2007). Additionally, let us emphasise that most of attempts to detect turbulence have been made at superKeplerian shear rates, 3/2 < q < 2. For example, Edlund & Ji (2014) and Shi et al. (2017) used q = 1.8 and q = 5/3, respectively^{1}. Indeed, there is a tendency to consider the shear rates starting from q = 2, which stems from the established opinion that the closer quasiKeplerian shear is to the Rayleigh line, R_{Ω} = −1 (q = 2), and the farther it is from the solidbody line, R_{Ω} → −∞ (q = 0), the less stable it becomes. However, an accompanying purpose of this work is to show that the most appropriate shear rate to test the transition to turbulence in the quasiKeplerian regime may be found in the vicinity of q = 1/2 (see Sect. 7).
Two decades ago it was established that shear flow turbulence is a result of the subtle interplay between nonnormal linear dynamics of perturbations and their nonlinear interaction (see for example the review by Grossmann 2000). Although the subcritical transition is in essence a loss of nonlinear stability, the phase of linear, nonmodal transient growth is of much importance in the general design of shear flow turbulence. Among others, this was elucidated by Henningson (1996). It was argued that (i) transient growth as linear process is the only source of energy for subcritical turbulence, (ii) transient growth must be large, that is the linear growth factor must be much greater than unity, for the transition to occur: for particular examples, the plane Couette and Poiseille flows may become turbulent at the Reynolds numbers as low as, respectively, R_{T} ≃ 350 and 1000, when according to Butler & Farrell (1992) and Reddy & Henningson (1993) the maximum growth factor of streamwise rolls is, respectively, ≃150 and 200. Also, a number of lowdimensional models of plane shear flow employing various representations of the nonlinear dynamics demonstrated qualitatively similar behaviour during the transition: substantial nonmodal growth followed by the “bootstrapping” effect (Trefethen et al. 1993; see Baggett & Trefethen 1997). These authors explained threshold exponents for transition to turbulence in each case using simple arguments based on the secondary linear instability of the transient perturbation attained the high amplitude. Thus, it was shown that transient dynamics can always be distinguished in the whole loop of the sustenance of turbulence (see caption for Fig. 3 by Baggett & Trefethen 1997)).
The reasoning just mentioned makes the transient growth a pivotal process responsible for the generation of the specific highamplitude inputs for the sequential nonlinear feedback. Such highamplitude inputs are the streamwise streaks in the case of plane shear flows, but what are they in the case of centrifugally stable homogeneous and boundless rotating shear flow? To answer this question in this paper, we (i) consider threedimensional (3D) perturbations, which are able to exhibit large transient growth, (ii) perform a critical test, which indicates that highamplitude perturbations generated during this transient growth have to do with the sustenance of turbulence. Issues (i) and (ii) are tackled solving, respectively, the linear problem (1)–(4) and the more general nonlinear problem (53)–(55).
It is known that shear flow allows for the existence of the socalled Kelvin modes, or alternatively, shearing vortices exhibiting the transient growth via the swing amplification^{2} (see Chagelishvili et al. 2003). The shearing vortices are spanwiseinvariant, nearly streamwise perturbations contracted by a background motion. In rotating flow, they take the form of spirals, whereas the instant of swing corresponds to the transformation of the leading spirals into trailing spirals or vice versa (see e.g. Razdoburdin & Zhuravlev 2015). Conservation of a spanwise vorticity perturbation is the physical reason for the growth of their velocity and pressure perturbations. The swing amplification is a twodimensional (2D) process, which takes place independently of whether the shear is plane or rotating: the Coriolis force does not affect the growth of shearing vortices in contrast to previously discussed rolls and streaks. However, shearing vortices grow considerably weaker rather than rolls and streaks, since they are much more suppressed by viscosity. The swing amplification gives the kinetic energy growth by a factor ∝ R^{2/3}, rather than ∝ R^{2} as in the case of lift and antiliftup. Mukhopadhyay et al. (2006) estimated that R ∼ 10^{6} is required to obtain the growth factor of shearing spirals comparable to the growth factor of rolls necessary to activate the SSP in the plane Couette flow. Remarkably, this claim is valid for both the plane and the rotating (no matter which cyclonicity) flows and does not change with the shear rate, since it can be shown that the maximum growth factor of shearing vortices G_{max} ≈ R^{2/3}exp ( − 2/3) provided that the time is measured in the inversed shear rate rather than in the inversed rotational frequency (cf. Eq. (83) by Afshordi et al. (2005) and Eq. (9) by Mukhopadhyay et al. 2005). Thus, the role of shearing vortices is negligible either at the Rayleigh line, or in absence of rotation, but as soon as R_{Ω} is either sufficiently smaller than −1, or sufficiently larger than 0, the growth of streaks and rolls vanishes and the shearing vortices become the kind of perturbations exhibiting the largest transient growth in the flow. This picture was confirmed for a variety of models by Ioannou & Kakouris (2001), Mukhopadhyay et al. (2005), Yecko (2004), and more recently by Maretzke et al. (2014), Zhuravlev & Razdoburdin (2014), and Razdoburdin & Zhuravlev (2017) who carried out studies of the transient dynamics of linear perturbations in the quasiKeplerian regime using rigorous optimisation methods, which allow us to obtain the optimal perturbations exhibiting the largest possible transient growth.
For very high Reynolds numbers typical in astrophysical boundary layers and the adjacent accretion discs, the growth of shearing spirals becomes quite significant. However, its role in the transition to turbulence remains unclear. On the one hand, it seems plausible that the shearing spirals are involved in 2D models of subcritical turbulence in shear flows simulated by Umurhan & Regev (2004), Johnson & Gammie (2005), and Horton et al. (2010). Furthermore, Lithwick (2007) suggested a 2D variant of positive nonlinear feedback sustaining the growth of shearing spirals. He showed that the coupling between shearing spirals at the instant of swing and smallamplitude axisymmetric vortices generates a new shearing spiral capable of transient growth. On the other hand, Lithwick (2009) revealed that the streamwise scale of a shearing spiral must exceed the disc scale height to make this process work in the 3D model, otherwise the shearing spiral is destroyed by resonant interaction with a 3D axisymmetric vortex (see also the simulations by Shen et al. 2006). Consequently, it seems reasonable that shearing spirals are only relevant to 2D turbulence, which may occur at scales above the accretion disc scale height: moreover, Razdoburdin & Zhuravlev (2017) have recently confirmed that shearing spirals keep the ability for significant transient growth even when their streamwise scale is comparable to the disc global radial scale. Also, the lack of the positive 3D feedback from shearing spirals at phase of swing is indirectly confirmed by the DNS of LL05, which show that R_{T} substantially depends on the shear rate for R_{Ω} ≲ 1 in the cyclonic case. In turn, this means that G_{max} strongly varies along the curve of R_{T}(q). Hence, positive nonlinear feedback would come into play at a different amplitude of working perturbation for different shear rates, which would require further clarification^{3}. Additionally, Maretzke et al. (2014) revealed a small difference between the optimal perturbations in cyclonic and quasiKeplerian regimes of the flow between the cylinders: in the cyclonic regime, the optimals are not shearing spirals but their slightly modified counterparts with a weak spanwise dependence. Coincidently or not, the nonlinear instability of the cyclonic flow between cylinders is already observed, rather than that of the quasiKeplerian flow. To summarise, we suppose that the shearing spirals considered as spanwise invariant perturbations seem unlikely to be an ingredient of the nonlinear feedback in the bypass scenario of subcritical transition to 3D turbulence on scales small compared to the disc scale height.
In this work we consider analytically the linear local dynamics of quasi2D shearing spirals in the rotating, homogeneous, boundless, and centrifugally stable shear flow. By a quasi2D shearing spiral, we mean a spatial Fourier harmonics (SFH) of small perturbation, which initially has a zero spanwise velocity perturbation, but a small nonzero spanwise wavenumber. Such perturbations are not optimal, meaning they do not attain the largest possible growth factor measured at the instant of swing. Instead, during the swing amplification they generate essentially 3D vortices of a new type, which we refer to as crossrolls, since those are aligned in the shearwise direction in contrast to ordinary streamwise rolls generated via the antiliftup mechanism on the Rayleigh line (see Fig. 9). We show that while the “planar” eddies constructed of shearwise and streamwise velocity perturbations decay as the shearing spiral is wound back by a background flow, the crossrolls conserve their amplitude and fade due to viscous dissipation only. The growth factor of crossrolls with respect to the growth factor of the shearing spiral at the instant of swing is less than unity and decreases when the shear rate alters from q = −∞ up to q = 2. Thus, the crossrolls amplitude is larger in the cyclonic regime, but smaller in the quasiKeplerian regime. Also, it vanishes both at the solidbody line and at the Rayleigh line. We note that the growth factor of crossrolls can be arbitrarily high provided that the initial quasi2D shearing spiral is sufficiently tightly wound. We suggest that the nonlinear instability and corresponding selfsustained turbulent solutions in the regime sufficiently far from both the centrifugally unstable interval and the nonrotating case are provided by the positive nonlinear feedback associated with either a secondary linear instability of the flow containing the finite amplitude crossrolls, or the interaction of crossrolls with each other. Either of these two processes may generate new small quasi2D shearing spirals capable of substantial swing amplification. Evidence for that comes from our simulations of subcritical transition to turbulence in the cyclonic regime. It turns out that in the shearing box model with periodic boundary conditions the transition is facilitated in the box stretched along the rotation axis, which indicates that perturbations with small but nonzero spanwise wavenumber are involved in that process. Indeed, the numerical domain of finite size does not affect the 2D shearing spirals while it works as an ideal highpass filter for harmonics with nonzero spanwise wavenumbers: it completely cuts off wavelengths larger than spanwise size of the numerical domain, whereas the first wavelength smaller than this size and satisfying the periodic boundary conditions is the best resolved by the numerical scheme. The inherent “flaw” of shearing box simulations becomes the useful tool in our case: the crossrolls with streamwise size comparable to the box streamwise size, so the least susceptible to viscous damping, are forbidden in a cubic box but allowed in a tall box. This way, in the tall box we manage to find the transition up to q = −3 (R_{Ω} ≈ 0.67). Furthermore, we find that R_{T} obtained in the tall box simulations corresponds approximately to the constant maximum growth factor of crossrolls in the range q ∼ −35 ÷ −3. Supposing that the transition to selfsustained turbulence in the quasiKeplerian regime corresponds to the same threshold growth factor of crossrolls, we predict the transition Reynolds number for 0 < q < 2 (see Fig. 8), which is R_{T} ≲ 10^{8} for the Keplerian shear rate. This is consistent with the results of LL05 concerning the asymmetry of R_{T}(q) in the vicinity of R_{Ω} = 0 and R_{Ω} = −1. Besides, it seems unrealistic to obtain the transition in the quasiKeplerian regime considering the superKeplerian shear rates. If our suggestion is true, the strategy of searching for hydrodynamical turbulence in the Keplerian flow should be changed: one should extend the turbulent solution from the cyclonic regime with high negative shear rate to the quasiKeplerian regime across the solidbody line, R_{Ω} = +∞ →R_{Ω} = −∞, rather than start from a turbulent solution obtained on the Rayleigh line, which has a different origin and properties. We predict that the most favourable value of the shear rate to test the existence of hydrodynamical turbulence in the quasiKeplerian regime is q = 1/2.
In the case of plane shear flow, R_{Ω} = 0, SFH with arbitrary nonzero spanwise wavenumber is described by an exact analytical solution obtained by Chagelishvili et al. (2016). At the limit of zero spanwise wavenumber, quasi2D shearing spirals turn into SFH, which are usually called shearing spirals. From now on, we omit prefix the “quasi2D”, so that by shearing spirals we mean all SFH with initially planar velocity perturbation, which may have a small nonzero spanwise wavenumber.
2. Linear perturbations in rotating shear flow
Local vortical 3D perturbations in viscous homogeneous rotating shear flow obey the following equations:
where u_{x}, u_{y}, and u_{z} are the Eulerian perturbations of velocity components, ρ_{0} is a constant background density, and p is the Eulerian perturbation of pressure. By definition,
and kinematic viscosity ν is assumed to be a constant. Variables x, y, and z are Cartesian coordinates, which locally correspond to radial, azimuthal, and axial directions in boundary layer (or disc), respectively, whereas Ω_{0} is angular velocity of fluid rotation at x = 0 corresponding to some radial distance r_{0} ≫ x. Also, q is a constant shear rate that defines the background azimuthal velocity as v_{y} = −qΩ_{0}x. Equations (1)–(4) are vertically unstratified linearised versions of small shearing box equations (see A.3 Umurhan & Regev 2004). Throughout the text, the variables x, y, z are alternatively referred to as, respectively, shearwise, streamwise, and spanwise coordinates (or directions). The shear rate is positive (negative) in the case of anticyclonic (cyclonic) flow. The anticyclonic case is usually divided into Rayleighstable regime 0 < q ≤ 2, also referred to as the quasiKeplerian regime with particular q = 3/2 corresponding to the Keplerian flow, and Rayleighunstable regime q > 2 subject to centrifugal linear instability. The solidbody line and the Rayleigh line are represented, respectively, by q = 0 and q = 2, whereas q → ±∞ are the limits corresponding to plane shear flow. In what follows, we are interested in the whole range −∞ < q < 2.
2.1. Dimensionless equations for shearing harmonics
We consider a boundless shear, which allows us to look for the partial solutions of Eqs. (1)–(4) in the form of SFH,
where f is any of perturbation quantities and f̂ is its Fourier amplitude specified by the dimensionless wavenumbers (k_{x}, k_{y}, k_{z}). By Eq. (5) it is implied that the problem is considered with respect to the dimensionless comoving Cartesian coordinates
where L is an auxiliary distance. Correspondingly, k_{x}, k_{y}, and k_{z} are expressed in units of L^{−1}. Let us assume from now on that k_{y} > 0 and L is chosen in such a way that k_{y} = 1. We omit the primes everywhere below and note that with respect to the original coordinates SFH acquires a varying radial wavenumber k̃_{x} ≡ k_{x} + qt. We arrive at the following dimensionless equations for SFH:
where velocities are expressed in units of Ω_{0}L and pressure normalised by a constant background density, W ≡ p/ρ_{0}, is expressed in units of L^{2}Ω_{0}^{2}. The nonzero kinematic viscosity yields a finite value of the Reynolds number
and
The normalisation used above is somewhat different from the accepted one (see e.g. LL05 and Mukhopadhyay et al. 2005). Specifically, the time is measured in units of Ω_{0}^{−1} rather than in units of (qΩ_{0})^{−1}. We choose the former variant, since it allows us to cross continuously the solidbody line.
Next, as soon as we are dealing with vortical dynamics, it is convenient to operate with equations for û_{x} and a shearwise component of vorticity SFH defined as ω̂ ≡ k × û:
3. Production of 3D vortices by shearing spirals: inviscid dynamics
We assume throughout this section that R → ∞.
3.1. Swing amplification
Here we revisit the wellknown mechanism of swing amplification of vortices, which takes place in the xyplane. Namely, in the case û_{z} = 0, k_{z} = 0, the problem (13)–(14) becomes onedimensional and the evolution of the corresponding 2D SFH is described by a single equation
provided that the velocity is divergencefree,
The solution of Eq. (15) is
where the normalisation k_{0} ≡ k(t = 0) is obtained from the condition of doubled kinetic energy density of 2D SFH equals to unity:
Accordingly, we have
The perturbation of pressure behaves in the following way:
The solution (17)–(18) conserves the spanwise component of the vorticity (see Chagelishvili et al. 2003):
which results in the transient growth of 2D SFH in the case of leading spirals, k_{x} < 0, for the quasiKeplerian regime, q > 0, and in the case of trailing spirals, k_{x} > 0, for the cyclonic regime, q < 0. Indeed, at the instant of swing, t_{s} ≡ −k_{x}/q, the shearwise wavenumber vanishes, k̃_{x} = 0, and their growth factor defined as the ratio of current energy density of 2D SFH to its initial energy density reads
For k_{x}≫1 g_{s} ≈ k_{x}^{2}, attaining an arbitrary high value for sufficiently wound spirals. The maximum possible g_{s} is limited by damping due to the microscopic viscosity of the fluid and will be evaluated in Sect. 4.
3.2. Dynamics of SFH with nonzero spanwise wavenumber
Let us consider SFH with small but nonzero k_{z}≪1^{4}. Without loss of generality, we assume that k_{z} > 0 from now on. We also assume that initially SFH has a planar velocity field, that is, û_{z}(t = 0) = 0. Combining (10) and (12) we formulate the initial condition in the form
It can be seen that Eq. (14) together with the initial condition Eq. (22) implies that ω̂_{x} ~ k_{z}, which makes the first term in the righthand side (RHS) of Eq. (13) as small as ∼O(k_{z}^{2}) comparing to the second term therein. Consequently, in order to obtain the solution for 3D SFH in the leading order in k_{z} ≪ 1, we drop the first term in the RHS of Eq. (13) going back to the 2D solution for û_{x} given by Eq. (15). Thus, for SFH weakly depending on the spanwise direction, û_{x}(t) is approximately given by Eq. (17).
Next, combining (10) and (12) we get an exact expression for û_{z} in terms of û_{x} and ω̂_{x}and :
Thereby, the term k_{z}û_{z} in the continuity Eq. (10) must be as small as ∼O(k_{z}^{2}) in comparison with the other two terms there. Thus, we assume that relation (16) holds for 3D SFH and û_{y} is approximately given by its 2D expression, (18). Of course, this implies that pressure also approximately obeys the 2D solution (19).
As we see, spanwise dynamics of 3D SFH is separated from the swing amplification dynamics described in Sect. 3.1 and ω̂_{x} is obtained from Eq. (14) by integrating the expression (17) over the time,
As t → ∞ and, respectively, k̃_{x} → ∞, the shearwise component of the vorticity perturbation tends to a constant nonzero value. We refer to this asymptote of SFH as a “plateau” of SFH. We conclude that along with the existing “planar” eddies, a set of other vortices associated with this nonzero limit of ω_{x} is generated as a “byproduct” of swing amplification. The latter will be referred to as “crossrolls” in this work. The crossrolls are directed across the stream in contrast to the streamwise rolls produced by the antiliftup mechanism in the flow on the Rayleigh line (see Fig. 9 and the details in Appendix A below). An important feature of the crossrolls is that the spanwise velocity perturbation gains a nonzero value at t → ∞, whereas the velocity perturbation components associated with the planar eddies vanish like û_{x} ∼ O(t^{−2}) and û_{y} ∼ O(t^{−1}). Equations (17) and (24) together with the relation (23) taken in the leading order in k_{z} ≪ 1 yield^{5}
Thus, Eq. (25) leads to the following asymptote of û_{z} at t → ∞,
As soon as the shearing spiral is wound up by the shear after the instant of swing, the ratio û_{z}/û_{y} increases, which implies that the crosssection of crossrolls becomes elongated in a spanwise direction. At the same time, crossrolls shrink along their axes, as far as the shearwise wavenumber increases. This, in turn, makes a streamwise component of the vorticity perturbation, ω_{y}, grow linearly with time. Since
Equation (26) demonstrates the basic features of crossrolls production. Clearly, their amplitude is proportional to pressure perturbation spanwise gradient at the instant of swing, which is ∝ k_{z}k_{x} (cf. Eq. (19)) at k̃_{x} = 0. Also, it decreases as the shear rate goes from q = −∞ to q = 2 and vanishes exactly at q = 2. Those features will be illustrated below in Fig. 1, where the growth factors of shearing spirals with nonzero spanwise number will be plotted for different shear rates.
Fig. 1.
Shortdashed, dotted, longdashed, solid, and dotdashed curves: g_{z}(k̃_{x}) (top panel) and g(k̃_{x}) (bottom panel) for transient growth of the particular SFH at q = −5, −1.5, 0.5, 1.5, and 1.6, respectively. In each case, k_{x} and k_{z} are chosen in such a way that SFH attains the constant value g_{z}(k̃_{x max}) = g_{z max} ≈ 180ϵ (see Sect. 3.3 for definition of ϵ). This value of g_{z max}/ϵ is obtained employing the procedure (ii) (see this section) in order to fit the transitional dependence R_{T}(q) found in the tall box simulations of turbulence in the cyclonic regime (see Sect. 6) for details and the dashed curve on the top panel in Figs. 6 and 8. Solid and dashed horizontal lines on the bottom panel mark the value of corresponding g_{z max} for ϵ = 1 and ϵ = 0.1, respectively. 

Open with DEXTER 
3.3. Restriction on value of spanwise wavenumber
Now we are a step away from estimating the maximum value of k_{z} = k_{z max}, which restricts the validity of the solution for 3D SFH obtained above in the leading order in k_{z} ≪ 1. At sufficiently large k_{z}, the swing amplification can no longer be separated from the spanwise dynamics, which starts to destruct the shearing spiral in the xyplane. Moreover, comparing the first term in the RHS of Eq. (13) to the second term therein, we find that these terms always become of the same order at a sufficiently large time for any small value of k_{z} provided that the first term contains ω̂_{x} given by Eq. (24), while the second term contains û_{x} given by Eq. (17). This takes place since û_{x} decreases as O(t^{−2}) long after the instant of swing. We must consequently provide a restriction on k_{z} estimating a higher order correction to the solution given in Sects. 3.1 and 3.2 at the selected period of time.
For that, let us assume that with the account of the first term in the RHS of Eq. (13), the 2D solution (17) acquires a small correction . Substituting
into Eq. (13), we derive the following equation for û′_{x}
Given the zero initial condition û′_{x}(t = 0) = 0, we find that
We estimate k_{z max} from the assumption that it corresponds to a constant ratio of the amplitude of and the amplitude of leading order û_{x} given by Eq. (17)
where k̃_{x max} is taken at some time after the instant of swing and corresponds to a certain location of SFH on its plateau. We note that the higher k̃_{x max} is, the smaller is k_{z max} and the longer is the duration of our analytical approximation. Most naturally, k̃_{x max} can be evaluated at the latest extremum of û_{z}(t), which emerges as one incorporates viscous damping into the problem.
4. Taking viscous damping into account
Primarily, viscous force confines the swing amplification of SFH as the initial winding of the spiral and so the duration of 2D transient growth become limited by viscous dissipation. Thus, the largest value of g_{s} corresponds to some k_{x} = k_{x max}. Additionally, viscous force causes damping of vertical motions and the asymptote (26) is replaced by extremum of û_{z} at some time after the instant of swing. Below we identify this time with k̃_{x} = k̃_{x max} used in Eq. (31) (see Sect. 5).
To determine both k_{x, max} and k̃_{x max} rigorously, let us suppose that û_{x}, ω̂_{x} is the solution of Eqs. (13)–(14) obtained in the inviscid limit R → ∞. Then, the solution of the same equations taken with the account of nonzero viscosity R < ∞ is represented by
where
For simplicity below we neglect by in comparison with unity in (32) since k_{z} ≲ 1 as far as q ϵ ≲ 1, which is assessed from the condition (31).
The spanwise velocity perturbation is modified in the same way:
where û_{z} is given by Eq. (25). This allows one to obtain k_{x, max} and k̃_{x max} jointly as the roots of the set of equations
provided that the initial SFH is assumed to be the leading spiral (k_{x max} < 0 and k̃_{x max} > 0) in the quasiKeplerian flow (q > 0) and, vice versa, the trailing spiral (k_{x max} > 0 and k̃_{x max} < 0) in cyclonic flow (q < 0).
4.1. The case of tightly wound SFH
In the limit of k_{x max}≫1 and k̃_{x max} ≫ 1, which is valid for sufficiently high R, the approximate relations
can be substituted into Eq. (25). Then, Eq. (34) simply recovers an approximate value of k_{x max}, which follows from the estimation of Razdoburdin & Zhuravlev (2017) of maximum duration of the swing amplification t_{max} (see their Eq. (B2)):
Thus, a very high transient growth of 2D SFH ∝ R^{2/3}, and so the amplitude of the crossrolls, can be reached taking into account that the Reynolds number may be as huge as 10^{10} in astrophysical boundary layers and discs.
Furthermore, Eq. (35) yields the estimation for k̃_{x max},
provided that R^{1/2} ≫ 1.
5. Maximum growth factor of crossrolls
The growth of the crossrolls can be measured by energy density stored in vertical motion. In this case, since the initial condition (22) sets SFH with the unit energy and the planar velocity field, the growth factor of the crossrolls is defined as
Then, the quantity
can be considered as an estimate of the largest possible amplification of crossrolls generated by shearing spirals. We note that g_{z max}/ϵ is the function of R and q. One may assume that violation of our analytical approximation discussed in Sect. 3.3 corresponds to the breakdown of swing amplification being a unique mechanism generating the crossrolls. Of course, an exact k_{z max} corresponding to the crossrolls of the largest amplitude can be obtained solving the full set of Eqs. (13)–(14) without an expansion over k_{z}. Nevertheless, it is known that k_{z max} should be at least k_{z max} ≲ q/(2 − q)^{6}, as the parameter β ≳ 1 introduced by Balbus & Hawley (2006) in their Eq. (39) substantially weakens the swing amplification and, consequently, the crossrolls amplitude. In this study, we intend to draw main conclusions assuming that Eq. (39) gives a reasonable estimate of the largest crossrolls generated by shearing spirals, whereas their exact identification from a strict solution of (13)–(14), as well as determination of the corresponding k_{z max}, is a subject for future studies with the use of an optimisation technique.
Below we determine the largest growth factor of crossrolls using accurate and approximate procedures.

i)
An accurate determination of g_{z max}(q, R) is done by solving Eqs. (34)–(35) numerically, which subsequently provides us with k_{z max} from the condition (31). Then, g_{z max} is obtained from Eq. (33) with k_{x} = k_{x max}, k̃_{x} = k̃_{x max}, and k_{z} = k_{z max}.

ii)
In the case of sufficiently high R_{gz max} can be estimated analytically. Indeed, as the conditions k_{x max} ≫ 1 and k̃_{x max} ≫ 1 are fulfilled, Eq. (30) simplifies to
and yields the following estimation for k_{z max}
Finally, in order to estimate the largest growth factor of crossrolls in this case, one must substitute Eqs. (36), (37), and (40) into
Now, setting g_{z max}/ϵ obtained according to the procedure (i) or (ii) to a constant value, we find R(q) as a root of Eq. (39). This also implies that k_{x max}(q), k̃_{x max}(q), and have been determined, so that we are ready to plot the growth factors of the particular SFH corresponding to constant g_{z max}/ϵ for several values of q (see Fig. 1). These growth factors have to be regarded as functions of time or, equivalently, k̃_{x}: so g(k̃_{x}) is defined by Eq. (21) and g_{z}(k̃_{x}) is defined by Eq. (38). As one can see on the top panel in Fig. 1, all curves of g_{z}(k̃_{x}) attain the same maximum value g_{z max}/ϵ, which occurs after the instant of swing at negative and positive k̃_{x} = k̃_{x max} for cyclonic and quasiKeplerian regimes, respectively. If it were not for viscous damping, g_{z} would tend to the plateau, which has been found previously for the inviscid solution (see Eq. (26)). While comparing two panels in Fig. 1, we see that the ratio (g_{s}/g_{z max})ϵ, where g_{s} is a factor of 2D swing amplification of SFH, dramatically (and monotonically) increases as one goes from large cyclonic shear rates across the solidbody line up to the marginal case of the Rayleigh line q = 2. For example, at the Keplerian shear rate this ratio attains almost 10^{2}. Thus, the crossrolls of the same intensity are generated by SFH with much larger k_{x} as one approaches the Rayleigh line. The larger is k_{x}, the tighter the winding of the initial shearing spiral is, the more it is susceptible to viscous dissipation at a typical time ∝ k_{x}^{−2}. Consequently, the closer one approaches the Rayleigh line, as starting from q = −∞, the higher the Reynolds number required in order to generate crossrolls of the same intensity. This is reflected in the fit of our simulations of the subcritical transition to turbulence in the cyclonic regime (see Sect. 6 and solid curve on the top panel in Fig. 6) as well as in our prediction of the subcritical transition to turbulence in the quasiKeplerian regime (see Sect. 7 and curve in Fig. 8).
Next, Fig. 2 illustrates the immediate cause of the suppression of the crossrolls generated in the lowshear (q∼1) cyclonic regime and even more in the quasiKeplerian regime in comparison with the highshear (q≫1) cyclonic regime. There are the profiles of pressure perturbation amplitude as a function of k̃_{x} plotted for shearing spirals shown in Fig. 1. For SFH these profiles also represent the behaviour of spanwise acceleration of fluid elements in the perturbed flow. At first, as the shear rate goes from −∞ up to 1/2, the net acceleration impulse in the spanwise direction decreases just because Ŵ in the vicinity of the instant of swing attains a smaller value. Next, at q = 1/2 the point Ŵ(k̃_{x} = 0) transforms from maximum to minimum and, further, for q > 1, Ŵ changes sign in the vicinity of the instant of swing. Eventually, as q → 2, the net acceleration of fluid elements in positive and negative spanwise directions balance each other, so that one is left with planar velocity perturbation as the shearing spiral is wound up by the shear after the instant of swing. Certainly, an additional dip in growth of the crossrolls, which must occur in the vicinity of the solidbody line q = 0, is imposed on the described picture. Therefore, there must exist a local maximum of growth of the crossrolls generated by the initial SFH with constant k_{x} and k_{z} in the quasiKeplerian regime. In Sect. 7 we show that this local maximum comes out through the local minimum of the transition Reynolds number occurring near the q = 1/2.
Fig. 2.
Shortdashed, dotted, longdashed, and solid curves: ℑ[Ŵ](k̃_{x}) according to Eq. (19) for q = −5, −1.5, 0.5, and 1.5, respectively. The growth factors of the corresponding SFHs are shown in Fig. 1. There is no curve for q = 1.6 since it virtually coincides with the curve for q = 1.5. 

Open with DEXTER 
5.1. Crossrolls in comparison with rolls and streaks
Let us estimate the growth factor of rolls generated by the antiliftup effect and streaks generated by the liftup effect in centrifugally stable rotating flow. This is to be done in order to evaluate the role of those known mechanisms relative to the generation of crossrolls discussed above.
Both liftup and antiliftup mechanisms describe the transient growth of axisymmetric perturbations. Consequently, we set k_{y} = 0 in Eqs. (7)–(10) and obtain the set of equations for û_{x} and ω̂_{x}
in the inviscid limit R → ∞, which yields the equation for û_{x}
where
with α ≡ k_{x}/k_{z}.
The solution of Eq. (44) reads
where amplitude A and phase φ are to be determined from initial conditions.
Further, ω̂_{x} is obtained with the help of Eq. (43)
whereas
since ω̂_{x} = k_{z}û_{y} for axisymmetric perturbations. Finally,
5.1.1. Growth of rolls
The liftup mechanism provides the transient growth of streamwise rolls, that is, initial vortices with û_{y} = 0. This implies that and for perturbations with initial . Thus, and φ = π/2.
Consequently, the growth factor of rolls defined as ḡ ≡ g + g_{z} reads
The corresponding maximum growth factor reads
and becomes unlimited as q → −∞, that is, in the nonrotating shear flow.
5.1.2. Growth of streaks
The antiliftup mechanism provides the transient growth of streamwise streaks, that is, initial vortices with û_{x} = 0, û_{y} = 1, û_{z} = 0. Thus, in this case , φ = 0.
The growth factor of streamwise streaks is
The corresponding maximum growth factor reads
and becomes unlimited as q → 2, on the Rayleigh line.
Hence, the maximum transient growth factor attained through either the liftup or antiliftup mechanism is independent of SFH wavenumbers. Since there is no low limit on their absolute values in the boundless shear, either the rolls or the streaks can be amplified up to ḡ_{max} given, respectively, by Eqs. (50) and (52) at any R. At the same time, the crossrolls cannot be amplified stronger than the value estimated by definition (39) with û_{z, v} substituted from Eq. (33). One can see that additionally to a constraint coming from the damping of the swing amplification, which is , the growth of the crossrolls is reduced by a factor of comparing to ḡ_{max}. We conclude that strictly in the absence of rotation or, alternatively, on the Rayleigh line, crossrolls are generally weaker than streaks or, respectively, rolls. However, when −∞ < q < 2, the transient growth of both rolls and streaks is highly suppressed. For example, the growth of rolls and streaks drops below ḡ_{max} ∼ 20, which is equal to g_{z max} represented in Fig. 1 at ϵ = 0.1, for q ≳ −40 and q ≲ 1.9, respectively. That is why in both the cyclonic and the quasiKeplerian regime of homogeneous boundless shear flow, crossrolls are the only high amplitude 3D perturbations that can be generated via the transient growth mechanisms.
6. Transition to turbulence: Cyclonic regime
In this section we consider 3D nonlinear dynamics of perturbations in the cyclonic regime of uniform and boundless rotating shear flow. Our goal is to recognize the subcritical transition to turbulence towards smaller q compared to the previous studies and to obtain the transition Reynolds number, R_{T}. This is done by performing hydrodynamical numerical simulations of the evolution of initial finite amplitude perturbations in the shearing box approximation. Importantly, we finally check how the crossrolls maximum growth factor estimated above behaves along the dependence R_{T}(q).
6.1. Model for nonlinear dynamics of perturbations and numerical details
Numerical simulations were carried out using the opensource code ATHENA for solving continuity and Euler equations (see Stone et al. (2008) and Stone & Gardiner (2010)),
supplemented by
where x, y, z, and t are the Cartesian coordinates and time used in Eqs. (1)–(4); v, P, and ρ are, respectively, velocity, pressure, and density of the perturbed flow. In Eqs. (53)–(55) kinematic viscosity ν and isothermal speed of sound c_{s} are constants and the full nonlinear dynamics of perturbations is considered in a box of size L_{x}, L_{y}, and L_{z} along x, y, and z directions, respectively. The usual periodic and shearingperiodic boundary conditions are imposed, respectively, along y, z, and x directions at each side of the box (see e.g. Hawley et al. 1995). We note that everywhere below L_{x} = L_{y} = 1 but L_{z} ≥ L_{x}, L_{y}.
In all our simulations
where H ≡ c_{s}/Ω_{0}. According to hydrodynamical equilibrium in accretion disc boundary layer, H is a typical scaleheight of the flow. The condition (56) ensures that we work with vortical dynamics of perturbations as long as vortical initial perturbations are imposed.
It is suitable to parametrize the problem by two dimensionless numbers: Mach number
and Reynolds number
We note that R_{nl} is different from R used above to describe the dynamics of the crossrolls in the boundless flow. The relation between both Reynolds numbers can be derived along the following line. Suppose that there is SFH with dimensional shearwise wavenumber K_{y}. In the box it takes values
where n is a natural number. On the other hand,
as far as in the analytical linear model of transient growth k_{y} = 1. Using the definitions (11) and (58), we find
for SFH with shearwise wavelength n times shorter than the shearwise box size.
In all simulations, the density of the unperturbed flow is set to ρ_{0} = 1 and the box rotation frequency Ω_{0} = 0.001. Setting q, M, and R_{nl} allows us to specify an unperturbed value of c_{s} and value of ν, as well as pressure in the unperturbed flow, P_{0}, through the equation of state (55). Numerical resolution is chosen to keep a cubic form of a gridcell, that is, the number of nodes in the x and y directions is equal to each other, N_{x} = N_{y} ≡ N, and the number of nodes in spanwise direction is N_{z} = l_{z}N, where the box ratio l_{z} ≡ L_{z}/L_{y}. We used N = 64 once in order to compare our results with those of LL05 for the case of q = −10. In all other simulations, N = 128.
6.2. Indication of transition to turbulence
In order to recognize turbulence, we employ a dimensionless angular momentum flux
arising due to velocity perturbation to the unperturbed value v_{y} = −qΩ_{0}x. In Eq. (61) v_{x} and v_{y}′ we introduce perturbations of shearwise and streamwise velocity components, respectively, and the brackets denote spatial averaging.
Angular momentum flux arising due to microscopic viscosity in unperturbed flow is
We compare α(t) with α_{ν} in the course of the nonlinear evolution of perturbations and interpret the instant enhanced transport of angular momentum α(t)> α_{ν} as the existence of turbulence at time t.
Further, in order to check the transition to turbulence in the flow for given R_{nl}, the timeaveraged variant of angular momentum flux is also used:
where qt_{1} = 500 Ω_{0}^{−1} and qt_{2} = 1000 Ω_{0}^{−1}. The transition to turbulence is defined by the condition
which implies that a selfsustained turbulent state exists for a sufficiently long time.
It was found that the evolution of initial perturbations crucially depends on R_{nl}. If R_{nl} exceeds some threshold value, which we refer to as transition Reynolds number, R_{T}, initial perturbations give rise to a selfsustained nonstationary solution. In the opposite case, R_{nl} < R_{T}, initial perturbations fade out (see Fig. 3). In order to determine R_{T}, we perform a set of simulations on a uniform grid of Reynolds number values. The transition Reynolds number is defined as R_{T} = (R_{turb} + R_{damp})/2, where R_{turb} equals the lowest R_{nl} at which the flow demonstrates the transition to turbulence and R_{damp} equals to the highest R_{nl} at which we do not observe a longlasting turbulent state of the flow and, consequently, the condition (64) is not satisfied. As we see in Fig. 4, the pass through R_{T} is associated with a significant jump of α̃.
Fig. 3.
Behaviour of a dimensionless angular momentum flux for initial vortex of perturbations defined by Eqs. (65)–(67) shown for q = −10 and M = 1. Solid and dashed lines correspond to R_{nl} = 1400 and R_{nl} = 1200, respectively. It was checked that turbulence for R_{nl} = 1400 does not decay at least during the three thousands shear times qt. The horizontal dotdashed line corresponds to the value of α_{ν} in this case. 

Open with DEXTER 
Fig. 4.
Timeaveraged angular momentum flux as function of Reynolds number for q = −10 and M = 1. Dotdashed line corresponds to value of α_{ν} in this case. The initial vortex of perturbations is determined by Eqs. (65)–(67). 

Open with DEXTER 
6.3. Initial conditions
At first, we check whether the subcritical transition to turbulence we observe in our numerical experiment is in accordance with the results of LL05. Performing simulations in a cubic box, they found R_{T} = 1200 at q = −10 for our definition of R_{nl} (LL05 measure time in units of the inversed shear rate, which makes their Reynolds number R_{nl} multiplied by q).
Generally, this is a feature of the subcritical turbulisation that R_{T} depends on the shape and amplitude of the initial perturbations. For example, this was shown experimentally by Darbyshire & Mullin (1995) and later numerically by Faisst & Eckhardt (2004) for pipe flow. One may assume that for any particular values of q and M there are optimal initial perturbations, which provide the lowest possible value of R_{T}. These optimals can be found by a variational approach to the nonlinear Cauchy problem (see Cherubini et al. 2010 and Pringle & Kerswell 2010 for such an example). Unfortunately, this method is not easy to implement and it is time consuming. In this work we restrict a comparison with LL05 to a single variant of initial perturbations. Also, we specify the same resolution as LL05, N_{x} = N_{y} = N_{z} = N ≡ 64 using the cubic box.
We choose an initial condition in the form of the vortex
where K_{y} and K_{z} are shearwise and spanwise dimensional wavenumbers of vortex, and the dimensionless constant A specifies the initial specific kinetic energy of perturbations, which equals A^{2}/2.
It was revealed that for l_{z} = 1, A = 1, and K_{y} = K_{z} = 2π the initial vortex (65)–(67) leads to selfsustained turbulence at R_{turb} = 1400 and damps at R_{damp} = 1200 (see Figs. 3 and 4. Thus, R_{T} = 1300, which is close to the value obtained by LL05. We additionally checked that both the decrease and the increase of the amplitude of initial vortex, respectively, up to A = 0.5 and A = 2 lead to the increase of the transition Reynolds number, respectively, up to R_{T} = 2600 and to R_{T} = 1500. Until the end of this section, we use the initial perturbations in the form (65)–(67) with A = 1, K_{y} = 2π, and K_{z} = 2π/l_{z}.
6.4. Results
The picture of the subcritical transition to turbulence we extract from our simulations is presented in Table 1. For the first time (with regards to the model of unbounded uniform hydrodynamical flow considered here), we use not only the cubic box, but also a tall box with the box ratio l_{z} > 1. The box acts as a filter with respect to perturbations of different wavelengths. It prevents the existence of both smallscale perturbations with a characteristic scale comparable and less than the size of a numerical gridcell, and largescale perturbations with a characteristic scale larger than the size of the box in any direction. The tall box allows perturbation harmonics with K_{z}/K_{y} = l_{z}^{−1} < 1 to take part in the nonlinear evolution of the perturbed flow and, probably, in the subcritical transition to turbulence. In particular, this applies to shearing spirals that effectively generate the crossrolls. Indeed, employing Eq. (40) one estimates K_{z}/K_{y} < k_{z max}/ϵ^{1/2} ≈ 0.8 for q = −10 and R_{T} = 1300. This ratio becomes smaller as one shifts to a smaller shear rate and higher Reynolds number. For example, at q = −6.67 and R_{T} = 10^{4}, which is the smallest shear rate examined by LL05 for the subcritical transition, its value is less than approximately 0.5. At the same time, the small size of the gridcell is always required in order for the tightly wound spirals exhibiting large transient growth and generating highamplitude crossrolls to be involved into simulations.
Summary of the numerical simulations for M = 1, N = 128.
The subcritical transition in the cubic box has been observed up to q = −4 (see Table 1 and Fig. 6). This is a higher value compared to q = −6.67 reached by LL05, which is due to a higher numerical resolution used in our simulations. Indeed, assuming that turbulisation, at least indirectly, is provided by the transient growth of shearing spirals and that the smallest shearwise scalelength of perturbations resolved in the numerical scheme comprises four gridcells, one can estimate the largest achievable R_{nl} provided that the shearing spirals subject to the largest swing amplification are resolved by the numerical scheme. The latter is done by Eq. (36),
For N = 128 and q = −4 Eq. (68) yields R_{nl} ≈ 330 000. This is close to the R_{T} obtained for q = −4 (see Table 1), which naturally explains why we did not see the subcritical transition at smaller q.
However, once we fix q = −4 and perform a set of simulations in the tall box, we find that R_{T} falls by a factor of ∼20 as l_{z} increases from 1.0 up to 1.5 ÷ 2 (see Fig. 5). We check that for l_{z} > 2.0, R_{T} tends to horizontal asymptote. Further, we managed to find the subcritical transition at even higher q = −3 in the box with the same resolution N = 128 but with the box ratio l_{z} = 2. We regard this result as evidence for the crucial role of harmonics of perturbations with certain 0 < (K_{z}/K_{y})_{T} < 1^{7}. It is plausible that R_{T}(l_{z}) acquires an asymptotic value as soon as l_{z}^{−1} ≤ (K_{z}/K_{y})_{T}. Such a harmonics already exists in the cubic box, however, it corresponds to n > 1 (see Eq. (59)), and, consequently, to the effective R smaller by a factor of n^{2} (see Eq. (60). Harmonics of a columnar shape corresponding to the highest effective R (i.e. to the largest shearwise wavelength equal to L_{y}), which implies the highest transient growth, emerges in the tall box only.
Fig. 5.
Transition Reynolds number R_{T} versus the box ratio l_{z} plotted for M = 1 and q = −4. 

Open with DEXTER 
The question arises, what can additionally testify in favour of crossrolls as moderators of subcritical transition? Trying to address this issue, we obtain R_{T} in the tall box with l_{z} = 2.0 for various q (see Table 1 and Fig. 6). It was additionally checked that the box ratio l_{z} sufficient for R_{T}(l_{z}) to approach the horizontal asymptote becomes smaller as one goes to higher q. This guarantees that subcritical transition observed in the tall box with l_{z} = 2.0 occurred in the presence of harmonics with (K_{z}/K_{y})_{T} in the range of the shear rate from q = −15 up to q = −3 covered by simulations. It was revealed that the difference between R_{T} obtained in cubic and tall boxes increases as q becomes smaller. The tall box result R_{T} = 16 250 at q = −4 allows us to evaluate the maximum growth factor of crossrolls according to the procedure (i) formulated in Sect. 5,
Fig. 6.
Transition Reynolds number R_{T} versus shear rate q in the cyclonic regime. Top panel: (1) squares represent DNS by LL05, see their Fig. 4; (2) triangles represent DNS performed in this work employing a cubic box; (3) circles represent DNS performed in this work employing a tall box with box ratio l_{z} = 2; (4) the solid curve shows R_{nl}(q) corresponding to constant maximum growth factor of the crossrolls given by equation (69) and obtained using the procedure (i) introduced in Sect. 5; (5) the dashed curve shows the same as described in (4), but with g_{z max}/ϵ obtained using the procedure (ii) introduced in Sect. 5; bottom panel: error bars for tall box DNS represented by circles on top panel are plotted in the range q = −10 ÷ −3 according to the values of R_{damp} and R_{turb} given in Table 1. The corresponding solid and dashed lines represent the same as on top panel but for g_{z max}/ϵ evaluated for R_{T}(q = −4)=15 000 (lower lines) and for R_{T}(q = −4)=17 500 (upper lines). 

Open with DEXTER 
Previous inspection of the maximum optimal growth factor, G_{max}, corresponding to transition to turbulence at R_{Ω} = 0 and R_{Ω} = −1 in the framework of the shearing sheet model (see e.g. Mukhopadhyay et al. 2005) has shown that G_{max} ∼ 150 in those cases, which coincides with Eq. (69) up to the smallness factor. Using procedure (i) further in order to look for the root of Eq. (39), we obtain the curve of constant (69) on the plane (R_{T}, q) (see the solid line in Fig. 6). Additionally, the dashed line in Fig. 6 represents the corresponding approximate R_{nl}(q) obtained according to procedure (ii) formulated in Sect. 5. A notable accordance is found between the growth of crossrolls and the subcritical transition: namely, the largest growth factor of crossrolls, as estimated in the main order in small k_{z}, remains almost constant along the transitional path, R_{T(q)}, deduced from our nonlinear stability analysis. Below we will refer to Eq. (69) as the threshold value of g_{z max}/ϵ. If there is a common nonlinear mechanism sustaining the subcritical turbulence in the uniform unbounded rotating shear flow of the cyclonic type, then it is natural to assume that it is triggered at the threshold value of the transient growth factor of relevant perturbations, at most weakly changing with the shear rate. Here we find such evidence in favour of the crossrolls and their maximum growth factor, g_{z max}, derived in the limit of small k_{z}. It can be suggested that the subcritical transition, which occurs at asymptotic value of R_{T(lz)} in the tall box simulations, takes place due to a positive nonlinear feedback from the harmonics of crossrolls with k_{z} ≃ (K_{z}/K_{y})_{T}. However, the ability of crossrolls with g_{z} ≳ g_{z max} to generate secondary shearing spirals subject to swing amplification sufficient to sustain turbulence remains to be shown.
We note that as q ≲ −5, the approximate analytical curve in Fig. 6 is no longer in accordance with the curve obtained accurately. The reason for that is shown in Fig. 7, where k_{x max}, k̃_{x max}, and along the corresponding curves of R_{nl}(q) are plotted. Clearly, the limit of tightly wound spiral, k_{x max}, k̃_{x max} ≫ 1, is valid as far as q∼1 and breaks when q ≲ −5. For q≫1 R_{nl}(q) corresponding to constant (69) becomes so small that spirals are no longer tightly wound (see Eqs. (36) and (37)). At the same time, we find that the first point of the subcritical transition obtained by LL05 in the presence of rotation, R_{Ω} > 0, which corresponds to q ≈ −33, is also in reasonable accordance with the threshold value (69). Along with the comparison of rolls and crossrolls growth factors with each other at q≫1 (see the end of Sect. 5.1), this argues in favour of the assumption that, at least, when q ≳ −35, SSP is already replaced by a different selfsustaining solution incorporating the growth of crossrolls.
Fig. 7.
Top panel, solid lines: roots of the set of Eqs. (34)–(35) along R_{nl}(q) plotted on the top panel in Fig. 6 by the solid line; dashed lines: Eqs. (36) and (37) along R_{nl}(q) plotted on the top panel in Fig. 6 by the dashed line. Bottom panel, solid line: obtained from the condition (31) with k_{x} = k_{x max} and k̃_{x max} along the corresponding solid lines on the top panel of this figure; dashed line: Eq. (40) along R_{nl}(q) plotted on the top panel in Fig. 6 by the dashed line. 

Open with DEXTER 
The dependence k_{z max}(q) obtained along the crossrolls threshold growth factor (69) (see Fig. 7) confirms that the latter is exhibited by shearing spirals with large spanwise length scale. As long as q becomes smaller, the destructive role of the Coriolis force in the evolution of SFH with nonzero k_{z} increases as compared to the pressure gradients enabling their swing amplification. That is why k_{z max} decreases as one approaches the line of rigid rotation, which is seen in Fig. 7 (see also Eq. (40)). In turn, this naturally explains why the difference between R_{T} obtained in cubic and tall boxes increases for smaller q: SFH with n = 1 (see Eq. (59)), corresponding to the growth factor (69) is progressively suppressed by periodic boundary conditions imposed in the cubic box in the spanwise direction. We note that SFH that generate the crossrolls of the largest amplitude will be highly columnar at q ≈1, as far as is estimated to become less than 0.1. If subcritical turbulence in the cyclonic regime is sustained by the nonlinear feedback from the crossrolls, it is worth checking its existence at q∼1 via the DNS in a very tall box.
7. Transition to turbulence: Towards the Keplerian shear across the solidbody line
Let us suppose that crossrolls are responsible for the subcritical transition in the whole range of centrifugally stable rotating shear flows. If so, one can assume that the corresponding threshold value of the crossrolls largest growth factor, g_{z max}, found in the cyclonic regime (see the previous section) remains independent of the shear rate further in the quasiKeplerian regime. Thus, we suppose that equality (69) holds for R_{nl} = R_{T} across the solid body line q = 0 up to the Rayleigh line q = 2. Since the shearing spirals, which produce the threshold growth factor, become increasingly wound up as one goes from q = −∞ to q = 2 (see the bottom panel in Fig. 1 and the top panel in Fig. 7), it is sufficient to use an approximate evaluation of g_{z max} provided by the procedure (ii) described in Sect. 5. So R_{nl}(q) corresponding to a constant threshold growth factor (69) extends to the quasiKeplerian regime as shown by the curve in Fig. 8. This curve may be considered as the prediction of the subcritical transition to turbulence in homogeneous boundless quasiKeplerian flow. Generally, the turbulence provided by the nonlinear feedback from the crossrolls is expected to occur at much higher Reynolds number ∼10^{7} ÷ 10^{9} in the quasiKeplerian regime when compared to the cyclonic regime. On the one hand, this is consistent with negative results from all existing studies of the nonlinear instability of the quasiKeplerian regime, whereas on the other, even such a huge Reynolds number is feasible in astrophysical discs (see e.g. Balbus 2003). Equations (40) and (26) show that such an asymmetry of growth of the crossrolls in the cyclonic and quasiKeplerian regimes is controlled by the factor
which stems from the suppression of fluid vertical motion in the course of the swing amplification of SFH (see comments to Fig. 2). Consequently, we predict an asymmetry of R_{T}(q) with respect to the line of rigid rotation: for the same absolute value of the shear rate, q, the subcritical transition in the quasiKeplerian regime should occur at the substantially higher R_{T} rather than in the cyclonic regime.
Fig. 8.
R_{nl}(q) along the constant maximum growth factor of the crossrolls given by Eq. (69), where g_{z max}/ϵ is obtained approximately using procedure (ii) introduced in Sect. 5. The circles represent DNS performed in this work using a tall box with box ratio l_{z} = 2 (see Sect. 6 for details). The curve is matched to the circle at q = −4, which corresponds to R_{T} = 16 250 (see Table 1). 

Open with DEXTER 
Let us note that the enhanced nonlinear stability of quasiKeplerian flow follows from our DNS in the cyclonic regime together with the results of LL05 and Rincon et al. (2007) discussed in Introduction. Indeed, the curves on top panel in Fig. 6 indicate that R_{T} ≈ 3 × 10^{5} at q = −1.5. At the same time, from extrapolation of DNS by LL05 carried out at R_{Ω} < −1, see expressions in their Fig. 7, it is expected that R_{T} ≳ 10^{9} at q = 1.5. This independently shows an asymmetry of R_{T}(q) with respect to change of sign of q. However, the swing amplification itself produces maximum growth factor , which would yield a symmetric R_{T} ∝ q^{−1} with respect to change from cyclonic to quasiKeplerian shear rate.
It is noteworthy that prediction of R_{T} in Fig. 8 combined with extrapolations of DNS suggested by LL05 for quasiKeplerian regime gives the most stable rotating flow located at superKeplerian shear rate q ∼ 1.7 ÷ 1.8.
Further, since in the vicinity of the solidbody line R_{T} naturally tends to infinity (cf. Eqs. (41) and (36), (40) at q→0), a minimal R_{T} is predicted for the quasiKeplerian regime. The corresponding most “favourable” value of q can be immediately estimated in the limit of tightly wound SFH. Using Eq. (41) in combination with Eqs. (36) and (40) and the condition R ≫ 1, we find that constant g_{z max} corresponds to
which yields minimal R_{T} at q = 1/2. So, R_{T}(q = 1/2) exceeds 10^{7}, while the transition at the Keplerian shear, q = 3/2, may occur at R approaching ∼10^{8}. Hence, the superKeplerian flow 3/2 < q < 2 is predicted to be the most nonlinearly stable throughout the whole range of centrifugally stable rotating shear flows. Notably, most of the efforts to detect turbulence in the quasiKeplerian regime in laboratory experiments and DNS have been made at superKeplerian shear rates. In contrast, our results suggest that one should extend the turbulent solutions obtained so far in the cyclonic regime to smaller q and then across the solidbody line, finally focusing the efforts on subKeplerian rotation. The effective numerical viscosity existing in DNS is approximately proportional to the distance between the nodes of the numerical grid. If so, the maximum Reynolds number accessible in DNS is ∝ N^{2}. Then, assuming that R_{T} ≈ 40 000 at q = −3, see Table 1, is the maximum accessible R_{nl} for N = 128 in tall box, we suppose that in order to be able to check hypothetical transition to hydrodynamic turbulence in quasiKeplerian flow, at least for “favourable” q = 1/2, N ≳ 2000 is required. This is quite a high DNS resolution, however, it is not unreachable for the most powerful modern supercomputers.
8. Summary
Angular momentum transport in boundary layers around the accreting weakly magnetised stars may be produced by turbulent fluctuations. Here we suggest that subcritical turbulence in boundless and homogeneous rotating shear flow of the cyclonic type, which is a local representative of the boundary layer, is supplied with energy by nearly optimal shearing spirals that generate essentially 3D vortices of high amplitude as a “byproduct” of swing amplification. We call those 3D vortices “crossrolls” in order to distinguish them from the known streamwise rolls generated via the antiliftup mechanism in a flow on the Rayleigh line. Such turbulence must be sustained by a novel type of positive nonlinear feedback from the crossrolls, which would recover the basin of growing shearing spirals. Furthermore, we suppose that this scenario also works in the quasiKeplerian regime. These suggestions are supported by the following evidence.

i)
As well as known spanwise invariant (columnar) shearing spirals, shearing spirals with small but nonzero spanwise wavenumber are subject to swing amplification, which causes their transient growth almost up to the highest possible values corresponding to a growth factor much greater that unity. Thus, these vortices satisfy a necessary condition for the subcritical transition in both the cyclonic and quasiKeplerian regimes. At the same time, it is known that all shearing harmonics with k_{z} ≳ 1 do not exhibit considerable transient growth since they drastically degenerate into epicyclic oscillations.

ii)
The substantial decrease of R_{T} and the existence of an asymptotic value of R_{T}(l_{z}) for large l_{z} found in the tall box simulations (see Fig. 5) indicates that these perturbations are crucial for the transition and sustenance of subcritical turbulence, at least in the cyclonic regime. This is the most striking evidence that the crossrolls could be ingredients of positive nonlinear feedback completing the bypass scenario in centrifugally stable rotating shear flows.

iii)
It is revealed that the maximum growth factor of crossrolls decreases as the shear rate alters from q = −∞ through the solidbody line q = 0 up to q = 2 (see Fig. 1). Moreover, its constant value accords with the curve of R_{T}(q < 0) obtained in tall box simulations (see Fig. 6). We regard this value as the threshold growth factor, which triggers the nonlinear sustenance of subcritical turbulence in the cyclonic regime. On the contrary, the spanwise invariant shearing spirals with k_{z} = 0 do not produce any 3D vortices and acquire a growth factor symmetric with respect to change between the cyclonic and quasiKeplerian regimes at the same q. However, DNS of turbulence in the cyclonic regime performed by LL05 and in this work far beyond the zone of the liftup effect show a strong variation of the 2D growth factor along the obtained R_{T}(q < 0).

iv)
The disappearance of the crossrolls in the vicinity of q = 2 naturally explains the asymmetry of the transition to turbulence revealed in DNS close to R_{Ω} = −1 versus R_{Ω} = 0: the bypass scenario based on the transient growth of shearing spirals with small but nonzero k_{z} for q < 0 replaces SSP, retaining the nonlinear instability of the cyclonic flow immediately after the liftup ceases to operate, whereas a similar situation with the antiliftup mechanism is not realised in the case of quasiKeplerian rotation.

v)
The crossrolls are the only highamplitude essentially 3D vortices that can be generated via the transient growth mechanism in the whole range of centrifugally stable shear rates, provided that the Reynolds number is sufficiently high. They create the vorticity (namely, its shearwise and streamwise components), which must be of crucial importance for 3D turbulence.

vi)
The transition to turbulence in the tall box is observed until smaller q than that in the cubic box at the same numerical resolution is obtained (see. Fig. 6). Generally, a sensitivity of the transition to the shape of a box indicates the anisotropic nature of subcritical turbulence. In future work, the properties of subcritical turbulence including details of the nonlinear sustaining process, as well as the precise role of the crossrolls in it, should be studied employing, in particular, 3D Fourier analysis in kspace (see Mamatsashvili et al. (2016) and Gogichaishvili et al. (2017) for an examples of such studies in application to nonrotating shear flow and Keplerian flow with an azimuthal magnetic field).

i)
The dynamics of perturbations in the plane orthogonal to the spanwise direction remain approximately the same as for spanwise invariant (columnar) shearing spirals, which is the 2D swing amplification process.

ii)
Spanwise acceleration of fluid elements is determined by the swing amplification.

iii)
Spanwise velocity perturbation increases while the shearing spiral swings from trailing to leading in the cyclonic flow and vice versa in the quasiKeplerian flow. After that, the spiral is wound by the shear, which causes the decay of streamwise and shearwise velocity perturbations, whereas the spanwise velocity perturbation tends to an asymptotic inviscid nonzero value proportional to the net impulse given to fluid elements in spanwise direction. The latter depends on the shear rate and vanishes as q → 2. Viscous damping results in a maximum value of the initial winding of shearing spiral and additionally in the maximum of spanwise velocity, which appears after the instant of swing (see Eqs. (36) and (37) giving the locations of those extrema in the limit of tightly wound spirals).

iv)
Assuming that the crossrolls of the largest amplitude sustain turbulence throughout the range of centrifugally stable rotating shear flows and the threshold value of their growth factor found in the cyclonic regime (see Eq. (69)) remains the same, we obtain a tentative prediction of R_{T}(q) for the quasiKeplerian regime, where hydrodynamical turbulence is not yet discovered. In particular, we find that the lowest R_{T} ∼ 10^{7} in the quasiKeplerian regime must be attained at a subKeplerian shear rate q ≈ 1/2.

v)
The threshold value of the crossrolls growth factor (69) exceeds the maximum growth factors of streaks and rolls generated via liftup and antiliftup up to the shear rates, very close, respectively, to nonrotating flow and to flow on the Rayleigh line (see Sect. 5.1).
The crossrolls’ appearance is shown in the top panels in Fig. 9 (see also Appendix A for the explanation of this representation). The growth factor of the corresponding shearing spiral is plotted by a solid curve in Fig. 1. For comparison, the rolls generated from axisymmetric streaks via antiliftup with the same k_{x} and k_{z} but in the flow close to the Rayleigh line are shown on the bottom panels in Fig. 9. The former and the latter ones are aligned in the shearwise and the streamwise directions, respectively. Furthermore, in contrast to rolls, the crossrolls are essentially 3D vortices, since they are of a finite length along their axes of rotation, being contracted by the background shear ∝ k̃_{x}^{−1}. In combination with the opposite rotation of the adjacent crossrolls, this causes monotonic growth of ω_{y} ∝ k̃_{x} in areas between the adjacent crossrolls, in the vicinity of for example point (x = λ_{x}/4, y = 0, z = λ_{z}/4). As one can see, the crossrolls coexist with “planar” eddies fading after the instant of swing. Such a complicated structure of the perturbation velocity field requires a new concept of nonlinear feedback, which may be completely different from both SSP and a scheme associated with streamwise rolls in the flow on the Rayleigh line. The latter itself remains unknown.
Fig. 9.
Left, middle, and right panels: crosssections of shearing spiral velocity pattern, respectively, in (y, z), (x, z), and (x, y)planes, which contain the point x = y = z = 0. The colour represents the velocity absolute value normalised by its maximum value separately for the top and bottom panels. Top panels: k_{x} = −141, k_{z} = 0.1, k̃_{x} = 43, and q = 1.5. The wavelengths are λ_{x} = 2π/k̃_{x}, λ_{y} = 2π, λ_{z} = 2π/k_{z}. Bottom panels: k_{x} = −141, k_{z} = 0.1, and t = π/(2σ), where σ is defined in Sect. 5.1; the shear rate is q = 1.99. The wavelengths are λ_{x} = 2π/k_{x}, λ_{y} = ∞, λ_{z} = 2π/k_{z}. The extraction of the shearing spiral velocity pattern is described in detail in Appendix A. 

Open with DEXTER 
The prediction of R_{T}(q) obtained in this work (see Fig. 8) makes it clear that, at least for boundless homogeneous flows, the issue of the transition to turbulence in the Keplerian regime should be attacked from the side of the cyclonic regime. This implies the crossing of the solidbody line, R_{Ω} = +∞ →R_{Ω} = −∞, which is not evident while looking at the axis of R_{Ω}. However, we suppose that R_{T} ∼ 10^{7} at the most “favourable” subKeplerian shear rate q ∼ 1/2, which presumably requires a resolution of N ≳ 2000 for the spectral Fourier code to be able to test the nonlinear instability of the quasiKeplerian regime. Still, such an extraordinary nonlinear stability of Keplerian flow may be defeated by the huge spatial scales of real astrophysical discs, leading to R = 10^{10} ÷ 10^{13}.
Especially, bearing in mind a result of Meseguer (2002) who revealed that subcritical transition observed in counterrotating flow between the cylinders occurs along the line of constant G_{max} on the plane of the Reynolds numbers corresponding to inner and outer cylinders (see Fig. 2 therein).
A rigorous restriction on the value of k_{z} including the possible smallness of q will be made after we have the corresponding modified solution at hand (see Sect. 3.3).
Cf. Eq. (40) below.
Acknowledgments
We kindly thank Professor K. Postnov and Professor P. Ivanov for genuine interest in the subject of this study and for a number of discussions during the preparation of the manuscript. The research is carried out using the equipment of the shared research facilities of HPC computing resources at Lomonosov Moscow State University (see paper Sadovnichy et al. 2013 for its detailed description). Equipment for the reported study was granted in part by the M. V. Lomonosov Moscow State University Programme of Development. DNR was supported by grant RSF 141200146 when writing Sect. 6 of this paper. VVZ was supported in part by RFBR grant 150208476 and by programme 7 of the Presidium of Russian Academy of Sciences.
References
 Abramowicz, M. A., Czerny, B., Lasota, J. P., & Szuszkiewicz, E. 1988, ApJ, 332, 646 [NASA ADS] [CrossRef] [Google Scholar]
 Afshordi, N., Mukhopadhyay, B., & Narayan, R. 2005, ApJ, 629, 373 [NASA ADS] [CrossRef] [Google Scholar]
 Baggett, J. S., & Trefethen, L. N. 1997, Phys. Fluids, 9, 1043 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Balbus, S. A. 2003, ARA&A, 41, 555 [NASA ADS] [CrossRef] [Google Scholar]
 Balbus, S. A., & Hawley, J. F. 2006, ApJ, 652, 1020 [NASA ADS] [CrossRef] [Google Scholar]
 Belyaev, M. A., Rafikov, R. R., & Stone, J. M. 2013, ApJ, 770, 67 [NASA ADS] [CrossRef] [Google Scholar]
 BisnovatyiKogan, G. S. 1994, MNRAS, 269, 557 [NASA ADS] [CrossRef] [Google Scholar]
 Burin, M. J., & Czarnocki, C. J. 2012, J. Fluid Mech., 709, 106 [NASA ADS] [CrossRef] [Google Scholar]
 Butler, K. M., & Farrell, B. F. 1992, Phys. Fluids A, 4, 1637 [NASA ADS] [CrossRef] [Google Scholar]
 Chagelishvili, G. D., Zahn, J.P., Tevzadze, A. G., & Lominadze, J. G. 2003, A&A, 402, 401 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Chagelishvili, G., Hau, J.N., Khujadze, G., & Oberlack, M. 2016, Phys. Rev. Fluids, 1, 043603 [NASA ADS] [CrossRef] [Google Scholar]
 Cherubini, S., de Palma, P., Robinet, J.C., & Bottaro, A. 2010, Phys. Rev. E, 82, 066302 [NASA ADS] [CrossRef] [Google Scholar]
 Darbyshire, A. G., & Mullin, T. 1995, J. Fluid Mech., 289, 83 [NASA ADS] [CrossRef] [Google Scholar]
 Edlund, E. M., & Ji, H. 2014, Phys. Rev. Lett., 89, 021004 [NASA ADS] [Google Scholar]
 Ellingsen, T., & Palm, E. 1975, Phys. Fluids, 18, 487 [NASA ADS] [CrossRef] [Google Scholar]
 Faisst, H., & Eckhardt, B. 2004, J. Fluid Mech., 504, 343 [NASA ADS] [CrossRef] [Google Scholar]
 Gogichaishvili, D., Mamatsashvili, G., Horton, W., Chagelishvili, G., & Bodo, G. 2017, ApJ, 845, 70 [NASA ADS] [CrossRef] [Google Scholar]
 Grossmann, S. 2000, Rev. Mod. Phys., 72, 603 [NASA ADS] [CrossRef] [Google Scholar]
 Grossmann, S., Lohse, D., & Sun, C. 2016, Ann. Rev. Fluid Mech., 48, 150724171740009 [CrossRef] [Google Scholar]
 Hamilton, J. M., Kim, J., & Waleffe, F. 1995, J. Fluid Mech., 287, 317 [NASA ADS] [CrossRef] [Google Scholar]
 Hawley, J. F., Gammie, C. F., & Balbus, S. A. 1995, ApJ, 440, 742 [NASA ADS] [CrossRef] [Google Scholar]
 Hawley, J. F., Balbus, S. A., & Winters, W. F. 1999, ApJ, 518, 394 [NASA ADS] [CrossRef] [Google Scholar]
 Henningson, D. 1996, Phys. Fluids, 8, 2257 [NASA ADS] [CrossRef] [Google Scholar]
 Horton, W., Kim, J.H., Chagelishvili, G. D., Bowman, J. C., & Lominadze, J. G. 2010, Phys. Rev. E, 81, 066304 [NASA ADS] [CrossRef] [Google Scholar]
 Inogamov, N. A., & Sunyaev, R. A. 1999, Astron. Lett., 25, 269 [NASA ADS] [Google Scholar]
 Ioannou, P. J., & Kakouris, A. 2001, ApJ, 550, 931 [NASA ADS] [CrossRef] [Google Scholar]
 Johnson, B. M., & Gammie, C. F. 2005, ApJ, 635, 149 [NASA ADS] [CrossRef] [Google Scholar]
 Lesur, G., & Longaretti, P.Y. 2005, A&A, 444, 25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lithwick, Y. 2007, ApJ, 670, 789 [NASA ADS] [CrossRef] [Google Scholar]
 Lithwick, Y. 2009, ApJ, 693, 85 [NASA ADS] [CrossRef] [Google Scholar]
 Mamatsashvili, G., Khujadze, G., Chagelishvili, G., et al. 2016, Phys. Rev. E, 94, 023111 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Maretzke, S., Hof, B., & Avila, M. 2014, J. Fluid Mech., 742, 254 [NASA ADS] [CrossRef] [Google Scholar]
 Meseguer, Á. 2002, Phys. Fluids, 14, 1655 [NASA ADS] [CrossRef] [Google Scholar]
 Mukhopadhyay, B., Afshordi, N., & Narayan, R. 2005, ApJ, 629, 383 [NASA ADS] [CrossRef] [Google Scholar]
 Mukhopadhyay, B., Afshordi, N., & Narayan, R. 2006, Adv. Space Res., 38, 2877 [NASA ADS] [CrossRef] [Google Scholar]
 Narayan, R., & Popham, R. 1993, Nature, 362, 820 [NASA ADS] [CrossRef] [Google Scholar]
 OstillaMónico, R., Verzicco, R., Grossmann, S., & Lohse, D. 2014, J. Fluid Mech., 748, R3 [CrossRef] [Google Scholar]
 OstillaMónico, R., Verzicco, R., & Lohse, D. 2016, J. Fluid Mech., 799, R1 [CrossRef] [Google Scholar]
 Philippov, A. A., Rafikov, R. R., & Stone, J. M. 2016, ApJ, 817, 62 [NASA ADS] [CrossRef] [Google Scholar]
 Popham, R., & Sunyaev, R. 2001, ApJ, 547, 355 [NASA ADS] [CrossRef] [Google Scholar]
 Popham, R., Narayan, R., Hartmann, L., & Kenyon, S. 1993, ApJ, 415, L127 [NASA ADS] [CrossRef] [Google Scholar]
 Pringle, C. C. T., & Kerswell, R. R. 2010, Phys. Rev. Lett., 105, 154502 [NASA ADS] [CrossRef] [Google Scholar]
 Razdoburdin, D. N., & Zhuravlev, V. V. 2015, Phys. Usp., 58, 1031 [NASA ADS] [CrossRef] [Google Scholar]
 Razdoburdin, D. N., & Zhuravlev, V. V. 2017, MNRAS, 467, 849 [NASA ADS] [CrossRef] [Google Scholar]
 Reddy, S. C., & Henningson, D. S. 1993, J. Fluid Mech., 252, 209 [NASA ADS] [CrossRef] [Google Scholar]
 Rincon, F., Ogilvie, G. I., & Cossu, C. 2007, A&A, 463, 817 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rincon, F., Ogilvie, G. I., Proctor, M. R. E., & Cossu, C. 2008, Astron. Nachr., 329, 750 [NASA ADS] [CrossRef] [Google Scholar]
 Sadovnichy, V., Tikhonravov, A., Voevodin, V., & Opanasenko, V. 2013, Contemporary High Performance Computing: From Petascale toward Exascale , (Boca Raton, USA: Chapman & Hall/CRC Computational Science), 283 [Google Scholar]
 Schartman, E., Ji, H., Burin, M. J., & Goodman, J. 2012, A&A, 543, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
 Shakura, N. I., & Sunyaev, R. A. 1988, Adv. Space Res., 8, 135 [NASA ADS] [CrossRef] [Google Scholar]
 Shen, Y., Stone, J. M., & Gardiner, T. A. 2006, ApJ, 653, 513 [NASA ADS] [CrossRef] [Google Scholar]
 Shi, L., Hof, B., Rampp, M., & Avila, M. 2017, Phys. Fluids, 29, 044107 [NASA ADS] [CrossRef] [Google Scholar]
 Stone, J. M., & Gardiner, T. A. 2010, ApJS, 189, 142 [NASA ADS] [CrossRef] [Google Scholar]
 Stone, J. M., Gardiner, T. A., Teuben, P., Hawley, J. F., & Simon, J. B. 2008, ApJS, 178, 137 [NASA ADS] [CrossRef] [Google Scholar]
 Trefethen, L. N., Trefethen, A. E., Reddy, S. C., & Driscoll, T. A. 1993, Science, 261, 578 [NASA ADS] [CrossRef] [MathSciNet] [PubMed] [Google Scholar]
 Umurhan, O. M., & Regev, O. 2004, ApJ, 427, 855 [NASA ADS] [Google Scholar]
 Waleffe, F. 1997, Phys. Fluids, 9, 883 [NASA ADS] [CrossRef] [Google Scholar]
 Yecko, P. A. 2004, A&A, 425, 385 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Zhuravlev, V. V., & Razdoburdin, D. N. 2014, MNRAS, 442, 870 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Representation of crossrolls
In order to plot the velocity perturbation in Fig. 9, we take the imaginary part of SFH given by Eq. (5), where the Fourier amplitudes of the velocity perturbations, û_{x,y,z}, are given by the analytic solution (17,18, 25). Explicitly, we have
where
Therefore, SFH is nothing but a combination of the single velocity pattern, C_{1}, with its duplicates translated along the coordinate axes on onefourth of the corresponding wavelength. This velocity pattern has the form
and is represented in Fig. 9 by the three mutually perpendicular crosssections having the common point {x = 0, y = 0, z = 0}.
Appendix B: Effect of finite numerical resolution
In order to reveal the influence of the numerical viscosity on the transition to turbulence at terminally high R_{nl}, we fix q = −4 and perform a set of additional simulations with exactly the same setup as described in Sect. 6. The results are shown in table B.1.
Summary of additional numerical simulations for M = 1 and q = −4.
The value of R_{T} obtained in cubic box substantially drops as we increase resolution up to N = 192. This is to be expected, since at N = 128 R_{T} is close to the restriction (68). However, an even higher resolution, N = 256, does not lead to further decrease of R_{T}: we conclude that convergence is achieved in this case. At the same time, R_{T} in tall box is only weakly affected by the increase of resolution. Thus, there is a difference between the transition Reynolds numbers in cubic and tall boxes, which cannot be attributed to effect of finite resolution.
All Tables
All Figures
Fig. 1.
Shortdashed, dotted, longdashed, solid, and dotdashed curves: g_{z}(k̃_{x}) (top panel) and g(k̃_{x}) (bottom panel) for transient growth of the particular SFH at q = −5, −1.5, 0.5, 1.5, and 1.6, respectively. In each case, k_{x} and k_{z} are chosen in such a way that SFH attains the constant value g_{z}(k̃_{x max}) = g_{z max} ≈ 180ϵ (see Sect. 3.3 for definition of ϵ). This value of g_{z max}/ϵ is obtained employing the procedure (ii) (see this section) in order to fit the transitional dependence R_{T}(q) found in the tall box simulations of turbulence in the cyclonic regime (see Sect. 6) for details and the dashed curve on the top panel in Figs. 6 and 8. Solid and dashed horizontal lines on the bottom panel mark the value of corresponding g_{z max} for ϵ = 1 and ϵ = 0.1, respectively. 

Open with DEXTER  
In the text 
Fig. 2.
Shortdashed, dotted, longdashed, and solid curves: ℑ[Ŵ](k̃_{x}) according to Eq. (19) for q = −5, −1.5, 0.5, and 1.5, respectively. The growth factors of the corresponding SFHs are shown in Fig. 1. There is no curve for q = 1.6 since it virtually coincides with the curve for q = 1.5. 

Open with DEXTER  
In the text 
Fig. 3.
Behaviour of a dimensionless angular momentum flux for initial vortex of perturbations defined by Eqs. (65)–(67) shown for q = −10 and M = 1. Solid and dashed lines correspond to R_{nl} = 1400 and R_{nl} = 1200, respectively. It was checked that turbulence for R_{nl} = 1400 does not decay at least during the three thousands shear times qt. The horizontal dotdashed line corresponds to the value of α_{ν} in this case. 

Open with DEXTER  
In the text 
Fig. 4.
Timeaveraged angular momentum flux as function of Reynolds number for q = −10 and M = 1. Dotdashed line corresponds to value of α_{ν} in this case. The initial vortex of perturbations is determined by Eqs. (65)–(67). 

Open with DEXTER  
In the text 
Fig. 5.
Transition Reynolds number R_{T} versus the box ratio l_{z} plotted for M = 1 and q = −4. 

Open with DEXTER  
In the text 
Fig. 6.
Transition Reynolds number R_{T} versus shear rate q in the cyclonic regime. Top panel: (1) squares represent DNS by LL05, see their Fig. 4; (2) triangles represent DNS performed in this work employing a cubic box; (3) circles represent DNS performed in this work employing a tall box with box ratio l_{z} = 2; (4) the solid curve shows R_{nl}(q) corresponding to constant maximum growth factor of the crossrolls given by equation (69) and obtained using the procedure (i) introduced in Sect. 5; (5) the dashed curve shows the same as described in (4), but with g_{z max}/ϵ obtained using the procedure (ii) introduced in Sect. 5; bottom panel: error bars for tall box DNS represented by circles on top panel are plotted in the range q = −10 ÷ −3 according to the values of R_{damp} and R_{turb} given in Table 1. The corresponding solid and dashed lines represent the same as on top panel but for g_{z max}/ϵ evaluated for R_{T}(q = −4)=15 000 (lower lines) and for R_{T}(q = −4)=17 500 (upper lines). 

Open with DEXTER  
In the text 
Fig. 7.
Top panel, solid lines: roots of the set of Eqs. (34)–(35) along R_{nl}(q) plotted on the top panel in Fig. 6 by the solid line; dashed lines: Eqs. (36) and (37) along R_{nl}(q) plotted on the top panel in Fig. 6 by the dashed line. Bottom panel, solid line: obtained from the condition (31) with k_{x} = k_{x max} and k̃_{x max} along the corresponding solid lines on the top panel of this figure; dashed line: Eq. (40) along R_{nl}(q) plotted on the top panel in Fig. 6 by the dashed line. 

Open with DEXTER  
In the text 
Fig. 8.
R_{nl}(q) along the constant maximum growth factor of the crossrolls given by Eq. (69), where g_{z max}/ϵ is obtained approximately using procedure (ii) introduced in Sect. 5. The circles represent DNS performed in this work using a tall box with box ratio l_{z} = 2 (see Sect. 6 for details). The curve is matched to the circle at q = −4, which corresponds to R_{T} = 16 250 (see Table 1). 

Open with DEXTER  
In the text 
Fig. 9.
Left, middle, and right panels: crosssections of shearing spiral velocity pattern, respectively, in (y, z), (x, z), and (x, y)planes, which contain the point x = y = z = 0. The colour represents the velocity absolute value normalised by its maximum value separately for the top and bottom panels. Top panels: k_{x} = −141, k_{z} = 0.1, k̃_{x} = 43, and q = 1.5. The wavelengths are λ_{x} = 2π/k̃_{x}, λ_{y} = 2π, λ_{z} = 2π/k_{z}. Bottom panels: k_{x} = −141, k_{z} = 0.1, and t = π/(2σ), where σ is defined in Sect. 5.1; the shear rate is q = 1.99. The wavelengths are λ_{x} = 2π/k_{x}, λ_{y} = ∞, λ_{z} = 2π/k_{z}. The extraction of the shearing spiral velocity pattern is described in detail in Appendix A. 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.