1. Introduction
Turbulence in wall-bounded shear flows in the transitional regime is characterised by coexisting turbulent and laminar regions, with the turbulent fraction increasing with Reynolds number. This phenomenon was first described by Coles & van Atta (Reference Coles and van Atta1966) and Andereck, Liu & Swinney (Reference Andereck, Liu and Swinney1986) in Taylor–Couette flow. Later, by constructing Taylor–Couette and plane Couette experiments with very large aspect ratios, Prigent et al. (Reference Prigent, Grégoire, Chaté, Dauchot and van Saarloos2002, Reference Prigent, Grégoire, Chaté and Dauchot2003) showed that these coexisting turbulent and laminar regions, called bands and gaps, respectively, spontaneously formed regular patterns with a selected wavelength and orientation that depend systematically on 
$Re$. These patterns have been simulated numerically and studied intensively in plane Couette flow (Barkley & Tuckerman Reference Barkley and Tuckerman2005, Reference Barkley and Tuckerman2007; Duguet, Schlatter & Henningson Reference Duguet, Schlatter and Henningson2010; Rolland & Manneville Reference Rolland and Manneville2011; Tuckerman & Barkley Reference Tuckerman and Barkley2011), plane Poiseuille flow (Tsukahara et al. Reference Tsukahara, Seki, Kawamura and Tochio2005; Tuckerman et al. Reference Tuckerman, Kreilos, Schrobsdorff, Schneider and Gibson2014; Shimizu & Manneville Reference Shimizu and Manneville2019; Kashyap Reference Kashyap2021) and Taylor–Couette flow (Dong Reference Dong2009; Meseguer et al. Reference Meseguer, Mellibovsky, Avila and Marques2009; Wang et al. Reference Wang, Ayats, Deguchi, Mellibovsky and Meseguer2022).
In pipe flow, the other canonical wall-bounded shear flow, only the streamwise direction is long, and transitional turbulence takes the form of puffs, also called flashes (Reynolds Reference Reynolds1883; Wygnanski & Champagne Reference Wygnanski and Champagne1973), which are the one-dimensional analogues of turbulent bands. In contrast to bands in planar shear flows, experiments and direct numerical simulations show that puffs do not spontaneously form spatially periodic patterns (Moxey & Barkley Reference Moxey and Barkley2010; Avila & Hof Reference Avila and Hof2013). Instead, the spacing between them is dictated by short-range interactions (Hof et al. Reference Hof, De Lozar, Avila, Tu and Schneider2010; Samanta, De Lozar & Hof Reference Samanta, De Lozar and Hof2011). Puffs have been studied extensively, especially in the context of the model derived by Barkley (Reference Barkley2011a,Reference Barkleyb, Reference Barkley2016) from the viewpoint of excitable media. In this framework, fluctuations from uniform turbulence trigger quasi-laminar gaps (i.e. low-turbulent-energy holes within the flow) at random instants and locations, as has been seen in direct numerical simulations of pipe flow. The bifurcation scenario giving rise to localised gaps has been investigated by Frishman & Grafke (Reference Frishman and Grafke2022), who called them anti-puffs. Interestingly, spatially periodic solutions like those observed in planar shear flows are produced in a centro-symmetric version of the Barkley model (Barkley Reference Barkley2011b), although the mechanism for their formation has not yet been clarified.
In this paper, we will show that in plane Couette flow, as in pipe flow, short-lived localised gaps emerge randomly from uniform turbulence at the highest Reynolds numbers in the transitional range, which we will see is 
$Re\simeq 470$ in the domain that we will study. The first purpose of this paper is to investigate these gaps. The emblematic regular oblique large-scale bands appear at slightly lower Reynolds numbers, which we will see is 
$Re\simeq 430$.
If the localised gaps are disregarded, then it is natural to associate the bands with a pattern-forming instability of the uniform turbulent flow. This was first suggested by Prigent et al. (Reference Prigent, Grégoire, Chaté and Dauchot2003) and later investigated by Rolland & Manneville (Reference Rolland and Manneville2011). Manneville (Reference Manneville2012) and Kashyap (Reference Kashyap2021) proposed a Turing mechanism to account for the appearance of patterns by constructing a reaction–diffusion model based on an extension of the Waleffe (Reference Waleffe1997) model of the streak–roll self-sustaining process. Reetz, Kreilos & Schneider (Reference Reetz, Kreilos and Schneider2019) discovered a sequence of bifurcations leading to a large-scale steady state that resembles a skeleton for the banded pattern, arising from tiled copies of the exact Nagata (Reference Nagata1990) solutions of plane Couette flow. The relationship between these pattern-forming frameworks and local nucleation of gaps is unclear.
The adaptation of classic stability concepts to turbulent flows is currently a major research topic. At the simplest level, it is always formally possible to carry out linear stability analysis of a mean flow, as was done by Barkley (Reference Barkley2006) for a limit cycle in the cylinder wake. The mean flow of uniformly turbulent plane Couette flow has been found to be linearly stable (Tuckerman, Barkley & Dauchot Reference Tuckerman, Barkley and Dauchot2010). However, this procedure makes the drastic simplification of neglecting the Reynolds stress entirely in the stability problem, hence its interpretation is uncertain (e.g. Bengana & Tuckerman Reference Bengana and Tuckerman2021). The next level of complexity and accuracy is to represent the Reynolds stress via a closure model. However, classic closure models for homogeneous turbulence (e.g. 
$(K,\varOmega )$) have yielded predictions that are completely incompatible with results from full numerical simulation or experiment (Tuckerman et al. Reference Tuckerman, Barkley and Dauchot2010). Another turbulent configuration in which large, spatially periodic scales emerge are zonal jets, characteristic of geophysical turbulence. For zonal jets, a closure model provided by a cumulant expansion (Srinivasan & Young Reference Srinivasan and Young2012; Tobias & Marston Reference Tobias and Marston2013) has led to a plausible stability analysis (Parker & Krommes Reference Parker and Krommes2013). Other strategies are possible for turbulent flows in general: Kashyap, Duguet & Dauchot (Reference Kashyap, Duguet and Dauchot2022) examined the averaged time-dependent response of uniform turbulence to large-wavelength perturbations, and provided evidence for a linear instability in plane channel flow. They computed a dispersion relation that is in good agreement with the natural spacing and angle of patterns.
Classic analyses for non-turbulent pattern-forming flows, such as Rayleigh–Bénard convection or Taylor–Couette flow, yield not only a threshold and a preferred wavelength, but also existence and stability ranges for other wavelengths through the Eckhaus instability (Busse Reference Busse1981; Ahlers et al. Reference Ahlers, Cannell, Dominguez-Lerma and Heinrichs1986; Riecke & Paap Reference Riecke and Paap1986; Tuckerman & Barkley Reference Tuckerman and Barkley1990; Cross & Greenside Reference Cross and Greenside2009). As the control parameter is varied, this instability causes spatially periodic states to make transitions to other periodic states whose wavelength is preferred. The Eckhaus instability is also invoked in turbulent zonal jets (Parker & Krommes Reference Parker and Krommes2013). The second goal of this paper is to study the regular patterns of transitional plane Couette flow and to determine the wavelengths at which they can exist and thrive. At low enough Reynolds numbers, patterns will be shown to destabilise and to acquire a different wavelength.
Pattern formation is sometimes associated with maximisation principles obeyed by the preferred wavelength, as in the canonical Rayleigh–Bénard convection. Such principles, like maximal dissipation, also have a long history for turbulent solutions. Malkus (Reference Malkus1954) and Busse (Reference Busse1981) proposed a principle of maximal heat transport, or equivalently maximal dissipation, obeyed by convective turbulent states. The maximal dissipation principle, as formulated by Malkus (Reference Malkus1956) in shear flows, occurs in other systems, such as von Kármán flow (Ozawa, Shimokawa & Sakuma Reference Ozawa, Shimokawa and Sakuma2001; Mihelich et al. Reference Mihelich, Faranda, Paillard and Dubrulle2017). (This principle has been somewhat controversial and was challenged by Reynolds & Tiederman (Reference Reynolds and Tiederman1967) within the context of stability theory. See a modern revisit of Malkus stability theory with statistical closures by Markeviciute & Kerswell (Reference Markeviciute and Kerswell2022).) Using the energy analysis formulated in our companion paper (Gomé, Tuckerman & Barkley Reference Gomé, Tuckerman and Barkley2023), we will associate the selected wavelength to a maximal dissipation observed for the large-scale flow along the bands.
2. Numerical set-up
Plane Couette flow consists of two parallel rigid plates moving at different velocities, here equal and opposite velocities 
$\pm U_{wall}$. Lengths are non-dimensionalised by the half-gap 
$h$ between the plates, and velocities by 
$U_{wall}$. The Reynolds number is defined to be 
$Re \equiv U_{wall}\,h/\nu$. We will require one further dimensional quantity that appears in the friction coefficient – the mean horizontal shear at the walls, which we denote by 
$U^\prime _{wall}$. We will use non-dimensional variables throughout except when specified. We simulate the incompressible Navier–Stokes equations
$$\begin{gather} \frac{\partial \boldsymbol{u}}{\partial t} + \left(\boldsymbol{u} \boldsymbol{\cdot}\boldsymbol{\nabla}\right) \boldsymbol{u} =-\boldsymbol{\nabla} p + \frac{1}{Re}\,\nabla^2 \boldsymbol{u}, \end{gather}$$
$$\begin{gather}\boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u} = 0, \end{gather}$$
using the pseudo-spectral parallel code Channelflow (Gibson Reference Frishman and Grafke2012). Since the bands are found to be oriented obliquely with respect to the streamwise direction, we use a doubly periodic numerical domain that is tilted with respect to the streamwise direction of the flow, shown as the oblique rectangle in figure 1. This choice was introduced by Barkley & Tuckerman (Reference Barkley and Tuckerman2005) and has become common in studying turbulent bands (Shi, Avila & Hof Reference Shi, Avila and Hof2013; Lemoult et al. Reference Lemoult, Shi, Avila, Jalikop, Avila and Hof2016; Paranjape, Duguet & Hof Reference Paranjape, Duguet and Hof2020; Tuckerman, Chantry & Barkley Reference Tuckerman, Chantry and Barkley2020). The 
$x$ direction is chosen to be aligned with a typical turbulent band, and the 
$z$ coordinate to be orthogonal to the band. The relationship between streamwise–spanwise coordinates and tilted band-oriented coordinates is
$$\begin{gather} \boldsymbol{e}_{strm} = \cos\theta \, \boldsymbol{e}_x + \sin\theta \, \boldsymbol{e}_z, \end{gather}$$
$$\begin{gather}\boldsymbol{e}_{span} =-\sin\theta \, \boldsymbol{e}_x + \cos\theta \, \boldsymbol{e}_z. \end{gather}$$
The usual wall-normal coordinate is denoted by 
$y$, and the corresponding velocity by 
$v$. Thus the boundary conditions are 
$\boldsymbol {u}(y=\pm 1)=\pm \boldsymbol {e}_{strm}$ in 
$y$, and periodic in 
$x$ and 
$z$, together with a zero-flux constraint on the flow in the 
$x$ and 
$z$ directions. The field visualised in figure 1 comes from an additional simulation that we carried out in a domain of size (
$L_{strm},L_y,L_{span}) = (200,2,100)$ aligned with the streamwise–spanwise directions. Exploiting the periodic boundary conditions of the simulation, the visualisation shows four copies of the instantaneous field.
Figure 1. Spatial visualisation of our numerical domains at 
$Re=360$. Colours show the wall-normal velocity 
$v$ at the mid-plane 
$y=0$ (blue 
$-0.2$, white 0, red 0.2) in a domain of size 
$L_{strm}=400$, 
$L_{span}=200$. Red and blue boxes respectively show a minimal band unit and a long slender box.
The tilted box effectively reduces the dimensionality of the system by disallowing large-scale variation along the short 
$x$ direction. The flow in this direction is considered to be statistically homogeneous as it is dictated only by small turbulent scales. In a large non-tilted domain, bands with opposite orientations coexist (Prigent et al. Reference Prigent, Grégoire, Chaté and Dauchot2003; Duguet et al. Reference Duguet, Schlatter and Henningson2010; Klotz et al. Reference Klotz, Lemoult, Avila and Hof2022), but only one orientation is permitted in the tilted box.
We will use two types of numerical domains, with different lengths 
$L_z$. Both have fixed resolution 
$\Delta z = L_z/N_z = 0.08$, along with fixed 
$L_x=10$ (
$N_x=120$), 
$L_y=2$ 
$(N_y=33)$ and 
$\theta = 24^{\circ }$. These domains are shown in figure 1.
(i) Minimal band units, an example of which is shown as the red box in figure 1, accommodate a single band–gap pair and so are used to study a strictly periodic pattern of imposed wavelength
$\lambda =L_z$. ($L_z$ must typically be below $\simeq 65$ to contain a unique band.)(ii) Long slender boxes have a large
$L_z$ direction that can accommodate a large and variable number of gaps and bands in the system. The blue box in figure 1 is an example of such a domain with $L_z = 240$, but larger sizes ($L_z=400$ or $L_z=800$) will be used in our study.
3. Nucleation of laminar gaps and pattern emergence
We have carried out simulations in a long slender box of size 
$L_z=800$ for various 
$Re$, with the uniform turbulent state from a simulation at 
$Re=500$ as an initial condition, a protocol called a quench. Figure 2, an extension of figure 2 of Gomé et al. (Reference Gomé, Tuckerman and Barkley2023), displays the resulting spatio-temporal dynamics at six Reynolds numbers. Plotted is the 
$(z,t)$ dependence of the cross-flow energy 
$(v^2 + u_{span}^2)/2$ at 
$x=L_x/2$, 
$y=0$. The cross-flow energy is a useful diagnostic because it is zero for laminar flow and is therefore a proxy for turbulent kinetic energy (TKE). The choice 
$x=L_x/2$ is arbitrary since there is no large-scale variation of the flow field in the short 
$x$ direction of the simulation.
Figure 2. Spatio-temporal visualisation of pattern formation with 
$L_z=800$, for 
$Re$ values (a) 
$500$, (b) 
$460$, (c) 
$440$, (d) 
$420$, (e) 
$400$, and (f) 
$380$. Flow at 
$t=0$ is initiated from uniform turbulence at 
$Re=500$. Colour shows cross-flow energy 
$(v^2 + u_{span}^2)/2$ at 
$x=L_x/2$, 
$y=0$ (white 0, red 0.02). At high 
$Re$, weak local gaps appear sparsely. When 
$Re$ is decreased, spatio-temporally intermittent patterns of finite spatial extent emerge. These consist of turbulent cores (dark red) and quasi-laminar gaps (white). For still lower 
$Re$, quasi-laminar regions live longer, and patterns are more regular and steady.
Figure 2 encapsulates the main message of this section: the emergence of patterns out of uniform turbulence is a gradual process involving spatio-temporal intermittency of turbulent and quasi-laminar flow. At 
$Re=500$, barely discernible low-turbulent-energy regions appear randomly within the turbulent background. At 
$Re=460$, these regions are more pronounced and begin to constitute localised, short-lived quasi-laminar gaps within the turbulent flow. As 
$Re$ is decreased further, these gaps are more probable and last for longer times. Eventually, the gaps self-organise into persistent, albeit fluctuating, patterns. The remainder of the section will quantify the evolution of states seen in figure 2.
3.1. Statistics of laminar and turbulent zones
We consider the 
$(x, y)$-averaged cross-flow energy
\begin{equation} e(z,t) \equiv \frac{1}{L_x L_y} \int_{-1}^1 \int_0^{L_x} \frac{1}{2}\,( v^2 + u_{span}^{2})(x,y,z,t) \, \mathrm{d}\kern0.7pt x \,\mathrm{d}y \end{equation}
as a useful diagnostic of quasi-laminar and turbulent zones. The probability density functions (p.d.f.s) of 
$e(z,t)$ are shown in figure 3(a) for various values of 
$Re$. The right tails, corresponding to high-energy events, are broad and exponential for all 
$Re$. The left, low-energy portions of the p.d.f.s vary qualitatively with 
$Re$, unsurprisingly since these portions correspond to the weak turbulent events and hence include the gaps. For large 
$Re$, the p.d.f.s are maximal at 
$e \simeq 0.007$. As 
$Re$ is decreased, a low-energy peak emerges at 
$e\simeq 0.002$, corresponding to the emergence of long-lived quasi-laminar gaps seen in figure 2. The peak at 
$e\simeq 0.007$ flattens and gradually disappears. An interesting feature is that the distributions broaden with decreasing 
$Re$, with both low-energy and high-energy events becoming more likely. This reflects a spatial redistribution of energy that accompanies the formation of gaps, with turbulent bands extracting energy from quasi-laminar regions and consequently becoming more intense (see figure 6 of Gomé et al. (Reference Gomé, Tuckerman and Barkley2023)).
Figure 3. (a) P.d.f.s of local cross-flow energy 
$e(z,t)$ defined in (3.1). Maximum at 
$e\simeq 0.002$ appears for 
$Re\leqslant 420$. (b) Illustration of the thresholding 
$e(z,t) < e_{turb}$ of a turbulent–laminar field at 
$Re=440$, with turbulent regions 
$e(z,t) > e_{ turb}$ in white, and quasi-laminar regions in blue. Definitions of 
$L_{lam}$ and 
$L_{turb}$, the lengths of quasi-laminar and turbulent regions, are illustrated. (c) P.d.f.s of laminar gap widths 
$L_{lam}$ showing plateaux near 15 appearing for 
$Re\leqslant 440$. (d) P.d.f.s of widths of turbulent regions 
$L_{turb}$ showing local increase near 20 for 
$Re\leqslant 420$.
An intuitive way to define turbulent and quasi-laminar regions is by thresholding the values of 
$e(z,t)$. In the following, a region will be called quasi-laminar if 
$e(z, t) < e_{turb}$ and turbulent if 
$e(z, t) \geqslant e_{turb}$. As the p.d.f. of 
$e(z,t)$ evolves with 
$Re$, we define a 
$Re$-dependent threshold as a fraction of its average value 
$e_{turb} = 0.75 \,\bar {e}$. The thresholding is illustrated in figure 3(b), which is an enlargement of the flow at 
$Re=440$ that shows turbulent and quasi-laminar zones as white and blue areas, respectively. Thresholding within a fluctuating turbulent environment can conflate long-lived gaps with tiny, short-lived regions in which the energy fluctuates below the threshold 
$e_{turb}$. These are seen as the numerous small blue spots in figure 3(b) that differ from the wider and longer-lived gaps. This deficiency is addressed by examining the statistics of the spatial and temporal sizes of quasi-laminar gaps.
We present the length distributions of laminar (
$L_{lam}$) and turbulent (
$L_{turb}$) zones in figures 3(c) and 3(d) at various Reynolds numbers. These distributions have their maxima at very small lengths, reflecting the large number of small-scale, low-turbulent-energy regions that arise due to thresholding the fluctuating turbulent field. As 
$Re$ is decreased, the p.d.f. for 
$L_{lam}$ begins to develop a plateau at 
$L_{lam}\simeq 15$, corresponding to the scale of the gaps visible in figure 2. The right tails of the distribution are exponential and shift upwards with decreasing 
$Re$. The p.d.f. of 
$L_{turb}$ also varies with 
$Re$, but in a somewhat different way. As 
$Re$ decreases, the likelihood of a turbulent length in the range 
$15 \lesssim L_{turb} \lesssim 35$ increases above the exponential background, but at least over the range of 
$Re$ considered, a maximum does not develop.
The laminar length distributions show the emergence of structure at 
$Re$ higher than the turbulent length distributions. This is visible at 
$Re=440$, where the distribution of 
$L_{turb}$ is indistinguishable from those at higher 
$Re$, while the distribution of 
$L_{lam}$ is substantially altered. This is entirely consistent with the impression from the visualisation in figure 2(c) that quasi-laminar gaps emerge from a uniform turbulent background. Although the distributions of 
$L_{lam}$ and 
$L_{turb}$ behave differently, the length scales emerging as 
$Re$ decreases are within a factor of two. This aspect is not present in the pipe flow results of Avila & Hof (Reference Avila and Hof2013). (See Appendix A for a more detailed comparison.)
3.2. Gap lifetimes and transition to patterns
Temporal measurements of the gaps are depicted in figure 4. Figure 4(a) shows the procedure by which we define the temporal extents 
$t_{gap}$ of quasi-laminar gaps. For each gap, i.e. a connected zone in 
$(z,t)$ satisfying 
$e(z,t)< e_{turb}$, we locate its latest and earliest times and define 
$t_{gap}$ as the distance between them. Here again, we fix the threshold at 
$e_{turb}=0.75\,\bar {e}$. Figure 4(b) shows the temporal distribution of gaps, via the survival function of their lifetimes. In a similar vein to the spatial gap lengths, two characteristic behaviours are observed: for small times, many points are distributed near zero (as a result of frequent fluctuations near the threshold 
$e_{turb}$), while for large enough times, an exponential regime is seen:
\begin{equation} P(t_{gap} > t) \propto \exp({-t/\tau_{gap}(Re)})\quad \text{for } t>t_0, \end{equation}
where 
$t_0=500$ has been used for all 
$Re$, although the exponential range begins slightly earlier for larger values of 
$Re$.
Figure 4. (a) Same as figure 3(b), but illustrating the definition of 
$t_{gap}$, the lifetime of a quasi-laminar gap. (b) Survival functions of 
$t_{gap}$. After initial steep portions, slopes yield the characteristic times. (c) Evolution with 
$Re$ of characteristic time 
$\tau _{gap}$ and of ratio of large- to small-scale energy 
$e_{L/S}$ defined by (3.5). Both of these quantities present two exponential regimes, with the same slopes and a common crossover at 
$Re_{gu}$. The horizontal dashed line delimits the region 
$e_{L/S} > 1$, defining 
$Re_{pg}$ below which regular patterns dominate. We estimate 
$Re_{pg} \simeq 430$ and 
$Re_{gu}\simeq 470$ (to two significant figures). (d) Evolution of friction coefficient 
$C_f$ with 
$Re$, with the three regimes delimited by 
$Re_{pg}$ and 
$Re_{gu}$, as defined in (c).
The slope of the exponential tail is extracted at each 
$Re$, and the resulting characteristic time scale 
$\tau _{gap}$ is shown in figure 4(c). The evolution of 
$\tau _{gap}$ with 
$Re$ displays two regimes, each with nearly exponential dependence on 
$Re$, but with very different slopes on the semi-log plot. For 
$Re\geqslant 470$, the characteristic lifetimes are 
$\tau _{gap} = O(10^2)$ and vary weakly with 
$Re$. These short time scales correspond to the small white events visible in figure 2(a), and are associated with low-energy values on the left tails of the p.d.f.s for 
$e(z,t)$ in figure 3(a). Discounting these events, we refer to such states as uniform turbulence. For 
$Re<470$, 
$\tau _{gap}$ varies rapidly with 
$Re$, increasing by two orders of magnitude between 
$Re = 470$ and 
$Re=380$. The abrupt change in slope seen in figure 4(c), which we denote by 
$Re_{gu}$, marks the transition between gaps and uniform turbulence; we estimate 
$Re_{gu} = 470$ (to two significant figures). We stress that as far as we have been able to determine, there is no critical phenomenon associated with this change of behaviour. That is, the transition is smooth and lacks a true critical point. It is nevertheless evident that the dynamics of quasi-laminar gaps changes significantly in the region of 
$Re = 470$, therefore it is useful to define a reference Reynolds number marking this change in behaviour.
Note that typical lifetimes of laminar gaps must become infinite by the threshold 
$Re \simeq 325$ below which turbulence is no longer sustained (Lemoult et al. Reference Lemoult, Shi, Avila, Jalikop, Avila and Hof2016). (We believe gap lifetimes to be infinite even for 
$Re\lesssim 380$ when the permanent banded regime is attained, although this is not shown here.) For this reason, we have restricted our study of gap lifetimes to 
$Re\gtrsim 380$, and we have limited our maximal simulation time to 
${\sim }10^4$.
To quantify the distinction between localised gaps and patterns, we introduce a variable 
$e_{L/S}$ as follows. Using the Fourier transform in 
$z$,
\begin{equation} \hat{\boldsymbol{u}}(x,y,k_z,t) \equiv \frac{1}{L_z}\int_{0}^{L_z} \boldsymbol{u}(x,y,z,t)\,{\rm e}^{- {\rm i} k_z z } \,\text{d}z , \end{equation}we compute the averaged spectral energy
\begin{equation} \hat{E}(y,k_z)\equiv\tfrac{1}{2}\,\overline{\hat{\boldsymbol{u}}\boldsymbol{\cdot}\hat{\boldsymbol{u}}^\ast}, \quad \hat{E}(k_z)\equiv \langle\hat{E}(y,k_z)\rangle_y, \end{equation}
where the overbar designates an average in 
$x$ and 
$t$. This spectral energy is described in figure 3(a) of our companion paper (Gomé et al. Reference Gomé, Tuckerman and Barkley2023). We are interested in the ratio of 
$\hat {E}(k_z)$ at large scales (pattern scale) to small scales (roll–streak scale), as it evolves with 
$Re$. For this purpose, we define the ratio of large-scale to small-scale maximal energy:
\begin{equation} e_{L/S} = \frac{\underset{ k_z < 0.5 }{ \max}\hat{E} (k_z) }{\underset{ k_z \geqslant 0.5 }{ \max}\hat{E} (k_z)}. \end{equation}
The choice of wavenumber 
$k_z=0.5$ to delimit large and small scales comes from the change in sign of nonlinear transfers, as established in Gomé et al. (Reference Gomé, Tuckerman and Barkley2023). This quantity is shown as blue squares in figure 4(c), and is highly correlated to 
$\tau _{gap}$. This correlation is in itself a surprising observation for which we have no explanation.
For 
$Re \gtrsim 430$, we have 
$e_{L/S} < 1$, signalling that the dominant peak in the energy spectrum is at the roll–streak scale, while for 
$Re\lesssim 430$, the large-scale pattern begins to dominate the streaks and rolls, as indicated by 
$e_{L/S} > 1$ (dashed blue line in figure 4c). Note that 
$Re=430$ is also the demarcation between unimodal and bimodal p.d.f.s of 
$e(z,t)$ in figure 3(a). The transition from gaps to patterns is smooth. In fact, we do not even observe a qualitative feature sharply distinguishing gaps and patterns. Nevertheless, we find it useful to define a reference Reynolds number associated with patterns starting to dominate the energy spectrum. This choice has the advantage of yielding a quantitative criterion, which we estimate as 
$Re_{pg} \simeq 430$ (to two significant figures). We find a similar estimation of the value of 
$Re$ below which patterns start to dominate via a wavelet-based measurement; see Appendix B.
In addition to the previous quantitative measures, we also extract the friction coefficient. This is defined as the ratio of the mean wall shear stress 
$\mu U^\prime _{wall}$ to the dynamic pressure 
$\rho U_{wall}^2/2$, which we write in physical units and then in non-dimensional variables as
\begin{equation} C_f \equiv \frac{\mu U^\prime_{wall}}{\frac{1}{2}\rho U_{wall}^2} = \frac{2\nu}{h U_{wall}}\,\frac{U^\prime_{wall}}{U_{wall}/h} = \frac{2}{Re} \left.\frac{\partial \langle u_{strm} \rangle_{x,z,t}}{ \partial y} \right\vert_{wall}. \end{equation}
In (3.6), the dimensional quantities 
$h$, 
$\rho$, 
$\mu$ and 
$\nu$ are the half-height, density, and dynamic and kinematic viscosities, and 
$U_{wall}$ and 
$U^\prime _{wall}$ are the velocity and mean velocity gradient at the wall. We note that the behaviour of 
$C_f$ in the transitional region has been investigated in plane channel flow by Shimizu & Manneville (Reference Shimizu and Manneville2019) and Kashyap, Duguet & Dauchot (Reference Kashyap, Duguet and Dauchot2020). Our measurements of 
$C_f$ are shown in figure 4(d). We distinguish different trends within each of the three regimes defined earlier in figure 4(c). In the uniform regime 
$Re>Re_{gu}=470$, 
$C_f$ increases with decreasing 
$Re$. In the patterned regime 
$Re < Re_{pg}=430$, 
$C_f$ decreases with decreasing 
$Re$. The localised-gap regime 
$Re_{pg} < Re < Re_{gu}$ connects these two tendencies, with 
$C_f$ reaching a maximum at 
$Re=450$.
The presence of a region of maximal 
$C_f$ (or equivalently maximal total dissipation) echoes the results on the energy balance presented in Gomé et al. (Reference Gomé, Tuckerman and Barkley2023): the uniform regime dissipates more energy as 
$Re$ decreases, up to a point where this is mitigated by the many laminar gaps nucleated. This is presumably due to the mean flow in the turbulent region requiring energy influx from gaps to compensate for its increasing dissipation.
3.3. Laminar–turbulent correlation function
The changes in regimes and the distinction between local gaps and patterns can be studied further by measuring the spatial correlation between quasi-laminar regions within the flow. We define
\begin{equation} \varTheta(z,t) = \begin{cases} 1 & \text{if } e(z,t) < e_{turb} \text{ (laminar)},\\ 0 & \text{otherwise (turbulent)} \end{cases} \end{equation}(this is the quantity shown in blue and white in figures 3b and 4a). We then compute its spatial correlation function
\begin{equation} C (\delta z) = \frac{\left\langle {\varTheta}(z)\,{\varTheta}(z+\delta z) \right\rangle_{z, t} - \left\langle {\varTheta}(z)\right\rangle^2_{z, t}}{\left\langle \varTheta(z)^2\right\rangle_{z, t} - \left\langle\varTheta(z)\right\rangle^2_{z, t}}. \end{equation}
Along with 
$(z,t)$ averaging, 
$C$ is also averaged over multiple realisations of quench experiments. As 
$\varTheta$ is a Heaviside function, 
$C$ can be understood as the probability of finding a gap at a distance 
$\delta z$ from a gap at position 
$z$. The results are presented in figure 5(a). The comparative behaviour of 
$C$ at near-zero values is enhanced by plotting 
$\tanh (10C)$ in figure 5(b). At long range, 
$C$ approaches zero with small fluctuations at 
$Re=480$, a noisy periodicity at 
$Re=460$, and a nearly periodic behaviour for 
$Re\leq 420$.
Figure 5. (a) Gap-to-gap correlation function 
$C (\delta z)$ defined by (3.8) for various values of 
$Re$. (b) For 
$Re\gtrsim 440$, the weak variation and short-ranged maxima are enhanced by plotting 
$\tanh (10 C (\delta z))$. The dots correspond to the first local maximum, indicating the selection of a finite distance between two local gaps, including at the highest 
$Re$. Large-scale modulations smoothly leave room to weak short-range interaction as 
$Re$ increases and the flow visits patterned, local-gap and uniform regimes.
In all cases, 
$C$ initially decreases from 1 and reaches a first minimum at 
$\delta z \simeq 20$, due to the minimal possible size of a turbulent zone that suppresses the creation of neighbouring laminar gaps. Also, 
$C$ has a prominent local maximum 
$\delta z_{max}$ right after its initial decrease, at 
$\delta z_{max} \simeq 32$ at 
$Re=480$, which increases to 
$\delta z_{max} \simeq 41$ at 
$Re=420$. These maxima, shown as coloured circles in figure 5(b), indicate that gap nucleation is preferred at distance 
$\delta z_{max}$ from an existing gap. The increase in 
$\delta z_{max}$ and in the subsequent extrema as 
$Re$ is lowered agrees with the trend of increasing wavelength of turbulent bands as 
$Re$ is decreased in the fully banded regime at lower 
$Re$ (Prigent et al. Reference Prigent, Grégoire, Chaté and Dauchot2003; Barkley & Tuckerman Reference Barkley and Tuckerman2005).
The smooth transition from patterns to uniform flow is confirmed in the behaviour of the correlation function. Large-scale modulations characteristic of the patterned regime gradually disappear with increasing 
$Re$, as gaps become more and more isolated. Only a weak, finite-length interaction subsists in the local-gap and uniform regimes, and will disappear with increasing 
$Re$. This is the selection of this finite gap spacing that we will investigate in §§ 4 and 5.
4. Wavelength selection for turbulent–laminar patterns
In this section, we investigate the existence of a preferred pattern wavelength by using as a control parameter the length 
$L_z$ of the minimal band unit. In a minimal band unit, the system is constrained and the distinction between local gaps and patterns is lost; see § 3 of our companion paper Gomé et al. (Reference Gomé, Tuckerman and Barkley2023). Here, 
$L_z$ is chosen such as to accommodate at most a single turbulent zone and a single quasi-laminar zone, which due to imposed periodicity can be viewed as one period of a perfectly periodic pattern. By varying 
$L_z$, we can verify whether a regular pattern of given wavelength 
$L_z$ can emerge from uniform turbulence, disregarding the effect of scales larger than 
$L_z$ or of competition with wavelengths close to 
$L_z$. We refer to these simulations in minimal band units as existence experiments. Indeed, one of the main advantages of the minimal band unit is the ability to create patterns of a given angle and wavelength that may not be stable in a larger domain.
In contrast, in a long slender box, 
$L_z$ is large enough to accommodate multiple bands and possibly even patterns of different wavelengths. An initial condition consisting of a regular pattern of wavelength 
$\lambda$ can be constructed by concatenating bands produced from a minimal band unit of size 
$\lambda$. The stability of such a pattern is studied by allowing this initial state to evolve via the nonlinear Navier–Stokes equations. Both existence and stability studies can be understood in the framework of the Eckhaus instability (Kramer & Zimmermann Reference Kramer and Zimmermann1985; Ahlers et al. Reference Ahlers, Cannell, Dominguez-Lerma and Heinrichs1986; Tuckerman & Barkley Reference Tuckerman and Barkley1990; Cross & Greenside Reference Cross and Greenside2009).
In previous work on transitional regimes, Barkley & Tuckerman (Reference Barkley and Tuckerman2005) studied the evolution of patterns as 
$L_z$ was increased. In § 4.1, we extend this approach to multiple sizes of the minimal band unit by comparing lifetimes of patterns that arise naturally in this constrained geometry. The stability of regular patterns of various wavelengths will be studied in long slender box domains (
$L_z=400$) in § 4.2.
4.1. Temporal intermittency of regular patterns in a short-$L_z$ box
Figure 6(a) shows the formation of a typical pattern in a minimal band unit of size 
$L_z=40$ at 
$Re=440$. While the system cannot exhibit the spatial intermittency seen in figure 2(c), temporal intermittency is possible and is seen as alternation between uniform turbulence and a pattern. We plot the spanwise velocity at 
$y=0$ and 
$x=L_x/2$. This is a particularly useful measure of the large-scale flow associated with patterns, seen as red and blue zones surrounding a white quasi-laminar region, i.e. a gap. The patterned state emerges spontaneously from uniform turbulence and remains from 
$t \simeq 1500$ to 
$t \simeq 3400$. At 
$t\simeq 500$, a short-lived gap appears at 
$z=10$, which can be seen as an attempt to form a pattern.
Figure 6. Pattern lifetimes. (a) Space-time visualisation of a metastable pattern in a minimal band unit with 
$L_z=40$ at 
$Re=440$. Colours show spanwise velocity (blue 
$-0.1$, white 0, red 0.1). (b) Values of the dominant wavelength 
$\hat {\lambda }_{max}$ (light blue curve) and of its short-time average 
$\langle \hat {\lambda }_{max}\rangle _{t_a}$ (dark blue curve); see (4.1). A state is defined to be patterned if 
$\hat {\lambda }_{max} = L_z$. (c) Survival function of lifetimes of turbulent–laminar patterns in a minimal band unit with 
$L_z=40$ for various 
$Re$. The pattern lifetimes 
$t_{pat}$ are the lengths of the time intervals during which 
$\hat {\lambda }_{max} = L_z$. (d) The top plot shows characteristic times 
$\tau _{pat}$ extracted from survival functions as a function of 
$L_z$ and 
$Re$. The bottom plot shows the intermittency factor for the patterned state 
$\gamma _{pat}$, which is the fraction of time spent in the patterned state. The Re values are the same as in (c) with the same colour coding.
We characterise the pattern quantitatively as follows. For each time 
$t$, we compute 
$|\langle \hat {\boldsymbol {u}}(y=0,k_z,t)\rangle _x|^2$, which is (twice) the instantaneous energy contained in wavenumber 
$k_z$ at the mid-plane. We then determine the wavenumber that maximises this energy, and compute the corresponding wavelength. That is, we define
\begin{equation} \hat{\lambda}_{max}(t) \equiv \frac{2{\rm \pi}}{\underset{k_z>0}{\text{argmax}}\, |{\langle\hat{\boldsymbol{u}}(y=0,k_z,t)}\rangle_x|^2}. \end{equation}
The possible values of 
$\hat {\lambda }_{max}$ are integer divisors of 
$L_z$, here 40, 20, 10, etc. Figure 6(b) presents 
$\hat {\lambda }_{max}$ and its short-time average 
$\langle \hat {\lambda }_{max}\rangle _{t_a}$ with 
$t_a=30$ as light and dark blue curves, respectively. When turbulence is uniform, 
$\hat {\lambda }_{max}$ varies rapidly between its discrete allowed values, while 
$\langle \hat {\lambda }_{max}\rangle _{t_a}$ fluctuates more gently around 10. The flow state is deemed to be patterned when its dominant mode is 
$\langle \hat {\lambda }_{max}\rangle _{t_a} =L_z$. The long-lived pattern occurring for 
$1500 \leqslant t \leqslant 3400$ in figure 6(a) is seen as a plateau of 
$\langle \hat {\lambda }_{max}\rangle _{t_a}$ in figure 6(b). There are other shorter-lived plateaus, notably at for 
$500\leqslant t \leqslant 750$. A similar analysis was carried out by Barkley & Tuckerman (Reference Barkley and Tuckerman2005) and Tuckerman & Barkley (Reference Tuckerman and Barkley2011) using the Fourier component corresponding to wavelength 
$L_z$ of the spanwise mid-gap velocity.
Figure 6(c) shows the survival function 
$t_{pat}$ of the pattern lifetimes obtained from 
$\langle \hat {\lambda }_{max}\rangle _{t_a}$ over long simulation times for various 
$Re$. This measurement differs from figure 4(b), which showed lifetimes of gaps in a long slender box and not regular patterns obtained in a minimal band unit. The results are, however, qualitatively similar, with two characteristic zones in the distribution: at short times, many patterns appear due to fluctuations, while after 
$t\simeq 200$, the survival functions enter an approximately exponential regime, from which we extract the characteristic times 
$\tau _{pat}$ by taking the inverse of the slope.
We then vary 
$L_z$, staying within the minimal box regime 
$L_z \lesssim 65$ in which only one band can fit. The top plot of figure 6(d) shows that 
$\tau _{pat}$ presents a broad maximum in 
$L_z$ whose strength and position depend on 
$Re$: 
$L_z\simeq 42$ at 
$Re=440$, and 
$L_z\simeq 44$ at 
$Re=400$. This wavelength corresponds approximately to the natural spacing observed in a large slender box (figure 2). The bottom plot of figure 6(d) presents the fraction of time that is spent in a patterned state, denoted 
$\gamma _{pat}$ to reflect that this should be thought of as the intermittency factor for the patterned state. The dependence of 
$\gamma _{pat}$ on 
$L_z$ follows the same trend as 
$\tau _{pat}$, but less strongly (the scale of 
$\gamma_{pat}$ is linear, while that for 
$\tau _{pat}$ is logarithmic).
The results shown in figure 6(d) complement the Ginzburg–Landau description proposed by Prigent et al. (Reference Prigent, Grégoire, Chaté and Dauchot2003) and Rolland & Manneville (Reference Rolland and Manneville2011). To quantify the bifurcation from featureless to pattern turbulence, these authors defined an order parameter and showed that it has a quadratic maximum at an optimal wavenumber. This is consistent with the approximate quadratic maxima that we observe in 
$\tau _{pat}$ and 
$\gamma _{pat}$ with regard to 
$L_z$. Note that the scale of the pattern can be set roughly from the force balance in the laminar flow regions (Barkley & Tuckerman Reference Barkley and Tuckerman2007), 
$\lambda \simeq Re\sin \theta /{\rm \pi}$, which yields wavelength 52 at 
$Re=400$ (not too distant from the value 44 found in figure 6d).
4.2. Pattern stability in a large domain
To study the stability of a pattern of wavelength 
$\lambda$, we prepare an initial condition for a long slender box by concatenating repetitions of a single band produced in a minimal band unit. We add small-amplitude noise to this initial pattern so that the repeated bands do not all evolve identically. Figures 7(a) and 7(b) show two examples of such simulations. Depending on the values of 
$Re$ and the initial wavelength 
$\lambda$, the pattern destabilises to either another periodic pattern (figure 7(a) for 
$Re=400$) or to localised patterns surrounded by patches of featureless turbulence (figure 7(b) for 
$Re=430$).
Figure 7. Simulation in a long slender box from a noise-perturbed periodic pattern with (a) initial 
$\lambda = 57$ at 
$Re=400$ and (b) initial 
$\lambda = 40$ at 
$Re=430$. Colours show spanwise velocity (red 0.1, white 0, blue 
$-0.1$). (c,d) Local dominant wavelength 
$\tilde {\lambda }_{max}(z,t)$ determined by wavelet analysis (see Appendix B) corresponding to the simulations shown in (a,b). Colour at 
$t=0$ shows the wavelength 
$\lambda$ of the initial condition. (e) Wavelet-defined 
$H_\lambda (t)$ given by (4.2), which quantifies the proportion of the domain that retains initial wavelength 
$\lambda$ as a function of time for cases in (a,b). Circles indicate the times for (a,b) after which 
$H_\lambda$ is below the threshold value 
$H_{stab}$ for a sufficiently long time. (f) Ensemble-averaged 
$\bar {t}_{stab}$ of the decay time of an imposed pattern of wavelength 
$\lambda$ for various values of 
$Re$. The relative stability of a wavelength, whether localised or not, is measured by 
$\bar {t}_{stab}$ via the wavelet analysis.
It can be seen that patterns often occupy only part of the domain. For this reason, we turn to the wavelet decomposition (Meneveau Reference Meneveau1991; Farge Reference Farge1992) to quantify patterns locally. In contrast to a Fourier decomposition, the wavelet decomposition quantifies the signal as a function of space and scale. From this, we are able to define a local dominant wavelength, 
$\tilde {\lambda }_{max}(z,t)$, similar in spirit to 
$\hat {\lambda }_{max}(t)$ in (4.1), but now at each space–time point (see Appendix B for details). Figures 7(c) and 7(d) show 
$\tilde {\lambda }_{max}(z,t)$ obtained from wavelet analysis of the simulations visualised in figures 7(a) and 7(b).
We now use the local wavelength 
$\tilde {\lambda }_{max}(z,t)$ to quantify the stability of an initial wavelength. We use a domain of length 
$L_z=400$, and we concatenate 
$n=7\unicode{x2013}13$ repetitions of a single band to produce a pattern with initial wavelength 
$\lambda (n)\equiv 400/n\simeq 57,50,44,\ldots,31$. (We have rounded 
$\lambda$ to the nearest integer value here and in what follows.) After adding low-amplitude noise, we run a simulation lasting 5000 time units, compute the wavelet transform, and calculate from it the local wavelengths 
$\tilde {\lambda }_{max}(z,t)$. We define 
$\epsilon _\lambda \equiv \min ((\lambda (n+1)-\lambda (n))/2, (\lambda (n)-\lambda (n-1))/2)$ such that 
$|\lambda - \tilde {\lambda }_{max}(z, t)|< \epsilon _\lambda$ if 
$\tilde {\lambda }_{max}$ is closer to 
$\lambda (n)$ than to its two neighbouring values. Finally, in order to measure the proportion of 
$L_z$ in the dominant mode 
$\tilde {\lambda }_{max}$ is 
$\lambda$, we compute
\begin{equation} H_\lambda(t) = \left\langle\frac{1}{L_z}\int_0^{L_z} \varTheta\left(\epsilon_\lambda - |\lambda - \tilde{\lambda}_{max}(z, t)|\right)\text{d}z \right\rangle_{t_a}, \end{equation}
where 
$\varTheta$ is the Heaviside function, and the short-time average 
$\langle {\cdot }\rangle _{t_a}$ is taken over time 
$t_a=30$ as before. In practice, because patterns in a long slender box still fluctuate in width, a steady pattern may have 
$H_\lambda$ somewhat less than 1. If 
$H_\lambda \ll 1$, then a pattern of wavelength 
$\lambda$ is present in at most a very small part of the flow.
Figure 7(e) shows how wavelet analysis via the Heaviside-like function 
$H_\lambda (t)$ quantifies the relative stability of the pattern in the examples shown in figures 7(a) and 7(b). The flow in figure 7(a) at 
$Re=400$ begins with 
$\lambda =57$, i.e. 7 bands. Figure 7(c) retains the red colour corresponding to 
$\lambda =57$ over the entire domain for 
$t\lesssim 1200$, and over most of it until 
$t\lesssim 2300$. The red curve in figure 7(e) shows 
$H_\lambda$ decaying quickly and roughly monotonically. One additional gap appears at approximately 
$t=2300$, and starting from then, 
$H_\lambda$ remains low. This corresponds to the initial wavelength 
$\lambda =57$ losing its dominance to 
$\lambda =40$, 44 and 50 in the visualisation of 
$\tilde {\lambda }_{max}(z,t)$ in figure 7(c). By 
$t=5000$, the flow shows 9 bands with local wavenumber 
$\lambda$ between 40 and 50.
The flow in figure 7(b) at 
$Re=430$ begins with 
$\lambda =40$, i.e. 10 bands. Figure 7(d) shows that the initial light green colour corresponding to 40 is retained until 
$t\lesssim 800$. The blue curve in figure 7(e) representing 
$H_\lambda$ initially decreases and drops precipitously at 
$t\simeq 1000$ as several gaps disappear in figure 7(b). Then 
$H_\lambda$ fluctuates around a finite value, which is correlated to the presence of gaps whose local wavelength is the same as the initial 
$\lambda$, visible as zones where 
$\tilde {\lambda }_{max} = 40$ in figure 7(d). The rest of the flow can be seen mostly as locally featureless turbulence, where the dominant wavelength is small (
$\tilde {\lambda }_{max} \leqslant 10$). The local patterns fluctuate in width and strength, and 
$H_\lambda$ evolves correspondingly after 
$t=1000$. The final state reached in figure 7(a) at 
$Re=430$ is characterised by the presence of intermittent local gaps.
The lifetime of an initially imposed pattern wavelength 
$\lambda$ is denoted 
$t_{stab}$ and is defined as follows. We first define a threshold 
$H_{stab} \equiv 0.2$ (marked by a horizontal dashed line in figure 7e). If 
$H_\lambda (t)$ is statistically below 
$H_{stab}$, then the imposed pattern will be considered as unstable. Following this principle, 
$t_{stab}$ is defined as the first time when 
$H_\lambda$ is below 
$H_{stab}$, with a further condition to dampen the effect of short-term fluctuations: 
$\left \langle H_\lambda (t)\right \rangle _{t\in [t_{stab}, t_{stab} + 2000]}$ must remain below 
$H_{stab}$. The corresponding times in the cases of figures 7(a) and 7(b) are marked respectively by red and blue circles in figure 7(e).
Repeating this experiment over multiple realisations of the initial pattern (i.e. different noise realisations) yields an ensemble-averaged 
$\bar {t}_{stab}$. This procedure estimates the time for an initially regular and dominant wavelength to disappear from the flow domain, regardless of the way in which it does so and of the final state approached. Figure 7(f) presents the dependence of 
$\bar {t}_{stab}$ on 
$\lambda$ for different values of 
$Re$. Although our procedure relies on highly fluctuating signals (like those presented in figure 7e) and on a number of arbitrary choices (
$H_{stab}$, 
$\epsilon _{\lambda }$, etc.) that alter the exact values of stability times, we find that the trends visualised in figure 7(f) are robust. (The sensitivity of 
$\bar {t}_{stab}$ with 
$H_{stab}$ is shown in figure 13(b) of Appendix B.)
A most-stable wavelength ranging between 40 and 44 dominates the stability times for all the values of 
$Re$ under study. This is similar to the results from the existence study in figure 6(d), which showed a preferred wavelength emerging from the uniform state at approximately 
$42$ at 
$Re=440$. Consistently with what was observed in minimal band units of different sizes, the most stable wavelength grows with decreasing 
$Re$.
4.3. Discussion
Our study of the existence and stability of large-scale modulations of the turbulent flow is summarised in figure 8. This figure resembles the existence and stability diagrams presented for usual (non-turbulent) hydrodynamic instabilities such as Rayleigh–Bénard convection and Taylor–Couette flow (Busse Reference Busse1981; Ahlers et al. Reference Ahlers, Cannell, Dominguez-Lerma and Heinrichs1986; Cross & Greenside Reference Cross and Greenside2009). In classic systems, instabilities appear with increasing control parameter, while here gaps and bands emerge from uniform turbulent flow as 
$Re$ is lowered. Therefore, we plot the vertical axis in figure 8 with decreasing Reynolds number upwards.
Figure 8. Visualisation of the pattern selection in the phase space 
$(\lambda, Re)$. Colours show the stability times 
$\bar {t}_{stab}$, while open circles are points 
$\gamma _{ pat}(\lambda, Re)=0.45$. The dashed line is an illustrative fit of these points. 
$Re_{pg}$ and 
$Re_{gu}$ delimit the patterned, local gap and uniform regimes defined in § 3 (see figure 4).
We recall that the existence study of § 4.1 culminated in the measurement of 
$\gamma _{pat}(\lambda, Re)$, the fraction of simulation time that is spent in a patterned state, plotted in figure 6(d). The parameter values at which 
$\gamma _{pat}(\lambda, Re) = 0.45$ (an arbitrary threshold that covers most of our data range) are shown as black circles in figure 8. The dashed curve is an interpolation of the iso-
$\gamma _{pat}$ points and separates two regions, with patterns more likely to exist above the curve than below. The minimum of this curve is estimated to be 
$\lambda \simeq 42$. This is a preferred wavelength at which patterns first emerge statistically as 
$Re$ is decreased from large values.
The final result of the stability study in § 4.2, shown in figure 7(f), was 
$\bar {t}_{stab} (Re,\lambda )$, a typical duration over which a pattern initialised with wavelength 
$\lambda$ would persist. The colours in figure 8 show 
$\bar {t}_{stab}$. The peak in 
$\bar {t}_{stab}$ is first discernible at 
$Re \simeq 440$ and occurs at 
$\lambda \simeq 40$. The pattern existence and stability zones are similar in shape and in their lack of symmetry with respect to line 
$\lambda =42$. The transition seen in figures 7(a) and 7(c) from 
$\lambda =57$ to 
$\lambda =44$ at 
$Re=400$ corresponds to motion from a light blue to a dark blue area in the top row of figure 8. This change in pattern wavelength resembles the Eckhaus instability, which in classic hydrodynamics leads to transitions from unstable wavelengths outside a stability band to stable wavelengths inside.
The presence of a most-probable wavelength confirms the initial results of Prigent et al. (Reference Prigent, Grégoire, Chaté and Dauchot2003) and those of Rolland & Manneville (Reference Rolland and Manneville2011). This is also consistent with the instability study of Kashyap et al. (Reference Kashyap, Duguet and Dauchot2022) in plane Poiseuille flow. However, contrary to classic pattern-forming instabilities, the turbulent–laminar pattern emerges not from an exactly uniform state, but instead from a state in which local gaps are intermittent, as established in § 3. In § 5, we will emphasise the importance of the mean flow in the wavelength selection that we have described.
5. Optimisation of the large-scale flow
This section is devoted to the dependence of various energetic features of the patterned flow on the domain length 
$L_z$ of a minimal band unit. We fix the Reynolds number at 
$Re=400$. In the existence study of § 4, the wavelength 
$\lambda \simeq 44$ was found to be selected. (Recall the uppermost curves corresponding to 
$Re=400$ in figure 6d.) We will show that this wavelength also extremises quantities in the energy balances of the flow.
5.1. Average energies in the patterned state
We first decompose the flow into a mean and fluctuations, 
$\boldsymbol {u} = \bar {\boldsymbol {u}} + \boldsymbol {u}^\prime$, where the mean (overbar) is taken over the statistically homogeneous directions 
$x$ and 
$t$. We compute energies of the total flow 
$\langle E\rangle \equiv \langle \boldsymbol {u} \boldsymbol {\cdot } \boldsymbol {u}\rangle /2$ and the fluctuations (TKE) 
$\langle K\rangle \equiv \langle \boldsymbol {u}^\prime \boldsymbol {\cdot } \boldsymbol {u}^\prime \rangle /2$, where 
$\langle {\cdot }\rangle$ is the 
$(x,y,z,t)$ average. Figure 9(a) shows these quantities as functions of 
$L_z$ for the patterned state at 
$Re=400$. At 
$L_z=44$, 
$\langle E\rangle$ is maximal and 
$\langle K\rangle$ is minimal. As a consequence, the mean-flow energy 
$\frac {1}{2} \langle \bar {\boldsymbol {u}} \boldsymbol {\cdot } \bar {\boldsymbol {u}} \rangle = \langle E\rangle - \langle K\rangle$ is also maximal at 
$L_z=44$. Figure 9(a) additionally shows average dissipation of the total flow 
$\langle D\rangle \equiv \left \langle |\boldsymbol {\nabla } \times \boldsymbol {u} |^2\right \rangle /Re$ and average dissipation of TKE 
$\langle \epsilon \rangle \equiv \left \langle |\boldsymbol {\nabla } \times \boldsymbol {u}^\prime |^2\right \rangle /Re$, both of which are minimal at 
$L_z=44$. Note that these total energy and dissipation terms change very weakly with 
$L_z$, with a variation of less than 
$6\,\%$.
Figure 9. Energy analysis for the patterned state at 
$Re=400$ as a function of the size 
$L_z$ of a minimal band unit (
$L_z$ is the wavelength of the imposed pattern). (a) Spatially-averaged total energy 
$\langle E\rangle$, TKE 
$\langle K\rangle$ (
$\times 5$), total dissipation 
$\langle D\rangle$, turbulent dissipation 
$\langle \epsilon \rangle$ (
$\times 3$), for the patterned state at 
$Re=400$ as a function of 
$L_z$. (b) Energy in each of the 
$z$ Fourier components of the mean flow (5.1) and (5.2a–c).
The mean flow is further analysed by computing the energy of each spectral component of the mean flow. For this, the 
$(x,t)$ averaged flow 
$\bar {\boldsymbol {u}}$ is decomposed into Fourier modes in 
$z$:
\begin{equation} \bar{\boldsymbol{u}}(y, z) = \bar{\boldsymbol{u}}_0(y) + 2\,\mathcal{R} \left(\bar{\boldsymbol{u}}_1(y) \exp({{\rm i} 2{\rm \pi} z/L_z }) \right) + \bar{\boldsymbol{u}}_{{>}1} (y,z), \end{equation}
where 
$\bar {\boldsymbol {u}}_0$ is the uniform component of the mean flow, 
$\bar {\boldsymbol {u}}_1$ is the trigonometric Fourier coefficient corresponding to 
$k_z= 2{\rm \pi} /L_z$, and 
$\bar {\boldsymbol {u}}_{>1}$ is the remainder of the decomposition, for 
$k_z> 2{\rm \pi} /L_z$. (We have omitted the hats on the 
$z$ Fourier components of 
$\bar {\boldsymbol {u}}$.) The energies of the spectral components relative to the total mean energy are
\begin{equation} e_0 = \frac{\left\langle\bar{\boldsymbol{u}}_0\boldsymbol{\cdot}\bar{\boldsymbol{u}}_0\right\rangle}{\left\langle\bar{\boldsymbol{u}}\boldsymbol{\cdot}\bar{\boldsymbol{u}}\right\rangle},\quad e_1 =\frac{\left\langle\bar{\boldsymbol{u}}_1\boldsymbol{\cdot}\bar{\boldsymbol{u}}_1\right\rangle}{\left\langle\bar{\boldsymbol{u}}\boldsymbol{\cdot}\bar{\boldsymbol{u}}\right\rangle},\quad e_{{>}1} =\frac{\left\langle\bar{\boldsymbol{u}}_{{>}1}\boldsymbol{\cdot}\bar{\boldsymbol{u}}_{{>}1}\right\rangle}{\left\langle\bar{\boldsymbol{u}}\boldsymbol{\cdot}\bar{\boldsymbol{u}}\right\rangle}. \end{equation}
These are presented in figure 9(b). It can be seen that 
$e_0 \gg e_1 > e_{>1}$, and also that all have an extremum at 
$L_z=44$. In particular, 
$L_z=44$ minimises 
$e_0$ (
$e_0=0.95$) while maximising the trigonometric component (
$e_1=0.025$) along with the remaining components (
$e_{>1} \simeq 0.011$). Note that for a banded state at 
$Re=350$, 
$L_z=40$, Barkley & Tuckerman (Reference Barkley and Tuckerman2007) found that 
$e_0\simeq 0.70$, 
$e_1\simeq 0.30$ and 
$e_{>1}\simeq 0.004$, consistent with a strengthening of the bands as 
$Re$ is decreased.
5.2. Mean-flow spectral balance
We now investigate the spectral contributions to the budget of the mean flow 
$\bar {\boldsymbol {u}}$, dominated by the mean flow's two main spectral components 
$\bar {\boldsymbol {u}}_0$ and 
$\bar {\boldsymbol {u}}_1$. The balances can be expressed as in Gomé et al. (Reference Gomé, Tuckerman and Barkley2023):
\begin{equation} \hat{\bar{A}}_0 -\hat{\bar{\varPi}}_0 - \hat{\bar{D}}_0 + I = 0\ \text{for } \bar{\boldsymbol{u}}_0 \quad \text{and} \quad \hat{\bar{A}}_1 - \hat{\bar{\varPi}}_1 - \hat{\bar{D}}_1 = 0\ \text{for } \bar{\boldsymbol{u}}_1, \end{equation}
where 
$I$ is the rate of energy injection by the viscous shear, and 
$\hat {\bar {\varPi }}_0$, 
$\hat {\bar {D}}_0$ and 
$\hat {\bar {A}}_0$ stand for, respectively, production, dissipation and advection (i.e. nonlinear interaction) contributions to the energy balance of mode 
$\bar {\boldsymbol {u}}_0$, and similarly for 
$\bar {\boldsymbol {u}}_1$. These are defined by
$$\begin{gather} I = \frac{2}{Re}\,\mathcal{R}\left\{ \int_{-1}^1 \frac{\partial}{\partial y} \left(\hat{\bar{u}}_j^* (k_z=0)\,\hat{\overline{s}}_{yj} (k_z=0) \right) \text{d} y \right\} = \frac{1}{Re}\left( \left.\frac{\partial \bar{u}_{strm}}{\partial y}\right\vert_{1} + \left.\frac{\partial \bar{u}_{strm}}{\partial y}\right\vert_{-1} \right), \end{gather}$$
$$\begin{gather}\hat{\bar{\varPi}}_0 = \mathcal{R}\left\{\int_{-1}^{1} \frac{\partial \hat{\bar{u}}_j^*}{\partial x_i} (k_z=0)\,\widehat{\overline{u_i^\prime u_j^\prime}} (k_z=0) \,\text{d} y \right\}, \end{gather}$$
$$\begin{gather}\hat{\bar{D}}_0 = \frac{2}{Re}\,\mathcal{R}\left\{\int_{-1}^{1} \hat{\bar{s}}_{ij} (k_z=0)\,\hat{\bar{s}}_{ij}^* (k_z=0) \,\text{d} y \right\}, \end{gather}$$
$$\begin{gather}\hat{\bar{A}}_0 =- \mathcal{R}\left\{\int_{-1}^{1}\hat{\bar{u}}_j^* (k_z=0)\, \widehat{\bar{u}_i\,\frac{\partial \bar{u}_j}{\partial x_i}}(k_z=0) \,\text{d} y \right\}, \end{gather}$$
where 
$\mathcal {R}$ denotes the real part. We define 
$\hat {\bar {\varPi }}_1$, 
$\hat {\bar {D}}_1$ and 
$\hat {\bar {A}}_1$ similarly by replacing 
$k_z=0$ by 
$k_z=2{\rm \pi} /L_z$ in (5.4a)–(5.4d).
We recall two main results from Gomé et al. (Reference Gomé, Tuckerman and Barkley2023). First, 
$\hat {\bar {A}}_1 \simeq - \hat {\bar {A}}_0$. This term represents the energetic transfer between modes 
$\bar {\boldsymbol {u}}_0$ and 
$\bar {\boldsymbol {u}}_1$ via the self-advection of the mean flow (the energetic spectral transfer due to 
$(\bar {\boldsymbol {u}} \boldsymbol {\cdot }\boldsymbol {\nabla }) \bar {\boldsymbol {u}}$). Second, 
$\hat {\bar {\varPi }}_1 <0$, and this term approximately balances the negative part of TKE production. This is an energy transfer from turbulent fluctuations to the component 
$\bar {\boldsymbol {u}}_1$ of the mean flow.
Each term contributing to the balance of 
$\bar {\boldsymbol {u}}_0$ and 
$\bar {\boldsymbol {u}}_1$ is shown as a function of 
$L_z$ in figures 10(a) and 10(b). We do not show 
$\hat {\bar {A}}_0$ because 
$\hat {\bar {A}}_0\simeq -\hat {\bar {A}}_1$.
Figure 10. Spectral energy balance of the mean flow components (a) 
$\bar {\boldsymbol {u}}_0$ (uniform component) and (b) 
$\bar {\boldsymbol {u}}_1$ (large-scale flow along the laminar–turbulent interface); see (5.3a,b). Advection and dissipation of the large-scale flow, 
$\hat {\bar {A}}_1$ and 
$\hat {\bar {D}}_1$, show the strongest variations with 
$L_z$ and are optimal at the preferred wavelength 
$L_z\simeq 44$.
We obtain the following results.
(i) Production
$\hat {\bar {\varPi }}_0$, dissipation $\hat {\bar {D}}_0$ and energy injection $I$ are nearly independent of $L_z$, varying by no more than 6 % over the range shown. These $k_z=0$ quantities correspond to uniform fields in $z$, hence it is unsurprising that they depend very little on $L_z$.(ii) The nonlinear term
$\hat {\bar {A}}_1 \simeq - \hat {\bar {A}}_0$, i.e. the transfer from $\bar {\boldsymbol {u}}_0$ to $\bar {\boldsymbol {u}}_1$ that is the principal source of energy of $\bar {\boldsymbol {u}}_1$, varies strongly with $L_z$ and has a maximum at $L_z\simeq 44$. This is the reason why $\bar {\boldsymbol {u}}_0$ is minimised by $L_z\simeq 44$ (see figure 9b): more energy is transferred from $\bar {\boldsymbol {u}}_0$ to $\bar {\boldsymbol {u}}_1$.(iii) Production
$\hat {\bar {\varPi }}_1$ increases with $L_z$ and does not show an extremum at $L_z\simeq 44$ (instead, it has a weak maximum at $L_z\simeq 50$). In all cases, $\hat {\bar {\varPi }}_1< \hat {\bar {A}}_1$: the TKE feedback on the mean flow, although present, is not dominant and not selective.(iv) Dissipation
$\hat {\bar {D}}_1$ accounts for the remaining budget, and its extremum at $L_z\simeq 44$ corresponds to maximal dissipation.
The TKE balance is also modified with changing 
$L_z$. This is presented in Appendix C. The impact of TKE is, however, secondary, because of the results established in item (iii).
6. Conclusion and discussion
We have explored the appearance of patterns from uniform turbulence in plane Couette flow at 
$Re\leq 500$. We used numerical domains of different sizes to quantify the competition between featureless (or uniform) turbulence and (quasi-)laminar gaps. In minimal band units, intermittency reduces to a random alternation between two states: uniform or patterned. In large slender domains, however, gaps nucleate randomly and locally in space, and the transition to patterns takes place continuously via the regimes presented in § 3: the uniform regime in which gaps are rare and short-lived (above 
$Re \simeq 470$), and another regime (
$Re<470$) in which gaps are more numerous and long-lived. Below 
$Re\simeq 430$, the large-scale spacing of these gaps starts to dominate the energy spectrum, which is a possible demarcation of the patterned regime. With further decrease in 
$Re$, gaps eventually fill the entire flow domain, forming regular patterns. The distinction between these regimes is observed in both gap lifetimes and in the friction factor.
Spatially isolated gaps were already observed by Prigent et al. (Reference Prigent, Grégoire, Chaté and Dauchot2003), Barkley & Tuckerman (Reference Barkley and Tuckerman2005) and Rolland & Manneville (Reference Rolland and Manneville2011) (see also Manneville (Reference Manneville2015, Reference Manneville2017) and references therein). Our results confirm that pattern emergence, mediated by randomly nucleated gaps, is necessarily more complex than the supercritical Ginzburg–Landau framework initially proposed by Prigent et al. (Reference Prigent, Grégoire, Chaté and Dauchot2003) and later developed by Rolland & Manneville (Reference Rolland and Manneville2011). However, this does not preclude linear processes in the appearance of patterns, such as those reported by Kashyap et al. (Reference Kashyap, Duguet and Dauchot2022) from an ensemble-averaged linear response analysis.
The intermittency between uniform turbulence and gaps that we quantify here in the range 
$380\lesssim Re\lesssim 500$ is not comparable to that between laminar flow and bands present for 
$325 \lesssim Re \lesssim 340$. The latter is a continuous phase transition in which the laminar flow is absorbing: laminar regions cannot spontaneously develop into turbulence and can become turbulent only by contamination from neighbouring turbulent flow. This is connected to the existence of a critical point at which the correlation length diverges with a power-law scaling with 
$Re$, as characterised by important past studies (Shi et al. Reference Shi, Avila and Hof2013; Lemoult et al. Reference Lemoult, Shi, Avila, Jalikop, Avila and Hof2016; Chantry, Tuckerman & Barkley Reference Chantry, Tuckerman and Barkley2017) that demonstrated a connection to directed percolation. The emergence of gaps from uniform turbulence is of a different nature. Neither uniform turbulence nor gaps are absorbing states, since gaps can always appear spontaneously and can also disappear, returning the flow locally to a turbulent state. While the lifetimes of quasi-laminar gaps do exhibit an abrupt change in behaviour at 
$Re = 470$ (figure 4c), we observe no evidence of critical phenomena associated with the emergence of gaps from uniform turbulence. Hence the change in behaviour appears to be smooth, in fact. This is also true in pipe flow where quasi-laminar gaps form, but not patterns (Avila & Hof Reference Avila and Hof2013; Frishman & Grafke Reference Frishman and Grafke2022).
We used the pattern wavelength as a control parameter, via either the domain size or the initial condition, to investigate the existence of a preferred pattern wavelength. We propose that the finite spacing between gaps, visible in both local gaps and patterned regimes, is selected by the preferred size of their associated large-scale flow. Once gaps are sufficiently numerous and patterns are established, their average wavelength increases with decreasing 
$Re$, with changes in wavelength in a similar vein to the Eckhaus picture (see, e.g., Busse (Reference Busse1981)).
The influence of the large-scale flow in wavelength selection is analysed in § 5, where we carried out a spectral analysis like that in Gomé et al. (Reference Gomé, Tuckerman and Barkley2023) for various sizes of the minimal band unit. In particular, we investigated the roles of the turbulent fluctuations and the mean flow, which is in turn decomposed into its uniform component 
$\bar {\boldsymbol {u}}_0$ and its trigonometric component 
$\bar {\boldsymbol {u}}_1$, associated with the large-scale flow along the laminar–turbulent interface. Our results demonstrate a maximisation of the energy and dissipation of 
$\bar {\boldsymbol {u}}_1$ by the wavelength naturally preferred by the flow, and this is associated primarily with an optimisation of the advective term 
$(\bar {\boldsymbol {u}} \boldsymbol {\cdot }\boldsymbol {\nabla }) \bar {\boldsymbol {u}}$ in the mean flow equation. This term redistributes energy between modes 
$\bar {\boldsymbol {u}}_0$ and 
$\bar {\boldsymbol {u}}_1$, and is mostly responsible for energising the large-scale along-band flow. Turbulent fluctuations are of secondary importance in driving the large-scale flow and do not play a significant role in the wavelength selection. Our results of maximal transport of momentum and dissipation of the large-scale flow are therefore analogous to the principles mentioned by Malkus (Reference Malkus1956) and Busse (Reference Busse1981). Explaining this observation from first principles remains a prodigious task.
It is essential to understand the creation of the large-scale flow around a randomly emerging laminar hole. The statistics obtained in our tilted configuration should be extended to large streamwise–spanwise oriented domains, where short-lived and randomly nucleated holes might align in the streamwise direction (Manneville Reference Manneville2017, figure 5). Presumably, this occurs at 
$Re$ above the long-lived-gap regime, in which the two gap orientations 
$\pm \theta$ compete. The selected pattern angle might also maximise the dissipation of the large-scale flow, similarly to what we found for the preferred wavelength. Furthermore, a more complete dynamical picture of gap creation is needed. The excitable model of Barkley (Reference Barkley2011a) might provide a proper framework, as it accounts for the emergence of both anti-puffs (Frishman & Grafke Reference Frishman and Grafke2022) and periodic solutions (Barkley Reference Barkley2011b). Connecting this model to the Navier–Stokes equations is, however, a formidable challenge. Our work emphasises the necessity of including the effects of the advective large-scale flow (Barkley & Tuckerman Reference Barkley and Tuckerman2007; Duguet & Schlatter Reference Duguet and Schlatter2013; Klotz, Pavlenko & Wesfreid Reference Klotz, Pavlenko and Wesfreid2021; Marensi, Yalnız & Hof Reference Marensi, Yalnız and Hof2023) to adapt this model to the establishment of the turbulent–laminar patterns of both preferred wavelength and angle observed in planar shear flows.
Acknowledgements
The authors wish to thank A. Frishman, T. Grafke and Y. Duguet for fruitful discussions, as well as the referees for their useful suggestions.
Funding
The calculations for this work were performed using high performance computing resources provided by the Grand Equipement National de Calcul Intensif at the Institut du Développement et des Ressources en Informatique Scientifique (IDRIS, CNRS) through grant no. A0102A01119. This work was supported by a grant from the Simons Foundation (grant no. 662985).
Declaration of interests
The authors report no conflict of interest.
Appendix A. Laminar and turbulent distributions in pipe versus Couette flows
From figures 3(c) and 3(d), both distributions of laminar and turbulent lengths, 
$L_{lam}$ and 
$L_{turb}$, are exponential for large enough lengths, similarly to pipe flow (Avila & Hof Reference Avila and Hof2013). It is, however, striking that the distributions of 
$L_{lam}$ and 
$L_{turb}$ have different shapes for 
$L_{lam}$ or 
$L_{turb}>10$ in plane Couette flow: 
$L_{lam}$ shows a sharper distribution, whereas 
$L_{turb}$ is more broadly distributed. We present in figures 11(a) and 11(b) the cumulative distributions of 
$L_{lam}$ and 
$L_{turb}$ for a complementary analysis. We focus on the characteristic length 
$l_{turb}^*$ or 
$l_{lam}^*$ for which 
$P(L_{lam}>l_{lam}^*) = P(L_{turb}>l_{turb}^*) = 10^{-2}$: for example, 
$l_{lam}^*=15.5$ and 
$l_{turb}^*=26.5$ at 
$Re=440$, and 
$l_{lam}^*=23.4$ and 
$l_{turb}^*=30.3$ at 
$Re=400$. We see that 
$l_{turb}^*$ and 
$l_{lam}^*$ are of the same order of magnitude. This differs from the same measurement in pipe flow, carried out by Avila & Hof (Reference Avila and Hof2013, figure 2): 
$l_{lam}^*=6$ and 
$l_{turb}^* \simeq 50$ at 
$Re=2800$, and 
$l_{lam}^*\simeq 17$ and 
$l_{turb}^* \simeq 160$ at 
$Re=2500$ (as extracted from their figure 2). This confirms that turbulent and laminar spacings are of the same order of magnitude in plane Couette flow, contrary to pipe flow.
Figure 11. Cumulative distribution of (a) laminar gaps and (b) turbulent zones, for various 
$Re$.
Appendix B. Wavelet transform
We introduce the one-dimensional continuous wavelet transform of the velocity 
$\boldsymbol {u}(z,t)$ taken along the line 
$(x,y)=(L_x/2, 0)$:
\begin{equation} \tilde{\boldsymbol{u}}(z, r, t) = C_\psi^{-1/2} r^{-1/2} \int_{0}^{L_z} \psi^* \left(\frac{z^\prime-z}{r}\right) \boldsymbol{u}(z^\prime, t) \,{\rm d}z^\prime. \end{equation}
Here, 
$\psi$ is the Morlet basis function, defined in Fourier space as 
$\hat {\psi }(k) = {\rm \pi}^{-1/4} \exp ({-(k - k_\psi )^2/2})$ for 
$k>0$. Its central wavenumber is 
$k_\psi = 6 /\Delta z$, where 
$\Delta z$ is the grid spacing. The scale factor 
$r$ is related to wavelength via 
$\lambda \simeq 2{\rm \pi} r/k_\psi$, and 
$C_\psi \equiv \int |k|^{-1} |\hat {\psi }(k)|^2 \,\text {d}k$ is a normalization constant. Tildes are used to designate wavelet transformed quantities. The inverse transform is
\begin{equation} \boldsymbol{u}(z, t) = C_\psi^{-1/2} \int_0^{\infty} \int_{-\infty}^{\infty} r^{-1/2}\,\psi\left(\frac{z-z^\prime}{r}\right)\tilde{\boldsymbol{u}}(z^\prime, r, t) \, \frac{{\rm d}z^\prime\,{\rm d}r}{r^2}. \end{equation}
The wavelet transform is related to the Fourier transform in 
$z$ by
\begin{equation} \tilde{\boldsymbol{u}}(z, r, t) = \frac{1}{2{\rm \pi}}\,C_\psi^{-1/2} r^{1/2}\int_{-\infty}^{\infty} \hat{\psi}(r k_z)\,\hat{\boldsymbol{u}} (k_z, t)\,{\rm e}^{{\rm i} k_z z} \,\text{d} k_z. \end{equation}We then define the most energetic instantaneous wavelength as
\begin{equation} \tilde{\lambda}_{max}(z,t)= \frac{2{\rm \pi}}{k_\psi}\ \underset{r}{\text{argmax}} \ |\tilde{\boldsymbol{u}}(z,r,t)|^2. \end{equation}
The characteristic evolution of 
$\tilde {\lambda }_{max}(z,t)$ is illustrated in figure 12(b) for the flow case corresponding to figure 12(a). Regions in which 
$\tilde {\lambda }_{max}$ is large 
$({>}10)$ and dominated by a single value correspond to the local patterns observed in figure 12(a). In contrast, in regions where 
$\tilde {\lambda }_{max}$ is small 
$({<}10)$ and fluctuating, the turbulence is locally uniform.
Figure 12. Space–time visualisation of a quench experiment at 
$Re=430$. (a) Spanwise velocity (blue 
$-0.2$, white 0, red 0.2). (b) Plot of 
$\tilde {\lambda }_{max}(z,t)$, defined by (B4), which quantifies the presence of local large-scale modulations within the flow. Dark blue zones where 
$\tilde {\lambda }_{max}(z,t)<10$ correspond to locally featureless turbulence in (a). Large-scale modulation of gaps at different wavelengths are visible as the green-to-red spots in (b).
This space–time intermittency of the patterns is quantified by measuring
\begin{equation} f_{L/S} = \left\langle \varTheta (\tilde{\lambda}_{max} (z,t) - 10) \right\rangle_{z,t} \end{equation}
and is shown in figure 13(a) as a function of 
$Re$ giving a behaviour similar to that of Fourier-based variable 
$e_{L/S}$ defined in the main text (3.5).
Figure 13. (a) Space–time fraction of large to small wavelengths obtained by wavelet transform; 
$f_{L/S}$ crosses 0.5 at 
$Re \simeq 427 \simeq Re_{pg}$. (b) Sensitivity of the stability analysis in 4.2 with regard to threshold 
$H_{stab}$, at 
$Re=430$.
Appendix C. Turbulent kinetic energy balance for various $L_z$
In this appendix, we address the balance of TKE 
$\hat {K}(k_z)$, written here in 
$y$-integrated form at a specific mode 
$k_z$ (see (5.3) of Gomé et al. (Reference Gomé, Tuckerman and Barkley2023) and the methodology in e.g. Bolotnov et al. Reference Bolotnov, Lahey, Drew, Jansen and Oberai2010; Lee & Moser Reference Lee and Moser2015; Mizuno Reference Mizuno2016; Cho, Hwang & Choi Reference Cho, Hwang and Choi2018):
\begin{equation} 0 = \hat{\varPi} - \hat{D} + \hat{A} + \hat{T}_{nl}, \end{equation}
where the variables in (C1) indicate 
$y$-integrated quantities
\begin{equation} \left.\begin{gathered} \hat{\varPi} (k_z)\equiv- \mathcal{R}\left\{ \int_{-1}^1 \overline{\widehat{u_j^\prime}^{*} \widehat{ \bar{u}_i\,\frac{\partial u_j^\prime}{\partial x_i} }} \, \text{d}y \right\},\quad \hat{D} (k_z)\equiv \frac{2}{Re} \int_{-1}^1 \overline{ \widehat{s_{ij}^\prime} \,\widehat{s_{ij}^\prime}^*}\,\text{d}y, \\ \hat{T}_{nl} (k_z)\equiv-\mathcal{R}\left\{ \int_{-1}^1 \overline{\widehat{u_j^\prime}^{*} \widehat{ u_i^\prime\,\frac{\partial u_j^\prime}{\partial x_i} }} \, \text{d}y \right\},\quad \hat{A} (k_z)\equiv-\mathcal{R}\left\{ \int_{-1}^1 \overline{\widehat{u_j^\prime}^{*} \widehat{ \bar{u}_i\,\frac{\partial u_j^\prime}{\partial x_i} }} \,\text{d}y \right\}, \end{gathered}\right\} \end{equation}
respectively standing for production, dissipation, triadic interaction and advection terms. We recall that 
$\overline {({\cdot })}$ is an average in 
$(x,t)$. The 
$y$ evolution of the energy balance was analysed in Gomé et al. (Reference Gomé, Tuckerman and Barkley2023).
Gomé et al. (Reference Gomé, Tuckerman and Barkley2023) reported robust negative production at large scales, along with inverse nonlinear transfers to large scales. If 
$k_{rolls}=1.41$ denotes the scale of rolls and streaks, then this inverse transfer occurs for 
$k_z < k_{LS} = 0.94$, while a downward transfer occurs for 
$k_z>k_{SS}= 3.6$. (We refer the reader to figure 5 of Gomé et al. Reference Gomé, Tuckerman and Barkley2023.) This spectral organisation of the energy balance will be quantified by the following transfer terms arising from (C2):
\begin{align} \hat{T}_{LS} \equiv \sum_{k_z=0}^{k_{LS}} \hat{T}_{nl} (k_z), \quad \hat{T}_{SS} \equiv \sum_{k_z=k_{SS}}^{\infty} \hat{T}_{nl}(k_z) ,\quad \hat{D}_{LS} \equiv \sum_{k_z=0}^{k_{LS}} \hat{D} (k_z), \quad \hat{A}_{LS} \equiv \sum_{k_z=0}^{k_{LS}} \hat{A} (k_z), \end{align}
where 
$\hat {T}_{LS}$ quantifies transfer to large scales, 
$\hat {T}_{SS}$ quantifies the transfer to small scales, 
$\hat {D}_{LS}$ quantifies the dissipation at large scales, and 
$\hat {A}_{LS}$ is a transfer of energy from the mean flow to the large fluctuating scales. Large-scale production is not shown here, as we presented in figure 10(b) a similar measurement of large-scale turbulent transfer to the mean flow, via 
$\hat {\bar {\varPi }}_1$.
The variables defined in (C3a–d) are displayed in figure 14 as functions of 
$L_z$. Here, 
$\hat {T}_{LS}$ is minimal at 
$L_z \simeq 44$, and 
$\hat {D}_{LS}$ is minimal at 
$L_z \simeq 40$. Contrary to 
$\hat {T}_{LS}$, 
$\hat {T}_{SS}$ is relatively constant with 
$L_z$ (green dashed line in figure 14), with variation approximately 
$6\,\%$. This demonstrates that transfers to small scales are unchanged with 
$L_z$. Large-scale TKE advection decays with increasing 
$L_z$, hence it does not play a role in the preference of a wavelength. Our results show that the balance at large scale is minimised at 
$L_z\simeq 44$, confirming the less important role played by turbulent fluctuations in the wavelength selection, compared to that of the mean-flow advection reported in the main text.
Figure 14. Evolution of the large-scale TKE balance with 
$L_z$ (see (C3a–d)).
Open access
