2.1. Large-Reynolds-number asymptotic analysis

Vortex-wave interaction theory describes the high-Re asymptotic nature of typical coherent structures in shear flows. The leading-order part of the theory naturally shows that the nonlinear interaction between two constituents of flows; the streamwise average flow (i.e. roll–streak field) and a streamwise-oscillatory wave which in streamwise Fourier space is represented by a complex-conjugate pair of monochromatic modes. In Hall & Sherwin (Reference Hall and Sherwin2010) plane Couette flow bounded in the normal (vertical, $y$) direction by counter-translating walls was considered in order to compare the asymptotic VWI approach with the finite-Reynolds-number results of Wang et al. (Reference Wang, Gibson and Waleffe2007). The VWI array theory of H18 removes the imposition of walls in the normal direction and replaces this with additional periodicity of both the roll and the streak, the latter becoming a periodic variation on a background linear shear flow. The wave is no longer required to be periodic in $y$ but instead is replaced with a periodic array of waves each allowed to be of large vertical extent, a relaxation allowed by the asymptotic structure of this flow constituent, whose interaction with the streamwise average flow is only significant in the critical layer of each wave.

Below we write down the VWI expansions and the interaction equations describing a periodic array of VWI states of unit wavelength in the $y$ direction. The basic scalings for this type of roll–streak–wave interaction remain those introduced by Hall & Smith (Reference Hall and Smith1991) and first used in the high-Reynolds-number asymptotic description of numerically exact Navier–Stokes solutions by Hall & Sherwin (Reference Hall and Sherwin2010). For later convenience, we define $(X,Z)=({\alpha ^{{\dagger} }} x, {\beta ^{{\dagger} }} z)$. In these new coordinates, the streamwise and spanwise periodic box dimensions are normalised to be $2{\rm \pi}$. In most parts of the flow, the velocity and pressure expansions of the periodic VWI array found in H18 are

(2.5a)\begin{gather} u=y+U(y,Z)+ \rho\frac{\mathrm{e}^{\mathrm{i} X}}{Re^{7/6}} \sum_{n={-}\infty}^{\infty} \mathrm{e}^{-\mathrm{i}{\alpha^{{{\dagger}}}} nt} \tilde{u}(y-n,Z)+\mathrm{c.c.}+\cdots, \end{gather}

(2.5b)\begin{gather}v =\frac{V(y,Z)}{Re}+ \rho\frac{\mathrm{e}^{\mathrm{i} X}}{Re^{7/6}} \sum_{n={-}\infty}^{\infty} \mathrm{e}^{-\mathrm{i}{\alpha^{{{\dagger}}}} n t} \tilde{v}(y-n,Z) +\mathrm{c.c.}+\cdots, \end{gather}

(2.5c)\begin{gather}w =\frac{W(y,Z)}{Re}+ \rho\frac{\mathrm{e}^{\mathrm{i} X}}{Re^{7/6}} \sum_{n={-}\infty}^{\infty} \mathrm{e}^{-\mathrm{i}{\alpha^{{{\dagger}}}} n t} \tilde{w}(y-n,Z)+\mathrm{c.c.}+\cdots, \end{gather}

(2.5d)\begin{gather}p =\frac{P(y,Z)}{Re^{2}}+ \rho\frac{\mathrm{e}^{\mathrm{i} X}}{Re^{7/6}} \sum_{n={-}\infty}^{\infty} \mathrm{e}^{-\mathrm{i}{\alpha^{{{\dagger}}}} n t}\tilde{p}(y-n,Z)+\mathrm{c.c.}+\cdots, \end{gather}

where $\mathrm {c.c.}$ stands for complex conjugate. The two-dimensional complex functions $\tilde {u},\tilde {v},\tilde {w},\tilde {p}$ are referred to as the wave part, while the roll–streak part $U,V,W,P$ are two-dimensional real functions. The ellipses $\cdots$ represent higher-order terms including higher harmonics of the wave. The real constant $\rho$ is the amplitude of the wave. The streak $U$ makes an $O(1)$ modulation to the mean flow so that it generates three-dimensional waves via the inviscid Rayleigh instability mechanism. The $1/Re$ prefactor of vertical–spanwise roll velocity $(V,W)$ is chosen to ensure the viscous-convective balance in the streamwise-averaged flow, noting that only the roll velocity contributes in the convective effect. The boundary conditions to be imposed on the roll–streak flow are periodicity in $y,Z$ so that

(2.6a)\begin{align} \boldsymbol{U}(y+1,Z) & = \boldsymbol{U}(y,Z+2{\rm \pi}) =\boldsymbol{U}(y,Z), \end{align}

(2.6b)\begin{align} P(y+1,Z) & = P(y,Z+2{\rm \pi}) =P(y,Z). \end{align}

The three-dimensional/three-component wave velocity field in (2.5) is, in the majority of the domain, smaller than that of the roll by a factor $Re^{1/6}$ and, thus, at first glance does not cause any feedback effect to the leading-order roll–streak. However, the wave becomes stronger within the thin critical layer of thickness $O(Re^{-1/3})$, where the feedback effect occurs (the amplification can be seen in figures 2 and 3). The scaling of the wave is chosen in order that the Reynolds stresses in the critical layer produced by the wave drive the roll flow; see Hall & Sherwin (Reference Hall and Sherwin2010) for a more detailed discussion of the scaling. The summations in (2.5) represent contributions from the infinite vertical array of waves. Individual waves are periodic in the streamwise and spanwise directions, and repetitive, but not assumed periodic in the vertical direction. As described in H18, the magnitude of each wave component takes on its maximum values within its associated critical layer and then decays to zero exponentially in neighbouring cells. For example, the wave pressure satisfies

(2.6c)\begin{gather} \tilde{p}(y,Z) =\tilde{p}(y,Z+2{\rm \pi}), \end{gather}

(2.6d)\begin{gather}\tilde{p}(y,Z) \to 0,~y \to \pm \infty. \end{gather}

Although an individual wave is not periodic in $y$, the entire superimposed array of the roll–streak–waves (2.5) does satisfy the spectral form of the shear-periodic conditions (2.4).

Figure 2. Solution for the mirror-symmetric mode at ${\alpha ^{{\dagger} }}=2$, ${\beta ^{{\dagger} }}=2{\rm \pi}$, $R_f=20\times 10^{3}$. Horizontal axis is $z$, vertical axis is $y$. The domain shown is $y\in [-5/2,5/2]$, $z\in [0,2{\rm \pi} /{\beta ^{{\dagger} }}]$. (a) Contours of roll stream function $\psi$ in the $(y,z)$ plane shown for three fundamental vertical modules. A slightly lightened colour is used except for the central module $y\in [-1/2,1/2]$. (b) Contours of streak flow $y+U$ on the domain used to solve for the wave $\tilde {\boldsymbol {u}}$. Critical layers for waves with speeds $0,\pm 1$ are shown as black lines. (c) Contours of kinetic energy $|\tilde {\boldsymbol {u}}|^{2}$ for the central stationary wave. (d) Three-dimensional structure of $\tilde {\boldsymbol {u}}$ for the central stationary wave shown as red/yellow isosurfaces of positive/negative streamwise vorticity. The critical layers for the waves with speeds $0,\pm 1$ are shown as translucent planes.

Figure 3. A lower-branch solution for the sinuous mode at ${\alpha ^{{\dagger} }}=2$, ${\beta ^{{\dagger} }}={\rm \pi}$, $R_f=20\times 10^{3}$. Panels as for figure 2.

The factors $\mathrm {e}^{-\mathrm {i}{\alpha ^{{\dagger} }} nt}$ in the summation account for the fact that the speeds of waves centred in vertically neighbouring cells differ by unity to match the change in mean flow speed between adjacent cells. We note that in the expansions above one could replace $nt$ in the exponential terms by $(n+c)t$ in order to describe the situation where each wave propagates at a phase speed $c$ different from the unperturbed flow speed at the centre of each cell. Here we describe the case where all the waves have an identical structure so that we could also then shift the origin in $y$ to reduce $c$ to zero, thereby recovering (2.5). However, a more detailed consideration of the phase speed would be essential in cases where waves of distinct properties participate in the array. For example, in § 3.2 we treat the flow maintained by two waves that move in opposite directions. In that case the $c$-dependence would need to be retained and we would split the summations in (2.5) into two sums representing the two sets of waves.

Substituting the expansions (2.5) into (2.1), we can find a leading-order system valid outside of the critical layer. As anticipated, the streamwise-invariant roll–streak part satisfies advection-diffusion equations

(2.7a)\begin{gather} V+VU_y+{\beta^{{{\dagger}}}} WU_{Z} =U_{yy}+{\beta^{{{\dagger}}}}^{2} U_{ZZ}, \end{gather}

(2.7b)\begin{gather}VV_{y}+{\beta^{{{\dagger}}}} WV_{Z} ={-}P_{y}+V_{yy}+{\beta^{{{\dagger}}}}^{2}V_{ZZ}, \end{gather}

(2.7c)\begin{gather}VW_{y}+{\beta^{{{\dagger}}}} WW_{Z} ={-}{\beta^{{{\dagger}}}} P_Z +W_{yy}+{\beta^{{{\dagger}}}}^{2} W_{ZZ}, \end{gather}

(2.7d)\begin{gather}V_{y}+{\beta^{{{\dagger}}}} W_{Z} =0. \end{gather}

