2.1. Computing coherent (filtered) motions

We consider the plane Couette flow of a viscous fluid of constant density ${\it\rho}$ and kinematic viscosity ${\it\nu}$ between two parallel plates located at $y=\pm h$ . The streamwise, wall-normal and spanwise coordinates are denoted by $x$ , $y$ and $z$ respectively. The two plates move in opposite directions with velocity $(\pm U_{w},0,0)$ and the Reynolds number is defined as $Re=hU_{w}/{\it\nu}$ . We are mainly interested in the fully developed turbulent regime where averaged quantities are of interest. The mean velocity profile $U(y)$ and the root-mean-square (r.m.s.) velocity fluctuations $u_{rms}(y)$ , $v_{rms}(y)$ , $w_{rms}(y)$ are obtained by averaging the instantaneous fields over horizontal $x{-}z$ planes and in time (after the extinction of transients). The friction Reynolds number $Re_{{\it\tau}}=u_{{\it\tau}}h/{\it\nu}$ is based on the friction velocity $u_{{\it\tau}}=\sqrt{{\it\tau}_{w}/{\it\rho}}$ , where ${\it\tau}_{w}/{\it\rho}={\it\nu}\,\text{d}U/\text{d}y|_{w}$ . We denote by a $+$ superscript variables expressed in wall units, i.e. made dimensionless with respect to ${\it\nu}/u_{{\it\tau}}$ for lengths and $u_{{\it\tau}}$ for velocities. Lengths expressed in wall units can be obtained by multiplication by $Re_{{\it\tau}}$ of those expressed in terms of the outer length $h$ , so that, for instance, $y^{+}=Re_{{\it\tau}}y/h$ .

In this investigation, large-eddy simulations are used to study the dynamics of large- and very large-scale motions. The equations for the filtered motions are the usual ones (see e.g. Deardorff Reference Deardorff1970; Pope Reference Pope2000):

where the overbar denotes the filtering action, $\overline{{\bf\tau}}^{r}=\overline{{\bf\tau}}^{R}-\text{tr}(\overline{{\bf\tau}}^{R})\unicode[STIX]{x1D644}/3$ , with $\overline{{\it\tau}}_{ij}^{R}=\overline{u_{i}u_{j}}-\overline{u}_{i}\overline{u}_{j}$ and $\overline{q}=\overline{p}+\text{tr}(\overline{{\bf\tau}}^{R})/3$ . The anisotropic residual stress tensor $\overline{{\it\tau}}_{ij}$ is modelled by choosing an appropriate subgrid model in terms of eddy viscosity ${\it\nu}_{t}$ as $\overline{{\it\tau}}_{ij}^{r}=-2{\it\nu}_{t}\overline{\unicode[STIX]{x1D61A}}_{ij}$ , where $\overline{\unicode[STIX]{x1D61A}}_{ij}$ is the rate of strain tensor associated with the filtered velocity field. For the eddy viscosity, we choose the static (Smagorinsky Reference Smagorinsky1963) model:

where $\overline{\mathscr{S}}\equiv (2\overline{\unicode[STIX]{x1D61A}}_{ij}\overline{\unicode[STIX]{x1D61A}}_{ij})^{1/2}$ , $\overline{{\it\Delta}}=(\overline{{\it\Delta}}_{x}\overline{{\it\Delta}}_{y}\overline{{\it\Delta}}_{z})^{1/3}$ is the average length scale of the filter based on the mean grid spacing and $C_{s}$ is the Smagorinsky constant. In the following we will use as a reference value $C_{s}=0.05$ , like Hwang & Cossu (Reference Hwang and Cossu2010c , Reference Hwang and Cossu2011). This value is known to provide the best turbulence statistics compared to those of direct numerical simulations (i.e. the best performance for a posteriori tests as discussed by Härtel & Kleiser Reference Härtel and Kleiser1998). To avoid non-zero residual velocity and shear stress at the wall we use the wall (damping) function $D=1-\exp (-(y^{+}/A^{+})^{2})$ proposed by Kim & Menon (Reference Kim and Menon1999) with $A^{+}=25$ .

