---
abstract: |
    A turbulent-laminar banded pattern in plane Couette flow is studied numerically. This pattern is statistically steady, is oriented obliquely to the streamwise direction, and has a very large wavelength relative to the gap. The mean flow, averaged in time and in the homogeneous direction, is analysed. The flow in the quasi-laminar region is not the linear Couette profile, but results from a non-trivial balance between advection and diffusion. This force balance yields a first approximation to the relationship between the Reynolds number, angle, and wavelength of the pattern. Remarkably, the variation of the mean flow along the pattern wavevector is found to be almost exactly harmonic: the flow can be represented via only three cross-channel profiles as ${{\bf U}}(x,y,z) \approx {{\bf U}}_0(y) +
    {{\bf U}}_c(y) \cos(kz) + {{\bf U}}_s(y)\sin(kz)$. A model is formulated which relates the cross-channel profiles of the mean flow and of the Reynolds stress. Regimes computed for a full range of angle and Reynolds number in a tilted rectangular periodic computational domain are presented. Observations of regular turbulent-laminar patterns in other shear flows – Taylor-Couette, rotor-stator, and plane Poiseuille – are compared.
bibliography:
- 'bands.bib'
date: '?? and in revised form ??'
title: |
     Mean Flow of Turbulent-Laminar Patterns\
    in Plane Couette Flow
---

Introduction
============

Pattern formation is associated with the spontaneous breaking of spatial symmetry. Many of the most famous and well-studied examples of pattern formation come from fluid dynamics. Among these are the convection rolls which spontaneously form in a uniform layer of fluid heated from below and the Taylor cells which form between concentric rotating cylinders. In these cases continuous translational symmetries are broken by the cellular flows beyond critical values of the control parameter – the Rayleigh number or Taylor number.

A fundamentally new type of pattern has been discovered in large-aspect-ratio shear flows in recent years by researchers at GIT-Saclay [@Prigent_arxiv; @Prigent_PRL; @Prigent_PhysD; @Prigent_IUTAM; @Bottin_EPF]. Figure \[fig:experiment\] shows an example from plane Couette experiments performed by these researchers. One sees a remarkable spatially-periodic pattern composed of distinct regions of turbulent and laminar flow. The pattern itself is essentially stationary. The pattern wavelength is large compared with the gap between the plates and its wavevector is oriented obliquely to the streamwise direction.

The pattern emerges spontaneously from featureless turbulence as the Reynolds number is [*decreased*]{}. This is illustrated in figure \[fig:R350isbanded\] with time series from our numerical simulations of plane Couette flow for decreasing Reynolds number (conventionally defined based on half the velocity difference between the plates and half the gap). At Reynolds number 500, the flow is uniformly turbulent. Following a decrease in the Reynolds number below 400 (specifically 350 in figure \[fig:R350isbanded\]) the flow organises into three regions of relatively laminar flow and three regions of more strongly turbulent flow. While the fluid in the turbulent regions is very dynamic, the pattern is essentially steady.

Shear flows exhibiting regular coexisting turbulent and laminar regions have a been known for many years. In the mid 1960’s, a state known as spiral turbulence was discovered [@Coles; @vanAtta; @ColesvanAtta] in counter-rotating Taylor-Couette flow. Consisting of a turbulent and a laminar region, each with a spiral shape, spiral turbulence was further studied in the 1980s [-@Andereck; -@Hegseth]. Experiments by the Saclay researchers [@Prigent_arxiv; @Prigent_PRL; @Prigent_PhysD; @Prigent_IUTAM] in a very large aspect-ratio Taylor-Couette system have shown that the turbulent and laminar regions in fact form a periodic pattern, of which the original observations of Coles and van Atta comprised only one wavelength. Analogues of these states occur in other shear flows as well. [@Cros] discovered large-scale turbulent spirals in the shear flow between a stationary and a rotating disk. [-@Tsukahara] observed oblique turbulent-laminar bands in plane Poiseuille flow. A unified Reynolds number based on the shear and the half-gap can be defined for these different flows [@Prigent_PhysD] and is described in the Appendix. When converted to comparable quantities in this way, the Reynolds-number thresholds, wavelengths, and angles are similar for all of these turbulent patterned flows. The patterns are always found near the minimum Reynolds numbers for which turbulence can exist in the flow.

In this paper we present a detailed analysis of these turbulent-laminar patterns. We will focus on a single case – the periodic pattern at Reynolds number 350. From computer simulations, we obtain the flow and identify the symmetries of the patterned state. We consider in detail the force balance responsible for maintaining the pattern. From the symmetries and harmonic content we are able to reduce the description to six ordinary-differential equations which very accurately describe the patterned mean flow.

![Photograph of a turbulent-laminar pattern in plane Couette flow from the Saclay experiment. Light regions correspond to turbulent flow and dark regions to laminar flow. The striped pattern of alternating laminar and turbulent flow forms with a wavevector ${\bf k}$ oblique to the streamwise direction. The wavelength is approximately 40 times the half-gap between the moving walls. The lateral dimensions are 770 by 340 half-gaps and the Reynolds number is $Re=385$. Figure reproduced with permission from Prigent [[*et al.*]{}]{}[]{data-label="fig:experiment"}](Figures/prigent_bands_1.ps "fig:"){width="8cm"} ![Photograph of a turbulent-laminar pattern in plane Couette flow from the Saclay experiment. Light regions correspond to turbulent flow and dark regions to laminar flow. The striped pattern of alternating laminar and turbulent flow forms with a wavevector ${\bf k}$ oblique to the streamwise direction. The wavelength is approximately 40 times the half-gap between the moving walls. The lateral dimensions are 770 by 340 half-gaps and the Reynolds number is $Re=385$. Figure reproduced with permission from Prigent [[*et al.*]{}]{}[]{data-label="fig:experiment"}](Figures/exp_labels.eps "fig:"){width="2.3cm"}

![Space-time plot from numerical simulations of plane Couette flow showing the spontaneous formation of a turbulent-laminar pattern at $Re=350$. The kinetic energy in the mid-plane is sampled at 32 equally spaced points along an oblique cut (in the direction of pattern wavevector) through three wavelengths of the pattern. At time zero, $Re=500$ and the flow is uniformly turbulent. Over about 3000 time units $Re$ is decreased in steps to 350, and then held constant.[]{data-label="fig:R350isbanded"}](Figures/space_time.eps){width="8cm"}

Preliminaries
=============

Geometry
--------

The unusual but key feature of our study of turbulent-laminar patterns is the use of simulation domains aligned with the pattern wavevector and thus tilted relative to the streamwise-spanwise directions of the flow. Figure \[fig:domains\] illustrates this and defines our coordinate system. In figure \[fig:domains\](a) a simulation domain is shown as it would appear relative to an experiment, figure \[fig:experiment\], in which the streamwise direction (defined by the direction of plate motion) is horizontal. The near (upper) plate moves to the right and the far (lower) plate to the left in the figure. As we have discussed in detail [@Barkley_PRL; @Barkley_IUTAM], simulating the flow in a tilted geometry has advantages in reducing computational expense and in facilitating the study of pattern orientation and wavelength selection. The important point for the present study is that the coordinates are aligned to the patterns. The $z$ direction is parallel to the pattern wavevector while the $x$ direction is perpendicular to the wavevector (compare figure \[fig:domains\](a) with figure \[fig:experiment\]).

Figures \[fig:domains\](b) and (c) show the simulation domain as it will be oriented in this paper. In this orientation the streamwise direction is tilted at angle $\theta$ (here $24^\circ$) to the $x$ direction. This choice of angle is guided by the experimental results and by our previous simulations. (In past publications [@Barkley_PRL; @Barkley_IUTAM] we have used un-primed $x-z$ coordinates for those aligned along spanwise-streamwise directions and primes for coordinates tilted with the simulation domain. Here we focus exclusively on coordinates fixed to the simulation domain and so for convenience denote them without primes.) In these tilted coordinates, the streamwise direction is $${{\bf e_x}}\cos\theta +{{\bf e_z}}\sin\theta \equiv \alpha {{\bf e_x}}+ \beta{{\bf e_z}}\label{eq:streamwise}$$ where $$\alpha \equiv  \cos\theta = \cos(24^\circ)=0.913,  ~~~~~~~ 
\beta \equiv \sin\theta = \sin(24^\circ)=0.407.
\label{eq:alphabeta}$$

We take $L_x=10$, for the reasons explained in [@Jimenez; @Hamilton; @Waleffe_03; @Barkley_PRL; @Barkley_IUTAM]. Essentially, $L_x \sin\theta$ must be near 4 in order to contain one pair of streaks or spanwise vortices, which are necessary to the maintenance of low Reynolds number wall-bounded turbulence. Although our simulations are in a three-dimensional domain, we will average the results in the homogeneous $x$ direction, as will be explained in section \[sec:average\]. For most purposes it is sufficient to view the flow in the $z-y$ coordinates illustrated in figure \[fig:domains\](c). The midplane between the plates corresponds to $y=0$.

The length $L_z$ of our computational domain is guided by the experimental results and by our previous simulations. One of the distinctive features of the turbulent-laminar patterns is their long wavelength relative to the gap between the plates. A standard choice for length units in plane Couette flow is the half-gap between the plates. In the simulation with $L_z=120$ and $\theta=24^\circ$ shown in figure \[fig:R350isbanded\], a pattern of wavelength 40 emerged spontaneously from uniform turbulence when the Reynolds number was lowered to $Re=350$. For this reason, the simulations we will describe below are conducted with $L_z=\lambda_z=40$. The corresponding wavenumber is $$k\equiv \frac{2\pi}{40}=0.157.$$ This large wavelength, or small wavenumber, expresses the fact that the pattern wavelength in $z$ is far greater than the cross-channel dimension.

![Computational domain oriented at angle $\theta$ to the streamwise-spanwise directions. The $z$ direction is aligned to the pattern wavevector while the $x$ direction is perpendicular to the pattern wavevector. The turbulent region is represented schematically by hatching. (a) Domain oriented with streamwise velocity horizontal, as in figure \[fig:experiment\]. (b) Domain oriented with $z$ horizontal, as it will be represented in this paper. In (a), (b) the near (upper) plate moves in the streamwise direction; the far (lower) plate in the opposite direction. (c) View between the plates. []{data-label="fig:domains"}](Figures/domains.eps){width="10cm"}