The wave part satisfies the inviscid instability problem for the base flow $y+U(y,Z)$, which after some manipulation becomes the single equation for the wave pressure function $\tilde {p}$,

(2.7e)\begin{equation} \left (\frac{\tilde{p}_{y}}{\{y+U(y,Z) \}^{2}} \right )_{y}+\left (\frac{\beta^{{{\dagger}} 2}\tilde{p}_Z}{\{y+U(y,Z) \}^{2}} \right )_Z-\frac{\alpha^{{{\dagger}} 2}\tilde{p}}{\{y+U(y,Z)\}^{2}}=0. \end{equation}

Given the base flow $y+U(y,Z)$, the wave equations constitute a linear eigenvalue problem. The wave pressure function $\tilde {p}(y,Z)$ is normalised such that the integral of $|\tilde {p}|^{2}$ over the streamwise-normal area of the periodic box is unity. In order to find the size of the wave measured by the amplitude $\rho$ we consider the feedback effect of the wave to the roll–streak.

The periodic roll–streak flow is supported by a periodic array of waves each of which interacts only significantly in the critical layer. The amplification mechanism of the wave there is due to the singularity of (2.7e) occurring when $y+U(y,Z)$ vanishes. The classical critical layer argument can be used to show that within a thin layer of thickness $O(Re^{-1/3})$ surrounding the singularity, one must alter the form of the asymptotic expansions. The viscous effect in the wave equations is no longer negligible there, and singularity observed in the outer scale is regularised in the inner scale. Remarkably, an analytical solution exists for the critical layer problem, and we can write the feedback effect down explicitly in terms of the jump of the roll–streak quantity across the critical layer. The jump conditions are most easily expressed in terms of $V^{s}$, the component of the roll flow locally parallel to the critical layer, using the new variables $s,N$ which measures distance along and normal to the layer in the $(y,z)$ coordinate. From Hall & Smith (Reference Hall and Smith1991) and Hall & Sherwin (Reference Hall and Sherwin2010), the jumps are found to be

(2.8a,b)\begin{equation} \left[ \frac{\partial {V^{s}}}{\partial N}\right]_-^{+} =J'(s), \quad [P]_-^{+}={-}\chi(s) J(s), \end{equation}

where $\chi (s)$ is the curvature of the critical layer and

(2.9a,b)\begin{equation} J(s)= \left.\frac{n_0\rho^{2}|\tilde{p}_s|^{2}}{\{{\alpha^{{{\dagger}}}}(y+U)_N\}^{5/3}}\right |_{N=0}, \quad n_{0}=2{\rm \pi}\left(\frac{2}{3}\right)^{\frac{2}{3}}\left (-\frac{2}{3}\right )!. \end{equation}

If the critical layer is flat the above jumps reduce to

(2.10a–c)\begin{equation} \left[ \frac{\partial W}{\partial N}\right] _{-}^{+} ={\beta^{{{\dagger}}}}\hat{J}'(Z), \quad [P] _{-}^{+}=0,\quad \hat{J}(Z)=\left. \frac{n_0\rho^{2}\beta^{{{\dagger}} 2}|\tilde{p}_Z|^{2}} {\{{\alpha^{{{\dagger}}}}(1+U_y)\}^{5/3}}\right |_{N=0}, \end{equation}

as used in H18. Once again we note that if there is more than a single critical layer in each cell, one must compute the jumps across each layer induced by the relevant wave.

The equations (2.7), the boundary conditions (2.6), the jump conditions (2.8a,b) and the wave pressure normalisation constitute a closed nonlinear eigenvalue problem. The solution of that canonical problem for given ${\alpha ^{{\dagger} }},{\beta ^{{\dagger} }}$ determines the leading-order structure $(\boldsymbol {U},\tilde {p})$ in the leading-order part of the expansion (2.5), together with the amplitude $\rho ({\alpha ^{{\dagger} }},{\beta ^{{\dagger} }})$.

In the high-Reynolds-number limit any neutral wave travelling in the flow will reinforce locally a flow in the $y$–$z$ plane predominantly by producing a stress jump associated with its own critical layer. The overall structure of the single wave interaction problem is unchanged from the usual VWI for parallel flows. This is advantageous for computations because the numerical solution of the canonical problem can be found using existing VWI codes for travelling wave solutions in plane Couette flow with some minor modifications (see § 2.2, following).

As remarked earlier, the boundary conditions used for the canonical problem are consistent with the shear-periodic conditions. This means that some time-periodic solutions of the finite-Reynolds-number shearing box problem as found, for example, by Sekimoto & Jiménez (Reference Sekimoto and Jiménez2017) might asymptote to the VWI structure described by the canonical problem even though the asymptotic solution is of travelling wave form. On first consideration, this might seem a rather paradoxical outcome. However, we recall that the travelling wave form of the leading-order structure is possible because in the asymptotic framework the fundamental mode is induced by the linear stability problem of the streak field. Thus, the interaction of fundamental waves centred in different cells induces higher streamwise harmonics of a non-travelling wave form in (2.5), and formally they are smaller than the fundamental mode. This means that an intrinsic property of periodic box solutions – that no pure travelling wave states are possible – is actually a finite-Reynolds-number effect.

A referee of this work has kindly pointed out similarities of H18 and the present paper to Montemuro et al. (Reference Montemuro, White, Klewicki and Chini2020). The latter was a development of an earlier paper, Chini et al. (Reference Chini, Montemuro, White and Klewicki2017), with significant modifications to the streamwise wavelength scales. The essential structure proposed in those papers is closely related to H18. However, the roll activity is taken to be sufficiently large that the streamwise velocity becomes uniform and the critical layer driving the interaction becomes embedded in an outer layer where the streamwise velocity adjusts between its different constant values in adjacent cells. In previous VWI and SSP problems it has been found that for a given periodic box size all flow quantities are determined by solving a nonlinear eigenvalue problem, while Montemuro et al. (Reference Montemuro, White, Klewicki and Chini2020) found that for a given periodic box size solutions exist for a continuum of values of the wave amplitude. This apparently surprising result perhaps arises because the streak profile in the new layer surrounding the critical layer becomes independent of the wave amplitude and roll driven by the wave. One might anticipate that similar structures could be found numerically from the hybrid formulation below but exhaustive searches did not find any such solutions.

2.2. Hybrid approach

As remarked earlier, the VWI formulation of H18 reduces the Navier–Stokes equations to a coupled pair of two- rather than three-dimensional momentum equations and takes the Reynolds number out of the problem. Nevertheless, the computational treatment of the asymptotically reduced problem (2.6)–(2.8a,b) is challenging because of the stress jump conditions to be applied across the wave's critical layer, whose position is a priori unknown. Here we adopt the ‘hybrid’ methodology first advanced in Blackburn et al. (Reference Blackburn, Hall and Sherwin2013) and used subsequently by Deguchi & Hall (Reference Deguchi and Hall2016). The key aspect of the methodology is to regularise the critical layer singularity by retaining some viscous effects in the reduced system, allowing us to use with only minor variations conventional numerical methods developed for solution of Navier–Stokes problems. The strength of the regularisation is measured by the ‘fictitious Reynolds number’ $R_f$, here based on the depth of the fundamental vertical module and the background shear, and as in Blackburn et al. (Reference Blackburn, Hall and Sherwin2013), we check that outcomes become independent of this parameter as it increases.

Here we remark that the reduced Navier–Stokes computations employed later by Thomas et al. (Reference Thomas, Lieu, Jovanović, Farrell, Ioannou and Gayme2014) and Beaume et al. (Reference Beaume, Chini, Julien and Knobloch2015, Reference Beaume, Knobloch, Chini and Julien2016) are essentially equivalent to the hybrid approach. All these methods rely on solving the finite-Reynolds-number Navier–Stokes equations reduced to just two streamwise Fourier modes: the streamwise average (the ‘vortex’, or roll–streak) and a single streamwise wave. The hybrid equations can be obtained by further assuming that there is no wave Reynolds stress driving the streak ($U$), and that the roll components ($V,W$) can be neglected in the wave equations. Both assumptions follow from the fact that the roll velocity is much smaller than that of the streak (see (2.5)). All reductions make explicit recognition of the crucial tripartite interaction between rolls, streaks and waves in the VWI theory.

Starting the hybrid reduction process from the Navier–Stokes equations (2.1) with the leading-order part of the VWI expansions shown in (2.5), and setting $Re=R_f$, we get the roll–streak equations from the $x$-averaged part

(2.11)\begin{align} &(V\partial_y +{\beta^{{{\dagger}}}} W\partial_Z) \begin{bmatrix} y+U\\ V\\ W \end{bmatrix} \nonumber\\ &\quad ={-} \begin{bmatrix} 0\\ \partial_y P\\ {\beta^{{{\dagger}}}}\partial_Z P \end{bmatrix}+ (\partial_{yy} + {\beta^{{{\dagger}}}}^{2}\partial_{ZZ}) \begin{bmatrix} U\\ V\\ W \end{bmatrix} -R_f^{{-}1/3}\rho^{2}\sum_{n={-}N}^{N} \begin{bmatrix} 0\\ F_1^{(n)}\\ F_2^{(n)} \end{bmatrix}, \end{align}

which, note, are solved at a unit effective Reynolds number owing to the adoption of the asymptotic scaling (Hall & Smith Reference Hall and Smith1991; Hall & Sherwin Reference Hall and Sherwin2010), while the $x$-fluctuation part yields the wave equations