The use of the static Smagorinsky model ensures that the residual motions cannot transfer energy to the filtered motions, i.e. there is no ‘backscatter’ of energy. This is essential in our approach of determining if large-scale motions can be sustained in the absence of forcing by smaller-scale motions. In particular, in the following we will use the same technique as used by Hwang & Cossu (Reference Hwang and Cossu2010c , Reference Hwang and Cossu2011) to quench small-scale motions, and investigate if the large-scale motions survive despite this quenching. The technique is simply based on increasing the Smagorinsky constant $C_{s}$ , which is equivalent to an increase of the ‘Smagorinsky mixing length’ $l_{0}=C_{s}{\it\Delta}$ , as shown by Mason & Callen (Reference Mason and Callen1986). As $C_{s}$ is increased, therefore, an increasing range of small-scale motions become inactive as they are modelled with an increasingly large (positive) eddy viscosity. In ‘over-damped’ simulations the value of $C_{s}$ will be increased to values in a typical range $C_{s}\approx 0.1{-}0.18$ . Two warnings must be issued here in order to avoid possible misunderstanding of the ‘over-damped LES’ technique. First, it must be recalled that the aim of using ‘over-damped LES’ is, of course, not to provide a quantitatively accurate representation of large-scale motions (which can be done using state-of-the-art models for the LES) but to show that large-scale motions can survive in the absence of the near-wall process. For e.g.  $C_{s}=0.1$ or $C_{s}=0.14$ , the near-wall region is (and must be) inaccurate because the buffer-layer cycle has been shut down. However, as already found by Hwang & Cossu (Reference Hwang and Cossu2010c , Reference Hwang and Cossu2011) and Hwang (Reference Hwang2015) and in the following, the over-damped flow displays key features of the ‘real’ large-scale motions despite the inaccuracy of the near-wall region. Secondly, it must be noted that increasing $C_{s}$ is not equivalent to a reduction of the ‘effective’ Reynolds number in the simulations because in the Smagorinsky model the viscosity depends on the local rate of strain (while it is constant in the Navier–Stokes equations). In practice it is observed that the friction Reynolds number of over-damped solutions is only weakly affected by the increased $C_{s}$ .

Table 1.

Numerical domains and discretization parameters used in the present study. The VLSM-Box  and the LSM-Box  are used to perform the numerical experiments respectively described in § 3 and appendix A.

Large-eddy simulations are performed with the code diablo (Bewley, Moin & Temam Reference Bewley, Moin and Temam2001) which implements the fractional-step method based on a semi-implicit time integration scheme and a mixed finite-difference and Fourier discretization in space. The computational domain, which extends from $0$ to $L_{x}$ and from $0$ to $L_{z}$ in the streamwise and spanwise directions respectively and from $-h$ to $h$ in the wall-normal direction, is discretized with $N_{x}\times N_{y}\times N_{z}$ points in respectively the streamwise ( $x$ ), wall-normal ( $y$ ) and spanwise ( $z$ ) directions. Grid stretching is applied in the wall-normal direction in order to refine the grid near the wall. No-slip boundary conditions are applied at the walls and periodic boundary conditions at other boundaries. The numerical domains and the associated discretization parameters used in the following are listed in table 1.

2.2. Computing three-dimensional nonlinear steady solutions

In the second part of this study (§ 4), we compute steady three-dimensional nonlinear solutions of the filtered equations by using a modified version of the code peanuts (Herault et al. Reference Herault, Rincon, Cossu, Lesur, Ogilvie and Longaretti2011; Riols et al. Reference Riols, Rincon, Cossu, Lesur, Longaretti, Ogilvie and Herault2013), which provides an implementation of a modified Newton–Krylov method through the PETSc toolkit (Balay et al. Reference Balay, Brown, Buschelman, Gropp, Kaushik, Knepley, McInnes, Smith and Zhang2011). An important feature of the algorithm is that it is matrix-free and therefore only relies on (many) calls to the large-eddy simulation solver diablo.

