On the relevance of Reynolds stresses in resolvent analyses of turbulent wall-bounded flows

TL;DR

This study demonstrates that including Reynolds stresses via eddy-viscosity in resolvent analysis significantly improves the estimation of turbulent coherent structures in wall-bounded flows, especially at buffer-layer scales.

Contribution

It shows the importance of incorporating Reynolds stresses in resolvent models to accurately predict turbulent flow structures and spectra.

Findings

Reynolds stresses improve resolvent spectral estimates.

Eddy-viscosity inclusion enhances peak frequency and location predictions.

Matching forcing spectra further improves model accuracy.

Abstract

The ability of linear stochastic response analysis to estimate coherent motions is investigated in turbulent channel flow at friction Reynolds number Re$_\tau$ = 1007. The analysis is performed for spatial scales characteristic of buffer-layer and large-scale motions by separating the contributions of different temporal frequencies. Good agreement between the measured spatio-temporal power spectral densities and those estimated by means of the resolvent is found when the effect of turbulent Reynolds stresses, modelled with an eddy-viscosity associated to the turbulent mean flow, is included in the resolvent operator. The agreement is further improved when the flat forcing power spectrum (white noise) is replaced with a power spectrum matching the measures. Such a good agreement is not observed when the eddy-viscosity terms are not included in the resolvent operator. In this case, the estimation based on the resolvent is unable to select the right peak frequency and wall-normal location of buffer-layer motions. Similar results are found when comparing truncated expansions of measured streamwise velocity power spectral densities based on a spectral proper orthogonal decomposition to those obtained with optimal resolvent modes.

Figures

Equations

$$\frac{\partial \mathbf{u}}{\partial t}+\nabla \mathbf{u}\cdot \mathbf{U}+\nabla \mathbf{U}\cdot \mathbf{u}=-\frac{1}{\rho }\nabla p+\nabla \cdot \left[{\nu }_T\left(\nabla \mathbf{u}+\nabla {\mathbf{u}}^T\right)\right]+\mathbf{f},$$

$$\nu_{t}={\nu\over 2}\big{\{}1+{{\kappa^{2}{\textit{Re}_{\tau}}^{2}}\over 9}(1-y^{2})^{2}(1+2y^{2})^{2}(1-e^{y^{+}/A})^{2}\big{\}}^{1/2}-{\nu\over 2},$$

$$\frac{\partial \hat{\mathbf{q}}}{\partial t}=\mathbf{A}\hat{\mathbf{q}}+\mathbf{B}\hat{\mathbf{f}},$$

$$\tilde{\mathbf{u}}$$

$$\mathbf{R}(\xi ,y,y^{\prime },\zeta ,\tau )=⟨\mathbf{u}(x,y,z,t){\mathbf{u}}^{\ast }(x+\xi ,y^{\prime },z+\zeta ,t+\tau )⟩,$$

$$\mathbf{S}(\alpha ,y,y^{\prime },\beta ,\omega )=\frac{1}{(2\pi {)}^3}{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }\mathbf{R}(\xi ,y,y^{\prime },\zeta ,\tau )e^{i(\alpha \xi +\beta \zeta -\omega \tau )}d\xi d\zeta d\tau .$$

$$\mathbf{S}=⟨\tilde{\mathbf{u}}(\alpha ,y,\beta ,\omega ){\tilde{\mathbf{u}}}^{\ast }(\alpha ,y^{\prime },\beta ,\omega )⟩;\mathbf{P}=⟨\tilde{\mathbf{f}}(\alpha ,y,\beta ,\omega ){\tilde{\mathbf{f}}}^{\ast }(\alpha ,y^{\prime },\beta ,\omega )⟩.$$

$${\mathbf{S}}^{(res)}=p\mathbf{H}{\mathbf{H}}^{\ast }$$

$$\displaystyle\mathbf{A}=\left[\begin{array}[]{cc}\Delta^{-1}\mathcal{L_{OS}}&0\\ -i\beta U^{\prime}&\mathcal{L_{SQ}}\\ \end{array}\right],~{}~{}~{}\mathbf{B}=\left[\begin{array}[]{ccc}-i\alpha\Delta^{-1}\mathcal{D}&-k^{2}\Delta^{-1}&-i\beta\Delta^{-1}\mathcal{D}\\ i\beta&0&-i\alpha\end{array}\right]$$

$${\mathcal{L}}_{\mathcal{O}\mathcal{S}}$$

$${\mathcal{L}}_{\mathcal{S}\mathcal{Q}}$$

$$\displaystyle\mathbf{C}={1\over k^{2}}\left[\begin{array}[]{ccc}i\alpha\mathcal{D}&-i\beta\\ k^{2}&0\\ i\beta\mathcal{D}&i\alpha\\ \end{array}\right]~{}~{}~{}~{}~{}~{};~{}~{}~{}~{}~{}~{}~{}\mathbf{D}=\left[\begin{array}[]{ccc}0&1&0\\ i\beta&0&-i\alpha\\ \end{array}\right].$$