(2.12)\begin{align} \mathrm{i}{\alpha^{{{\dagger}}}}(y+U-c) \begin{bmatrix} \tilde{u}\\ \tilde{v}\\ \tilde{w} \end{bmatrix} + \begin{bmatrix} \tilde{v}\partial_y (y+U)+{\beta^{{{\dagger}}}}\tilde{w}\partial_Z U\\ 0\\ 0 \end{bmatrix} ={-} \begin{bmatrix} \mathrm{i}{\alpha^{{{\dagger}}}} \tilde{p}\\ \partial_y \tilde{p}\\ {\beta^{{{\dagger}}}}\partial_Z \tilde{p} \end{bmatrix}+R_f^{{-}1}\nabla^{2} \begin{bmatrix} \tilde{u}\\ \tilde{v}\\ \tilde{w} \end{bmatrix}. \end{align}

Owing to the presence of the viscous term, (2.12) are not singular. We have here allowed for a non-zero wave speed $c$, allowing us to obtain the type of results to be described in § 3.2 below, though, as discussed in § 2.1, we could equally well have used a slight elaboration to the form of (2.5) (where $c=0$). The final terms in (2.11) represent the Reynolds stresses from the waves

\begin{align*}F_1^{(n)}&=\partial_y\{\tilde{v}(y-n,Z)\tilde{v}^{{\ast}}(y-n,Z)\}+ {\beta^{{{\dagger}}}}\partial_Z\{\tilde{v}(y-n,Z)\tilde{w}^{{\ast}}(y-n,Z)\}+\mathrm{c.c.},\\ F_2^{(n)}&=\partial_y\{\tilde{v}(y-n,Z)\tilde{w}^{{\ast}}(y-n,Z)\}+ {\beta^{{{\dagger}}}}\partial_Z\{\tilde{w}(y-n,Z)\tilde{w}^{{\ast}}(y-n,Z)\}+\mathrm{c.c.}, \end{align*}

where $*$ denotes complex conjugation and $N$ should formally be taken to be large in order to recover an infinite array system. These stress terms are inserted in body force form because, owing to the viscous term in the wave equations, the thickness of the critical layer becomes finite. Note that these body force terms were absent in the formal asymptotic reduction (2.7) since the stress terms at the asymptotic limit are concentrated within the infinitely thin critical layer and, hence, become discontinuous stress jumps across the layer as given in (2.8a,b).

Now consider the fundamental vertical cell $y\in [-1/2,1/2]$ for the numerical computations. We enforce periodic conditions of period unity for the roll–streak. In general, the Reynolds stress term from waves in cells above and below the module must be accounted for. However, for many of the numerical solutions which follow, the critical layer is trapped within the fundamental module. In that case we can set $N=0$ such that each roll is influenced only by the wave inside its local module, since the stress is significant only around the critical layer (it is possible to relax this restriction if required). Each wave decays exponentially with increasing distance from its corresponding module. For the wave component, it is typically sufficient to consider a domain of depth three to five units with no-slip or no-stress conditions imposed at its upper and lower edges.

The preceding discussion assumes that all waves have an identical structure. Next we consider the case when two distinct waves coexist in one period of the roll–streak flow. In order to compute this second mode, we exploit the property that the two waves possess certain symmetries, allowing us to compute over only half of the periodic structure, and, hence, it is convenient to choose the depth of the fundamental module to be half of the roll–streak period, such that it encloses only one wave critical layer.

We solve the hybrid system numerically at finite but large enough $R_f$, checking that the solution is in the asymptotic range where the wave amplitude $\rho$ required to produce an equilibrium state becomes independent of $R_f$. We have employed two alternative, and equivalent, means to finding solutions to the system ((2.11), (2.12)). The first method couples a large-time $\tau$-integration of the unsteady version of (2.11) (i.e. $\partial _{\tau } [U,V,W]^{\text {T}}$ is added in the left-hand side) with the eigenvalue analysis of (2.12), using spectral elements to discretise both equations. The former time integration problem can be solved if the wave forcing is given, and the latter eigenvalue problem returns the complex phase speed $c$ for a given streak field. The wave amplitude $\rho$ is tuned to neutralise the growth rate of the wave in a relaxation iteration based on the successive updates of the wave and roll–streak fields until reaching a fixed point; this was the approach adopted in Blackburn et al. (Reference Blackburn, Hall and Sherwin2013). The convergence of a time-stepper-based eigensolution method used in the scheme depends on the stability of the solution. For lower-branch solutions, the neutrally stable Rayleigh mode is the least stable one and there is no difficulty in obtaining convergence. However, for upper-branch solutions, unstable Rayleigh modes typically exist in addition to a neutrally stable one (Gibson, Halcrow & Cvitanović Reference Gibson, Halcrow and Cvitanović2009; Deguchi & Hall Reference Deguchi and Hall2016), and, hence, a shift-invert solution method outlined in Gómez et al. (Reference Gómez, Pérez, Blackburn and Theofilis2015) is used to isolate neutrally stable modes. Using this slight modification to our original methodology we were here able to find both upper- and lower-branch solutions with the fixed-point technique, something which did not prove possible in the earlier work, which contained only lower-branch solutions. The second method uses Chebyshev–Fourier/Fourier–Fourier global expansions for discretising the waves/roll–streaks and directly solves the resultant algebraic equations by the Newton–Raphson method. The Navier–Stokes solver developed in Deguchi, Hall & Walton (Reference Deguchi, Hall and Walton2013) can be modified to find solutions of ((2.11), (2.12)); this was the approach adopted in, for example, Deguchi & Hall (Reference Deguchi and Hall2016). The fixed-point iteration code is useful in finding solutions from fairly arbitrary initial conditions, while the Newton–Raphson code provides an efficient continuation method to trace the solution branch without regarding its stability. We confirmed that the two methods produce the same outcomes at matching values of ${\alpha ^{{\dagger} }},{\beta ^{{\dagger} }},R_f$.

3.1. Mirror-symmetric and sinuous modes

The numerical results presented here were motivated by the previous findings of two typical steady solutions in plane Couette flow. The first is the mirror-symmetric mode found by Itano & Generalis (Reference Itano and Generalis2009) and Gibson et al. (Reference Gibson, Halcrow and Cvitanović2009) where in the latter, the labels EQ7 and EQ8 were used. This mode has some similarity to the VWI array solution found in H18. The second is the sinuous-mode solution found by Nagata (Reference Nagata1990), Clever & Busse (Reference Clever and Busse1992) and Waleffe (Reference Waleffe1998) and labelled EQ1 and EQ2 in Gibson et al. (Reference Gibson, Halcrow and Cvitanović2009). This mode has been the most intensively studied among other numerous solutions, because it is the first to appear when the Reynolds number is increased, and the lower branch created by that saddle-node bifurcation defines the edge of the laminar and turbulent attractors (Wang et al. Reference Wang, Gibson and Waleffe2007). The large-Reynolds-number asymptotic development of the branch has been explained by applying the VWI theory to that branch (Hall & Sherwin Reference Hall and Sherwin2010; Blackburn et al. Reference Blackburn, Hall and Sherwin2013; Deguchi & Hall Reference Deguchi and Hall2014b, Reference Deguchi and Hall2016). The sinuous mode possesses so-called shift-reflection and shift-rotation symmetries

(3.1a)\begin{gather}{} [u,v,w](x,y,z)=[u,v,-w](x+{\rm \pi}/{\alpha^{{{\dagger}}}},y,-z), \end{gather}

(3.1b)\begin{gather}{} [u,v,w](x,y,z)=[{-}u,-v,w]({-}x,-y,z+{\rm \pi}/{\beta^{{{\dagger}}}}), \end{gather}

while the mirror-symmetric mode satisfies those two symmetries and

(3.1c)\begin{equation} [u,v,w](x,y,z)=[{-}u,-v,w]({-}x+{\rm \pi}/{\alpha^{{{\dagger}}}},-y,z).\end{equation}

Indeed, Deguchi & Hall (Reference Deguchi and Hall2014a) found a bifurcation point at which the sinuous-mode branch bifurcates from the mirror-symmetric mode branch by a symmetry breaking of (3.1c).

Here we consider how those symmetries could be carried over to our triply periodic problem. It is convenient to begin by first demonstrating that the mirror symmetry (3.1) imposed for the roll–streak field can reproduce features of the VWI array solutions presented in H18. The roll–streak field version of (3.1) can be found as

(3.2a)\begin{gather}{} [U,V,W](y,z)=[U,V,-W](y,-z), \end{gather}

(3.2b)\begin{gather}{} [U,V,W](y,z)=[{-}U,-V,W]({-}y,z+{\rm \pi}/{\beta^{{{\dagger}}}}), \end{gather}

(3.2c)\begin{gather}{} [U,V,W](y,z)=[{-}U,-V,W]({-}y,z). \end{gather}

From the third condition we see that $U=0$ at $y=0$. Such a flat critical layer is special in the sense that the wave pressure admits a regular expansion there and, hence, presents less numerical difficulty (see Deguchi & Hall Reference Deguchi and Hall2016, H18 and Deguchi Reference Deguchi2019). In fact, the numerical technique used by H18 can be applied only for this case. As we shall see shortly, for the sinuous mode, the critical layer is curved, owing to the absence of the symmetry (3.2c), and, hence, in order to resolve the singular critical layer structure the hybrid method must be employed.