For a Newton-iteration-based method to converge, it is crucial to have a good initial guess for the solution. One possibility to obtain such an initial guess is to compute the edge state of the system, i.e. the relative attractor on the surface which (in phase space) is the boundary of the basin of attraction of self-sustained motions. In cases where the dimension of this boundary is $N-1$ , $N$ being the dimension of the phase space, it is possible to constrain the solution to remain in a neighbourhood of the boundary by a one-parameter bisection technique (see e.g. Itano & Toh Reference Itano and Toh2001; Toh & Itano Reference Toh and Itano2003), sometimes labelled ‘edge tracking’. Edge tracking has indeed been successfully used by e.g. Viswanath (Reference Viswanath2007) and Schneider et al. (Reference Schneider, Gibson, Lagha, De Lillo and Eckhardt2008) to compute lower-branch non-trivial steady solutions of the Navier–Stokes equations (unfiltered motions) in plane Couette flow. In the present study, edge tracking is implemented by selecting as initial condition of the LES the mean flow profile $\boldsymbol{U}$ to which we add a coherent, three-dimensional large-scale perturbation $\boldsymbol{u}_{0}$ of amplitude $A_{0}$ : $\boldsymbol{u}_{0}=\boldsymbol{U}+A_{0}\boldsymbol{u}_{0}$ . A standard bisection is then performed to compute the threshold value of $A_{0}$ which separates solutions evolving to active large-scale motions from the laminar solution, as will be further shown in § 4.

3.1. Reference large-eddy simulations

As a first step, we assess the ability of our large-eddy simulations to reproduce the main features of near-wall, large-scale and very large-scale motions in turbulent Couette flow by performing a reference simulation with the Smagorinsky constant fixed to its reference value $C_{s}=0.05$ . The simulations are performed in a very large domain $L_{x}\times L_{z}=132h\times 18.85h$ (labelled VLSM-Box  in table 1), the size of which is comparable to the largest ones used by Tsukahara et al. (Reference Tsukahara, Kawamura and Shingai2006) in their direct numerical simulations (DNS) of turbulent Couette flow. Given the typical streamwise and spanwise size of the very large-scale structures ( ${\it\lambda}_{x}\simeq 42h{-}64h$ and ${\it\lambda}_{z}\simeq 4h{-}5h$ , according to Tsukahara et al. Reference Tsukahara, Kawamura and Shingai2006), the domain allows for two or three very large-scale structures in the streamwise direction and three or four in the spanwise direction. The grid spacings in the streamwise and the spanwise direction are set to ${\rm\Delta}x^{+}=26.3$ and ${\rm\Delta}z^{+}=18.2$ , respectively (before dealiasing), while the wall-normal grid spacing is chosen such that the minimum and maximum spacings respectively become ${\rm\Delta}y_{min}^{+}=1.27$ and ${\rm\Delta}y_{max}^{+}=9.7$ . We note that while grid spacing is small enough to reasonably resolve the near-wall structures (Härtel & Kleiser Reference Härtel and Kleiser1998), it is much coarser than the one use by Tsukahara et al. (Reference Tsukahara, Kawamura and Shingai2006) (who use ${\rm\Delta}x^{+}=8$ and ${\rm\Delta}z^{+}=6$ and ${\rm\Delta}y^{+}=0.2{-}5.7$ ).

Figure 1.

Comparison of the reference LES ( $C_{s}=0.05$ , lines) results obtained in the VLSM-Box  to the data (symbols) of Tsukahara et al. (Reference Tsukahara, Kawamura and Shingai2006) obtained by DNS in domains of similar extent at the same Reynolds number $Re=2150$ ( $Re_{{\it\tau}}=128$ ). The mean flow profiles $U^{+}$ are compared in (a), while in (b) the r.m.s. velocity profiles $u_{rms}^{+}$ , $v_{rms}^{+}$ and $w_{rms}^{+}$ are reported.

The mean velocity and r.m.s. velocity profiles of the reference LES are in reasonably good agreement with the DNS data of Tsukahara et al. (Reference Tsukahara, Kawamura and Shingai2006), as shown in figure 1. In the instantaneous field visualized in figure 2(a), very long streaky structures with the spanwise spacing roughly at ${\it\lambda}_{z}\approx 4h{-}5h$ are clearly visible throughout the entire domain on top of small-scale background turbulence, in agreement with previous simulations (see e.g. Komminaho et al. Reference Komminaho, Lundbladh and Johansson1996; Tsukahara et al. Reference Tsukahara, Kawamura and Shingai2006). The large-scale streaks also show the characteristic ‘ramp’ structure of the low-speed regions with angles ${\approx}10^{\circ }{-}12^{\circ }$ to the wall (see figure 3) already observed e.g. in turbulent boundary layers (see e.g. figure 8 of Dennis & Nickels Reference Dennis and Nickels2011).