Full text

1 Introduction

Most of the fluctuating energy in wall-bounded turbulent shear flows resides in highly correlated streamwise velocity fluctuations associated with coherent streamwise streaks, i.e. spanwise alternated high- and low-velocity regions elongated in the streamwise direction. The existence of these structures is known since the visualisations of Kline et al. (1967) which revealed that the near-wall region of turbulent boundary layers is populated by streaks with average spanwise spacing $\lambda_{z}^{+}\approx 100$ in the buffer layer and higher farther from the wall. Streaky motions with larger scales have been observed in the logarithmic and the outer region where ‘large-scale’ (LSM) and ‘very large scale’ motions (VLSM) have typical spanwise spacings $\lambda_{z}\approx\delta-1.5\delta$ (where $\delta$ is the outer length scale) and streamwise size of $\lambda_{x}\approx 2\lambda_{z}$ (Corrsin & Kistler, 1954; Kovasznay et al., 1970) and $\lambda_{x}\simeq O(10\lambda_{z})$ (Komminaho et al., 1996; Kim & Adrian, 1999; Hutchins & Marusic, 2007), respectively.

The ubiquity of streaks in transitional and turbulent flows has been related to the ‘lift-up’ effect where high-energy streamwise streaks are induced by low-energy quasi-streamwise vortices immersed in a shear flow (Moffatt, 1967; Ellingsen & Palm, 1975; Landahl, 1980) via algebraic growth (Ellingsen & Palm, 1975; Gustavsson, 1991). The large energy amplifications associated with this mechanism have been related to the non-normality of the linearised Navier-Stokes operator (e.g Böberg & Brosa, 1988; Reddy & Henningson, 1993; Trefethen et al., 1993) and much attention has been given to the computation of the largest (optimal) energy amplifications supported by laminar basic flows (Gustavsson, 1991; Butler & Farrell, 1992; Reddy & Henningson, 1993; Schmid & Henningson, 2001). Böberg & Brosa (1988) suggested that the subcritical onset of turbulence in wall-bounded flows can be attributed to a mechanism where the low-energy quasi-streamwise vortices leading to the amplification of high-energy streamwise streaks are then regenerated by nonlinear effects related to the breakdown of the streaks. A similar ‘self-sustained mechanism’ has been invoked by Hamilton et al. (1995) to explain the dynamics of buffer-layer streaks in the turbulent regime with Waleffe (1995) and Reddy et al. (1998) attributing the breakdown of the streaks to a modal secondary instability and Schoppa & Hussain (2002) to a secondary non-modal energy amplification. Hwang & Cossu (2010c, 2011) have found evidence that similar coherent self-sustained processes sustain all active streaky motions, with scales ranging from those of buffer-layer streaks to those of large-scale motions.

While a relative consensus begins to emerge about the existence of self-sustained processes, there are different interpretations of the mechanisms involved in these processes in turbulent flows, and, in particular, of the nature of the linear operator intervening in the non-normal amplification mechanism. In a first approach, pursued by Malkus (1956); Butler & Farrell (1993); Farrell & Ioannou (1993a); McKeon & Sharma (2010) and many others, the Navier-Stokes equations are rewritten in terms of perturbations to the turbulent mean velocity; the instantaneous and averaged (Reynolds stresses) perturbation nonlinear terms are accounted for as an external input which forces the response via the linear operator. In the second approach, initiated by Reynolds & Tiederman (1967); Reynolds & Hussain (1972) and pursued by Bottaro et al. (2006); del Álamo & Jiménez (2006); Cossu et al. (2009); Pujals et al. (2009); Hwang & Cossu (2010a, b) among others, the ‘incoherent’ part of turbulent Reynolds stresses is included in the linear operator by means of a $\nu_{t}$ eddy-viscosity modelling. We will refer to the latter approach as ‘$\nu_{t}$-model’ and to the former as ‘$\nu$-model’.

An important number of investigations have dealt with the computation of optimal energy amplifications and of the associated optimal inputs and outputs for the linear initial-value problem and the response to harmonic and stochastic forcing. While a review of this research effort and of the respective merits attributed to the $\nu$ and $\nu_{t}$ formulations is beyond the scope of this paper (we refer the reader to McKeon, 2017; Cossu & Hwang, 2017, for a summary of recent results), it appears that most of the comparisons between linear models predictions and real turbulent flows are either qualitative or concern wall-normal integrated energy densities and the reproduction of spatial spectra (the Fourier transform of the second-order velocity spatial correlations). Also, with very few exceptions (e.g. Illingworth et al., 2018), these analyses very often focus only on the most amplified modes without resulting in a detailed quantitative comparison of the performance of the $\nu$ and $\nu_{t}$ models.