Imposing (3.2) and periodicity in $y$ with unit wavelength for the roll–streak field, we are able to generate the mirror-symmetric mode solution shown in figure 2. In figure 2(a), for presentation of roll velocity field $(V,W)$, we use contours of stream function $\psi$, where ${\beta ^{{\dagger} }}\partial _Z\psi =-V$ and $\partial _y\psi =W$. Three modules of depth unity are shown; except for the centre module $y\in [-1/2,1/2]$, the colour map is slightly lightened (and owing to the colour map chosen, variations in streak speed around the central critical layer may be difficult to perceive, but are similar in nature to the more visible variations near the other two critical layers shown). The rolls redistribute the streamwise velocity to create the streaky pattern shown in (b), which in turn generates the wave field shown by figure 2(c,d) via the wave equations (2.12) that are solved on domains of large but finite vertical extent. We note that within the domain shown in figure 2(c), three critical layers for waves of speeds $0,{\pm }1$ may coexist, as shown by the black lines. Our focus in the computation is on resolving the central wave that is both neutral and stationary; the nature of other neutrally stable waves is identical to the centre wave and, hence, are not shown in figure 2(c,d). The generated centre wave $\tilde {\boldsymbol {u}}$ indeed satisfies the symmetries (3.1). In order to solve (2.12) for the centre wave, here a domain three times greater than the roll–streak module (i.e. $y\in [-3/2,3/2]$) was used. As may be seen in figure 2(c,d) each wave localises around its critical layer, and falls rapidly in strength with distance from the layer, consistent with the asymptotic theory. Only one critical layer is present in each roll–streak module and, thus, we can set $N=0$ in (2.11) to calculate the feedback effect from the waves to the rolls.

Now we turn our attention to the sinuous mode. From the plane Couette flow results we anticipate that a pair of counter-rotating rolls form within a periodic box. If such rolls are generated in the fundamental module, it would be more natural to combine a half-shift in $z$ with the periodic boundary conditions in $y$, in order to tile the infinite domain. Thus, we impose

(3.3)\begin{equation} [U,V,W](y,z)=[U,V,W](y+1,z+{\rm \pi}/{\beta^{{{\dagger}}}}),\end{equation}

rather than the usual periodicity. A solution generated by imposing symmetries (3.2a), (3.2b), (3.3) is shown in figure 3. From figure 3(a) it is evident that the shift-periodic symmetry (3.3) generates a feature characteristic of sinuous-mode arrays: rolls in vertically and horizontally adjacent cells counter-rotate, chequerboard fashion. In other words, the period of the roll–streak field is twice unity in $y$ rather than unity; thus, except for this minor rescaling for numerical convenience, the structure here is fully covered by the asymptotic theory of § 2. Combining (3.2b) and (3.3), we see that the roll–streak ($U,V,W$) has reflectional symmetry with respect to the edge of the fundamental box. As seen in figure 3(b), the critical layers for the waves of speeds ${\pm }1$ are placed symmetrically to that for the centre wave. The kinetic energy of the centre wave is trapped within the fundamental module (see figure 3c), and, thus, we can use the same numerical technique as adopted for the mirror-symmetric mode, again with $N=0$. As should be apparent from an examination of figure 3(d), the centre wave $\tilde {\boldsymbol {u}}$ satisfies symmetries (3.1a) and (3.1b).

3.2. Spatial intermittency

Interestingly, when a wider domain is used, our periodic array computation is able to capture spanwise-isolated solutions similar to those found in the plane Couette flow by Schneider, Gibson & Burke (Reference Schneider, Gibson and Burke2010a), Schneider, Marinc & Eckhardt (Reference Schneider, Marinc and Eckhardt2010b), Deguchi et al. (Reference Deguchi, Hall and Walton2013) and Gibson & Brand (Reference Gibson and Brand2014). The solution shown in figure 4 is an example, which may be considered a counterpart of one of the plane Couette flow solutions in Schneider et al. (Reference Schneider, Gibson and Burke2010a,Reference Schneider, Marinc and Eckhardtb), originating in a symmetry breaking bifurcation of the shift-rotation symmetry (3.1b) from the sinuous mode. By analogy, only (3.2a) and (3.3) are imposed for the roll–streak to generate the solution shown in figure 4. The structure of each roll cell is similar to that observed in the sinuous mode (see figure 3a), but these now only appear in a subinterval of $z$ smaller than the box dimension. The central wave is solved in a domain five times deeper than the fundamental module using the total streak field shown in figure 4(b). The wave satisfies (3.1a) but not (3.1b), and, hence, it is no longer stationary. Nevertheless, the reflectional symmetry with respect to the edge of the modules survives, meaning that the wave in the module adjacent to the central module propagates in the opposite direction with the same absolute speed, if the unit increment of the wave speed is ignored. In figure 4 (b) critical layers exist at $U=2+c, 1-c,c,-1-c,-2+c$, from top to bottom, where $c=-0.24$ is the speed of the most central wave.

Figure 4. A spanwise-isolated solution at ${\alpha ^{{\dagger} }}=0.8$, ${\beta ^{{\dagger} }}=0.35$, $R_f=40\times 10^{3}$. Domain shown is $y\in [-5/2,5/2],z\in [0,2{\rm \pi} /{\beta ^{{\dagger} }}]$. (a) Contours of roll stream function $\psi$. (b) Contours of total streak flow $y+U$. Critical layers for waves with speeds $2+c, 1-c,c,-1-c,-2+c$ with $c=-0.24$ are shown as black curves.

It is clear from the flow visualisation of figure 4 that the typical scale of the flow variation in the spanwise direction is comparable to the vertical module, rather than the spanwise dimension of the computational box size. In fact, the perturbation field decays exponentially in large $z$, and, hence, for sufficiently small ${\beta ^{{\dagger} }}$, the coherent structure becomes insensitive to changes of the width of the box. This strongly suggests that there is a limiting solution in the limit of ${\beta ^{{\dagger} }}\to 0$. The existence of such a spanwise-isolated canonical array solution is the basis of the discussion in § 4, where the multiple-scale version of the VWI array theory will be used to derive the logarithmic law of the wall.

3.3. Momentum transport and emergence of uniform momentum zones

Now let us examine the typical parameter dependence of the solutions drawing the solution branch by varying ${\alpha ^{{\dagger} }}$ for fixed ${\beta ^{{\dagger} }}$. Here we choose the sinuous mode as the representative solution.

In order to illustrate the bifurcation diagram, we need to select a typical scalar flow quantity. The shear of the mean flow on the wall is commonly used as a measure of nonlinear states in the computational study of wall-bounded shear flows, because it can be related to the momentum transport. For our triply periodic VWI problem, by integrating the streak (2.7a) over one period in $Z$ we can deduce the fact that the momentum transport $\sigma$ defined by

(3.4)\begin{equation} \sigma=\frac{\mathrm{d} U_m}{\mathrm{d} y}-\frac{1}{2{\rm \pi}}\int^{2{\rm \pi}}_0UV\,\mathrm{d} Z \end{equation}

must be a constant for all $y$. Here we have defined the mean flow distortion by

(3.5)\begin{equation} U_m(y)=\frac{1}{2{\rm \pi}}\int^{2{\rm \pi}}_0U\,\mathrm{d} Z. \end{equation}

For the sinuous and mirror-symmetric modes, there are axes of symmetry (e.g. $y=\pm 1/2$) on which $V$ vanishes. Thus, for those modes the momentum transport and the shear of the mean flow are related by $\sigma =(\mathrm {d} U_m/\mathrm {d} y)|_{y=\pm 1/2}$, similar to wall-bounded flows. We also note that by further integrating (3.4) over one period in $y$ we can establish the simple relationship

(3.6)\begin{equation} \sigma={-}\langle UV\rangle, \end{equation}

where $\langle \,\rangle$ denotes an average over one periodic box in the $y,Z$ plane. As shown by H18 and to be discussed in § 4, the quantity $\sigma ({\alpha ^{{\dagger} }},{\beta ^{{\dagger} }})$ plays an important role when the VWI array theory for the homogeneous shear flow is extended by the method of multiple scales.

Figure 5 illustrates the variation of $\sigma$ along the sinuous-mode solution branch. Here we use ${\beta ^{{\dagger} }}=4$, i.e. the same box aspect ratio used to draw the bifurcation diagram of plane Couette flow study in figure 1 of Deguchi & Hall (Reference Deguchi and Hall2014c) (hybrid systems for both flows become identical except for boundary conditions). The overall bifurcation diagram of the array problem (solid curve) is similar to that for the plane Couette flow computed in Deguchi & Hall (Reference Deguchi and Hall2014c) (dashed curve). The branch turns back around ${\alpha ^{{\dagger} }}\approx 2.86$. The existence of the turning point suggests that the flow may have a minimal flow unit, consistent with the observation in finite-Reynolds-number shearing box simulations by Sekimoto et al. (Reference Sekimoto, Dong and Jiménez2016). According to the empirical observation of that paper, the shearing box turbulence generates typical properties of near-wall turbulence when ${\beta ^{{\dagger} }}/{\alpha ^{{\dagger} }} \gtrsim 2$ and ${\beta ^{{\dagger} }} \gtrsim {\rm \pi}$. The computational boxes we selected in figure 6(a,b) satisfy those conditions.