Figure 2.

Iso-surface ( $u^{+}=-2.5$ ) of the instantaneous negative streamwise velocity fluctuations (low-speed streaks) on the lower wall at $Re=2150$ : (a)  $C_{s}=0.05$ (reference case); (b) $C_{s}=0.10$ ; (c) $C_{s}=0.14$ ; (d) $C_{s}=0.18$ . The visualization illustrates that the progressive quenching of small-scale motions for increasing values of $C_{s}$ does not affect the self-sustainment of large- and very large-scale streaky motions.

One-dimensional spectra of the streamwise velocity shown in figure 4 are also in good agreement with the DNS ones reported by Tsukahara et al. (Reference Tsukahara, Kawamura and Shingai2006). The spanwise spectra reveal two well-separated peaks depending on the wall-normal location (figure 4 a). The first peak appears at ${\it\lambda}_{z}^{+}\simeq 100$ in the near-wall region ( $y^{+}\leqslant 30$ ), while the second peak is at ${\it\lambda}_{z}\approx 4h{-}5h$ in the outer region ( $y\geqslant 0.2h$ ); the former corresponds to near-wall streaks whereas the latter corresponds to large- and very large-scale streaky structures. The streamwise spectra have a peak at ${\it\lambda}_{x}^{+}\simeq 700$ in the near-wall region, which is representative of near-wall streaks. In the outer region, the spectra exhibit two peaks, one at ${\it\lambda}_{x}\approx 10h$ and the other at ${\it\lambda}_{x}\approx 42h{-}64h$ . The two peaks at ${\it\lambda}_{x}\approx 10h$ and ${\it\lambda}_{x}\approx 42h{-}64h$ , which are also in good agreement with the DNS results, correspond respectively to large- (LSM) and very large-scale motions (VLSM).

A terminological note is necessary here. We choose to use the LSM and VLSM terminology in Couette flow also, where it does not appear to be currently standard, because we think that the mechanisms underlying their dynamics are common to other wall-bounded turbulent shear flows. Also, we will denote by VLSM only structures that have a very large (e.g.  ${\approx}40{-}60h$ ) but finite streamwise wavelength. These motions can be detected only in very long domains. Structures with infinite wavelength in the streamwise direction will not be referred to as VLSM because streamwise-uniform motions cannot self-sustain (see e.g. Waleffe Reference Waleffe1995) and therefore are not of interest here.

3.2. Over-damped simulations: self-sustained large-scale motions

Having validated the ability of the reference LES to quantitatively reproduce the main features of the turbulent Couette flow obtained by DNS in previous studies, we now examine if the motions at large scales persist even in the absence of smaller-scale active motions in the near-wall region. As in Hwang & Cossu (Reference Hwang and Cossu2010b ), we therefore gradually increase the Smagorinsky constant $C_{s}$ , which determines the filter width of the LES, to remove the small-scale active structures from the near-wall region. It should be noted that the increase of $C_{s}$ does not significantly affect the friction Reynolds number $Re_{{\it\tau}}=128$ of the reference simulation ( $C_{s}=0.05$ ): it is indeed found that $Re_{{\it\tau}}=119$ for $C_{s}=0.10$ , $Re_{{\it\tau}}=117$ for $C_{s}=0.14$ , $Re_{{\it\tau}}=119$ for $C_{s}=0.18$ .

The flow fields visualized in figure 2 show that the small-scale turbulent motions are gradually damped out as $C_{s}$ is increased. For $C_{s}=0.14$ (figure 2 c), only the motions at large and very large scale are left in the flow field, suggesting the outer motions are likely to be sustained by themselves. The survival of the large-scale motions persists at $C_{s}=0.18$ (see figure 2 d), although they are eventually quenched for $C_{s}=0.25$ (not shown). An examination of figure 3 further confirms that when small-scale motions are smoothed out the surviving coherent large- and very large-scale structures are reasonably similar to the ‘natural’ ones obtained in the reference case. In particular, they show the typical ‘ramp’ structure of the low-speed regions with angles ${\approx}10^{\circ }{-}12^{\circ }$ to the wall.

Figure 3.