A similar situation has emerged in the evaluation of the performance of Karhunen-Loève decompositions (proper orthogonal decompositions, POD) in the approximation of the power spectral density. In this context it has been recognised that the aggregation of all frequencies in the analysis, typically used when only spatial correlations are estimated, often blurs the interpretation of the results. It has hence been suggested to analyse separately each contributing frequency by making use of the Fourier transform of the spatio-temporal correlation tensor $\mathbf{R}(\mathbf{x},\mathbf{x}+\Delta\mathbf{x},t,t+\tau)$ and the associated ‘spectral’ POD modes (Lumley, 1970; Picard & Delville, 2000). It has then been highlighted that a strong relation exists between these spectral POD (SPOD) modes and those issued from the Schmidt (singular value, SVD) decomposition of the resolvent operator (Semeraro et al., 2016; Towne et al., 2018). However, in those investigations, only the performance of the resolvent analysis based on the $\nu$-model was evaluated.

The scope of the present study, inspired by the approach recently revived in the context of POD analyses, is therefore to evaluate the respective performance of the $\nu$-model and the $\nu_{t}$-model in the estimation of the velocity spatio-temporal power spectral density  and power cross-spectral density. The turbulent channel flow at ${\textit{Re}_{\tau}}=1007$ is used as a testbed for this analysis.

The paper is organized as follows: the mathematical formulation of the problem addressed is briefly summarized in §2. The turbulent channel flow obtained by direct numerical simulation and the relevant spectral information are described in §3 and compared to the estimations obtained via the resolvent-based models in §4. The results are summarized and discussed in §5. Further details on the operators involved in the resolvent-based models and on the data analysis of direct numerical simulations data are provided in App. A and B, respectively.

2 Background

2.1 Linear models for the dynamics of turbulent fluctuations

We are interested in the dynamics of coherent perturbations in a pressure-driven turbulent channel flow of incompressible viscous fluid with density $\rho$, kinematic viscosity $\nu$ between two infinite parallel walls located at $y=\pm h$. Here, we denote by $x$, $y$ and $z$ the streamwise, wall-normal and spanwise coordinates, respectively. In the following we will consider dimensionless variables based on the reference length $h$ and the reference velocity $(3/2)U_{bulk}$, where $U_{bulk}$ is the constant mass-averaged streamwise velocity.

Following Reynolds & Tiederman (1967); Reynolds & Hussain (1972); del Álamo & Jiménez (2006); Pujals et al. (2009) and many others, the evolution of small coherent perturbations to the turbulent mean flow can be modelled with linearised equations which include the effect of the turbulent Reynolds stresses in the linear model by means of an eddy viscosity $\nu_{t}$ corresponding to the turbulent mean flow profile $\mathbf{U}=(U(y),0,0)$:

[TABLE]

where $\mathbf{u}=(u,v,w)$ and $p$ are the coherent perturbation velocity and pressure, $\nu_{T}(y)=\nu+\nu_{t}(y)$ is the total effective viscosity and $\mathbf{f}$ is the forcing term. These equations are completed by the incompressibility condition $\nabla\cdot\mathbf{u}=0$.

For the eddy viscosity, we adopt the semi-empirical expression proposed by Cess (1958), as reported by Reynolds & Tiederman (1967):

[TABLE]

where $y\in[-1,1]$ is the centre-channel based dimensionless wall-normal coordinate, $y^{+}={\textit{Re}_{\tau}}(1-|y|)$ and ${\textit{Re}_{\tau}}=u_{\tau}h/\nu$ is the Reynolds number based on the friction velocity $u_{\tau}$. We set the von Kármán constant $\kappa=0.426$ and the constant $A=25.4$ as in Pujals et al. (2009); Hwang & Cossu (2010b). In the following we will refer to this model as the ‘$\nu_{t}$-model’.

In an alternative approach, used e.g. by Butler & Farrell (1993); Farrell & Ioannou (1993b); McKeon & Sharma (2010) and many others, turbulent Reynolds stresses and nonlinear terms are both included in the forcing term $\mathbf{f}$ (which is therefore different from the forcing term of the $\nu_{t}$-model) so that the linear model reduces to Eq. (1) but now only including the molecular kinematic viscosity $\nu_{T}=\nu$. We will refer to this model as ‘$\nu$-model’ (also known as ‘quasi-laminar model’).

2.2 Analyses of harmonic and stochastic forcing