Equations and Numerics {#sec:eqnum}
----------------------

The flow is governed by the incompressible Navier–Stokes equations

&=& -([[**u**]{}]{})[[**u**]{}]{} - p + [\^2]{}[[**u**]{}]{},\
0&=&[[**u**]{}]{}, \[eq:nse\]

where ${{\bf u}}({{\bf x}},t)$ is the velocity field and $p({{\bf x}},t)$ is the static pressure. Without loss of generality the density is taken to be one. The equations have been nondimensionalized by the plate speed and the half gap between the plates. $\Omega$ is the tilted computational domain discussed in the previous section.

No-slip boundary conditions are imposed at the plates and periodic boundary conditions are imposed in the lateral directions. In our coordinates the conditions are

[[**u**]{}]{}(x,y=1,z) & = & ([[**e\_x**]{}]{}+ [[**e\_z**]{}]{})\
[[**u**]{}]{}(x + L\_x,y,z) & = & [[**u**]{}]{}(x,y,z)\
[[**u**]{}]{}(x,y,z + L\_z) & = & [[**u**]{}]{}(x,y,z). \[eq:bcs\]

Linear Couette flow ${{\bf u}}^L$ is a solution to –, which is stable for all ${\mbox{\textit{Re}}}$ and satisfies $${{\boldsymbol{\nabla}}^2}{{\bf u}}^L = ({{\bf u}}^L\cdot{{\boldsymbol{\nabla}}}){{\bf u}}^L = 0$$ In our tilted coordinate system, $${{\bf u}}^L= y({{\bf e_x}}\cos\theta +{{\bf e_z}}\sin\theta) = y(\alpha {{\bf e_x}}+ \beta{{\bf e_z}}) 
= u^L{{\bf e_x}}+w^L{{\bf e_z}}\label{eq:uL}$$

The Navier-Stokes equations  with boundary conditions  are simulated using the spectral-element ($x$-$y$) – Fourier ($z$) code [Prism]{} [@Henderson]. We use a spatial resolution consistent with previous studies [@Hamilton; @Waleffe_03]. Specifically, for a domain with dimensions $L_x=10$ and $L_y=2$, we use a computational grid with 10 elements in the $x$ direction and 5 elements in the $y$ direction. Within each element, we use $8$th order polynomial expansions for the primitive variables. In the $z$ direction, a Fourier representation is used and the code is parallelized over the Fourier modes. Our domain with $L_z=40$ is discretized with 512 Fourier modes or gridpoints. Thus the total spatial resolution we use for the $L_x\times L_y
\times L_z = 10\times 2\times 40$ domain can be expressed as $N_x \times
N_{y} \times N_z = 81 \times 41 \times 512 = 1.7\times 10^6$ modes or gridpoints.

Dataset and averaging {#sec:average}
---------------------

The focus of this paper is the mean field calculated from the simulation illustrated by the spatio-temporal diagram in figure \[fig:average\](a). The velocity field in the portion of the domain shows high-frequency and high-amplitude fluctuations, while the flow in the right portion is basically quiescent. We will call the flow on the left turbulent, even though it could be argued that it is not fully developped turbulence. We will call the flow on the right laminar, even though occasional small fluctuations can be seen in this region.

The turbulent-laminar pattern subsists during the entire simulation of $14
\times 10^3$ time units. However the pattern undergoes short-scale “jiggling”, seen particularly at the edges of the turbulent regions, and longer-scale drifting or wandering in the periodic $z$ direction. We seek to describe the field which results from smoothing the turbulent fluctuations, but for which drifting is minimal, by averaging over an appropriate time interval. The desired averaging time interval represents a compromise between the short and long timescales. We have chosen to average the flow in figure \[fig:average\](a) over the shaded time interval $[t,t+T]=[6000,8000]$, during which the pattern is approximately stationary.

The time-averaged flow is homogeneous in $x$-direction. This is illustrated in figure \[fig:average\](b) where we plot one of the velocity components time-averaged flow over the interval $[6000,8000]$. Cuts at different $x$ locations show that there is essentially no variation in the $x$-direction. All other quantities are similarly independent of $x$. It is therefore appropriate to consider mean flows as averages over the $x$ direction as well as over the time.

We define mean flows as

(y,z) & & \_t\^[t+T]{} \_0\^[L\_x]{} [[**u**]{}]{}(x,y,z,t) dx dt\
[p ]{}(y,z) && \_t\^[t+T]{} \_0\^[L\_x]{} p(x,y,z,t) dx dt.

The mean fields obey the averaged Navier-Stokes equations

\[eq:NSavg\] 0&=&-() - ( ) -[p ]{}+ [\^2]{}\
0&=&,

where $${\tilde{{{\bf u}}}}\equiv {{\bf u}}- {\langle {{\bf u}}\rangle}$$ is the fluctuating field and ${\langle}{\rangle}$ denotes $x$-$t$ average. The mean fields are subject to the same boundary conditions as equations . We denote the Reynolds-stress force from the fluctuating field in equations  by ${\bf F}$: $${\bf F}\equiv - {\langle}\left({\tilde{{{\bf u}}}}\cdot {{\boldsymbol{\nabla}}}\right) {\tilde{{{\bf u}}}}{\rangle}= - {{\boldsymbol{\nabla}}\cdot}{\langle}{\tilde{{{\bf u}}}}{\tilde{{{\bf u}}}}{\rangle},$$ We shall focus almost exclusively on the difference between the mean flow and linear Couette flow, for which we introduce the notation $${{\bf U}}\equiv {\langle {{\bf u}}\rangle}- {{\bf u}}^L,
\label{eq:UminusC}$$ as well as $P\equiv {\langle}p {\rangle}$.

Letting the components of ${{\bf U}}$ be denoted by $(U,V,W)$ and the components of ${\bf F}$ be denoted by $({F}^U,{F}^V,{F}^W)$, then the averaged Navier-Stokes equations for the deviation from linear Couette flow in component form become

\[eq:RANS\] $$\begin{aligned}
0 &=& -\left(V {\ensuremath{\partial}}_y + (W + \beta y)~ {\ensuremath{\partial}}_z \right) (U + \alpha y) 
~~~~~~~~~~+ {\frac{1}{{\mbox{\textit{Re}}}}}({\ensuremath{\partial}}^2_y + {\ensuremath{\partial}}^2_z) U + {F}^U \label{eq:NSavgU}  \\
0 &=& -\left(V {\ensuremath{\partial}}_y + (W + \beta y)~ {\ensuremath{\partial}}_z  \right) V ~~~~~~~~~~
 -{\ensuremath{\partial}}_y P + {\frac{1}{{\mbox{\textit{Re}}}}}({\ensuremath{\partial}}^2_y + {\ensuremath{\partial}}^2_z) V + {F}^V \label{eq:NSavgV}\\
0 &=& -\left(V {\ensuremath{\partial}}_y + (W + \beta y)~ {\ensuremath{\partial}}_z \right) (W + \beta y) 
-{\ensuremath{\partial}}_z P + {\frac{1}{{\mbox{\textit{Re}}}}}({\ensuremath{\partial}}^2_y + {\ensuremath{\partial}}^2_z) W + {F}^W \label{eq:NSavgW}\end{aligned}$$ $$\begin{aligned}
0 &=& {\ensuremath{\partial}}_y V + {\ensuremath{\partial}}_z W. \label{eq:RANSdiv}\end{aligned}$$

${{\bf U}}$ is required to satisfy homogeneous boundary conditions at the plates $${{\bf U}}(y=\pm1,z) = {\bf 0}
\label{eq:homogbc}$$ and periodic boundary conditions in $z$.

A system of this type, with three components depending on two coordinates, is sometimes called 2.5 dimensional. The transverse, or out-of-plane flow $U(y,z)$ appears only in the first equation and is effectively a passive scalar advected by the in-plane flow $(V,W)$ and driven by the Reynolds-stress force ${F}^U$. The in-plane flow can be expressed in terms of a streamfunction $\Psi$ where $$V{{\bf e_y}}+W{{\bf e_z}}={{\bf e_x}}\times \nabla \Psi = 
-{\ensuremath{\partial}}_z \Psi {{\bf e_y}}+ {\ensuremath{\partial}}_y \Psi {{\bf e_z}}.
\label{eq:definepsi}$$ We shall use both $(U,V,W)(y,z)$ and $(U,\Psi)(y,z)$ to describe the mean flows.

[ ]{}

[ ]{}

![ $U(y,z)$: transverse component of mean flow. $\Psi(y,z)$: streamfunction of in-plane mean flow. A long cell extends from one laminar-turbulent boundary to the other. Gradients of $\Psi$ are much larger in $y$ than in $z$, [i.e.]{} $|W|\gg|V|$. In the laminar region at the center, $W, V\approx 0$. $E_{\rm turb}(y,z)$: mean turbulent kinetic energy $\langle
{\tilde{{{\bf u}}}}\cdot{\tilde{{{\bf u}}}}\rangle/2$. There is a phase difference of $\lambda_z/4=10$ between extrema of $E_{\rm turb}$ and of $U$. $P(y,z)$: mean pressure field. Pressure gradients are primarily in the $y$ direction and within the turbulent region. Color ranges for each field from blue to red: $U$ \[–0.4, 0.4\], $\Psi$ \[0, 0.09\], $E_{\rm turb}$ \[0, 0.4\], P \[0, 0.007\] .[]{data-label="fig:upsikep"}](Figures/average.ps){width="14cm"}

![ $U(y,z)$: transverse component of mean flow. $\Psi(y,z)$: streamfunction of in-plane mean flow. A long cell extends from one laminar-turbulent boundary to the other. Gradients of $\Psi$ are much larger in $y$ than in $z$, [i.e.]{} $|W|\gg|V|$. In the laminar region at the center, $W, V\approx 0$. $E_{\rm turb}(y,z)$: mean turbulent kinetic energy $\langle
{\tilde{{{\bf u}}}}\cdot{\tilde{{{\bf u}}}}\rangle/2$. There is a phase difference of $\lambda_z/4=10$ between extrema of $E_{\rm turb}$ and of $U$. $P(y,z)$: mean pressure field. Pressure gradients are primarily in the $y$ direction and within the turbulent region. Color ranges for each field from blue to red: $U$ \[–0.4, 0.4\], $\Psi$ \[0, 0.09\], $E_{\rm turb}$ \[0, 0.4\], P \[0, 0.007\] .[]{data-label="fig:upsikep"}](Figures/upsikepmake.ps){width="14cm"}

Results
=======

We present a characterisation of the turbulent-laminar pattern at ${\mbox{\textit{Re}}}=350$. We describe in detail the mean flow, its symmetries, and the dominant force balances within the flow. Our goal here is not to consider closures for averaged Navier-Stokes equations (\[eq:RANS\]). We will make no attempt to model the turbulence, [i.e.]{} to relate the Reynolds-stress tensor ${\langle}{\tilde{{{\bf u}}}}{\tilde{{{\bf u}}}}{\rangle}$ to the mean flow ${{\bf U}}$. Instead we use fully resolved (three-dimensional, time-dependent) numerical simulations of the turbulent flow to measure both the mean field ${{\bf U}}$ and Reynolds-stress force ${\bf F}$. From these we extract the structure of these fields and the dominant force balances at play in sustaining turbulent-laminar patterns.

Mean flow {#sec:meanflow}
---------

The mean flow is visualised in figure \[fig:upsikep\] via the transverse, out-of-plane flow $U(y,z)$ and the in-plane streamfunction $\Psi(y,z)$. Recall \[equation \] that these fields are the deviations of the mean flow from linear Couette flow ${{\bf u}}^L$. The mean turbulent kinetic energy $$E_{\rm turb}\equiv\frac{1}{2}\langle {\tilde{{{\bf u}}}}\cdot{\tilde{{{\bf u}}}}\rangle$$ serves to clearly identify the turbulent region. In these and subsequent plots, the middle of the laminar region is positioned at the centre of the figure and the turbulent region at the periodic boundaries of the computational domain. In figure \[fig:upsikep\] (but not in subsequent figures), plots are extended in the $z$-direction one quarter-period beyond each periodic boundary to help visualise the flow in the turbulent region. The pattern wavelength is $\lambda_z=40$, so that $z=30$ and $z=-10$ describe the same point, as do $z=-30$ and $z=10$.

The mean flow can be described as follows. $U$ is strongest in the turbulent-laminar transition regions. In the transition region to the left of centre ($z=-10$) in figure \[fig:upsikep\], $U$ is negative and primarily in the upper half of the channel. To the right of centre ($z=10$), $U$ is positive and is seen primarily in the lower half of the channel. Comparison with turbulent kinetic energy shows that the transverse mean flow $U$ is out of phase with respect to the fluctuating field ${\tilde{{{\bf u}}}}$ by $\lambda_z/4$. This has been seen experimentally by [@ColesvanAtta] and Prigent [[*et al.*]{}]{} ([-@Prigent_PRL], [-@Prigent_PhysD], [-@Prigent_IUTAM]).

The in-plane flow $\Psi$ in figure \[fig:upsikep\] has a large-aspect ratio cellular structure consisting of alternating elliptical and hyperbolic points. The flow around the elliptical points, located in the centre of the turbulent regions, rotates in a counter-clockwise sense, opposing linear Couette flow. In the vicinity of the hyperbolic points, centred in the laminar regions, the in-plane deviation from linear Couette flow is very weak ($W$ and $V$ nearly zero).

Figure \[fig:meanflow\] shows $y$-profiles at four key points equally spaced along the pattern: centre of the laminar region, turbulent-laminar transition region, centre of the turbulent region, and the other turbulent-laminar transition region. While the $V$ profile is plotted, its variation is very small on the scale of $U$ and $W$ and can essentially be used to indicate the axis. Figure \[fig:meanflowwithCou\] shows profiles for the full mean flow ${\langle {{\bf u}}\rangle}= {{\bf U}}+ {{\bf u}}^L$ containing the linear Couette profile.

![ Same as in figure \[fig:meanflow\], but with laminar Couette flow ${{\bf u}}^L$ included.[]{data-label="fig:meanflowwithCou"}](Figures/Meanuvwmc.ps){width="14cm"}

![ Same as in figure \[fig:meanflow\], but with laminar Couette flow ${{\bf u}}^L$ included.[]{data-label="fig:meanflowwithCou"}](Figures/Meanuvw.ps){width="14cm"}

![ Mean velocity components seen in three planes with standard orientation for Couette flow. The turbulent regions are shaded. Top: velocity components in the streamwise-spanwise plane at $y=0.725$ (upper part of the channel). Middle: same except $y=-0.725$ (lower part of the channel). Bottom: flow in a constant spanwise cut. The mean velocity is shown in the enlarged region. []{data-label="fig:coles"}](Figures/coles.ps){width="10cm"}

The $U$ profiles in figure \[fig:meanflowwithCou\] are S-shaped, of the type found in turbulent Couette flow. This is to be expected in the turbulent region, even at these low Reynolds numbers. However, it is very surprising that the $U$ profile in the laminar region is also of this form. In the laminar region, local Reynolds stresses are absent (see figure \[fig:upsikep\]) and so cannot be responsible for maintaining the S-shaped velocity profile in the laminar regions. The other prominent features in figures \[fig:meanflow\] and \[fig:meanflowwithCou\] are the asymmetric profiles at the transition regions.

The relationship between the mean flow field and the regions of turbulence can be seen in figure \[fig:coles\]. Here the flow is shown in the standard orientations. In each view, greyscale indicates the size of the turbulent energy and the arrows show the mean flow within the plane. In the top two views, the flow is shown in the streamwise-spanwise planes located at $y=0.725$ and at $y=-0.725$. The next view shows the flow between the plates, [i.e.]{}in a streamwise-cross-channel plane, and the last shows an enlargement of one of the laminar-turbulent transition regions. Note that the length $L_z=40$ of our tilted computational domain corresponds to a streamwise length of $L_z/\sin\theta = 40/.407=98.3\approx
100$ and to a spanwise length of $L_z/\cos\theta= 43.78$.

The flow in figure \[fig:coles\] can be compared with the mean flow reported by [@ColesvanAtta] in experiments on turbulent spirals in Taylor-Coutte flow. Coles and van Atta measured the mean flow near the midgap between the rotating cylinders and noted an asymmetry between the mean flow into and out of turbulent regions. They found that the mean flow into turbulent regions was almost perpendicular to the turbulent-laminar interface whereas flow out of the turbulent region was almost parallel to the turbulent-laminar interface. We also observe a striking asymmetry between the mean flow into and out of the turbulent regions. The orientation of our mean flow does not agree in detail with that of Coles and van Atta, but this is most likely due to the fact that Coles and van Atta considered circular Taylor-Couette flow and measured the flow near the mid-gap. Referring to figures \[fig:upsikep\] and \[fig:meanflow\] one sees that the mid-plane ($y=0$) is not the ideal plane on which to obverve the mean flow since its structure is most pronounced between the midplane and the upper or lower walls.

Before considering the symmetries and force balances in detail, it is instructive to consider the dominant force balance just at the centre of the laminar region. Recall that one of the more interesting features of the mean flow is that the $U$ profile appears very similar to a turbulent profile, even in the absence of turbulence in the laminar region. Here the balance is dominated by advection and viscous diffusion, as shown in figure \[fig:forcepreface\]. Equation for flow in the $x$-direction is $$0 = -\left(V {\ensuremath{\partial}}_y + (W  + \beta y) {\ensuremath{\partial}}_z \right) (U + \alpha y)
 + {\frac{1}{{\mbox{\textit{Re}}}}}({\ensuremath{\partial}}^2_y + {\ensuremath{\partial}}^2_z) U + {F}^U.$$

Variations in $y$ dominate variations in $z$, [i.e.]{} the usual boundary-layer approximation $({\ensuremath{\partial}}^2_y + {\ensuremath{\partial}}^2_z) U \simeq {\ensuremath{\partial}}^2_y U$ holds; see, [e.g.]{} [@Pope]. Indeed, approximating the $y$ dependence of $U$ by the functional form $\sin(\pi y)$ suggested by figure \[fig:meanflow\], we have $$O\left(\frac{{\ensuremath{\partial}}^2_y U}{{\ensuremath{\partial}}^2_z U}\right) = 
\frac{\pi^2}{k^2} = \frac{\pi^2}{(2\pi/40)^2}=400
\label{eq:d2ytod2z}$$ This is confirmed by the second panel of figure \[fig:forcepreface\]. In the centre of the laminar region ${F}^U$, $V$, and $W$ are all negligible, so that $-\beta y {\ensuremath{\partial}}_z U$ dominates the advective terms, as shown in the third panel of figure \[fig:forcepreface\]. Thus the balance between advection and viscosity in the laminar region is $$\beta\: y \: {\ensuremath{\partial}}_z U \approx {\frac{1}{{\mbox{\textit{Re}}}}}{\ensuremath{\partial}}^2_y U.
\label{eq:firstapprox}$$

This equation is appealingly simple and yet leads immediately to some interesting conclusions. The first is that a non-zero tilt angle $\theta$ is necessary to maintain the S-shaped $U$ profile in the laminar region, since otherwise $\beta = \sin \theta = 0$ and $U$ could be at most linear in $y$ and would in fact be zero, due to the homogeneous boundary conditions (\[eq:homogbc\]). The second conclusion follows from consideration of $y$ parity. The multiplication by $y$ on the left-hand-side reverses $y$-parity, while the second derivative operator on the right-hand-side preserves $y$ parity. The conclusion is that $U$ should be decomposed into odd and even components in $y$ and equation  should actually be understood as two equations coupling the two components. Specifically, as can be seen in figure \[fig:meanflow\], $U$ is odd in $y$ in the centre of the laminar region, yet ${\ensuremath{\partial}}_z U$ must be even for equation to hold.

The remainder of the paper is devoted to formalising, demonstrating and extending this basic idea.

![Balance of forces in the center of the laminar region. Left: Forces in the $U$ direction. Advective terms (blue, solid), viscous terms (red, dashed), Reynolds-stress terms (green, dotted). Middle: Viscous terms in the $U$ direction. $(1/Re)\,{\ensuremath{\partial}}^2_y U$ (red, dashed) dominates dominates $(1/Re)\,{\ensuremath{\partial}}^2_z U$ (green, dotted). Right: Advective terms in the $U$ direction. Curves show $-\beta\,y\,{\ensuremath{\partial}}_z U$ (blue, solid) $-W{\ensuremath{\partial}}_z U$ (black, dash-dot), $-\alpha V$ (green, dotted), $-V{\ensuremath{\partial}}_y U$ (red, dashed). Right-most: Advective terms in the $W$ direction (for later reference). Curves show $-\beta\, y\, {\ensuremath{\partial}}_z W$ (blue, solid) $-W{\ensuremath{\partial}}_z W$ (black, dash-dot), $-\beta V$ (green, dotted), $-V{\ensuremath{\partial}}_y W$ (red, dashed). []{data-label="fig:forcepreface"}](Figures/forcepreface.ps){width="13cm"}

[ ]{}

Symmetry and Fourier modes {#sec:Symmetry}
--------------------------

We now consider in depth the symmetry properties of the flow. We start with the symmetries of the system before averaging, that is, the Navier-Stokes equations and boundary conditions . The system has translation symmetry in $x$ and $z$ as well as centrosymmetry under combined reflection in $x$, $y$ and $z$: $${\kappa_{xyz}}(u,v,w)(x,y,z) \equiv (-u,-v,-w)(x_0-x,-y,z_0-z)
\label{eq:Cfdef}$$ where the origin $x_0$, $z_0$ can be chosen arbitrarily. Linear Couette flow ${{\bf u}}^L$ possesses all the system symmetries, as does the mean flow at Reynolds numbers for which turbulence is statistically homogeneous in $x$ and $z$.

Note that in the absence of tilt ($\theta=0$), the system possesses two reflection symmetries. These can be taken to be ${\kappa_{xyz}}$ and reflection in the spanwise direction. For the tilted domain (at angles other than multiples of $90^\circ$), the only reflection symmetry is ${\kappa_{xyz}}$. This can be seen in figure \[fig:domains\](a): for general tilt angles $\theta$, spanwise reflection does not preserve the domain, [i.e.]{} does not leave the periodic boundaries in place. The experimental system shown in figure \[fig:experiment\] possesses spanwise reflection symmetry and hence bands can be observed in the either of the two symmetrically related angles, the choice is dictated by factor such as initial conditions. By design, our tiled computational domain precludes the symmetry-related pattern given by spanwise reflection.

The transition to the turbulent-laminar patterned state breaks symmetry. Specifically, both the mean flow and the Reynolds-stress force break $z$-translation symmetry but break neither $x$-translation symmetry nor centrosymmetry. The spatial phase of the pattern in $z$ is arbitrary, but given a phase there are two values of $z_0$, separated by half a period, for which the flow is invariant under ${\kappa_{xyz}}$, as is typical for a circle pitchfork bifurcation [-@Crawford]. As can be seen in figure \[fig:upsikep\], the values of $z_0$ about which the patterns are centrosymmetric are the centres of the laminar ($z_0=0$) and of the turbulent ($z_0=\pm 20$) regions.

The centrosymmetry operator for our averaged fields ${{\bf U}}$, which depend only on $y$ and $z$, is $${\kappa_{yz}}(U,V,W)(y,z) \equiv (-U,-V,-W)(-y,z_0-z)
\label{eq:Cdef}$$ Since the Reynolds-stress force $({F}^U, {F}^V, {F}^W)$ is centrosymmetric in the case we consider, then the averaged equations  for the mean field have centrosymmetry.

We formalise this further as follows. Any $x$-independent field ${{\bf g}}$ can be decomposed into even and odd functions of $y$ and $z$ as $${{\bf g}}(y,z)={{\bf g}}_{oo}(y,z) + {{\bf g}}_{oe}(y,z) + {{\bf g}}_{eo}(y,z) + {{\bf g}}_{ee}(y,z)
\label{eq:4fcns}$$ where, for example, ${{\bf g}}_{oe}$ is odd in $y$ and even in $z-z_0$. Applying the operator in (\[eq:Cdef\]) to (\[eq:4fcns\]), we obtain $$\begin{aligned}
{\kappa_{yz}}{{\bf g}}(y,z)
&=&-{{\bf g}}_{oo}(-y,z_0-z) - {{\bf g}}_{oe}(-y,z_0-z) 
- {{\bf g}}_{eo}(-y,z_0-z) - {{\bf g}}_{ee}(-y,z_0-z)
\nonumber\\
&=&-{{\bf g}}_{oo}(y,z) + {{\bf g}}_{oe}(y,z) + {{\bf g}}_{eo}(y,z) - {{\bf g}}_{ee}(y,z)\end{aligned}$$ For the field ${{\bf g}}$ to be centrosymmetric requires ${\kappa_{yz}}{{\bf g}}= {{\bf g}}$, so that in fact $${{\bf g}}(y,z)={{\bf g}}_{oe}(y,z) + {{\bf g}}_{eo}(y,z)
\label{eq:2fcns}$$ Table \[tab:symmetry\], as well as figure \[fig:upsikep\], shows that this is indeed the case for ${{\bf U}}$; it holds for ${\bf F}$ as well.

[ ]{}

$$\begin{aligned}
\begin{array}{c||c|c|} 
&  {\mbox{$z$ even}} & {\mbox{$z$ odd}} \\ 
\cline{1-3}
{\mbox{$y$ even}} & {0.03\%} & {25.48\%} \\ 
\cline{1-3}
{\mbox{$y$ odd}} & {74.48\%} & {0.01\%} \\ 
\cline{1-3}
\end{array}
\nonumber\end{aligned}$$

$$\begin{aligned}
\begin{array}{c||ccc|cc|} 
& \multicolumn{3}{|c|} {\mbox{$z$ even (cosine)}} & 
\multicolumn{2}{|c|} {\mbox{$z$ odd (sine)}}\\ 
\mbox{$z$ wavenumber}& 0 & k & \geq 2k & k & \geq 2k \\ \cline{1-6}
\mbox{$y$ even} &&&&\fbox{25.2\%} & {0.3\%} \\ \cline{1-6}
\mbox{$y$ odd} & \fbox{69.7\%} & \fbox{4.7\%} & 0.1\% & & \\ \cline{1-6}
\end{array}
\nonumber\end{aligned}$$

We now Fourier transform in $z$ to further decompose the mean velocity and the Reynolds-stress force. We find that the $z$-wavenumbers $0$ and $\pm k$ have contributions to ${{\bf U}}$ which are an order of magnitude higher than the remaining wavenumber combinations. See table \[tab:fourier\]. The deviation from the $z$ average is thus almost exactly trigonometric, with almost no higher harmonic content. The dominance of these terms in the Fourier series means that ${{\bf U}}$ and ${\bf F}$ can be represented by only three functions of $y$, namely: $${{\bf g}}(x,y,z)={{\bf g}}_0(y) + {{\bf g}}_c(y) \cos(kz) + {{\bf g}}_s(y) \sin(kz)
\label{eq:3fcns}$$ which is a special case of , with the first two terms of coinciding with ${{\bf g}}_{oe}(y,z)$ and the last to ${{\bf g}}_{eo}$. Thus, ${{\bf g}}_0$ and ${{\bf g}}_c$ are odd functions of $y$, while ${{\bf g}}_s$ is even. The fields thus consist of a $z$-independent component ${{\bf g}}_0$ and two components which vary trigonometrically and out of phase with one another, ${{\bf g}}_c$ dominating in the laminar and turbulent regions and ${{\bf g}}_s$ dominating in the boundaries between them. Moreover, ${{\bf g}}_s$ dominates in the bulk, since ${{\bf g}}_0$ and ${{\bf g}}_c$ are odd in $y$ and thus zero in the channel centre.

\[eq:trigexpand\] $$\begin{aligned}
&{{\bf g}}={{\bf g}}_0(y)+{{\bf g}}_c(y) && z=0:& 
\mbox{Centre of laminar region} \label{eq:g0}\\
&{{\bf g}}={{\bf g}}_0(y)+{{\bf g}}_s(y) && z=\lambda_z/4=10:& 
\mbox{Laminar-turbulent boundary} \label{eq:g10}\\
&{{\bf g}}={{\bf g}}_0(y)-{{\bf g}}_c(y) && z=\lambda_z/2=20:& 
\mbox{Centre of turbulent region} \label{eq:g20}\\
&{{\bf g}}={{\bf g}}_0(y)-{{\bf g}}_s(y) && z=3\lambda_z/4=30:& 
\mbox{Turbulent-laminar boundary} \label{eq:g30}\end{aligned}$$

Figure \[fig:fourmean\] shows the three trigonometric components, each a function of $y$, obtained by Fourier transforming $U$, $V$, and $W$. Figure \[fig:profileZuvw\] shows $U$, $V$ and $W$ as functions of $z$ at locations in the upper and lower channel and compares them with the values obtained from the trigonometric formula using the functions shown in figure \[fig:fourmean\]. Figures \[fig:udecomp\], \[fig:psidecomp\] and \[fig:Reyudecomp\] depict $U(y,z)$, $\Psi(y,z)$ and ${F}^U(y,z)$ with their trigonometric decompositions. Each of these figures uses only the three scalar functions of $y$, figure \[fig:fourmean\], to reproduce the corresponding two-dimensional field. As shown by equation , the streamfunction $\Psi$ of a centrosymmetric field has symmetry opposite to that of the velocity components, [i.e.]{} it is composed of functions of the same parity in $y$ and $z$.

![Fourier decomposition of mean velocity. $U$ component (blue, solid), $V$ component (green, dotted), $W$ component (red, dashed). $W_c \approx
  -W_0$, corresponding to the fact that $W$ shows no deviation from the linear in the laminar region. []{data-label="fig:fourmean"}](Figures/fourmean.ps){width="14cm"}

![Mean flow as a function of $z$ at $y=0.725$ (lower curves) and $y=-0.725$ (upper curves). $U$ (blue, solid), $V$ (green, dotted), $W$ (red, dashed). Dots show values calculated from trigonometric formula .[]{data-label="fig:profileZuvw"}](Figures/profileZuvw.ps){width="8cm"}

![Reynolds-stress force ${F}^U$ and its trigonometric decomposition. Color scale is ${F}^U~(-0.017,0.017)$, ${F}^U_0 (-0.0085,0.0085)$, ${F}^U_c (-0.0085,0.0085)$, ${F}^U_s (-0.0085,0.0085)$.[]{data-label="fig:Reyudecomp"}](Figures/udecompmake.ps){width="14cm"}

![Reynolds-stress force ${F}^U$ and its trigonometric decomposition. Color scale is ${F}^U~(-0.017,0.017)$, ${F}^U_0 (-0.0085,0.0085)$, ${F}^U_c (-0.0085,0.0085)$, ${F}^U_s (-0.0085,0.0085)$.[]{data-label="fig:Reyudecomp"}](Figures/psidecompmake.ps){width="14cm"}

![Reynolds-stress force ${F}^U$ and its trigonometric decomposition. Color scale is ${F}^U~(-0.017,0.017)$, ${F}^U_0 (-0.0085,0.0085)$, ${F}^U_c (-0.0085,0.0085)$, ${F}^U_s (-0.0085,0.0085)$.[]{data-label="fig:Reyudecomp"}](Figures/Reyudecompmake.ps){width="14cm"}

Figures \[fig:Reynolds\] and \[fig:fourrey\] show the three Reynolds-stress forces and their Fourier decompositions. Each component obeys ${\bf F}_c\approx-{\bf F}_0$, a necessary condition for ${\bf F}$ to vanish at the centre of the laminar region, as shown by equation and also illustrated in figure \[fig:Reyudecomp\]. More precisely, $${F}^U_c = -1.09 {F}^U_0, ~~~~~ {F}^V_c = -1.22 {F}^V_0, ~~~~~ {F}^W_c = -1.16 {F}^W_0
\label{eq:falmostzero}$$ In addition, $${\bf F}\approx -{\ensuremath{\partial}}_y {\langle}{\tilde{{{\bf u}}}}{\tilde{v}}{\rangle}.$$ as is typical for turbulent channel flows; see, [e.g.]{} [@Pope]. [ ]{}

![Reynolds-stress force ${\bf F}=-{\langle}({\tilde{{{\bf u}}}}\cdot\nabla){\tilde{{{\bf u}}}}{\rangle}$ as a function of $y$. Curves show ${F}^U$ (blue, solid), ${F}^V$ (green, dotted) and ${F}^W$ (red, dashed).[]{data-label="fig:Reynolds"}](Figures/Reynolds.ps){width="14cm"}

![Fourier decomposition of Reynolds-stress force. ${F}^U$ components (blue, solid), ${F}^V$ components (green, dotted), ${F}^W$ components (red, dashed). ${\bf F}_c\approx-{\bf F}_0$, as required for the vanishing of ${\bf F}$ at $z=0$.[]{data-label="fig:fourrey"}](Figures/fourrey.ps){width="14cm"}

Force balance for $U$ {#sec:ForceU}
---------------------

We now turn to understanding the balance of forces responsible for maintaining the mean flow profiles. We focus primarily on $U$, both because it is the component of largest amplitude and also because it appears only in equation : $U$ is subject to Reynolds-stress and viscous forces, and is advected by $(V,W)$ but is not self-advected. We begin by showing the balance of forces in the $U$ direction as a function of $z$ at locations in the upper and lower channel in figure \[fig:balance\]. One can again see the centrosymmetry of each of the forces, [i.e.]{} invariance under the combined operations of reflection in $y$ and $z$ and change of sign. The Reynolds-stress force disappears at the center of the laminar region and the advective and viscous forces exactly counterbalance, as emphasized in the figures on the right. Figure \[fig:Forceu\] shows another view of this balance, displaying the forces as a function of $y$ at four locations in $z$. As previously stated, ${\nabla^2}U$ is dominated by ${\ensuremath{\partial}}^2_y U$ and ${F}^U$ by $-{\ensuremath{\partial}}_y {\langle}{\tilde{u}}{\tilde{v}}{\rangle}$. In figure \[fig:FourForceu\], we show the Fourier-space analogue of figure \[fig:Forceu\].

![Mean forces in $U$ direction as a function of $z$ at $y=\pm 0.725$ for turbulent-laminar pattern at $Re=350$. Advective $-({{\bf U}}\cdot\nabla) U$ (blue, solid), viscous ${\nabla^2}U$ (red, dashed), and turbulent $-{\langle}({\tilde{{{\bf u}}}}\cdot\nabla) {\tilde{u}}{\rangle}$ (green, dotted) forces. In the laminar region ($z\approx 0$), the Reynolds-stress force vanishes and the viscous and advective forces are equal and opposite to one another. In figures on right, enlarged around the laminar region, ${\nabla^2}U$ and $+({{\bf U}}\cdot\nabla) U$ are shown to emphasize equality between viscous and advective forces.[]{data-label="fig:balance"}](Figures/newbalance.ps){width="17cm"}

![Advective terms in the $U$ direction, decomposed into modes. Curves show $-\beta y {\ensuremath{\partial}}_z U$ (black, solid) $-W{\ensuremath{\partial}}_z U$ (blue, dash-dot), $-\alpha V$ (green, dotted), $-V{\ensuremath{\partial}}_y U$ (red, dashed). The 1 mode is generated by the product $-k W_c U_s/2$; a second harmonic of the same small size is also generated. The $\cos(kz)$ mode is dominated by $-\beta\,y\,k\,U_s$. The $\sin(kz)$ term is dominated by $\beta\,y\,k\,U_c$ near the boundaries and by $-V_s{\ensuremath{\partial}}_y (U_0+\alpha y)$ in the bulk.[]{data-label="fig:FourNonu"}](Figures/Forceu.ps){width="14cm"}

![Advective terms in the $U$ direction, decomposed into modes. Curves show $-\beta y {\ensuremath{\partial}}_z U$ (black, solid) $-W{\ensuremath{\partial}}_z U$ (blue, dash-dot), $-\alpha V$ (green, dotted), $-V{\ensuremath{\partial}}_y U$ (red, dashed). The 1 mode is generated by the product $-k W_c U_s/2$; a second harmonic of the same small size is also generated. The $\cos(kz)$ mode is dominated by $-\beta\,y\,k\,U_s$. The $\sin(kz)$ term is dominated by $\beta\,y\,k\,U_c$ near the boundaries and by $-V_s{\ensuremath{\partial}}_y (U_0+\alpha y)$ in the bulk.[]{data-label="fig:FourNonu"}](Figures/fourtotu.ps){width="14cm"}

![Advective terms in the $U$ direction, decomposed into modes. Curves show $-\beta y {\ensuremath{\partial}}_z U$ (black, solid) $-W{\ensuremath{\partial}}_z U$ (blue, dash-dot), $-\alpha V$ (green, dotted), $-V{\ensuremath{\partial}}_y U$ (red, dashed). The 1 mode is generated by the product $-k W_c U_s/2$; a second harmonic of the same small size is also generated. The $\cos(kz)$ mode is dominated by $-\beta\,y\,k\,U_s$. The $\sin(kz)$ term is dominated by $\beta\,y\,k\,U_c$ near the boundaries and by $-V_s{\ensuremath{\partial}}_y (U_0+\alpha y)$ in the bulk.[]{data-label="fig:FourNonu"}](Figures/fournonudetail.ps){width="14cm"}

We now turn to the more complex advective forces, whose Fourier decompositions are shown in figure \[fig:FourNonu\]. The $\cos(0z)$ component of the advective force is small but non-zero. Because this term results from the product of trigonometric functions, it also provides a measure of the generation of higher harmonics, a point which we will explore further in section \[sec:Model\]. The advective $\cos(kz)$ term is well approximated by the contribution from advection by $w^L=\beta y$. The advective $\sin(kz)$ term is dominated near the walls by advection by $w^L$, but in the bulk by advection by $V$. Properties of the $\cos(kz)$ and $\sin(kz)$ modes echo their physical space counterparts: the advective term is well approximated by advection by $w^L$ in the laminar region, as was shown in figure \[fig:forcepreface\], while the advective forces in the laminar-turbulent boundaries combine advection by $w^L$ near the walls and by $V$ in the bulk.

We illustrate these conclusions via schematic visualizations of the dynamics of $U$. Figure \[fig:schematic\_cos\] illustrates the dynamics in the laminar and turbulent regions. The dynamics in the laminar region are essentially described by the simple balance between viscous diffusion of $U$ profiles and advection by linear Couette flow in $z$, given by equation . Viscous diffusion tends to reduce curvature, but the profiles have greater curvature upstream (to the left for the upper channel, to the right for the lower channel). Hence advection replenishes the curvature damped by viscosity. However, this trend towards greater curvature upstream cannot continue indefinitely, since the pattern is periodic in $z$. Hence eventually a maximum is reached (at a turbulent-laminar boundary), beyond which the curvature decreases upstream. Thus, in the turbulent region, advection and diffusion act together to decrease curvature and must both be counter-balanced by turbulent forcing. These features are essentially described by the $\cos(0z)=1$ and $\cos(kz)$ modes. Figure \[fig:schematic\_sin\] illustrates the dynamics in the turbulent-laminar boundaries. These dynamics include advection by $V$ in the bulk, leading to the $U>0$ ($U<0$) patch in the lower right (upper left) of figure \[fig:udecomp\] and are described by the $\sin(kz)$ mode. [ ]{}

![Schematic depiction of the dynamics of $U$ near the centres of the laminar and of the turbulent regions. The cross-channel direction is exaggerated. Shown are $U$ (profiles) and $W+\beta y$ (arrows). Viscous diffusion tends to diminish both peaks of the profile. In laminar region surrounding $z=0$, the peaks in the upper half-channel increase in amplitude with decreasing $z$; advection towards positive $z$ (upper arrow) replenishes these peaks, maintaining $U$. Conversely, the peaks in the lower half-channel increase with $z$; advection towards negative $z$ (lower arrow) replenishes these peaks. That is, the sign of $-(W+\beta y)\,{\ensuremath{\partial}}_z U$ is opposite to that of ${\ensuremath{\partial}}^2_y U$ in both the upper and lower parts of the laminar region. In the turbulent region around $z=\pm 20$, the size of the upper (lower) peak decreases to the left (right) and so advection, like viscous diffusion, acts instead to diminish the peaks. $U$ is maintained by the Reynolds-stress force, which counterbalances both. The effect is to modulate the amplitude of the $U$ profiles periodically in $z$.[]{data-label="fig:schematic_cos"}](Figures/schematic_cos.ps){width="13cm"}

![Schematic depiction of the dynamics of $U$ near the turbulent-laminar boundaries. The cross-channel direction is exaggerated. Shown are $U+\alpha
  y$ (profiles), and $(V, W+\beta y)$ (arrows). Near the upper and lower walls, the $U+\alpha y$ profiles are advected towards increasing/decreasing $z$ by $W + \beta y$. In the bulk, advection by $V$ is significant. At $z\approx 10$, $V$ advects downwards the right-moving fluid in the upper portion of the channel. At $z\approx -10$, $V$ advects upwards the left-moving fluid in the lower portion of the channel. The effect is to tilt the $U=0$ boundary periodically in $z$.[]{data-label="fig:schematic_sin"}](Figures/schematic_sin.ps){width="13cm"}

Force balance for $W$ and $V$ {#sec:ForceWV}
-----------------------------

![Advective terms in the $W$ direction, decomposed into modes. Curves show $-\beta y {\ensuremath{\partial}}_z W$ (black, solid) $-W{\ensuremath{\partial}}_z W$ (blue, dash-dot), $-\beta V$ (green, dotted), $-V{\ensuremath{\partial}}_y W$ (red, dashed). The $\cos(kz)$ mode is dominated by $-\beta\,y\,k\,W_s$. The $\sin(kz)$ term is dominated by $\beta\,y\,k\,W_c$ near the boundaries and by $-V_s{\ensuremath{\partial}}_y (W_0+\beta y)$ in the bulk.[]{data-label="fig:FourNonw"}](Figures/Forcew.ps){width="14cm"}

[ ]{}

![Advective terms in the $W$ direction, decomposed into modes. Curves show $-\beta y {\ensuremath{\partial}}_z W$ (black, solid) $-W{\ensuremath{\partial}}_z W$ (blue, dash-dot), $-\beta V$ (green, dotted), $-V{\ensuremath{\partial}}_y W$ (red, dashed). The $\cos(kz)$ mode is dominated by $-\beta\,y\,k\,W_s$. The $\sin(kz)$ term is dominated by $\beta\,y\,k\,W_c$ near the boundaries and by $-V_s{\ensuremath{\partial}}_y (W_0+\beta y)$ in the bulk.[]{data-label="fig:FourNonw"}](Figures/fourtotw.ps){width="14cm"}

![Advective terms in the $W$ direction, decomposed into modes. Curves show $-\beta y {\ensuremath{\partial}}_z W$ (black, solid) $-W{\ensuremath{\partial}}_z W$ (blue, dash-dot), $-\beta V$ (green, dotted), $-V{\ensuremath{\partial}}_y W$ (red, dashed). The $\cos(kz)$ mode is dominated by $-\beta\,y\,k\,W_s$. The $\sin(kz)$ term is dominated by $\beta\,y\,k\,W_c$ near the boundaries and by $-V_s{\ensuremath{\partial}}_y (W_0+\beta y)$ in the bulk.[]{data-label="fig:FourNonw"}](Figures/fournonwdetail.ps){width="14cm"}

Figure \[fig:Forcew\] shows the balance of forces in the $W$ direction and figure \[fig:FourForcew\] its analogue in Fourier space. This balance resembles that in the $U$ direction shown in figures \[fig:Forceu\] and \[fig:FourForceu\]. In physical space (compare the leftmost panels of figures \[fig:Forcew\] and \[fig:Forceu\]), the main difference is that the advective and viscous forces are both small in the laminar region, in keeping with the fact that $W\approx 0$. The pressure gradient ${\ensuremath{\partial}}_z P$ is far smaller than the other forces throughout (see below). In Fourier space (compare the middle panels of figures \[fig:FourForcew\] and \[fig:FourForceu\]), the main difference with the $U$ balance is that the relative importance of the viscous and advective forces in the $\cos(kz)$ balance is reversed from that in the case of $U$: for $W$, the viscous component is larger than the advective component, which is especially small in the bulk. The decomposition of the advective terms (figure \[fig:FourNonw\]) shows that, as is the case for $U$, the advective $\cos(kz)$ term is well approximated by the contribution from advection by $w^L=\beta y$, whereas all four advective components contribute to the $\sin(kz)$ term.

The balance of forces in the $V$ direction is entirely different. The dominant balance in this equation is: $$0 = -{\ensuremath{\partial}}_y P + {F}^V
\label{eq:Vdombal}$$ as shown in figure \[fig:forcev\]. This is typical for turbulent channel flows; see, [e.g.]{} [@Pope]. This balance between the mean pressure gradient $P$ and the Reynolds-stress force ${F}^V$ does not constrain or provide information about any of the velocity components. Since $${F}^V=-{{\boldsymbol{\nabla}}\cdot}{\langle}{\tilde{{{\bf u}}}}{\tilde{v}}{\rangle}\approx -{\ensuremath{\partial}}_y{\langle}{\tilde{v}}^2{\rangle}$$ [ ]{} we in fact have $$P\approx-{\langle}{\tilde{v}}^2{\rangle}\label{eq:v2fluc}$$ up to a small $z$-dependent correction. Figure \[fig:upsikep\] shows the pressure field $P$ calculated from and suggests that its $y$ dependence can be approximated by the functional form $\cos(\pi y/2)$. This leads to an estimate of the relative importance of the pressure gradients in the $y$ and $z$ directions: $$O\left(\frac{{\ensuremath{\partial}}_y P}{{\ensuremath{\partial}}_z P}\right)= \frac{\pi/2}{k} \approx 10
\label{eq:dyPtodzP}$$ while our data shows $$\frac{({\ensuremath{\partial}}_y P)_{\rm max}}{({\ensuremath{\partial}}_z P)_{\rm max}} = \frac{0.012}{0.0017} = 7.05$$ The same estimate applies to the relative magnitudes of $V$ and $W$, using the streamfunction shown in figure \[fig:upsikep\]: $$O\left(\frac{W}{V}\right) = 
O\left(\frac{{\ensuremath{\partial}}_y \Psi}{{\ensuremath{\partial}}_z \Psi}\right) = \frac{\pi/2}{2\pi/40} = 10
\label{eq:WtoV}$$ while the actual ratio of maximum values is $$\frac{W_{\rm max}}{V_{\rm max}} = \frac{0.15}{0.013} = 11.$$

[ ]{}

![Forces in the $V$ direction. Reynolds-stress force ${F}^V$(red, dashed) is counterbalanced by pressure gradient $-{\ensuremath{\partial}}_y P$ (blue, solid). Both are zero in the laminar region. Advective and viscous forces (green, dotted) are negligible throughout.[]{data-label="fig:forcev"}](Figures/Forcev.ps){width="14cm"}

Model equations {#sec:Model}
---------------

We now derive a system of ordinary differential equations by substituting the trigonometric form into the Reynolds-averaged Navier-Stokes equations . The drawback in this procedure is the usual one, namely that this form is not preserved by multiplication. However, Table \[tab:fourier\] shows that higher harmonics contribute very little to ${{\bf U}}$.

We expand the advective term as:

\[eq:nonnon\] $$\begin{aligned}
(({{\bf U}}+{{\bf u}}^L)\cdot{{\boldsymbol{\nabla}}})&&({{\bf U}}+{{\bf u}}^L) 
= (V\:{\ensuremath{\partial}}_y + (W+\beta y)\:{\ensuremath{\partial}}_z) ({{\bf U}}+\alpha y {{\bf e_x}}+\beta y {{\bf e_z}})\nonumber\\
&=& (V_c\cos(kz)+V_s\sin(kz))
({{\bf U}}{^\prime}_0+{{\bf U}}{^\prime}_c\cos(kz) + {{\bf U}}{^\prime}_s\sin(kz) +\alpha {{\bf e_x}}+\beta {{\bf e_z}})\nonumber\\
&+&(W_0+\beta y +W_c\cos(kz)+W_s\sin(kz))
(-k\:{{\bf U}}_c\sin(kz)+k\:{{\bf U}}_s\cos(kz))\nonumber\\
&=&{\frac{1}{2}}\left(V_c\:{{\bf U}}{^\prime}_c+V_s\:{{\bf U}}{^\prime}_s+k\:(W_c\:{{\bf U}}_s-W_s\:{{\bf U}}_c)\right)
\label{eq:non0}\\
&+&(V_c ({{\bf U}}{^\prime}_0+\alpha{{\bf e_x}}+\beta{{\bf e_z}})  +k\:(W_0+\beta y){{\bf U}}_s)\cos(kz) 
\label{eq:nonc1}\\
&+&(V_s ({{\bf U}}{^\prime}_0 +\alpha{{\bf e_x}}+ \beta{{\bf e_z}}) -k\:(W_0+\beta y){{\bf U}}_c)\sin(kz)
\label{eq:nons1}\\
&+&{\frac{1}{2}}\left(V_c\:{{\bf U}}{^\prime}_c-V_s\:{{\bf U}}{^\prime}_s+k\:(W_c\:{{\bf U}}_s+W_s\:{{\bf U}}_c)\right)
\cos(2kz),\label{eq:non2}\end{aligned}$$

where primes denote $y$ differentiation. We neglect the second harmonic term , and will discuss the accuracy of this approximation below. We now rewrite the $U$ and $W$ components of the averaged momentum equations, neglecting the $z$-derivatives ${\ensuremath{\partial}}^2_z U$, ${\ensuremath{\partial}}^2_z W$ and ${\ensuremath{\partial}}_z P$, as justified by equations and :

\[eq:RANSnew\] $$\begin{aligned}
0 &=& -\left(V {\ensuremath{\partial}}_y + (W + \beta y)~ {\ensuremath{\partial}}_z \right) (U + \alpha y) 
+ {\frac{1}{{\mbox{\textit{Re}}}}}{\ensuremath{\partial}}^2_y U + {F}^U \\
0 &=& -\left(V {\ensuremath{\partial}}_y + (W + \beta y)~ {\ensuremath{\partial}}_z \right) (W + \beta y) 
+ {\frac{1}{{\mbox{\textit{Re}}}}}{\ensuremath{\partial}}^2_y W + {F}^W \end{aligned}$$

Substituting – in and separating terms in $\cos(0z)=1$, $\cos(kz)$ and $\sin(kz)$, we obtain

$$\begin{aligned}
0 &=& -{\frac{1}{2}}\left[V_c U{^\prime}_c + V_s U{^\prime}_s + k(W_c U_s  - W_s U_c) \right]
+{\frac{1}{{\mbox{\textit{Re}}}}}U{^{\prime\prime}}_0 + {F}^U_0 \label{eq:u0}\\
0 &=& -V_c (U{^\prime}_0+\alpha) - k\:(W_0+\beta y) U_s + 
{\frac{1}{{\mbox{\textit{Re}}}}}U{^{\prime\prime}}_c + {F}^U_c \label{eq:uc}\\
0 &=& -V_s (U{^\prime}_0+\alpha) + k\:(W_0 +\beta y) U_c + 
{\frac{1}{{\mbox{\textit{Re}}}}}U{^{\prime\prime}}_s + {F}^U_s \label{eq:us}\\
0 &=& -{\frac{1}{2}}\left[V_c W{^\prime}_c+ V_s W{^\prime}_s \right]
+{\frac{1}{{\mbox{\textit{Re}}}}}W{^{\prime\prime}}_0 + {F}^W_0 \label{eq:w0}\\
0 &=& -V_c (W{^\prime}_0+\beta) - k\:(W_0+\beta y) W_s +
{\frac{1}{{\mbox{\textit{Re}}}}}W{^{\prime\prime}}_c + {F}^W_c \label{eq:wc}\\
0 &=& -V_s (W{^\prime}_0+\beta)  +k\:(W_0 +\beta y) W_c +
{\frac{1}{{\mbox{\textit{Re}}}}}W{^{\prime\prime}}_s + {F}^W_s \label{eq:ws}\end{aligned}$$

\[eq:UWeqs\]

where the Fourier modes of $V$ and $W$ are related via those of the streamfunction $\Psi$ of :

\[eq:Psicomps\] $$\begin{aligned}
&V_0 = 0 & W_0 = \Psi{^\prime}_0 \\
&V_c = -k\Psi_s & W_c = \Psi{^\prime}_c \\
&V_s = k\Psi_c & W_s = \Psi{^\prime}_s\end{aligned}$$

and where homogeneous boundary conditions are imposed:

0&=& U\_0 = U\_c = U\_s y=1\
0&=& W\_0 = W\_c = W\_s y=1 \[eq:UWbcs\]

[ ]{}

System with boundary conditions is composed of six ordinary differential equations coupling the six scalar functions $U_0,U_c,U_s,\Psi_0,\Psi_c,\Psi_s$ of $y$, with six turbulent forces ${F}^U_0,{F}^U_c,{F}^U_s,{F}^W_0,{F}^W_c,{F}^W_s$.

We have solved – numerically, using as inputs ${F}^U$ and ${F}^W$ obtained from our full simulations, [i.e.]{} the ${\bf F}$ modes shown in figure \[fig:fourrey\]. The resulting solutions are shown in figure \[fig:fourdb\]. For comparison, we reproduce from figure \[fig:fourmean\] the mean velocity fields, in Fourier representation, from our full simulations (DNS). The ODE solutions are virtually indistinguishable from the mean fields from DNS. Only in the sine component of $U$ can the ODE solutions be distiguished (and only very slightly) from the DNS results. From the profiles in figure \[fig:fourdb\], the full mean fields could be constructed as in figures \[fig:udecomp\] and \[fig:psidecomp\]. Thus, while the ODE model requires input of the Reynolds-stress force terms, ${F}^U$ and ${F}^W$, it demonstrates the simplicity of the force balance responsible for generating the patterned flow when viewed in the Fourier representation. Considering higher harmonics would be straightforward, but would serve little purpose.

![Comparison between mean velocities (in Fourier representation) from full DNS and ODE models. Top row: Curves show $U$ (blue, solid) and $W$ (red, dashed) from DNS. Dots show solution to the full ODE model , essentially indistinguishable from the solid curves. Bottom row: Curves again show $U$ (blue, solid) and $W$ (red, dashed) from DNS. Dots show solution to simplified ODE model . The agreement with DNS is very good, though there are small differences particularly in the $U$ component.[]{data-label="fig:fourdb"}](Figures/modeldb.ps){width="14cm"}

We can go in the other direction and attempt to simplify system . The approximate equalities ${F}^U_c\approx -{F}^U_0$, ${F}^W_c\approx -{F}^W_0$ \[see equation \], necessary for ${\bf F}$ to vanish at the center of the laminar region, can be imposed exactly, reducing the number of turbulent forcing input functions to four. The terms arising from the advective forces can be reduced by making approximations justified from figures \[fig:FourNonu\] and \[fig:FourNonw\]. The nonlinear terms in and can be neglected. The advective terms in and can be approximated by $-k\,\beta \,y \,U_s$ and $-k\,\beta\, y \,W_s$. Making these approximations, we obtain:

$$\begin{aligned}
0 &=& {\frac{1}{{\mbox{\textit{Re}}}}}U{^{\prime\prime}}_0 + {F}^U_0 \label{eq:u0simp}\\
0 &=& - k \,\beta \,y \,U_s + {\frac{1}{{\mbox{\textit{Re}}}}}U{^{\prime\prime}}_c - {F}^U_0 \label{eq:ucsimp}\\
0 &=& - V_s (U{^\prime}_0+\alpha) + k (W_0+\beta y) U_c 
+ {\frac{1}{{\mbox{\textit{Re}}}}}U{^{\prime\prime}}_s +{F}^U_s \label{eq:ussimp}\\
0 &=&  {\frac{1}{{\mbox{\textit{Re}}}}}W{^{\prime\prime}}_0 + {F}^W_0\label{eq:w0simp}\\
0 &=& - k \,\beta \,y \,W_s + {\frac{1}{{\mbox{\textit{Re}}}}}W{^{\prime\prime}}_c - {F}^W_0\label{eq:wcsimp}\\
0 &=& -V_s (W{^\prime}_0 +\beta) + k (W_0+\beta y) W_c 
+ {\frac{1}{{\mbox{\textit{Re}}}}}W{^{\prime\prime}}_s + {F}^W_s \label{eq:wssimp}\end{aligned}$$

\[eq:UWeqssimp\]

The solutions to this simplified ODE model are also presented in figure \[fig:fourdb\]. There is quite good agreement with full DNS results, thus demonstrating that the dominant force balance is captured by this very simple system of ODEs. We stress that the only nonlinearities in this model are in equations  and . This reflects the complexity of the dynamics in the turbulent-laminar boundaries regions (and the simplicity of the dynamics in the centre of the turbulent and laminar regions.)

From the simplified ODE model we can obtain the approximate equation satisfied at the centre of the laminar region by adding equations (\[eq:u0simp\]) and (\[eq:ucsimp\]): $$k\,\beta\, y \, U_s= {\frac{1}{{\mbox{\textit{Re}}}}}(U_0+U_c){^{\prime\prime}}.
\label{eq:lam}$$ This is a restatement in terms of Fourier components of the balance described by equation and figure \[fig:forcepreface\].

Discussion
==========

We have presented an analysis of a particular turbulent-laminar pattern obtained in simulations of large-aspect-ratio plane Couette flow. We have focused on a single example so as to understand in quantitative detail the structure of these unusual flows. The key findings obtained in our study are as follows. First we find that in the (quasi-) laminar flow region the velocity profiles are not simply those of linear Couette flow. Instead a non-trivial flow is maintained in the laminar regions by a balance between viscous diffusion and nonlinear advection. Next we have considered the symmetries of the flow. When the pattern forms, the time-averaged flow breaks the translation symmetry but not centrosymmetry. The patterned state is centrosymmetric about the centre of the laminar region and about the center of the turbulent region. Next we have considered a spatial Fourier decomposition of the mean flow in the direction of the pattern wavevector. From this we find that the lateral structure of the pattern is almost completely harmonic, [i.e.]{} composed of a constant and single harmonic. Thus the pattern description can be reduced to just three cross-channel functions for each field, in that ${{\bf U}}(x,y,z) \approx {{\bf U}}_0(y) + {{\bf U}}_c(y) \cos(kz) + {{\bf U}}_s(y)\sin(kz)$. The absence of higher harmonics suggests that the pattern is near the threshold, in some sense, of a linear instability of a uniform turbulent state. Such an instability would be governed by a linear equation with coefficients which are constant in $z$, whose solutions are necessarily trigonometric in $z$.

From our analysis of the turbulent-laminar pattern, in particular its Fourier decomposition, we derive a model which reproduces the patterned flow. The model is derived from the averaged Navier-Stokes equations with the following assumptions. The crucial assumption, which is strongly supported by our numerical computations, is that the mean flow can be expressed in terms of just three horizontal modes. Effectively the generation of higher harmonics via nonlinear terms in the Navier-Stokes equations is negligible in the mean flow. The model is then further simplified because viscous diffusion is dominated by cross-channel diffusion – the standard boundary-layer approximation – and because pressure variation is negligible along the pattern wavevector. We take as input to the model the Reynolds-stress forces measured from computations. Assuming that the Reynolds stresses exactly vanish in the centre of the laminar regions, the number of inputs to the model is just four cross-channel functions. The result is a system of six simple ordinary differential equations which depend on four forcing functions. The model equations accurately reproduce the mean flow from full direct numerical simulations.

A number of other researchers have attempted to reduce the description of turbulent or transitional plane Couette flow by various means. At these low Reynolds numbers, there is no doubt that fully resolved direct numerical simulation is feasible and gives accurate results. The purpose of formulating a reduced description is therefore to yield understanding. We now comment on the differences between the approaches used by other authors and our reduction.

In parallel with their experiments, Prigent [[*et al.*]{}]{}([-@Prigent_PRL; -@Prigent_PhysD]) considered a pair of coupled Ginzburg-Landau (GL) equations with additive noise as a model for the transition from uniform turbulence to turbulent-laminar banded patterns via noisy (intermittent) patterns. These equations describe the variation in time and spanwise coordinate of the amplitudes $A^{\pm}$ of two sets of laminar bands at opposite tilt angles. These laminar bands modulate the uniform turbulence in competition with one another. Each equation separately has one reflection symmetry which corresponds physically to the centrosymmetry ${\kappa_{yz}}$ (equation ) of a banded pattern. The coupled GL equations possess a second reflection symmetry, corresponding physically to a spanwise reflection, which takes the amplitude $A^+$ to $A^-$ and vice versa. By design, this symmetry is not present in our numerical computations. Prigent [[*et al.*]{}]{} used their experimental results to fit the parameters in the GL equation and then compared simulations of the equations with experimental results. Steady patterns in the resulting GL equations have only one non-zero amplitude and this amplitude possesses the reflection symmetry corresponding to ${\kappa_{yz}}$. Hence, the steady patterns in these simulations have exactly the symmetries of the patterns we have considered.

An important class of models aims at reproducing dynamics of streamwise vortices and streaks in plane Couette turbulence by using a small number of ordinary differential equations (ODEs). These equations describe the time-evolution of amplitudes of modes with fixed spatial dependence. [@Waleffe_97], guided by the discovery of the self-sustaining process (SSP) in direct numerical simulations [@Hamilton], derived a system of eight ODEs, whose variables represent amplitudes of the key ingredients of the SSP, namely longitudinal vortices, streaks, and streak waviness. This model was later also studied and extended by [@Dauchot_Vioujard] and by [-@Moehlis_NJP].

Two other Galerkin projection procedures have been used to derive ODE models. The most energetic streamwise-independent modes in a principle orthogonal composition has been used as a basis for a 13-equation model [-@Moehlis_PF] exhibiting heteroclinic cycles; when streamwise-dependent modes are added, the resulting 31-equation model [-@Moehlis_JFM] reproduces elements of the SSP cycle. Eckhardt and co-workers [@Schmiegel; @Mersmann] have proposed a Fourier space truncation of the Navier-Stokes equations in all three spatial directions leading to a 19-equation model. They calculated turbulent lifetimes and saddle-node bifurcations giving rise to new steady states in this model.

Manneville and co-workers [@Manneville_Locher; @Lagha_Manneville] have proposed a drastic Galerkin truncation in the cross-channel direction $y$, retaining one or two trigonometric (for free-slip boundary conditions) or polynomial (for rigid boundary conditions) basis functions, but fully resolving both lateral directions. Simulating the resulting PDEs using a Fourier basis, they have been able to study phenomena such as the statistics of lifetimes of turbulent spots in domains with very large lateral dimensions.

The reduction we have presented differs from the aforementioned studies in several respects. Most importantly, we do not describe any time-dependent behaviour. We consider here neither turbulent-laminar patterns which are themselves dynamic (as in Prigent [[*et al.*]{}]{}), nor do we consider the dynamics of streaks and vortices within the turbulence, nor do we consider the transient dynamics of turbulence. Instead we focus on the spatially periodic mean flow of steady turbulent-laminar patterns. While the turbulent portions of patterns are dynamic, containing streaks and streamwise vortices, these are on a fine scale relative to spatial scales of interest here. Our model description follows directly from an analysis of full numerical simulations (not from any [*a priori*]{} assumptions, physical or phenomenological), that show that all averaged velocity components and forces, including the Reynolds stress force, are almost exactly trigonometric in the direction of the pattern wavevector. It follows directly that the steady Reynolds-averaged Navier-Stokes equations can be reduced to 6 ODEs for cross-channel profiles of the Fourier modes.

One of the more significant aspects of this work is the consideration of the force balance in just the laminar region. This balance is expressed by simple equations either in physical space, equation , or in Fourier space, equation . These equations are particularly interesting because they do not contain the Reynolds stresses, as these are negligable in the laminar region, and hence their implications can be understood without the need for closure assumptions.

As noted in §\[sec:meanflow\], equation implies that a non-zero tilt angle is necessary to maintain the S-shaped $U$ profile in the laminar region. If the patterns were not tilted, the flow would necessarily be laminar Couette flow in the centre of the laminar regions where the turbulence vanishes. We can also derive implications for the relationship between Reynolds number, tilt angle and wavelength of the patterns from equation , which we rewrite as: $$\frac{Re\,\sin\theta}{\lambda} = \frac{(U_0+U_c){^{\prime\prime}}}{2\pi \, y \, U_s}
\label{eq:firstapproxrewrite}$$ Except where $y\,U_s\approx 0$, the function on the right-hand-side is indeed approximately constant across the channel, between about 2.8 and 3.6. The value of $Re\,\sin\theta/\lambda$ used in our simulations is $350 \,\sin(24^\circ)/40=3.56$.

We may obtain a qualitative understanding of this constant as follows; see figure \[fig:Re\_beta\_lambda\] (left). Observe that in the center of the laminar region, the functional form of $U = U_0+U_c$ is like $\sin(\pi y)$. Hence its second $y$ derivative can be approximated by multiplication by $-\pi^2$, or equivalently $(U_0+U_c){^{\prime\prime}}/(-\pi^2) \approx U_0+U_c$. We also find that the odd function $-2 \,y \,U_s$ is close to $U_0+U_c$ and is in fact almost indistinguishable from $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$. This implies that the right-hand-side of equation is nearly constant across the channel and equal to $\pi$, leading to: $$\frac{{\mbox{\textit{Re}}}\, \sin\theta}{\lambda} \approx \pi
\label{eq:firstapprox2}$$

![Left: Comparison of profiles of $U_0+U_c$ (solid), $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$ (dashed) and $-2\,y\,U_s$ (points). The profile of $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$ is close to $U_0+U_c$, in accordance with the approximate functional form $(U_0+U_c)\sim \sin(\pi y)$. Note that $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$ is almost indistinguishable from $-2\,y\,U_s$, showing that $(U_0+U_c){^{\prime\prime}}/(2\pi \, y \, U_s)$ is near $\pi$ over the entire channel. Right: Plot of ${\mbox{\textit{Re}}}\,\sin\theta/\lambda $ as a function of ${\mbox{\textit{Re}}}$ for the experimentally observed patterns of Prigent [[*et al.*]{}]{} ([-@Prigent_PhysD], [-@Prigent_IUTAM]). The open triangle shows ${\mbox{\textit{Re}}}\,\sin\theta/\lambda=3.56$ for the case studied numerically in this paper.[]{data-label="fig:Re_beta_lambda"}](Figures/wavelength.ps){width="5cm"}

![Left: Comparison of profiles of $U_0+U_c$ (solid), $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$ (dashed) and $-2\,y\,U_s$ (points). The profile of $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$ is close to $U_0+U_c$, in accordance with the approximate functional form $(U_0+U_c)\sim \sin(\pi y)$. Note that $(U_0+U_c){^{\prime\prime}}/(-\pi^2)$ is almost indistinguishable from $-2\,y\,U_s$, showing that $(U_0+U_c){^{\prime\prime}}/(2\pi \, y \, U_s)$ is near $\pi$ over the entire channel. Right: Plot of ${\mbox{\textit{Re}}}\,\sin\theta/\lambda $ as a function of ${\mbox{\textit{Re}}}$ for the experimentally observed patterns of Prigent [[*et al.*]{}]{} ([-@Prigent_PhysD], [-@Prigent_IUTAM]). The open triangle shows ${\mbox{\textit{Re}}}\,\sin\theta/\lambda=3.56$ for the case studied numerically in this paper.[]{data-label="fig:Re_beta_lambda"}](Figures/Re_beta_lambda.ps){width="6.2cm"}

We believe that equation provides a good first approximation for the relationship between ${\mbox{\textit{Re}}}$, $\lambda$, and $\theta$. Figure \[fig:Re\_beta\_lambda\] (right) shows a plot of ${\mbox{\textit{Re}}}\,
\sin\theta/\lambda$ as a function of ${\mbox{\textit{Re}}}$ from the experimental data of Prigent [[*et al.*]{}]{} ([-@Prigent_PhysD]). It can be seen that this combination of quantities is approximately constant with a value near $\pi$. The range of values of the individual factors ${\mbox{\textit{Re}}}$, $\theta$, and $\lambda$ can be seen in table \[tab:exp\_obs\]. In prior studies [@Barkley_PRL; @Barkley_IUTAM], we have studied a large range of Reynolds numbers and tilt angles in a domain of length $L_z=120$. In this domain, the wavelength of a periodic pattern is less constrained, though it must be a divisor of 120. Figure \[fig:angles\] shows the observed states as a function of ${\mbox{\textit{Re}}}$ and $\theta$. Equation captures the correct order of magnitude of ${\mbox{\textit{Re}}}\, \sin\theta/\lambda$; specifically $1.8 \lesssim {\mbox{\textit{Re}}}\, \sin\theta/\lambda \lesssim 5$. Moreover, in figure \[fig:angles\] one sees that for fixed $Re$, $\lambda$ increases with increasing $\theta$, as $\eqref{eq:firstapprox2}$ predicts.

Equation does not hold in detail, however. Most notably, figure \[fig:angles\] shows that when $Re$ is decreased at fixed $\theta$, the wavelength $\lambda$ increases rather than decreases as one would expect from $\eqref{eq:firstapprox2}$. We believe that the force balance holds for all patterns which possess a laminar region free of turbulence, but that the additional approximations made in deriving the simple relationship do not hold over the full range of conditions considered in figure \[fig:angles\]. In particular, the right-hand-side of depends implicitly on $Re$, $\theta$, and $\lambda$ via the dependence of $U_0$, $U_c$, and $U_s$ on these quantities. The approximate functional relationships between $U_0$, $U_c$ and $U_s$ that we have observed in our simulations and on which we have relied in deriving may not hold for other parameter values. Finer adjustments must come from another mechanism.

![Patterns as a function of Reynolds number $Re$ and $\theta$ in a computational domain of size $L_x \times L_y \times L_z= (4/\sin\theta)\times
2 \times 120$. Turbulent-laminar patterns with wavelength $\lambda=40$ ($\times$), $\lambda=60$ ($\bullet$), $\lambda=120$ ($\ast$). Uniform turbulence ($\blacksquare$), intermittent turbulence ($\boxtimes$), laminar Couette flow ($\square$). Wavelengths in computations are constrained to be divisors of 120. Numbers are wavelengths of experimentally observed patterns of Prigent [[*et al.*]{}]{} ([-@Prigent_PhysD], [-@Prigent_IUTAM]).[]{data-label="fig:angles"}](Figures/scat_th_Re.ps){width="11cm"}

The main issue not addressed in our study is closure. We have not attempted to relate the forcing of the mean flow due to Reynolds stresses back to the mean flow itself. In the future we will report on studies employing closure models.

We thank F. Daviaud, O. Dauchot, P. Le Gal, P. Manneville and A. Prigent for helpful comments. The simulations analyzed in this work were performed on the IBM Power 4 of the IDRIS-CNRS supercomputing center as part of project 1119. This work was supported in part by a CNRS-Royal Society grant.

Turbulent-laminar bands in other shear flows
============================================

[|c|cc|cc|cc|c|]{} & & & & PP\
$Re$ & 340 & 395 & 340 & 415 & 303 & 438 & 357\
$\lambda_\text{stream}$ & 110 & 110 & 145 & 95 & 71 & 106 & 103\
$\lambda_\text{span}$ & 83 & 52 & 70 & 35 & 24 & 36 & 45\
$\lambda_z$ & 60 & 46 & 63 & 33 & 23 & 34 & 41\
$\theta$ & $\:37^\circ$ & $\:25^\circ$ & $\:26^\circ$ & $\:20^\circ$ & $\:19^\circ$ & $\:19^\circ$ & $\:24^\circ$\

Turbulent-laminar banded patterns have been observed in a number of shear flows: plane Couette (PC) flow, Taylor-Couette (TC) flow, rotor-stator (RS) flow (torsional Couette flow; the flow between differentially rotating disks) and plane Poiseuille (PP) flow (channel flow). Comparisons between these flows are impeded by the fact that different conventions are used to non-dimensionalise each of them.

In order to compare their observations in Taylor-Couette flow with those in plane Couette flow, [-@Prigent_PhysD] generalise the Reynolds number used in plane Couette flow $U=y/h$ by considering it as based on the shear and the half-gap: $$Re^{PC} = \frac{(\text{Shear}^\text{PC})\;\text{(half-gap)}^2}{\nu} = 
\frac{\frac{U}{h} h^2}{\nu} = \frac{U h}{\nu}$$ For flows whose shear is not constant, the average shear is used. We also convert streamwise and spanwise wavelengths to total wavelength and angle of the pattern wavevector via $$\tan(\theta) = \frac{\lambda_\text{span}}{\lambda_\text{stream}} \qquad
\lambda_z = \lambda_\text{span} \cos(\theta)
\label{eq:convert}$$ Table \[tab:exp\_obs\] presents the Reynolds numbers, wavelengths, and angles for which turbulent-laminar patterns have been observed experimentally or numerically. The subsections which follow explain how Table \[tab:exp\_obs\] was obtained from the data in [@Prigent_PhysD; @Cros; @Tsukahara].

Taylor-Couette flow
-------------------

For Taylor-Couette flow between differentially rotating cylinders, the azimuthal and axial directions correspond to the streamwise and spanwise directions of plane Couette flow. For cylinders of radius $r_i$ and $r_o$, rotating at angular velocities $\omega_i$ and $\omega_o$ with $2h\equiv r_o-r_i$ and $\eta\equiv r_i/r_o$, the shear averaged over the gap is $$\langle \text{Shear}^\text{TC}\rangle = 
\frac{r_i\omega_i - \eta r_o\omega_o}{(1+\eta)h}$$ leading to the Reynolds number: $$Re^\text{TC} \equiv \frac{r_i \omega_i - \eta r_o \omega_o}{(1+\eta)h} \;
\frac{h^2}{\nu} \approx \frac{Re_i-Re_o}{4\nu}$$ where the last approximate equality corresponds to exact counter-rotation ($\omega_o=-\omega_i$) and the narrow gap limit ($\eta \rightarrow 1$), and $R_i$, $R_o$ are the conventionally defined inner and outer Reynolds numbers, [e.g.]{} $Re_i\equiv 2h r_i\omega_i/\nu$. The wavelengths and Reynolds numbers observed in Taylor-Couette and plane Couette flow are compared in Figure 5 of [@Prigent_PhysD].

Torsional Couette flow
----------------------

The laminar profile for torsional Couette flow between a rotating and a stationary disk (rotor-stator flow) is $${{\bf u}}={{\bf e}}_\theta \;\frac{\omega \,r \,z}{h}$$ and the Reynolds number based on axial shear and half-gap is $$Re^\text{RS}=\frac{\omega r}{h} \;\frac{h^2}{4\nu} = \frac{\omega r h}{4\nu}
\label{eq:Re_RS}$$ For $m$ spirals, the azimuthal wavelength in units of the half-gap is $$\lambda^\text{RS}_\text{stream} = \frac{2\pi r}{m h/2}=\frac{4\pi r}{m h}
\label{eq:lambda_RS}$$ Turbulent spiral patterns which are rather regular occur for a range of angular velocities and radii. In their figures 12, 16 and 18, [@Cros] focus particularly on the radius and gap: $$r = 0.8\times 140\,\text{mm}=11.2\,\text{cm} \qquad
h=0.22\,\text{cm}
\label{eq:rh_RS}$$ The highest and lowest rotation rates for which turbulent spirals are seen are

= 68 = 7.12 && m=6\
\[eq:omegahi\_RS\] = 47 = 4.92 && m=9 \[eq:omegalo\_RS\]

Substituting - and the viscosity $\nu=10^{-2}\,$cm$^2$/sec of water into - leads to the values shown in Table \[tab:exp\_obs\]. The pitch angle of the spirals remains approximately constant at $19^\circ$. We use to calculate $\lambda_\text{span}$ and $\lambda_z$, neglecting the variation in radius.

Plane Poiseuille flow
---------------------

Figure 14 of [-@Tsukahara] shows a visualisation from a direct numerical simulation of plane Poiseuille (PP) flow in a channel with domain and Reynolds number $$L_\text{stream}\times L_y\times L_\text{span} = 
51.2\;\delta \times 2\;\delta \times 22.5\;\delta
\qquad
Re_c \equiv \frac{u_c\delta}{\nu} = 1430$$ where $u_c$ is the centerline velocity. The domain contains a single wavelength of an oblique turbulent-laminar banded pattern oriented at $\theta=24^\circ$ to the streamwise direction. (Both the wavelength and the angle are dictated by the computational domain.) Following [@Waleffe_03], we view the Poiseuille profile in the half-channel $[-\delta,0]$, over which the shear has one sign, as comparable to the Couette profile in the channel $[-h,h]$, and thus take $\delta/2$ as the unit of length, rather than $\delta$. The shear is obtained by averaging over $[-\delta,0]$: $$\langle\text{Shear}^\text{PP}\rangle 
= \left\langle \frac{du}{dy} \right\rangle 
= \frac{u_c}{\delta}$$ For the Reynolds number based on the average shear and half-gap, we obtain $$Re^\text{PP}=\frac{\langle\text{Shear}^\text{PP}\rangle 
\;(\text{half-gap})^2}{\nu} = 
\frac{\frac{u_c}{\delta}\frac{\delta^2}{4}}{\nu} = 
\frac{u_c\delta}{4\nu} = \frac{Re_c}{4\nu} = \frac{1430}{4}=357.5$$