Side view of the low-velocity fluctuations in the VLSM-Box  for the reference case $C_{s}=0.05$ (a) and the over-damped $C_{s}=0.14$ (b) at $Re=2150$ . The data are the same as in figure 2(a,c). In order to highlight the ramp structure of large-scale motions, only a limited streamwise part of the computational domain is represented and the levels $u^{+}=-3.9$ for $C_{s}=0.05$ and $u^{+}=-5.3$ for $C_{s}=0.14$ are used. The $12^{\circ }$ angle is measured from the wall.

Figure 4.

Spanwise premultiplied power spectra $k_{z}E_{uu}({\it\lambda}_{z})$ (a,c,e,g) and streamwise premultiplied power spectra $k_{x}E_{uu}({\it\lambda}_{x})$ (b,d,f,h) for respectively the reference simulation with $C_{s}=0.05$ (a,b) and for the cases $C_{s}=0.10$ (c,d), $C_{s}=0.14$ (e,f) and $C_{s}=0.18$ (g,h). The premultiplied spectra are extracted in the inner layer ( $y^{+}=11,17,28$ based on $Re_{{\it\tau}}$ of the reference simulation, dashed lines) and in the outer layer (at distances from the bottom wall of $y+h=0.40h,0.57h,0.78h,h$ solid lines).

The premultiplied spectra reported in figure 4 provide a quantitative confirmation of this scenario. The original ( $C_{s}=0.05$ ) peak location of the near-wall motions at ${\it\lambda}_{z}^{+}\simeq 100$ in the spanwise spectra (figure 4 a) is initially shifted to a larger spanwise wavelength for $C_{s}=0.10$ (figure 4 c) and completely suppressed for $C_{s}=0.14$ (figure 4 e) and higher (figure 4 g), where only the footprint of large-scale motions appears in the near-wall region. Similar features are observed in the streamwise spectra where the original peak location at ${\it\lambda}_{x}^{+}\simeq 700$ in the near-wall region (figure 4 b) is shifted to higher wavelengths up to its complete suppression (figure 4 d,f,h).

In contrast, large-scale motions (LSM) appear to be fairly robust to the increase of $C_{s}$ . The peak at ${\it\lambda}_{z}\simeq 4.5h$ in the spanwise spectra robustly remains at the same location for all increasing values of $C_{s}$ (figure 4 a,c–e). In the streamwise spectra, for $C_{s}=0.10$ the LSM peak at ${\it\lambda}_{x}\simeq 10h$ remains unchanged, while the VLSM peak at ${\it\lambda}_{x}\simeq 42h{-}64h$ is blurred. For $C_{s}=0.14$ , the original location of the first peak ( ${\it\lambda}_{x}\simeq 10h$ ) has drifted to a slightly longer streamwise wavelength ( ${\it\lambda}_{x}\simeq 14h$ ) (figure 4 f), and it appears that this feature makes the presence of the VLSM structures even more blurred. Similar drifts of the streamwise wavelength of LSM for relatively large $C_{s}$ have also been observed in the case of the pressure-driven channel flow (Hwang & Cossu Reference Hwang and Cossu2010b ). When $C_{s}=0.18$ , the artificial over-damping begins to have a non-negligible effect on the dynamics of large-scale motions in the near-wall region where energy accumulates at the longest available streamwise wavelengths. This is why in figure 3 and in the following we only discuss the dynamics of structures corresponding to $C_{s}\lesssim 0.14$ .

These results suggest that motions with typical spanwise spacings ${\it\lambda}_{z}\simeq 4h{-}5h$ are sustained by a process which is not forced by motions at smaller scales. The relation between large-scale (LSM) and very large-scale motions (VLSM) is left partially unsolved by these simulations because they do not exclude the self-sustained process being associated with the VSLM very long wavelengths. The reference and over-damped simulations have been therefore repeated in a smaller domain (the LSM-Box  listed in table 1) having the typical dimensions of large-scale motions. The results of this additional series of LES, discussed in appendix A, are very similar to the ones discussed above and show that large-scale motions are also self-sustaining in the absence of potentially active motions of even larger scale when near-wall active processes are quenched. There is, therefore, an active process at precisely this scale ( ${\it\lambda}_{z}\approx 4h{-}5h$ , ${\it\lambda}_{x}\approx 10h{-}12h$ ). The nature of this process is further explored in the next section.