Figure 5. Bifurcation diagram for the sinuous mode at ${\beta ^{{\dagger} }}=4$ showing variation of momentum-transport parameter $\sigma$ (see text) as a function of streamwise wavenumber ${\alpha ^{{\dagger} }}$. Solid line: $R_f=280\times 10^{3}$; open circles: $R_f=200\times 10^{3}$; dashed line: equivalent result for plane Couette flow geometry by Deguchi & Hall (Reference Deguchi and Hall2014c).

Figure 6. Upper-branch sinuous-mode solutions at ${\beta ^{{\dagger} }}=4$, showing the effect of streamwise wavenumber ${\alpha ^{{\dagger} }}$. (a–c) Contours of total streak $(y+U)$ for ${\alpha ^{{\dagger} }}=0.8$, 1.6, 2.4, respectively, with critical layers for waves with speeds $-1,0,1$. (d–f) Corresponding contours of kinetic energy $|\tilde {\boldsymbol {u}}|^{2}$ for the stationary wave with its critical layer.

Upon passing the turning point, the upper-branch solutions are harder to resolve in the channel flow studies, and this is indeed the case here. The difficulty is due to the eventual emergence of a homogenised zone developing in the streak field. The upper panels of figure 6 show that the homogenisation develops with decreasing ${\alpha ^{{\dagger} }}$, associated with progressively more contorted critical layers. As the lower panels of figure 6 illustrate, when the critical layer straddles the homogenised zone, the local shear around the layer must necessarily be weak, meaning that we need very large $R_f$ in order to ensure the thinness of the critical layer. In figure 5 the asymptotic convergence of the solution was checked by increasing $R_f$ up to $280\times 10^{3}$. For such very large values of $R_f$, the critical layer near the non-homogenised zone becomes extremely thin, as seen in figure 6(a–c). This multi-scale feature of the flow along the critical layer requires very high resolution in the computation (to resolve the wave structure, here we have employed 640 Chebyshev modes in the vertical direction, while up to 42 Fourier harmonics are considered in $Z$). These features of the upper-branch solutions are reminiscent of the results of Deguchi & Hall (Reference Deguchi and Hall2014c), where uniform momentum zones were found for the first time for equilibrium solutions in plane Couette flow.

4.1. The slow-scale evolution problem

Equations (4.1) are parameter free and so we base our analysis on the assumption that the structures we seek are distributed over a large range of values of $y=O(\epsilon ^{-1})\gg 1$. The quantity $\epsilon$ is however not an arbitrarily small parameter but a quantity to be found numerically by full numerical solutions of the wall-layer equations. We will look for an array structure with the local cell depth a function of the slow variable spatial variable $Y=\epsilon y$. Thus, at the fast scale we define the phase function $\varPhi$ for the locally periodic flow structure by

(4.2)\begin{equation} \varPhi= \varTheta(Y)/\delta. \end{equation}

The above choice of scale means that the vertical derivative is transformed as

(4.3)\begin{equation} \partial_y\to \frac{\epsilon \varTheta_Y}{\delta}\partial_\varPhi+\cdots \end{equation}

so that the cell depth at a large height $y$ in wall units is $O(\delta /(\epsilon \varTheta _Y))$. For that reason, we will refer to $\varTheta _Y$ as the vertical wavenumber. We shall see from a higher-order analysis that $\delta =\epsilon ^{2/5}$ is the appropriate distinguished limit.

We now consider a slowly varying distributed VWI array driven by an infinite set of waves of different streamwise wavenumbers moving downstream with the local mean flow speed at the location $Y=Y_n$. The mean part of the flow is predominantly in the $x$ direction. In terms of the phase variable $\varPhi$, the waves have critical layers located at $\varPhi =\varPhi _n=\varTheta (Y_n)/\delta$ where $\varPhi _{n+1}-\varPhi _n=1$. We further assume that the $n$th wave has streamwise wavenumber $\alpha (Y_n)$. Note here that we are not considering a single wave with its wavenumber varying in the $Y$ direction but an infinite set of distinct waves with wavenumber assigned using the function $\alpha (Y)$. Since the local array structure has a roll–streak–wave flow with variations in all three directions on the same scale we write $(X,Z)= (\epsilon /\delta ){(x,\beta z)}$, and on the basis of the above discussion we look for a roll–streak flow of the form (see H18 also)

(4.4a)\begin{gather} u = \frac{\overline{u}(Y)}{\epsilon} + \frac{\delta \overline{u}_Y }{\epsilon \varTheta_Y} \{ U _{1}(\varPhi,Y,Z)+ \delta U _{2}(\varPhi,Y,Z)+ \cdots\}+\sum_{n={-}\infty}^{\infty}u_w^{(n)}, \end{gather}

(4.4b)\begin{gather}v = \frac{\epsilon \varTheta_Y}{\delta} \{V _{1}(\varPhi,Y,Z)+ \delta V _{2}(\varPhi,Y,Z)+\cdots\}+\sum_{n={-}\infty}^{\infty}v_w^{(n)}, \end{gather}

(4.4c)\begin{gather}w =\frac{\epsilon\varTheta_Y}{\delta} \{W _{1}(\varPhi,Y,Z)+ \delta W _{2}(\varPhi,Y,Z)+\cdots\}+\sum_{n={-}\infty}^{\infty}w_w^{(n)}, \end{gather}

(4.4d)\begin{gather}p ={Kx}+ \left(\frac{\epsilon \varTheta_Y}{\delta}\right )^{2} \{P _{1}(\varPhi,Y,Z)+\delta P _{2}(\varPhi,Y,Z)+\cdots\}+ \sum_{n={-}\infty}^{\infty}p_w^{(n)}, \end{gather}

where $K$ represents a streamwise pressure gradient driving the laminar flow, assumed to be smaller than $O(\epsilon )$ when measured in wall units. We shall shortly see that the solvability condition is needed for the terms with subscript 2. Components $(u_w^{(n)},v_w^{(n)},w_w^{(n)},p_w^{(n)})$ represent the $n$th wave. In particular, the wave pressure takes the form

(4.4e)\begin{equation} p_w^{(n)}= \overline{u}_Y^{5/6} \left (\frac{\epsilon\varTheta_Y}{\delta}\right )^{1/3}\rho_1(Y) E^{(n)}p_1(\varPhi-\varPhi_n,Y,Z)+\mathrm{c.c.}+\cdots, \end{equation}

where

(4.5)\begin{equation} E^{(n)}=\exp\left \{\mathrm{i} X-\mathrm{i}\alpha(Y_n) \left (\frac{\overline{u}(Y_n)}{\delta} \right )t \right \}. \end{equation}

Here $\rho _1$ is the amplitude of the wave, which is chosen so that the shear associated with the roll field is modified at leading order when crossing the critical layer; the procedure is essentially identical to that described in detail in Hall & Sherwin (Reference Hall and Sherwin2010).

The scalings in (4.4) were motivated as follows. First, the scaling of the mean flow $\bar {u}$ is chosen so that it is consistent to the base flow $\bar {u}=y$. The factor $(\epsilon \varTheta _Y)/\delta$ in the roll expansions are introduced so that, in view of (4.3), the roll advection term balances to the viscous term. The need for the factor $(\delta \overline {u}_Y)/(\epsilon \varTheta _Y)$ in the streak expansion can be seen by considering the mean flow equation, assuming that to leading order it is driven by the nonlinear interaction of the roll and the streak. The origin of the wave scaling is less obvious as it involves the critical layer analysis; however, the essential mechanism that is determined so that the wave drives the roll–streak is unchanged from the usual VWI. The derivation of the scaling (4.4e) is explained in detail in the appendix. For simplicity in the above expansion, we have assumed that each cell is maintained by a single wave. If, say, two waves are involved we must, as previously noted in § 2.1, split the contributions into two sums and allow for the wave speeds of each wave to differ from the mean flow speed at the middle of each cell. The more general case is easily accounted for by replacing the summation terms in (4.4) into two sums and inserting a wave speed into (4.5).

The factors dependent on $\varTheta _Y, \overline {u}_Y$ in (4.4) have been introduced in order to facilitate comparison with the uniform shear case described in § 2. Owing to those factors, at leading order the roll–streak–wave flow satisfies the triply periodic canonical array problem (2.6)–(2.8a,b) but with $y$ replaced by $\varPhi$, as shown in H18. In order to do so the local effective wavenumbers and the global wavenumbers must be related as

(4.6)\begin{equation} ({\alpha^{{{\dagger}}}},{\beta^{{{\dagger}}}})= \left (\frac{\alpha}{\varTheta_Y},\frac{\beta}{\varTheta_Y} \right), \end{equation}

in which case

\begin{gather}\boldsymbol{U}_{1}(\varPhi,Y,Z;\alpha,\beta,\varTheta_Y,\overline{u}_Y) = \boldsymbol{ {U}} (\varPhi,Z;{\alpha^{{{\dagger}}}},{\beta^{{{\dagger}}}}),\notag\\ p_1(\varPhi-\varPhi_n,Y,Z;\alpha,\beta,\varTheta_Y,\overline{u}_Y) = \tilde{p}(\varPhi-n,Z;{\alpha^{{{\dagger}}}},{\beta^{{{\dagger}}}})\equiv \tilde{p}^{(n)}(\varPhi,Z;{\alpha^{{{\dagger}}}},{\beta^{{{\dagger}}}}),\notag\\ \rho_1(Y;\alpha,\beta,\varTheta_Y,\overline{u}_Y) = \rho({\alpha^{{{\dagger}}}},{\beta^{{{\dagger}}}}), \notag \end{gather}