Since the mean flow is homogeneous in the streamwise and spanwise directions, we consider the Fourier-modes ${\hat{\mathbf{u}}}(\alpha,y,\beta,t)e^{i(\alpha x+\beta z)}$ and ${\hat{\mathbf{f}}}(\alpha,y,\beta,t)e^{i(\alpha x+\beta z)}$ of streamwise and spanwise wavenumbers $\alpha=2\pi/\lambda_{x}$ and $\beta=2\pi/\lambda_{z}$. Eq. (1) can be reduced to the following linear system, expressed in terms of the state vector $\hat{\mathbf{q}}=\left[\hat{v},\hat{\omega}_{y}\right]^{T}$ formed with the wall-normal components of the velocity and vorticity Fourier modes:

[TABLE]

where $\hat{\mathbf{u}}=\mathbf{C}\hat{\mathbf{q}}$ and $\hat{\mathbf{q}}=\mathbf{D}\hat{\mathbf{u}}$ (the explicit expressions of the operators $\mathbf{A}$, $\mathbf{B}$, $\mathbf{C}$ and $\mathbf{D}$ are given in App. A). As the system (3) is linearly stable (Reynolds & Tiederman, 1967), the response to deterministic finite-power forcing can be analysed by Fourier-transforming the linear system in time, hence considering the harmonic forcing $\hat{\mathbf{f}}=\widetilde{\mathbf{f}}e^{-i\omega t}$ with the harmonic response $\widetilde{\mathbf{u}}e^{-i\omega t}$ related by:

[TABLE]

where $\mathbf{H}$ is the resolvent operator (or transfer function). When system  (3) is driven by stochastic forcing, the response is also stochastic. For the considered horizontally- and time-invariant problem the velocity second-order spatio-temporal correlation tensor is defined as

[TABLE]

where $\langle{~{}}\rangle$ denotes ensemble averaging and ∗ denotes complex conjugate transpose. The (spatio-temporal) power cross-spectral density  tensor is obtained through Fourier transform in the horizontal plane and time:

[TABLE]

The velocity power cross-spectral density  tensor $\mathbf{S}$ can also be obtained directly as the average of the Fourier transform $\widetilde{\mathbf{u}}$ of the velocity (e.g. Bendat & Piersol, 1986) and the stochastic forcing power cross-spectral density  tensor $\mathbf{P}$ can be defined similarly:

[TABLE]

An estimation $\mathbf{S}^{(est)}=\mathbf{H}\mathbf{P}\mathbf{H}^{*}$ of the velocity power cross-spectral density  tensor is obtained by making use of the linear models by replacing the linear expression of Eq. (4) in Eq. (7). If it is assumed that $\widetilde{\mathbf{f}}$ is $\delta$-correlated in the wall-normal direction $y$, with decorrelated components (see e.g. Farrell & Ioannou, 1993a; Hwang & Cossu, 2010a, b), then $\mathbf{P}=p(\alpha,\beta,\omega)\,\mathbf{I}$ (where $\mathbf{I}$ is the identity operator). This assumption, which has been made in a number of previous studies, underlies the so-called resolvent analysis because in this case $\mathbf{S}^{(est)}$ reduces to:

[TABLE]

The scope of this study is to explore the validity of the estimation provided by Eq. (8).

3 Direct numerical simulations in a LSM flow unit

The turbulent channel flow has been simulated at ${\textit{Re}_{\tau}}=1007$ in a domain $L_{x}=3$ long and $L_{z}=1.5$ wide in order to focus on the physics of the near-wall and large-scale processes only and to keep the data processing manageable. This domain size also corresponds to that of the most energetic large-scale motions (LSM) in the channel (del Álamo et al., 2004) and represents the minimum flow unit for the self-sustainment of coherent LSM (Hwang & Cossu, 2010c).

We have verified that the computed mean flow and the $rms$ fluctuation profiles are in reasonable agreement with those of Lee & Moser (2015) obtained in the much larger domain $L_{x}=25,L_{z}=9.5$ as shown in Fig. 1. We have further verified that Cess’s analytic fit works well by comparing the mean flow profile and the corresponding eddy viscosity to the ones issued by DNS data. It should be remarked that the values of the eddy viscosity are small in the near-wall region but very large in the bulk of the flow where $\nu_{T}/\nu\sim O(10-100)$ at the considered ${\textit{Re}_{\tau}}$.

The premultiplied streamwise kinetic energy spectral density $\alpha\beta E_{uu}(\alpha,y,\beta)=\alpha\beta\int_{-\infty}^{\infty}{S_{uu}}(\alpha,y,y,\beta,\omega)d\omega$, evaluated in the $y^{+}=15$ and $y=0.5$ planes, is shown in Fig. 2. As expected, the most energetic structures at $y=0.5$ have spatial scales corresponding to the dimensions of the LSM flow unit $\lambda_{x}=L_{x}=3$, $\lambda_{z}=L_{z}=1.5$ (corresponding to $\alpha\approx 2$, $\beta\approx 4$). The most energetic structures at $y^{+}=15$ have spatial scales $\lambda_{x}^{+}\approx 450$ and $\lambda_{z}^{+}\approx 100$ (corresponding to $\alpha\approx 14$, $\beta\approx 63$ at ${\textit{Re}_{\tau}}=1007$) typical of the near-wall self-sustained process. In the following we will therefore concentrate on these two sets of wavenumbers $(\alpha,\beta)$, corresponding to $(\lambda_{x}=3,\lambda_{z}=1.5)$ and $(\lambda_{x}^{+}=450,\lambda_{z}^{+}=100)$, which represent the dynamics of large-scale and near-wall self-sustained motions, respectively.