4.1. Coherent steady solutions at $Re=750$

We begin by looking for exact solutions near the relatively low Reynolds number $Re=750$ (corresponding to $Re_{{\it\tau}}=52$ ), which roughly corresponds to twice the value of the transitional Reynolds number. The box discretization details are reported in table 1. It has been verified (Rawat Reference Rawat2014) that reference large-eddy simulations at the considered values $Re=750$ with $C_{s}=0.05$ in the LSM-Box compare well to the DNS results obtained in the same box and are relatively similar to those obtained in larger domains at the same Reynolds number (see also Komminaho et al. Reference Komminaho, Lundbladh and Johansson1996; Tsukahara et al. Reference Tsukahara, Kawamura and Shingai2006, for a discussion of the box size effects). The premultiplied spectra also show a reasonable agreement with the DNS results for $C_{s}=0.05$ and only large-scale motions survive when the near-wall cycle is artificially quenched by increasing the Smagorinsky constant $C_{s}$ (Rawat Reference Rawat2014).

Figure 5.

Edge tracking in the over-damped LES at $Re=750$ ( $Re_{{\it\tau}}=52$ ) and $C_{s}=0.14$ . Two initial conditions (circle symbol) lying almost on the edge of chaos surface initially remain near the surface while being attracted to the steady solution edge state (triangle symbol). The solution lying on the upper side of the edge of chaos (green line, dotted) finally evolves toward a state with irregular large-scale active motions, while the other one (solid, red line) finally relaxes to a state with no large-scale active motions. The trajectories are shown in both the $t{-}\Vert v\Vert$ (a) and in the $\Vert u\Vert {-}\Vert v\Vert$ plane (b), displaying the streaks–vortices dynamics.

Figure 6.

Continuation diagram of the steady solutions in the Smagorinsky constant $C_{s}$ at $Re=750$ ( $Re_{{\it\tau}}=52$ ).  $S={\it\tau}_{w}/({\it\mu}U_{w}/h)$ is the wall shear rate of the solution normalized by its laminar value (where ${\it\mu}$ is the viscosity of the fluid). The symbols denote the solutions found at $C_{s}=0$ (Navier–Stokes, squares, with $S=1.37$ on the lower branch and $S=4.59$ on the upper branch), $C_{s}=0.05$ (reference LES, circles) and $C_{s}=0.14$ (over-damped LES where the edge tracking has been performed, triangles). Empty symbols denote lower-branch solutions, filled symbols denote upper-branch solutions.

As discussed in § 2.2, an edge-tracking technique can be effective for computing lower-branch steady solutions in plane Couette flow. We therefore apply the edge-tracking technique to the over-damped LES simulations with $C_{s}=0.14$ at $Re=750$ in the LSM-Box where the near-wall cycle is quenched and only large-scale motions survive. The initial condition given to the LES consists of the mean turbulent flow, to which is added a pair of streamwise-uniform counter-rotating rolls along with a sinuous perturbation of the spanwise velocity, which is similar to that used by Toh & Itano (Reference Toh and Itano2003) and Cossu et al. (Reference Cossu, Brandt, Bagheri and Henningson2011): $\boldsymbol{u}_{0}=\{U(y),0,0\}+A_{0}\{0,\partial {\it\psi}_{0}/\partial z,-\partial {\it\psi}_{0}/\partial y\}+0.1\,A_{0}\{0,0,w_{sin}\}$ where ${\it\psi}_{0}(y,z)=(1-y^{2})\sin (2{\rm\pi}z/L_{z})$ and $w_{sin}(x,y)=(1-y^{2})\sin (2{\rm\pi}x/L_{x}).$ The bisection is performed by adjusting the amplitude $A_{0}$ of the perturbations. The results of the edge tracking, displayed in figure 5, show that the edge-state solution has constant energy. The initial transient clearly shows the coherent transient growth of the initial condition (slow decrease of the wall-normal velocity perturbations and large increase of the streamwise-streak norm $\Vert u\Vert$ ) leading to the edge state.