together with similar expressions for the wave velocity field. The left-hand sides of these equations are the quantities appearing in the expansions (4.4), while on the right-hand side the solution of the canonical equations (2.7) are used. The implication is that hereafter we can remove the subscripts 1 which appear in the expansions (4.4).

We note at this stage that $\overline {u}_Y(Y)$ and the normal and streamwise wavenumbers $\varTheta _Y(Y),\alpha (Y)$ are all unknown and it is of interest to know how many of these functions can be specified a priori. If we think of the array being formed by the bringing together of the vertically localised canonical states of Blackburn et al. (Reference Blackburn, Hall and Sherwin2013) then it might be anticipated that, certainly for large separations, the normal wavenumber $\varTheta _Y(Y)$ should be at our disposal. We might anticipate that the mean flow function is then fixed to be the function supporting that solution leaving $\alpha$ also to be chosen. However, we shall shortly see that when the flow in each cell strongly interacts, there are two constraints on the slow dynamics. It will turn out that here, unlike in H18, only one of $\varTheta _Y(Y), \alpha (Y), \overline {u}_Y(Y)$ can be chosen.

Hall (Reference Hall2018) obtained the first constraint by considering the development of $\overline {u}$. Integrating the streamwise momentum equation over one periodic cell $(X,\varPhi ,Z)\in [0,2{\rm \pi} ]\times [0,1]\times [0,2{\rm \pi} ]$, we obtain

(4.7)\begin{equation} {\overline{u}_{YY}}+{(\overline{u}}_Y\sigma)_Y=0. \end{equation}

The function $\sigma =-\langle UV\rangle$, where $\langle \,\rangle$ denotes the average over $0<\varPhi <1$, $0< Z<2{\rm \pi}$, can be obtained by the leading-order problem (see (3.6)) and depends on the slow variable $Y$ through ${\alpha ^{{\dagger} }},{\beta ^{{\dagger} }}$. As pointed out in H18, (4.7) is essentially the RANS equation for the mean flow but with the only difference that here the mean is taken over $X,\varPhi ,Z$ rather than $X,Z$. At this stage, in the absence of slow time dependence, one of the functions $\overline {u}_Y(Y), {\varTheta _Y}(Y),\alpha (Y)$ is fixed, and in H18 no further constraint was enforced.

We shall see below that in fact the second constraint is available by consideration of the higher-order equations for wave, roll and streak. The slow dependence on $Y$ of the leading-order roll–streak flow is through the different slow amplitude functions multiplying each quantity and through the local wavenumbers ${\alpha ^{{\dagger} }},{\beta ^{{\dagger} }}$. As is usual with the phase equation approach, the slow variation will be constrained at next order through a solvability condition.

As observed earlier, the nature of the wave field sustaining the locally periodic roll–streak structure differs according to the symmetries to be imposed. For simplicity, here we assume that (i) just one dominant sustaining wave is present in each cell and (ii) the local roll–streak field possesses the mirror symmetry mentioned in § 3.1. With this symmetry, we can assume that the critical layer is flat and located at $Y=Y_n$. Furthermore, this allows us to significantly contract the algebra that is needed in the higher-order critical layer analysis (see the appendix) to derive the solvability condition.

At next order the slow variation of the roll–streak with respect to $Y$ comes into play, and it is $O(\delta )$ smaller than the leading-order effect. We also remark that even when there is no slow effect in the flow there are some higher-order terms in the VWI. The size of the second type of higher-order terms are comparable to the critical layer thickness, which is the largest perturbative parameter in the asymptotic analysis. The critical layer analysis shown in the appendix reveals that the layer thickness is $O((\epsilon /\delta )^{2/3})$, and by balancing the two perturbative effects, we have the distinguished limit $\delta =\epsilon ^{2/5}$. Collecting the terms of this size, the second-lowest-order roll–streak equations are written as

(4.8)\begin{equation} (\mathbb{A}+\mathbb{B}+\mathbb{C}+\mathbb{D}+\mathbb{E})\boldsymbol{X}= \boldsymbol{F},\end{equation}

where $\boldsymbol {X},\mathbb {A},\mathbb {B},\mathbb {C},\mathbb {D},\mathbb {E}$ are defined by

\begin{gather} {{{\boldsymbol{X}}}} =\begin{bmatrix} U _2 \\ V_2 \\ W_2\\ P_2 \end{bmatrix},\quad {\mathbb{A}}=\begin{bmatrix} {\dfrac{\partial^{2}}{\partial \varPhi^{2}}+ {{\beta^{{{\dagger}}}}^{2}}\dfrac{\partial^{2}}{\partial Z^{2}}} & 0 & 0 & 0 \\ 0 & {\dfrac{\partial^{2}}{\partial \varPhi^{2}}+{{\beta^{{{\dagger}}}}^{2}}\dfrac{\partial^{2}}{\partial Z^{2}}} & 0 & 0\\ 0 & 0 & {\dfrac{\partial^{2}}{\partial \varPhi^{2}}+{{\beta^{{{\dagger}}}}^{2}}\dfrac{\partial^{2}}{\partial Z^{2}}} & 0\\ 0 & 0 & 0 & 0 \end{bmatrix}, \notag\end{gather}

\begin{gather} {\mathbb{B}}=\begin{bmatrix} -{{ {V}}\dfrac{\partial }{\partial \varPhi }} & 0 & 0 & 0 \\ 0 & -{{ {V}}\dfrac{\partial }{\partial \varPhi }} & 0 & 0\\ 0 & 0 & -{{ {V}}\dfrac{\partial }{\partial \varPhi }} & 0\\ 0 & { \dfrac{\partial }{\partial \varPhi }} & 0 & 0 \end{bmatrix}, \notag\end{gather}

\begin{gather} {\mathbb{C}}=\begin{bmatrix} -{{\beta^{{{\dagger}}} {W}}\dfrac{\partial }{\partial Z }} & 0 & 0 & 0 \\ 0 & -{{\beta^{{{\dagger}}} {W}}\dfrac{\partial }{\partial Z }} & 0 & 0\\ 0 & 0 & -{{\beta^{{{\dagger}}} {W}}\dfrac{\partial }{\partial Z }} & 0\\ 0 & 0 & {{\beta^{{{\dagger}}} }\dfrac{\partial }{\partial Z }} & 0 \end{bmatrix},\notag\\ {\mathbb{D}}=\begin{bmatrix} 0 & -{ {U}}_\varPhi-1 & 0 & 0 \\ 0 & -{ {V}}_\varPhi & 0 & 0\\ 0 & -{ {W}}_\varPhi & 0 & 0\\ 0 & 0 & 0 & 0 \end{bmatrix},\quad {\mathbb{E}}=\begin{bmatrix} 0 & 0 & -{\beta^{{{\dagger}}}}{ {U}}_Z & 0 \\ 0 & 0 & -{\beta^{{{\dagger}}}}{ {V}}_Z & 0 \\ 0 & 0 & -{\beta^{{{\dagger}}}}{ {W}}_Z & 0 \\ 0 & 0 & 0 & 0 \end{bmatrix}. \notag\end{gather}

The term $\boldsymbol {F}=[F_1,F_2,F_3,F_4]^{\text {T}}$ on the right-hand side of (4.8) arises from the slow variation of the coherent structure with $Y$ and is given by

(4.9a)\begin{gather} F_1 =\frac{{\overline{u}}_{YY}}{\varTheta_Y {\overline{u}}_Y}({ {U}}{ {V}}-2{{ {U}}}_\varPhi)+\frac{\varTheta_{YY}}{\varTheta_Y^{2}}({ {U}}_\varPhi-{ {U}}{ {V}}) +\frac{1}{\varTheta_Y }({ {V}}{ {U}}_Y-2{ {U}}_{Y \varPhi}), \end{gather}

(4.9b)\begin{gather}F_2 =\frac{\varTheta_{YY}}{\varTheta_Y^{2}}({ {V}}^{2}-3{ {V}}_\varPhi+{ 2{P}})+\frac{1}{\varTheta_Y}({ {P}}_Y-2{ {V}}_{Y \varPhi}+{ {V}}{ {V}}_Y) +f_2, \end{gather}

(4.9c)\begin{gather}F_3 =\frac{\varTheta_{YY}}{\varTheta_Y^{2}}({ {V}}{ {W}}-3{ {W}}_\phi)+\frac{1}{\varTheta_Y}({ {V}}{ {W}}_Y-2{ {W}}_{Y \varPhi})+ f_3, \end{gather}

(4.9d)\begin{gather}F_4 ={-}\frac{\varTheta_{YY}}{\varTheta_Y^{2}}{ {V}}-\frac{1}{\varTheta_Y}{ {V}}_Y. \end{gather}

Here

\begin{gather} f_2=\sum_{n={-}\infty}^{\infty}\frac{\varTheta_Y^{2/3}}{\overline{u}_Y^{1/3}\alpha^{{{\dagger}} 2}} \left \{ \left (\frac{|\tilde{p}^{(n)}_{\varPhi}|^{2}}{(\varPhi+U )^{2}}\right )_{\varPhi} +{\beta^{{{\dagger}}}}\left ( \frac{{\beta^{{{\dagger}}}}\tilde{p}^{(n)}_{\varPhi}\tilde{p}_{Z}^{(n)*}}{(\varPhi+U )^{2}} \right )_{Z}\right \}+\mathrm{c.c.},\notag\\ f_3=\sum_{n={-}\infty}^{\infty}\frac{\varTheta_Y^{2/3}}{\overline{u}_Y^{1/3}\alpha^{{{\dagger}} 2}} \left \{ \left ( \frac{{\beta^{{{\dagger}}}}\tilde{p}^{(n)}_{\varPhi}\tilde{p}_{Z}^{(n)*}}{(\varPhi+U )^{2}} \right )_{\varPhi} +{\beta^{{{\dagger}}}}\left ( \frac{\beta^{{{\dagger}} 2}|\tilde{p}^{(n)}_{Z}|^{2}}{(\varPhi+U )^{2}} \right )_{Z}\right \}+\mathrm{c.c.}, \notag \end{gather}