The power cross-spectral density  tensor has been computed for the two considered $(\alpha,\beta)$ pairs by means of Welch’s method (see e.g. Bendat & Piersol, 1986), as detailed in App. B. In Fig. 3 we report the streamwise velocity power spectral density  profiles ${S_{uu}}(\alpha,y,y,\beta,\omega)=\langle{\widetilde{u}(\alpha,y,\beta,\omega)\widetilde{u}^{*}(\alpha,y,\beta,\omega)}\rangle$ as a function of the phase speed $c=\omega/\alpha$, expressed in wall units, and of the wall normal coordinate for the considered $(\alpha,\beta)$ pairs. The peaks of these distributions are found at $\omega_{max}=4.3$ and $\omega_{max}=1.4$ for the near-wall and large-scale peaks, respectively. These peaks correspond to the two phase velocities $c^{+}_{max}\approx 10$ and $c^{+}_{max}\approx 20$ and their wall-normal locations ($y^{+}_{max}\approx 15$ and $y_{max}\approx 0.4-0.5$) approximatively correspond to the wall-normal position where $U^{+}=c^{+}$.

The ${S_{uu}}(\alpha,y,y^{\prime},\beta,\omega)=\langle{\widetilde{u}(\alpha,y,\beta,\omega)\widetilde{u}^{*}(\alpha,y^{\prime},\beta,\omega)}\rangle$ streamwise velocity power cross-spectral densities, are reported, for later comparison, in Fig. 4 as a function of $y$ and $y^{\prime}$ for the considered $(\alpha,\beta)$ at the peak $\omega_{max}$ corresponding to each pair.

4 Estimations of power cross-spectral density  based on linear models

4.1 Estimations of the power cross-spectral density

We now evaluate the capacity of the linear models of the coherent structures dynamics to reproduce statistics of the turbulent flow via the estimation $\mathbf{S}^{(res)}=p\mathbf{H}\mathbf{H}^{*}$ given by Eq. (8) using both the $\nu_{t}$ and the $\nu$ models to build the resolvent $\mathbf{H}$.

As a first case, we assume that the temporal power spectrum of the forcing is flat i.e. that $p$ does not depend on $\omega$ (white noise assumption) and $p=\overline{p}$ is chosen such that the estimated total spectral power of the streamwise velocity matches the one issued from direct numerical simulations $\int_{-\infty}^{\infty}\int_{-1}^{1}{S_{uu}}^{(est)}(\alpha,y,y,\beta,\omega)dyd\omega=\int_{-\infty}^{\infty}\int_{-1}^{1}{S_{uu}}^{(dns)}(\alpha,y,y,\beta,\omega)dyd\omega$. As a second case, we remove the white-noise assumption and determine $p$ (which depends now on $\omega$) such that the estimated power of the streamwise velocity matches that issued from direct numerical simulations at each selected frequency $\omega$: $\int_{-1}^{1}{S_{uu}}^{(est)}(\alpha,y,y,\beta,\omega)dy=\int_{-1}^{1}{S_{uu}}^{(dns)}(\alpha,y,y,\beta,\omega)dy$. We will refer to this second case as ‘coloured’ noise (see e.g. Zare et al., 2017).

In Figs. 5 and 6 we show, analogously to Fig. 3, the $y$-profiles of the estimated streamwise velocity power spectral density  versus the phase speed $c^{+}$.

The cases where the estimation is based on the $\nu$-model to compute the resolvent are reported in Fig. 5. For the case with white-noise (flat-spectrum) stochastic forcing (panels $a$ and $b$), the $\nu$-model does not select the correct values and locations of the power spectral density  peaks, which are predicted to be near the channel center, and therefore are predicted at too large $c^{+}$ values (and therefore too large $\omega$ values). For both near-wall and large-scale structures the power spectral density  appears also to be too narrowly concentrated near the critical layer (where $U^{+}=c^{+}$) when compared to the DNS data of Fig. 3. The estimation improves for the case of coloured noise (panels $c$ and $d$), where the selective power spectrum of the forcing is able to drive the response peaks nearer their DNS values. However, even in this case the power spectral density  spatial distributions remain too narrowly concentrated and therefore too large near the critical layer curves (the uniform dark red regions in panels $b$, $c$ and $d$ which are strongly offscale).