The edge-tracking solution is then used as an initial guess for Newton–Krylov iterations, performed with the peanuts  code. The solver converges to a steady solution (i.e. with zero phase speed, up to the precision of the Newton solver computation) of the (LES) equations for the filtered motions. The converged solution can be continued in $C_{s}$ using a Newton-based continuation. When continued to higher $C_{s}$ , the solution encounters a saddle–node bifurcation from which the lower branch and an associated upper branch originate, as shown in figure 6. As the solutions obtained for $C_{s}\gtrsim 0.18$ lose much of their physical relevance because the LSM become themselves over-damped, the high- $C_{s}$ continuation should be considered only as a method (just as an homotopy) to access upper-branch solutions. Also, it must be noted that the high- $C_{s}$ bifurcation is not simply reproducing the low- $Re$ bifurcation of Navier–Stokes solutions. Here the stress tensor depends nonlinearly on the rate of strain, which models turbulent dissipation. As a consequence, for example, the value of the shear parameter $S$ at the saddle–node bifurcation is more than double the value of $S$ at the low-Reynolds number saddle–node bifurcation of Navier–Stokes solutions.

Even more interestingly, both upper- and lower-branch solutions can be continued to lower values of $C_{s}$ and in particular to $C_{s}=0.05$ , the ‘reference’ value used to match DNS solutions, and even to $C_{s}=0$ , which corresponds to solutions of the (unfiltered) Navier–Stokes equations. The mean and r.m.s. velocity profiles corresponding to $C_{s}=0.14$ , $C_{s}=0.05$ and $C_{s}=0$ , reported in figure 7, are very similar to one another and are reminiscent of the NBCW solutions of the Navier–Stokes equations. For the lower-branch solutions, most of the energy is concentrated in the quasi-streamwise streaky motions (mainly the $u_{rms}$ component) which are forced by low-amplitude quasi-streamwise vortices (mainly the $v_{rms}$ component), while for the upper-branch solutions, quasi-streamwise streaks and vortices have comparable energy.

Figure 7.

Lower- (a–c) and upper- (d–f) branch steady solutions in the LSM-Box at $Re=750$ ( $Re_{{\it\tau}}=52$ ): mean flow (a,d), streamwise (b,e) and cross-stream (c,f) r.m.s. velocity profiles corresponding to $C_{s}=0.14$ (over-damped LES), $C_{s}=0.05$ (reference LES) and $C_{s}=0$ (Navier–Stokes steady solution).

Figure 8.

Continuation diagram of steady solutions in the $C_{s}{-}Re$ parameters plane. The solid (red) line continuation in $C_{s}$ is that shown in figure 6 (see its caption for the meaning of the symbols). The continuation in Reynolds number at $C_{s}=0$ (black line, long-dashed-dotted) shows that the solutions continued from the turbulent case (square symbols) belong to the NBCW branch of Navier–Stokes solutions originating in a saddle–node bifurcation at $Re=127$ (in this box with $L_{x}=10.89$ and $L_{z}=5.46$ ). It has been possible to directly continue the $C_{s}=0.05$ solution from $Re=750$ to $Re=1164$ (green line connecting to the solution denoted by a filled triangle). Nonlinear steady solutions at higher Reynolds numbers are found by first continuing the solutions obtained for $C_{s}=0.1$ to higher $Re$ (blue line, short-dashed-dotted) and then reducing $C_{s}$ to $C_{s}=0.05$ (filled diamond symbols), except for the $Re=2187$ solution which could only be continued to $C_{s}=0.06$ .

Figure 9.

Mean flow (a) and r.m.s. velocity profiles (b) corresponding to the upper-branch coherent steady solutions obtained for increasing Reynolds numbers, LSM-Box  ( $C_{s}=0.06$ for $Re=2187$ and $C_{s}=0.05$ otherwise).

Figure 10.

Visualization of the upper-branch large-scale coherent steady solutions of the LES equations obtained with the reference value $C_{s}=0.05$ and for increasing values of the Reynolds number: (a,b) $Re=750$ (corresponding to $Re_{{\it\tau}}=52$ ), (c,d) $Re=1600$ ( $Re_{{\it\tau}}=99$ ), (e,f) $Re=2187$ ( $Re_{{\it\tau}}=127$ , with $C_{s}=0.06$ ). (a,c,e) The streaks and quasi-streamwise vortices: the green surface corresponds to a streamwise velocity value $50\,\%$ of its maximum value while the blue and green surfaces correspond streamwise vorticity values equal to $\pm 70\,\%$ of the maximum. (b,d,e) The streaks and the relative eddy viscosity associated with the filtered small-scale motions. The green surface is the same as in (a,c,e), while the violet surfaces correspond to ${\it\nu}_{t}/{\it\nu}=40\,\%$ and the yellow one to ${\it\nu}_{t}/{\it\nu}=10\,\%$ .