are associated with the wave Reynolds stress outside of the critical layer. The above perturbed roll–streak equations must once again be solved subject to periodicity in $Z$. We do not write down the perturbed wave pressure equation because it is readily shown that, whatever the symmetries, the perturbed wave pressure equation does not enter the perturbed jump conditions and so we need not consider it here.

We see that the normal and spanwise momentum equations are independent of the streamwise momentum equation in (4.8), just as was the case for (2.7). The link between the roll and streak equations must therefore come from the jump conditions. For the mirror-symmetric mode, we can take the critical layers to be at $\varPhi =0, \pm 1, \pm 2,\ldots$ and so the boundary conditions in $\varPhi$ can be most conveniently written in the form

\begin{gather} U_2(1, Z) ={U_2} (0, Z) ,\quad {U_{2\varPhi}} (1, Z) ={U_{2 \varPhi}}(0, Z) ,\notag\\ V_2(1, Z) ={V_2} (0, Z) ,\quad {V_{2\varPhi}} (1, Z) ={V_{2 \varPhi}}(0, Z),\notag\\ W_{2\varPhi}(1, Z) =W_{2 \varPhi}(0, Z), \notag \end{gather}

and

(4.10)\begin{equation} W_2(0, Z)-W_2(1, Z)={\beta^{{{\dagger}}}}\left \{ \frac {{U_2} (0, Z) \hat{J}(Z)}{1+U_{\varPhi}(0,Z)}\right \}_Z, \end{equation}

where the jump is derived in (A13) and $\hat {J}(Z)$ is the function defined in (2.10a–c).

The crucial point here is that the homogeneous form of (4.8) and the associated boundary conditions now has a solution $(U_2,V_2,W_2,P_2)=\partial _\varPhi (U,V,W,P)$ representing the translational invariance of the roll–streak flow. The outcome is that the inhomogeneous problem specified by (4.8) and the associated boundary conditions require a solvability condition to be satisfied; in the usual manner this is done by constructing the system adjoint to the homogeneous problem. The adjoint is formed by multiplying the homogeneous form of (4.8) by the adjoint function ${\boldsymbol {S}}=[S_1,S_2,S_3,S_4]^{\text {T}}$ and integrating over $0<\varPhi < 1$, $0< Z<2 {\rm \pi}$ and integrating by parts. We find that the adjoint function must satisfy the equation

(4.11)\begin{equation} (\mathbb{A}-\mathbb{B}-\mathbb{C}+\mathbb{D}^{\text{T}}+\mathbb{E}^{\text{T}}) \boldsymbol{S} =\boldsymbol{0}, \end{equation}

together with the boundary conditions

\begin{gather} S_1(1, Z)=S_1(0, Z),\quad S_{1 \varPhi}(1, Z)=S_{1 \varPhi}(0, Z),\notag\\ S_2(1, Z)=S_2(0, Z),\quad S_{2 \varPhi}(1, Z)=S_{2 \varPhi}(0, Z),\notag\\ S_3(1, Z)=S_3(0, Z),\quad S_{3 \varPhi}(1, Z)=S_{3 \varPhi}(0, Z). \notag \end{gather}

Moreover, $S_1,S_2$ are odd, $S_3,S_4$ are even with respect to the $\varPhi =0$ axis.

The solvability condition then takes the form

(4.12)\begin{equation} 0=\sum_{k=1}^{k=4} \langle S_k F_k\rangle. \end{equation}

Here we note that there is no energy input for $S_1$ in (4.11) and, thus, $S_1\equiv 0$. Further using (4.9) and the fact that $\beta$ is constant, (4.12) may be written in the form

(4.13)\begin{equation} e_1\left (\frac{\alpha}{\varTheta_Y}\right )_Y+ e_2\left (\frac{\beta}{\varTheta_Y}\right )_Y+ e_3\frac{\varTheta_{YY}}{\varTheta_Y}=0, \end{equation}

which is the second constraint of our interest. Here $e_1,e_2,e_3$ are defined by

(4.14a)\begin{gather} e_1=\langle S_2(VV_{{\alpha^{{{\dagger}}}}}-2V_{{\alpha^{{{\dagger}}}}}+ P_{{\alpha^{{{\dagger}}}}})+S_3(VW_{{\alpha^{{{\dagger}}}}}-2W_{\varPhi {\alpha^{{{\dagger}}}}})+ S_4V_{{\alpha^{{{\dagger}}}}} \rangle, \end{gather}

(4.14b)\begin{gather}e_2=\langle S_2(VV_{{\beta^{{{\dagger}}}}}-2V_{{\beta^{{{\dagger}}}}}+ P_{{\beta^{{{\dagger}}}}})+ S_3(VW_{{\beta^{{{\dagger}}}}}-2W_{\varPhi {\beta^{{{\dagger}}}}})+S_4V_{{\beta^{{{\dagger}}}}} \rangle, \end{gather}

(4.14c)\begin{gather}e_3=\langle S_2(V^{2}-3V_{\varPhi}+2P)+ S_3(VW-3W_{\varPhi})+S_4V \rangle, \end{gather}

and are functions of the slow variable $Y$ through its dependence on ${\alpha ^{{\dagger} }},{\beta ^{{\dagger} }}$. Note that by symmetry $f_2$ and $f_3$ do not contribute to the solvability. In the above set of equations, the subscripts ${\alpha ^{{\dagger} }},{\beta ^{{\dagger} }}$ denote the derivative with respect to the local wavenumbers and, hence, in principle, the integrand can be computed by a parameter search of the local solutions and the solution of the adjoint problem (4.11).

Although the evolution (4.13) has been derived for the mirror-symmetric mode, the essential structure of the equation is unchanged for other symmetries. In any case, the evaluation of the integrals appearing in the definitions above requires the solution of the adjoint partial differential equations. That rather daunting calculation is not discussed here but in the next section, where we will see that important conclusions regarding the possibility of a totally self-similar state consistent with the attached eddy hypothesis can be drawn without doing so.

We now have sufficient information to see how the canonical vertically localised states may be moved close enough together for them to interact. We recall that we have three unknown functions ${u_Y},{\varTheta _Y},\alpha$ and the two evolution equations (4.7) and (4.13) linking those functions. If we think of each VWI cell as being initially far away from its neighbours to interact it seems most natural to determine how the mean flow $\overline {u}$ and wavenumber distribution change in order to support the mutual interaction. In other words, it is perhaps more natural to specify $\varTheta _Y$ and determine $\overline {u}_Y, \alpha$ from the two evolution equations.

4.2. Self-similar coherent structures with logarithmic mean flow profile

Here we shall explore the role of the VWI arrays in the logarithmic layer by seeking self-similar distributed solutions consistent with the attached eddy hypothesis. Firstly we note that the attached eddy hypothesis requires that the dominant streamwise and spanwise wavelengths vary linearly with distance from the wall. The apparently insurmountable constraint not permitting self-similarity is that the spanwise wavenumber $\beta$ is a constant so the spanwise variation cannot be self-similar with respect to $Y$. So at first consideration it appears that the slowly varying distributed structures we have constructed are not consistent with the hypothesis; the key to the resolution is the consideration of the local VWI array solutions fully isolated in the spanwise direction. Thus, we will see that self-similarity is achieved when the spanwise variations of the roll–streak structure are no longer influenced by the distant periodic constraint in the spanwise direction.

We have already seen such isolated solutions indeed exist in § 3.2, and confirmed that their flow structure becomes insensitive to the change of ${\beta ^{{\dagger} }}$ for any choice of sufficiently small ${\beta ^{{\dagger} }}$. First we note that from the invariance of the isolated flow structure against ${\beta ^{{\dagger} }}$, we can derive the asymptotic properties

(4.15a–c)\begin{equation} \sigma\rightarrow {\beta^{{{\dagger}}}}\sigma_0({\alpha^{{{\dagger}}}}),\quad \frac{e_1}{e_3}\rightarrow B_1({\alpha^{{{\dagger}}}}),\quad \frac{e_2}{e_3}\rightarrow \frac{B_2({\alpha^{{{\dagger}}}})}{{\beta^{{{\dagger}}}}}\end{equation}

as ${\beta ^{{\dagger} }}\to 0$. The scaling law for $\sigma$ is essentially the mean flow scaling found in Deguchi et al. (Reference Deguchi, Hall and Walton2013), where the asymptotic behaviour was observed for ${\beta ^{{\dagger} }}$ numerically smaller than unity. It can be deduced by recalling that $\sigma =-\langle UV\rangle$, where $\langle \,\rangle$ denotes the average over the periodic box. If $UV$ is only non-zero for an $O(1)$ spanwise region of the domain and the rest of $O(1/{\beta ^{{\dagger} }})$ region does not contribute to the average, obviously the scaling above is expected. The other two scaling laws can be similarly deduced from (4.14).