As shown in Fig. 6, the use of the $\nu_{t}$ model leads to a significant improvement because the effective eddy diffusivity dampens the response near the channel centre (where the eddy viscosity is high) and smoothens the critical layer peaks. Even with (non-selective, flat spectrum) white-noise stochastic forcing (panels $a$ and $b$), the $\nu_{t}$-model is able to select a reasonable $y^{+}_{max}$-$c^{+}_{max}$ location of the cross-spectral peak for buffer-layer structures. For large-scale structures (panel $b$), however, the amplitude of the response is still over-predicted near the channel centre and near the wall.

The estimation provided by the $\nu_{t}$ model further improves when the coloured-spectrum forcing is used (panels $c$ and $d$). A significant improvement is obtained, in particular, for the large-scale case (panel $d$) where the streamwise velocity power spectral density is not too large near the wall and the channel center and has a well defined peak located at reasonable $y_{max}$ and $c^{+}_{max}$ values.

These observations are further supported by the examination of streamwise velocity power cross-spectral density  ${S_{uu}}(\alpha,y,y^{\prime},\beta,\omega)=\langle{\widetilde{u}(\alpha,y,\beta,\omega)\widetilde{u}^{*}(\alpha,y^{\prime},\beta,\omega)}\rangle$ which are shown in Fig. 7 as a function of $y$ and $y^{\prime}$ for the same $(\alpha,\beta)$ and $\omega_{max}$ values already considered in Fig. 4. Also in this case, the $\nu$-model predicts too ‘narrow’ and too large power cross-spectral density  distributions which decay too fast from their peak location while the $\nu_{t}$ model is able to estimate reasonably well the power cross-spectral density  issued from the DNS data.

4.2 Expansion on SPOD and resolvent modes

Great attention has been given recently to the use of highly truncated expansions of the velocity power spectral densitybased on resolvent and spectral POD (SPOD) modes. Even if the main focus of this paper is not on these expansions, we find it useful to briefly examine the relevance of using the $\nu_{t}$ or the $\nu$ models in their computation.

Considering, for simplicity, only the streamwise velocity correlations, the spectral POD modes are the set of eigenfunctions $\psi_{j}^{(SPOD)}(y)$ of the self-adjoint operator ${S_{uu}}^{(dns)}$ ordered by decreasing corresponding real eigenvalues $e_{j}^{(SPOD)}$ (Lumley, 1970; Picard & Delville, 2000; Semeraro et al., 2016; Towne et al., 2018). These modes are orthonormal in the inner product associated to the power norm and therefore each SPOD eigenvalue, when normalized with respect to the sum of all eigenvalues, represents the fraction of the power associated with the corresponding SPOD mode. Thus, a salient property of the SPOD modes is that they provide the orthonormal basis of functions which is optimal (in terms of spectral power) for the expansion of the streamwise velocity coherent fluctuations.

Analogously, the resolvent modes $\psi_{j}^{(res)}(y)$ are the set of (orthonormal) eigenfunctions of ${S_{uu}}^{(res)}$ ordered by decreasing respective real eigenvalues $e_{j}^{(res)}$ and, in the present contest, they represent an estimation of the SPOD modes, under the assumption that $\mathbf{P}=p\,\mathbf{I}$. The leading resolvent eigenvalue $e_{1}^{(res)}=p\|\mathbf{H}\|^{2}$ also represents the optimal power amplification (modulo $p$) of the forcing supported by the resolvent operator $\mathbf{H}$ (see e.g. Trefethen et al., 1993; Schmid & Henningson, 2001).

We therefore examine the two already considered $(\alpha,\beta)$ pairs corresponding to buffer-layer and large-scale structures and the corresponding peak frequencies $\omega_{max}$ (the same as in Fig. 4 and Fig. 7); for these values, we compute the SPOD and resolvent modes and the corresponding eigenvalues based on the $\nu$ and the $\nu_{t}$ models, respectively including appropriate integration weights in the discretized operators (see e.g. Zare et al., 2017; Schmidt et al., 2018; Towne et al., 2018).

In what concerns the eigenvalues, a common trend is observed for the SPOD and the resolvent modes, as shown in Fig. 8: the two leading modes capture 70-90% of the total power, depending on the model; six modes are enough to capture more than 95% of the power of large-scale structures while they capture $\approx$85% of the power of near-wall structures (these values are similar to those found for linearly estimated space-only PODs by Farrell & Ioannou, 1993b; Hwang & Cossu, 2010a, b).