4.2. Continuation to higher Reynolds numbers

The steady solutions obtained at $Re=750$ are very similar to the NBCW solutions. We have, therefore, continued the solution corresponding to $C_{s}=0$ (i.e. the unfiltered Navier–Stokes solution) from $Re=750$ to $Re=400$ , where it appears to be identical to the one computed by Waleffe (Reference Waleffe2003) in the same box at the same $Re$ . Continuation to lower $Re$ does indeed show that the upper and the lower branches are connected by the saddle–node bifurcation at $Re=127$ (Waleffe Reference Waleffe2003) as shown in figure 8.

At $Re=750$ ( $Re_{{\it\tau}}\approx 52$ ), the near-wall and large-scale motions scales are not yet separated enough to allow definitive conclusions on the nature of the steady solutions of the filtered equations. It is, therefore, potentially unclear that the computed steady solutions are associated with the self-sustained large-scale motions discussed in § 3 and not with near-wall motions. In order to clarify this issue the coherent steady solutions of the over-damped LES have been continued to higher Reynolds numbers, keeping the box dimensions constants in outer units ( $L_{z}=5.5h$ and $L_{x}=10.9h$ ).

We concentrate on upper-branch solutions because they are most relevant to the turbulent dynamics (lower-branch solutions are related to the transition problem which is probably less relevant in the turbulent regime). A continuation to higher $Re$ fails to converge for Reynolds numbers larger than $Re\approx 1164$ ( $Re_{{\it\tau}}\approx 74$ ) for the reference $C_{s}=0.05$ (but also for $C_{s}=0$ ). To circumvent this difficulty, an alternative path has been followed in the parameter space by first computing solutions at $C_{s}=0.1$ , as shown in figure 8. The $C_{s}=0.1$ steady solutions can be continued up to $Re\approx 2500$ . At selected Reynolds numbers ranging up to $Re=2187$ , the $C_{s}=0.1$ solutions can then be tracked down to $C_{s}=0.05$ , except for $Re=2187$ , corresponding to $Re_{{\it\tau}}\approx 127$ , where the solution cannot be continued below $C_{s}=0.06$ .

The mean and r.m.s. velocity profiles of the ‘reference’ ( $C_{s}=0.05$ and $C_{s}=0.06$ for the highest $Re$ ) solutions obtained at these larger Reynolds numbers are reported in figure 9. The change of the wall-normal and spanwise r.m.s. velocity profiles with the Reynolds number is minor. The peak streamwise r.m.s. velocity, associated with the large-scale streaks, increases with the Reynolds number and is attained nearer to the wall. All the r.m.s. velocity components remain relatively large in the bulk of the flow. These features are consistent with the observed behaviour of coherent large-scale motions (see e.g. Avsarkisov et al. Reference Avsarkisov, Hoyas, Oberlack and Garcia-Galache2014; Pirozzoli, Bernardini & Orlandi Reference Pirozzoli, Bernardini and Orlandi2014, for the most recent DNS results). This similarity is confirmed by the analysis of the flow fields associated with the large-scale coherent steady solutions reported in figure 10. This figure also illustrates how the eddy viscosity associated with the residual (small-scale) motions increases with the Reynolds number and is located mainly on the flanks of the low-speed streak. The absolute levels of eddy viscosity are not huge even at $Re=2187$ , where the maximum of ${\it\nu}_{t}$ does not exceed $90\,\%$ of the molecular viscosity for those solutions but, contrary to the Navier–Stokes case, the eddy viscosity is not spatially uniform and depends on the local properties of the large-scale flow.

Having observed a clear relation between NBCW and coherent large-scale motions, it remains to verify if these solutions are also related to buffer-layer structures. This question is addressed in appendix B where it is shown that NBCW solutions do not converge to buffer-layer structures if continued to higher $Re$ in minimal flow units with $L_{x}^{+}\times L_{z}^{+}=250\times 100$ .