Noting that complete self-similarity in the global problem is only achieved if the local streamwise wavelength is everywhere constrained such that ${\alpha ^{{\dagger} }}$ is a constant, we find that in the limit (4.15a–c) the evolution (4.7) and (4.13) reduce to

(4.16)\begin{gather} \overline{u}_{YY}+\sigma_0\beta\left(\frac{\overline{u}_Y}{\varTheta_Y}\right)_Y=0, \end{gather}

(4.17)\begin{gather}(1-B_2)\frac{\varTheta_{YY}}{\varTheta_Y}=0, \end{gather}

respectively. From the second evolution (4.17) we see that the self-similar roll–streak field is possible only if

(4.18)\begin{equation} B_2({\alpha^{{{\dagger}}}})=1, \end{equation}

which fixes the values of ${\alpha ^{{\dagger} }}$ permitting self-similarity. If the flow is to be consistent with the attached eddy hypothesis, then the cell depth must locally scale linearly in $Y$. Thus, we choose the vertical wavenumber

(4.19)\begin{equation} \varTheta_Y=1/(Y-Y_0), \end{equation}

assuming that the self-similarity begins at $Y=Y_0$. For spanwise-isolated local solutions, whatever the value of ${\beta ^{{\dagger} }}$ is, the typical streamwise and spanwise scales are comparable to the depth of the vertical domain. Thus, (4.19) implies that in the global problem the cell depth and the typical spanwise and streamwise length scales of the equilibrium structure all vary like $Y-Y_0$, and, hence, the size variation of the VWI array is the same as envisaged in the attached eddy hypothesis.

Given a choice of $\varTheta (Y)$ and a spanwise-isolated solution of the canonical VWI array problem, we are able to construct a self-similar solution in terms of the wall-unit variables $y,z$. As an illustrative example, in figure 7 we show the global flow structure inferred by the local leading-order solution seen in figure 4 and the distribution $\varTheta =\ln Y$. In order to explain how that picture was generated, let us write $(y,z)$ used in §§ 2 and 3 as $(\hat {y},\hat {z})$ to avoid confusion. The choice of $\varTheta$ means that $\varTheta _Y=Y^{-1}$, i.e. $Y_0=0$ in (4.19) and the self-similarity begins at $y=0$. From (4.6), the local wavenumber varies like ${\beta ^{{\dagger} }}=\beta Y$. Further using the two definitions of $Z$ introduced in §§ 2 and 4, we have the relationship $Z=(\epsilon ^{\frac {3}{5}})\beta z={\beta ^{{\dagger} }}\hat {z}=\beta \epsilon y\hat {z}$, while for the vertical coordinate, we can simply use the definition of $\varPhi$ (which equals $\hat {y}$) to deduce

\[ \hat{y}=\frac{\ln y}{\epsilon^{\frac{2}{5}}},\quad \hat{z}=\frac{z}{\epsilon^{\frac{2}{5} y}}, \]

from which we can find the flow at any point in the $y$–$z$ plane through the local array solution at $(\hat {y}, \hat {z})$. As may be seen in figure 4, self-similarity consistent with the attached eddy hypothesis can indeed be reached with the isolated solutions.

Figure 7. The roll stream function of the self-similar solution in the wall-unit coordinate $(y,z)$. The isolated solution shown in figure 4 is used together with $\varTheta =\ln Y$ and $\epsilon =0.001$. The horizontal coordinate is $z\in [0,1200]$, the vertical coordinate is $y\in [50,500]$.

Note that we have not solved for the adjoint problem and so we can only assume here that the condition (4.18) may be achieved at some ${\alpha ^{{\dagger} }}$. On the assumption that ${\alpha ^{{\dagger} }}$ satisfying (4.18) exists, the mean flow (4.16) can be integrated to yield

(4.20)\begin{equation} \overline{u}(Y)=C_1\ln\left (Y-Y_0+\frac{1}{\sigma_0\beta}\right )+C_2, \end{equation}

where $C_1, C_2$ are integration constants.

Therefore, by writing the mean part of $u$ as $u_m(y)=\epsilon ^{-1}\overline {u}(Y)$, we can see that the self-similar VWI array leads naturally to the Kármán law of the wall

(4.21)\begin{equation} u_m(y)=\frac{1}{\kappa}\ln y+C, \end{equation}

in which we have chosen

(4.22)\begin{equation} Y_0=\frac{1}{\sigma_0\beta}. \end{equation}

We remark here that this choice and the linear wall-normal variation of the eddy length scale are assumptions to yield the log-law, although they seem physically reasonable. The use of those assumptions may conversely be justified by certain special scale invariances of the logarithmic profile, e.g. as noted in Moarref et al. (Reference Moarref, Sharma, Tropp and McKeon2013) and H18, but we do not go into detail here.

The Kármán constant $\kappa$ appearing in (4.21) may be found by solving the matching problem which connects the solution in what we believe to be the logarithmic layer to the adjacent regions such as the buffer and/or wake layers. The conclusion is therefore that our self-similar VWI structure does lead to the law of the wall, but the associated value of the Kármán constant can only be fixed by a more complete matched asymptotic expansion of the entire flow. This outcome is slightly disappointing but the structure we have constructed is arguably the first nonlinear equilibrium state built from slowly varying SSP or VWI states shown to be consistent with the law of the wall and the attached eddy hypothesis.

Although such a matching problem is out of the scope of the present study, it is nevertheless interesting to see how the self-similar structure breaks down on exiting the logarithmic layer. As remarked earlier, the self-similarity starts at $Y=Y_0$ meaning that the roll cell size formally shrink to zero there; see figure 8. However, the VWI theory breaks down before reaching this position. When the roll cell size becomes $O(1)$ in the wall unit, the flow reaches the Kolmogorov microscale and, hence, the assumption used for the VWI array theory that the cell size is asymptotically large in the wall unit is not satisfied. Since the size of the roll cell varies like $(\delta /\epsilon )(Y-Y_0)$, the VWI array is terminated at slightly ($O(\epsilon ^{3/4}$) in $Y$) above the location $Y=Y_0$. The small-scale flows may persist all the way down to the viscous sublayer, and here the localised solutions at the Kolmogorov microscale found by Deguchi (Reference Deguchi2015), Eckhardt & Zammert (Reference Eckhardt and Zammert2018) and Yang et al. (Reference Yang, Willis and Hwang2019), or their distributed form, (see Yang, Willis & Hwang Reference Yang, Willis and Hwang2018) would be relevant.

Figure 8. Schematic of coherent structures in near-wall turbulent boundary layers, with distances given in wall units. The self-similar VWI eddies of size $O(\epsilon ^{-3/5})$ produce a logarithmic mean flow profile in a layer of depth $O(\epsilon ^{-1})$. At the lower bound of the logarithmic layer, the eddy size reaches the Kolmogorov microscale $O(1)$. At the other bound $y=Y_c/\epsilon$, isolated eddies begin to interact with each other and, thus, self-similarity must be destroyed. The global flow has a spanwise period $2{\rm \pi} /\beta \epsilon ^{3/5}$.

The opposite end of the logarithmic region can be found by recalling that in the canonical VWI array problem localisation is possible only when ${\beta ^{{\dagger} }}$ is sufficiently small (i.e. the spanwise dimension of the periodic box is large in the local problem). Thus, when ${\beta ^{{\dagger} }}$ becomes sufficiently large and reaches a critical value, $\beta _c^{{\dagger} }$ say, the spanwise-isolated coherent structures originally placed at a large distance apart begin to interact with each other. If this saturation occurs in the global problem the flow should deviate from the self-similar state and, hence, the wake region may appear. Since the local wavenumber varies like ${\beta ^{{\dagger} }}=\beta (Y-Y_0)$, the vertical position in $Y$ at which the self-similarity halts can be estimated as

(4.23)\begin{equation} Y_c=Y_0+\frac{\beta_c^{{{\dagger}}}}{\beta}= (1+\sigma_0\beta_c^{{{\dagger}}})Y_0. \end{equation}

Here we recall that the slow evolution equations imply that the vertical mean flow variation is determined by the wavenumber dependence of the local solution, and we used the asymptotic behaviour of the isolated solution in the limit of ${\beta ^{{\dagger} }}\rightarrow 0$ to deduce the logarithmic mean flow profile. Thus, if the isolated structure were only occurring in an extremely wide box, so that $\beta _c^{{\dagger} }$ is extremely small, the depth of the logarithmic layer would have been very thin, as expected from (4.23). This is not the case at least for the isolated solution shown in figure 4, where $\sigma _0\approx 0.437$ and $\beta _c^{{\dagger} }\approx 1.5$ yield the ratio $Y_c/Y_0\approx 1.66$. This number is still not too large, but we note that the isolated solution on which we have based the above discussion would not be the only possible isolated state of the canonical VWI array problem. The ratio also depends on $\sigma _0$. The main reason why that value is not so large for our solution is that the coherent structure possesses the character of the typical lower-branch solutions (see figure 4). However, in view of the substantial evidence of the existence of uniform momentum zones in near-wall turbulence (Meinhart & Adrian Reference Meinhart and Adrian1995; Adrian, Meinhart & Tomkins Reference Adrian, Meinhart and Tomkins2000), there may also be isolated solutions of the upper-branch type, and they would play more major roles in the flow dynamics. For such upper-branch-type solutions, as seen in figure 5, the momentum transfer $\sigma$ is much larger, and, hence, the logarithmic law will be valid for a wider range of $Y$.