Finally, we expand the measured streamwise velocity root-mean-square ($rms$) profile $\widetilde{u}_{rms}^{(dns)}(y)=[{S_{uu}}^{(dns)}(y,y)]^{1/2}$ for the considered $\alpha,\beta,\omega_{max}$ on the leading six SPOD and resolvent modes getting the estimates $\widetilde{u}_{rms}^{(SPOD)}(y)=\sum_{j=1}^{6}a_{j}^{(SPOD)}\psi_{j}^{(SPOD)}(y)$, $\widetilde{u}_{rms}^{(res)}(y)=\sum_{j=1}^{6}a_{j}^{(res)}\psi_{j}^{(res)}(y)$.

These estimates are compared to the original profile in Fig. 9; from this figure we see that the SPOD expansion and the resolvent modes expansion based on the $\nu_{t}$ model provide a reasonable approximation of $\widetilde{u}_{rms}^{(dns)}$ for buffer-layer structures (Fig. 9$a$) and a quite good approximation for large-scale structures (Fig. 9$b$) while this is clearly not the case for the expansion based on the $\nu$ model resolvent modes, especially for large-scale structures.

5 Summary and conclusions

In this study we have compared the ‘measured’ spatio-temporal power cross-spectral density  $\mathbf{S}^{(dns)}$ obtained by direct numerical simulation of plane channel flow at ${\textit{Re}_{\tau}}=1007$ to the estimated one $\mathbf{S}^{(res)}=p\mathbf{H}\mathbf{H}^{*}$ based on the resolvent operator (linear transfer function) $\mathbf{H}$ under the customary assumption that the forcing spatio-temporal power cross-spectral density  is of the type $\mathbf{P}=p\,\mathbf{I}$. The main advantage of this type of analysis is that the contributions of different temporal frequencies to the spatial energy spectral density can be analysed separately, as emphasized in the recent studies of Semeraro et al. (2016); Schmidt et al. (2018) and Towne et al. (2018) who analyse the relation between proper orthogonal decomposition and resolvent analysis.

Two cases have been considered for the forcing power spectrum $p(\omega)$: a flat (white-noise) forcing power spectrum $p=\overline{p}$ and that of a (‘coloured’ noise) forcing power spectrum where $p(\omega)$ is chosen so as to match the power spectrum of the estimated streamwise velocity power spectral density  to the measured one. Two separate linear models have been used to evaluate the resolvent: the first, labelled ‘$\nu_{t}$ model’, includes the effect of the turbulent Reynolds stresses modelled with an eddy viscosity while the second, labelled ‘$\nu$ model’, does not include them.

The two linear models and the two types of forcing power spectrum, have been tested for structures with spatial scales typical of buffer-layer structures ($\lambda_{z}^{+}=100,\lambda_{x}^{+}=450$) and of large-scale motions ($\lambda_{z}=1.5,\lambda_{x}=3$). We find that estimations based the $\nu_{t}$ model resolvent lead to qualitatively correct estimations of power cross-spectral density  distributions and of the dependence of the power spectral density  wall-normal distributions on the frequency (or, equivalently, the phase speed). On the contrary, estimations based on the $\nu$ model lead to generally overestimated values of the streamwise velocity power spectral density  which are (too) narrowly concentrated near the critical layer. The $\nu$ model is also unable to even qualitatively reproduce the $c^{+}_{max}-y^{+}_{max}$ values of the peak amplitude of the power spectral density  for buffer-layer structures. We also show that in all cases, as expected, the estimation is improved by using the appropriately coloured input power spectrum $p(\omega)$ instead of white noise as already observed by Zare et al. (2017) using a different methodology.

Similar results are found when comparing highly truncated expansions of the measured streamwise velocity power spectral density  on the basis of the leading spectral POD (SPOD) modes to the expansions on leading resolvent modes. The six-modes expansion based on SPOD modes approximates well the original power spectral density  and a comparable accuracy is obtained by making use of the same number of resolvent modes based on the $\nu_{t}$ model, while the use of resolvent modes based on the $\nu$ model leads to a less accurate approximations, especially for large-scale structures.

These results can be understood by recalling that in estimates based on the $\nu_{t}$ model the effect of turbulent Reynolds stresses is included in the resolvent $\mathbf{H}$ while in the $\nu$ model the effect of turbulent Reynolds stresses resides only in the forcing and the associated power cross-spectral density  $\mathbf{P}$. As a consequence, when using linear models in the estimation $\mathbf{S}^{(est)}=\mathbf{H}\mathbf{P}\mathbf{H}^{*}$, when no accurate information on $\mathbf{P}$ is available (which is often the case, especially at high Reynolds numbers, and one therefore assumes the standard $\mathbf{P}=p\mathbf{I}$), the $\nu_{t}$ model performs relatively well because the effect of the turbulent Reynolds stresses is embedded in $\mathbf{H}$. It is likely that the $\nu$ model could perform equally well in the case where relevant information on turbulent Reynolds stresses is available leading to a realistic modelling of the forcing cross-spectral tensor $\mathbf{P}$; in this second case the scale selection would be dictated by $\mathbf{P}$ instead of $\mathbf{H}$.

These findings are consistent with those of Illingworth et al. (2018) who find that the use of the $\nu_{t}$ model leads to a better performance in the estimation of streamwise velocity fluctuations in the turbulent channel at the same friction Reynolds number ${\textit{Re}_{\tau}}=1007$ considered here (however, they assume a white-noise spectrum and do not separate the different frequency contributions in their analysis). It is also important to note that there is no contradiction between our findings and those of Semeraro et al. (2016); Schmidt et al. (2018) and Towne et al. (2018) who show a reasonable agreement of resolvent modes computed with the $\nu$ model to experimental power cross-spectral density  data in turbulent jets because, contrary to turbulent wall-bounded flows, in the case of jets the turbulent eddy viscosity has only a weak dependence on the (cross-stream) radial coordinate and therefore the $\nu$ and $\nu_{t}$ formulations almost coincide except for a Reynolds number rescaling by $\nu/\nu_{T}$.

In the present investigation turbulent Reynolds stresses have been modelled with eddy viscosity, as in a number of previous studies, but an improvement of quantitative predictions of the power cross-spectral density  tensor could probably come from a more advanced modelling of the Reynolds stress tensor in the linear operator. Also, as the linear operators that we have discussed do mainly model the non-normal streaks amplification mechanism, important progress could come from an improved, necessarily nonlinear, modelling of the coherent forcing terms which are related to the regeneration mechanism of the vortices in the self-sustained processes discussed by Waleffe (1995); Hwang & Cossu (2010c, 2011) and Cossu & Hwang (2017) among others. These issues are under current intensive scrutiny.

Appendix A Linear model operators

The linear system of Eq. (3) is derived (see e.g. Schmid & Henningson, 2001) from equations (1) with the operators $\mathbf{A}$ and $\mathbf{B}$ defined as

[TABLE]

where the generalised Orr-Sommerfeld and Squire operators Cossu et al. (2009); Pujals et al. (2009) are:

[TABLE]

with $\mathcal{D}$ and ′ denoting $d/dy$, $k^{2}=\alpha^{2}+\beta^{2}$ and $\Delta=\mathcal{D}^{2}-k^{2}$. Homogeneous boundary conditions are enforced on both walls: $\hat{v}(\pm 1)=\mathcal{D}\hat{v}(\pm h)=\hat{\omega}_{y}(\pm 1)=0$. The linear operators relating $\hat{\mathbf{u}}$ and $\hat{\mathbf{q}}$ are

[TABLE]

Appendix B Direct numerical simulations and data analysis

The simulations have been performed using the SIMSON code, which is pseudo-spectral and uses Fourier expansions in the streamwise and spanwise directions, where periodicity is enforced, and Chebyshev expansions in the wall-normal direction to solve the three-dimensional time-dependent incompressible Navier-Stokes equations in the channel (see Chevalier et al., 2007, for further details).

$N_{x}=128,N_{y}=129,N_{z}=128$ points (before 3/2-rule dealiasing of Fourier modes) have been used for the ${\textit{Re}_{\tau}}=1007$ simulations in the domain $L_{x}=3,L_{y}=2,L_{z}=1.5$, with a uniformly spaced grid in $x$ and $z$ and Gauss-Lobatto points in $y$. At ${\textit{Re}_{\tau}}=1007$ this corresponds to a grid spacing $\Delta x^{+}=23,\Delta z^{+}=12$ and to $\Delta y^{+}_{min}=0.3$ near the wall and $\Delta y^{+}_{max}=24$ at the channel centre. This chosen grid is coarser than those typically used in DNS at this Reynolds number (e.g. Lee & Moser, 2015, who use $\Delta x^{+}=11,\Delta z^{+}=5$) to keep the data analysis manageable but it has been verified that this has no significant impact on the first ans second-order flow statistics (see Fig. 1).

The initial transient of the simulation is discarded and statistics and samples are accumulated starting from $t=95000$. Welch’s method with Hamming windowing and $75$% overlap has been used to compute the velocity power cross-spectral density  $\mathbf{S}$ using a total of $N_{s}=8001$ snapshots of the DNS solutions with sampling interval $(\Delta t)_{s}=0.25$ for a total acquisition time $T_{max}=2000$ for the small scales, and using $N_{s}=20001$ snapshots sampled every $(\Delta t)_{s}=0.5$ for a total acquisition time $T_{max}=10000$ for the large scales. The average temporal step of the DNS during the acquisition time is $\overline{\Delta t}=0.0106$, i.e. $\overline{\Delta t}^{\,+}=0.35$. Data have also been averaged between the two walls.

Acknowledgements.

Direct numerical simulations were performed on resources provided by the Swedish National Infrastructure for Computing (SNIC) at NSC, HPC2N and PDC