1. Introduction
The transition of a boundary layer from the laminar regime to fully developed turbulence is a central problem in an immense range of technological applications because turbulent wall friction can be several times larger than that exerted by a laminar boundary layer. Frictional losses in the boundary layer are responsible for the performance degradation of engineering flow systems, such as turbomachinery and jet engines, for the enhanced aerodynamic drag of transport vehicles, and, in turn, for wasted fuel consumption, unwanted noise production and environmental pollution. For design purposes, it is therefore paramount to be able to predict under which conditions boundary-layer transition occurs. Free-stream turbulence acts as a triggering factor for transition, and it has been shown that the transition Reynolds number decreases as the free-stream turbulence level increases (Mayle Reference Mayle1991).
Dryden (Reference Dryden1936) and Taylor (Reference Taylor1939) were probably the first to study the effects of free-stream turbulence on a flat-plate boundary layer. They showed that the dominant streamwise velocity fluctuations generated by free-stream turbulence in the boundary layer are of very low frequency and reach amplitudes that can be several times larger than those in the free stream.
The Dryden–Taylor observations did not receive much attention until Klebanoff (Reference Klebanoff1971) carried out experiments in which he reproduced the earlier findings of Dryden and Taylor. Klebanoff demonstrated that the disturbances grow more or less linearly with the boundary-layer thickness, and they are quite narrow in the spanwise direction. Klebanoff referred to these disturbances as ‘breathing modes’ because, as noted earlier by Taylor (Reference Taylor1939), they appeared to correspond to a thickening and thinning of the boundary layer. Kendall (Reference Kendall1991) renamed them Klebanoff modes, and that name has taken hold even though these disturbances are not modes in the strict mathematical sense, i.e. they are not homogeneous solutions of differential equations.
The early transition experiments were conducted at very low free-stream turbulence levels (
$Tu<0.1\,\%$), but more recent experiments, such as those by Westin et al. (Reference Westin, Boiko, Klingmann, Kozlov and Alfredsson1994), Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001), Fransson, Matsubara & Alfredsson (Reference Fransson, Matsubara and Alfredsson2005), Fransson & Shahinfar (Reference Fransson and Shahinfar2020) and Mamidala, Weingärtner & Fransson (Reference Mamidala, Weingärtner and Fransson2022) were carried out at higher turbulence levels. However, the results are invariably the same. The dominant streamwise velocity fluctuations are always of the Klebanoff type, i.e. the boundary layer acts as a low-frequency-pass filter on the free-stream perturbation spectrum, and amplifies streamwise stretched streaky vortical structures. The spanwise wavelength of the Klebanoff modes is constant along the streamwise direction, and the peak amplitude occurs at the same Blasius-similarity wall-normal location. Direct numerical simulations have also been employed to study the development of low-frequency streaks and the induced bypass transition (Jacobs & Durbin Reference Jacobs and Durbin2001; Ovchinnikov, Choudhari & Piomelli Reference Ovchinnikov, Choudhari and Piomelli2008; Yao, Mollicone & Papadakis Reference Yao, Mollicone and Papadakis2022).
The mathematical framework describing the incompressible Klebanoff modes was developed by Leib, Wundrow & Goldstein (Reference Leib, Wundrow and Goldstein1999) (LWG99). They proved that these disturbances, near the leading edge, are well represented by forced solutions of the linearized unsteady boundary-layer equations for which the spanwise viscous effects are negligible. As the mean boundary layer grows downstream, these equations lose their validity because the spanwise length scale of the Klebanoff modes becomes comparable with the boundary-layer thickness. Their dynamics is then ruled by the unsteady boundary-region equations, i.e. the Navier–Stokes equations where the spanwise viscous terms are retained, while the streamwise pressure gradient and the viscous effects can be neglected because the perturbations are of low frequency and streamwise elongated. The boundary-region equations, and their terminology, were first used by Kemp (Reference Kemp1951) to study the corner boundary-layer problem. A crucial ingredient in the LWG99 formulation is the continuous action of the free-stream perturbations that are responsible for the generation and evolution of the Klebanoff modes. LWG99 utilized matched asymptotic expansions to obtain the initial and outer boundary conditions that synthesize the interaction between the free-stream flow and the boundary-layer flow. Wundrow & Goldstein (Reference Wundrow and Goldstein2001) and Ricco, Luo & Wu (Reference Ricco, Luo and Wu2011) extended the linearized study of LWG99 to include nonlinear effects, focusing on the steady and unsteady cases, respectively. Ricco et al. (Reference Ricco, Luo and Wu2011) also explained the occurrence of nonlinear effects in the results by Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001), and studied the secondary instability of the saturated Klebanoff modes, thereby describing the mechanism at the heart of bypass transition induced by free-stream turbulence. Extensions to the compressible regime include the investigations by Ricco & Wu (Reference Ricco and Wu2007), Ricco, Tran & Ye (Reference Ricco, Tran and Ye2009), Ricco, Shah & Hicks (Reference Ricco, Shah and Hicks2013) and Marensi, Ricco & Wu (Reference Marensi, Ricco and Wu2017).
Other theories describing the Klebanoff modes have been proposed. The non-modal growth theory (Schmid & Henningson Reference Schmid and Henningson2001) and the optimal growth theory (Andersson, Berggren & Henningson Reference Andersson, Berggren and Henningson1999; Luchini Reference Luchini2000) model the growth of streaky disturbances already present in the boundary layer, while allowing the disturbances to vanish in the free stream. Continuous Orr–Sommerfeld modes have also been used extensively since Jacobs & Durbin (Reference Jacobs and Durbin2001) to synthesize the penetration of free-stream disturbances into a boundary layer. Reviews of this approach are found in Dong & Wu (Reference Dong and Wu2013), Ricco et al. (Reference Ricco, Walsh, Brighenti and McEligot2016) and Durbin (Reference Durbin2017).
In the present study, we develop the theoretical background of previously unexplained experimental results of a transitional boundary layer exposed to free-stream turbulence, reported by Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001) (MA01). These findings are remarkable because the energy spectra at different streamwise locations were found to collapse on one another when scaled properly. MA01 described their discovery as ‘an unexpected new finding’ and their energy spectra showing ‘an astonishing similarity’ for which ‘there is no theoretical explanation’.
In § 2, the experimental findings of MA01 on the scaling of the Klebanoff modes are discussed. In § 3, we present the key features of the mathematical framework describing the Klebanoff modes, while the theoretical results behind the experimental findings of MA01 are found in § 4. Section 5 contains the conclusions.
2. Discussion of the experimental results of Matsubara & Alfredsson
MA01studied experimentally an incompressible flow of uniform velocity 
$U_\infty ^*$ past a thin flat plate located in a low-speed wind tunnel. Rigid grids were placed upstream of the leading edge of the plate to generate free-stream vortical disturbances. A thin laminar boundary layer developed over the flat plate and transitioned to a fully-developed turbulent boundary layer because of the perturbative action of the free-stream disturbances. The objective of the MA01 study was to fully characterize the transitional boundary layer. In our discussion of the MA01 results and in the theoretical analysis, the flow is described through a Cartesian coordinate system, i.e. 
$\boldsymbol {x}^*$ = 
$x^* \hat {\boldsymbol {i}} + y^* \hat {\boldsymbol {j}} + z^* \hat {\boldsymbol {k}}$, where 
$x^*,y^*,z^*$ define the streamwise, wall-normal and spanwise directions, respectively, and the superscript 
$^*$ indicates a dimensional quantity. The flat plate is located at 
$y^*=0$, and its leading edge is at 
$x^*=0$. Lengths are scaled by 
$\varLambda _z^*$, the integral spanwise length scale of the free-stream vortical disturbances, velocities are scaled by 
$U_\infty ^*$, pressure is scaled by 
$\rho ^* U_\infty ^{*2}$, where 
$\rho ^*$ is the density, and time is scaled by 
$\varLambda _z^*/U_\infty ^*$. The kinematic viscosity is denoted by 
$\nu ^*$. Non-dimensional quantities are not marked by any symbol.
2.1. Validity of linearized dynamics
As our theoretical framework hinges on the assumption that the boundary-layer disturbances are described by a linearized dynamics, we first examine the MA01 findings to support our hypothesis. Figure 1 shows the mean boundary-layer streamwise velocity profiles measured by MA01 at different streamwise locations. The data displayed by the black circles correspond to the three streamwise stations that are closest to the leading edge, i.e. 
$x^*=100, 300, 500$ mm. The data represented by the thin lines were acquired at 
$x^*>500$ mm. The black-circle data show excellent agreement with the numerical solution of the Blasius laminar boundary-layer flow, represented by the thick red line, while the thin-line data deviate progressively more and more from the laminar solution as 
$x^*$ increases. For 
$x^*>500$ mm, nonlinear effects become important as the boundary-layer perturbations grow in amplitude, and the wall-shear stress is enhanced as fully-developed turbulence ensues. These results are evidence of the perturbed flow obeying a linearized dynamics at the locations closest to the leading edge because the mean-flow profiles follow the laminar solution. Figure 3 in MA01 further reveals that the boundary-layer thickness and the shape factor match the laminar values for 
$x^*\leq 700$ mm. Additional support for these results is given by profiles of the root-mean-square (r.m.s.) of the streamwise velocity fluctuations, shown in figure 2(c) of MA01, which denote clear signs of nonlinear effects for 
$x^*\geq 1100$ mm, such as the disturbances growing in the outer part of the boundary layer, and the perturbation peak moving closer to the wall. The theoretical and numerical results that match quantitatively the nonlinear MA01 data are discussed in Ricco et al. (Reference Ricco, Luo and Wu2011). We conclude that a linearized dynamics can be utilized to study the perturbed flow for 
$x^*\leq 500$ mm, despite the free-stream turbulence intensity not being vanishingly small for these experiments, i.e. 
$Tu=2.2\,\%$ (refer to grid A in table 1 in MA01).
Figure 1. Mean boundary-layer streamwise velocity profiles reported by MA01 for 
$x^*\leq 500$ mm (black circles) and 
$700\,{\rm mm}\leq x^*\leq 1900$ mm (thin lines). The red thick line denotes the numerical solution of the Blasius laminar boundary-layer flow, and 
$\delta _d^*$ indicates the displacement thickness.
2.2. Scaling of experimental turbulence spectra
Figures 2(a,b), a reproduction of figure 13 in MA01, depict streamwise velocity energy spectra at 
$y^*/\delta _d^*=1.2$, where 
$\delta_d^*$ is the displacement thickness. For this experimental dataset, 
$U^{*}_\infty =5$ m s
$^{-1}$ and 
$\varLambda _z^*=7$ mm, computed from the autocorrelation of the streamwise velocity shown in figure 7 on p. 161 of MA01. The Reynolds number based on 
$\varLambda _z^*$ is 
$R_\lambda =U_\infty ^* \varLambda ^*_z/\nu ^*=2232$. The spectrum 
$E_\alpha$ is shown as a function of the streamwise wavenumber 
$k_x^*=2 {\rm \pi}f^*/U_\infty ^*$, where 
$f^*$ is the frequency (figure 2a), and the spectrum 
$E_\beta$ is shown as a function of the spanwise wavenumber 
$k_z^*=2 {\rm \pi}/\lambda _z^*$, where 
$\lambda _z^*$ is the spanwise wavelength (figure 2b).
Figure 2. (a,b) Reproduction of figure 13 in MA01. Energy spectra as functions of (a) the streamwise wavenumber and (b) the spanwise wavenumber. (c,d) Reproduction of figure 14 in MA01. Rescaled energy spectra as functions of (c) the scaled streamwise wavenumber and (d) the dimensional spanwise wavenumber. Data were acquired at 
$y^*/\delta _d^*=1.2$. The solid lines are for 
$x^*=120, 150, 200, 250, 300, 400, 500$ mm. Labels in the original graphs have been changed to conform to the present notation.
The wavenumbers in figures 2(a,b) are dimensional, while in our theoretical analysis they are scaled by 
$\varLambda _z^*$, that is, 
$k_x=k_x^* \varLambda _z^*$ and 
$k_z=k_z^* \varLambda _z^*$. The spectra 
$E_\alpha$ and 
$E_\beta$ are linked to the variance of the streamwise velocity fluctuations,
\begin{equation} \epsilon^2 \langle {u^{\prime 2}} \rangle_{zt}(x,y) = C_\alpha \int_0^\infty E_\alpha(k_x)\,\mathrm{d}k_x = C_\beta \int_0^\infty E_\beta(k_z)\,\mathrm{d}k_z, \end{equation}
where 
$C_\alpha$ and 
$C_\beta$ are constants, computed in § 2.3, and 
$\left \langle {\cdot } \right \rangle _{zt}$ indicates averaging along 
$z$ and over 
$t$. In figure 2, the dash-dotted lines refer to locations upstream of the solid lines, while the dashed lines correspond to locations downstream of the solid lines.
In figure 2(a), for 
$x^*\leq 500$ mm, the dash-dotted and solid lines show that the low-wavenumber portion of the spectrum grows downstream, while the high-wavenumber portion is unchanged. This behaviour confirms that the boundary layer acts as a low-frequency-pass filter (Durbin Reference Durbin2017), consistently with the algebraic growth of the streamwise-elongated, low-frequency Klebanoff modes. The high-frequency free-stream disturbances do not penetrate sufficiently into the boundary layer to reach these wall-normal locations. Nonlinear effects becomes predominant further downstream, where the high-wavenumber fluctuations grow more significantly than the low-wavenumber ones (dashed lines). Figure 2(b) shows that the spanwise energy spectrum grows uniformly for all the spanwise wavenumbers.
Figures 2(c,d) are a reproduction of figure 14 in MA01. The spectra 
$E_\alpha$ and 
$E_\beta$, shown in figures 2(a,b), are scaled as (the symbol 
$\hat {\cdot }$ is used here in lieu of 
$^*$ in MA01)
\begin{equation} \hat E_\alpha = \frac{E_\alpha}{C_e\,Re_x^{3/2}},\quad \hat E_\beta = \frac{E_\beta}{C_e\,Re_x}, \end{equation}
where 
$Re_x=U_\infty ^* x^*/\nu ^*$, and the constant 
$C_e=16$ is the same for the two spectra. The scaling of 
$E_\beta$ with 
$Re_x$ is expected because the integral of 
$E_\beta$ along 
$k_z$, given by the last equation of (2.1), is equal to the variance of the streamwise velocity fluctuations, which grows linearly with 
$Re_x$, as shown by the experimental results in figure 2(d) of MA01. On the abscissas of figures 2(c,d), the streamwise wavenumber is scaled by the displacement thickness 
$\delta _d^*$, while the spanwise wavenumber is dimensional. Both sets of profiles represented by the solid lines show excellent collapse when rescaled. The objective of our study is explain the scaling of those solid lines in figures 2(c,d).
This scaling demonstrates that the streamwise spectrum 
$E_\alpha$ grows downstream at a faster rate (proportional to 
$Re_x^{3/2}$) than its integral across the streamwise wavenumbers 
$\epsilon ^2 \langle {u^{\prime 2}} \rangle _{zt}$, which grows linearly with 
$Re_x$, as shown in figure 2(d) of MA01. The different growth rates are caused by the low-frequency fluctuations becoming larger more rapidly than the high-frequency ones, as shown in figure 2(a).
It is worth mentioning that Zhigulev, Uspenskii & Ustinov (Reference Zhigulev, Uspenskii and Ustinov2009), in their figures 7 and 8, reported similar scaling of streamwise spectra, in their case by 
$Re_x^2$ and 
$Re_x^{3/2}$, for different boundary-layer datasets collected in their low-turbulence wind tunnel (
$Re_x^2$ and 
$Re_x^{3/2}$ were written as 
$\epsilon ^2 \langle {u^{*\prime 2}} \rangle _{zt}$ and 
$\epsilon ^2 \langle {u^{*\prime 2}} \rangle _{zt} \delta _d^*$, respectively, in their formulas (2.8) and (2.9)). They attributed the scaling by 
$Re_x^{3/2}$ to nonlinear effects. We show in the following that the scaling of the MA01 spectra can be explained by asymptotic results emerging from the linearized theory of LWG99, although our form of free-stream spectrum does model nonlinear effects through its streamwise dependency.
2.3. Computation of 
${C_\alpha }$ and ${C_\beta }$
The constants 
${C_\alpha }$ and 
${C_\beta }$ in (2.1) are found as follows. The integrals in (2.1) are first computed by using the spectral data in figures 2(a,b) at different streamwise locations 
$Re_x$. For the experimental data of figure 2, MA01 do not report the values of 
$\epsilon ^2 \langle {u^{\prime 2}} \rangle _{zt}$ at different 
$Re_x$. The data shown in figure 2(d) on p. 156 of MA01 for a similar set of flow conditions can, however, be used for our purpose because that graph shows that the r.m.s. of the streamwise velocity starts to deviate from the linear behaviour when it reaches a value of about 
$9 \times 10^{-3}$. The constants 
${C_\alpha }$ and 
${C_\beta }$ can thus be found by linear fitting of the integrated experimental data in order to obtain 
$\epsilon ^2 \langle {u^{\prime 2}} \rangle _{zt}=9 \times 10^{-3}$ at 
$Re_x=159\,438$, which is the most downstream location where the data of figure 2 obey the scaling discussed in § 2.2 (denoted by solid lines). Data downstream of this location, displayed by dashed lines in figures 2(a,b), are affected by nonlinear effects, similarly to the r.m.s. data larger than 
$9 \times 10^{-3}$ in figure 2(d) on p. 156 of MA01. The computed values are 
${C_\alpha }=1.62 \times 10^{-10}$ and 
${C_\beta }=4 \times 10^{-12}$. Figure 3 shows that the r.m.s. values, obtained by integrating 
$E_\alpha$ and 
$E_\beta$, agree well with each other and grow linearly with 
$Re_x$ as expected. MA01 give the free-stream turbulence level for this experimental dataset, 
$Tu(\%)=0.022$, and we thus take 
$\epsilon =0.022$.
Figure 3. Growth of r.m.s. of streamwise velocity fluctuations as a function of 
$Re_x$, computed by integrating the spectra 
$E_\alpha$ (red circles) and 
$E_\beta$ (blue squares) shown by solid lines in figures 2(a,b).
2.4. Power-law dependence of scaled turbulence spectra
The data in figures 2(c,d) are replotted in figure 4, which reveals that the experimental data of the energy spectra by MA01 are well approximated by the power laws
$$\begin{gather} \hat E_\alpha = \frac{1.91\times 10^{{-}5}}{(k_x \delta_d)^{\tilde \alpha}},\quad \mbox{where}\ \tilde \alpha=2.82, \end{gather}$$
$$\begin{gather}\hat E_\beta = \frac{8.3\times 10^2}{k_z^{\tilde \beta}}, \quad \mbox{where}\ \tilde \beta=1.55. \end{gather}$$The power laws (2.3) and (2.4) are useful in our theoretical analysis of § 4.
Figure 4. (a) Scaled energy spectrum 
${\hat E}_\alpha$ as a function of the streamwise wavenumber 
$k_x \delta _d$, and (b) scaled energy spectrum 
${\hat E}_\beta$ as a function of the spanwise wavenumber 
$k_z$. The experimental data by MA01, also shown in figures 2(c,d), are represented by the red circles, and the algebraic best fitting lines in solid blue represent relations (a) (2.3) and (b) (2.4).
3. Theoretical framework for the Klebanoff modes
The theory of the Klebanoff modes is found in LWG99. Here, we report the main points that are useful for our analysis of the wind-tunnel flow studied by MA01.
3.1. The free-stream disturbance flow at short streamwise distances
A uniform flow of velocity 
$U_\infty ^*$ past an infinitely thin flat plate transports homogeneous, statistically stationary vortical fluctuations of the gust type, i.e. disturbances that are convected passively by the mean flow. These free-stream perturbations are assumed to be of small amplitude with respect to 
$U_\infty ^*$, so that the free-stream flow is represented as the sum of the mean uniform flow and the free-stream vortical disturbances, as
\begin{align} \boldsymbol{u}_\infty &= \hat{\boldsymbol{i}} + \epsilon\,\boldsymbol{u}'_\infty (x - t, y , z) = \hat{\boldsymbol{i}} + \epsilon \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \hat{\boldsymbol{u}}'_\infty(k_x,k_y,k_z) \nonumber\\ &\quad \times\exp({\mathrm{i}({\boldsymbol{k}\boldsymbol{\cdot}\boldsymbol{x}} - k_x t)}) \,\mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z, \end{align}
where 
$\epsilon \ll 1$, 
$\hat {\boldsymbol {u}}'_\infty =\{\hat u^\infty,\hat v^\infty,\hat w^\infty \}={O}(1)$, 
$\boldsymbol {k}=\{k_x,k_y,k_z\}$, and the streamwise wavenumber 
$k_x$ and the frequency 
$-k_x$ are related because of Taylor's hypothesis (Taylor Reference Taylor1938; Hunt Reference Hunt1973). In the experiments of MA01, the turbulence is generated by a grid located upstream of the leading edge of the plate, but we consider 
$x=0$ as the streamwise location where the free-stream turbulence starts influencing the system because that is where the turbulence intensity was measured by MA01, as explained in the second paragraph of p. 154 in MA01. The representation (3.1) is valid at wall-normal distances that are sufficiently large for the flow not to be influenced by the presence of the boundary layer and the flat plate. The free-stream perturbation (3.1) is not influenced by viscous dissipation while being transported downstream by the free-stream potential flow because it is valid only at sufficiently small 
$x$ location. The streamwise evolution of the free-stream flow is nevertheless taken into account at larger streamwise locations by the model of the free-stream spectrum studied in §§ 4.1 and 4.2, and by the numerical solution of the free-stream disturbance flow that includes the viscous dissipation and the inviscid displacement of the mean-flow streamlines due to the boundary layer, as discussed in § 3.2. Furthermore, expansion (3.1) is not valid for amplitudes of free-stream disturbances comparable with that of the mean flow and for a non-uniform free-stream mean flow because Taylor's hypothesis does not apply in those cases (Lundell & Alfredsson Reference Lundell and Alfredsson2004).
3.2. The Klebanoff modes
In the limit of large Reynolds number, 
$R_\lambda \gg 1$, the mean laminar boundary layer that develops over the flat plate is described by the steady boundary-layer equations (Schlichting & Gersten Reference Schlichting and Gersten2000). The mean-flow streamwise and wall-normal velocity components are 
$U(x,y)$ and 
$V(x,y)$, and the wall-normal similarity coordinate is 
$\eta =y/\delta =y \sqrt {R_\lambda /2 x}$, where 
$\delta =\sqrt {2 x/R_\lambda }=\sqrt {2}\,\delta _d/1.72$ is the boundary-layer thickness used in LWG99.
The free-stream vortical flow encounters the boundary layer and generates the Klebanoff modes, as documented by the experimental data of MA01 discussed in § 2. We consider the limit 
$k_x={O}(R_\lambda ^{-1})\ll k_y,k_z$ because the Klebanoff modes are of low frequency. The boundary layer indeed acts as a low-frequency-pass filter, thus only the low-frequency disturbances penetrate into the boundary layer, as evidenced in figure 9(b) on p. 162 of MA01. We study the flow at downstream locations where 
$\delta ^*={O}(\varLambda _z^*)$, and we scale the streamwise coordinate as 
$\bar {x}=k_xx={O}(1)$. As explained in LWG99, the condition for linearization in the boundary layer is 
$\epsilon /k_x\ll 1$. The boundary-layer flow is expressed as the sum of the mean boundary-layer flow 
$\boldsymbol {U}$ and the disturbance flow 
$\epsilon \boldsymbol {u}'$, as follows (LWG99; Hunt Reference Hunt1973; Hunt & Carruthers Reference Hunt and Carruthers1990):
\begin{align} \boldsymbol{u} &=
\boldsymbol{U}(x,y) +
\epsilon\,\boldsymbol{u}'(x,y,z,t) =
\boldsymbol{U}(x,y) + \epsilon \int_{-\infty}^\infty
\int_{-\infty}^\infty \hat{\boldsymbol{u}}'(x,y,k_x,k_z)
\nonumber\\ &\quad \times \exp({\mathrm{i}(k_z z -
k_xt)})\,\mathrm{d}k_x \,\mathrm{d}k_z \nonumber\\ &=
\{U,V,0\} + \epsilon\int_{-\infty}^\infty
\int_{-\infty}^\infty \left\{ \bar{u}_0 (\bar{x},\eta),
\left(\frac{2 \bar{x} k_x}{R_\lambda}\right)^{1/2}
\bar{v}_0 (\bar{x}, \eta), \bar{w}_0 (\bar{x},\eta)
\right\} \nonumber\\ &\quad \times \exp({\mathrm{i}(k_z z -
k_xt)}) \,\mathrm{d}k_x \,\mathrm{d}k_z +
{O}(\epsilon^2),
\end{align}
where the leading-order velocity components with respect to 
$k_x\ll 1$ are retained, i.e. 
$\{\bar {u}_0,\bar {v}_0,\bar {w}_0\}=[\hat w^\infty +\mathrm {i} k_z \hat v^\infty /(k_x^2+k_z^2)^{1/2}]$ 
$\{(\mathrm {i} k_z/k_x) \bar {u},(\mathrm {i} k_z/k_x) \bar {v},\bar {w}\}$. The components 
$\{\bar {u},\bar {v},\bar {w}\}$ satisfy the linearized unsteady boundary-region equations, complemented by initial and boundary conditions, all found in LWG99. Homogeneous boundary conditions at the wall represent the no-slip condition, while mixed boundary conditions in the free stream account for the boundary-layer inviscid displacement and the perturbation decay due to viscous dissipation. The system is solved by a second-order implicit finite-difference scheme and a standard block-elimination algorithm (Ricco & Wu Reference Ricco and Wu2007), described in Appendix A.
The scaled wavenumber 
$\kappa _z=k_z/(k_x R_\lambda )^{1/2}={O}(1)$ represents the relative importance between spanwise and wall-normal viscous effects at 
$\bar {x}={O}(1)$. In the limit 
$\kappa _z\ll 1$, the spanwise viscous diffusivity becomes negligible and the dynamics is ruled by the boundary-layer equations.
We now discuss an asymptotic result, based on the parameter 
$\kappa _z$, which is central in the analysis developed in § 4. LWG99 showed that an asymptotic solution exists in the low-frequency, large-spanwise-wavenumber limit 
$\kappa _z\gg 1$ with 
$\tilde \kappa =\kappa _y/\left |\kappa _z\right |={O}(1)$, where 
$\kappa _y=k_y/(k_x R_\lambda )^{1/2}$. In this limit, the leading-order velocity components 
$\{\bar {u},\bar {v},\bar {w}\}$ are rescaled and expressed as a function of the new streamwise coordinate 
$\tilde x=\kappa _z^2 \bar {x}={O}(1)$, i.e. 
$\tilde u(\tilde x,\eta,\tilde \kappa )=\kappa _z^2 \bar {u}={O}(1)$, 
$\{\tilde v,\tilde w\}(\tilde x,\eta,\tilde \kappa )=\{\bar {v},\bar {w}\}={O}(1)$. The rescaled velocity components 
$\{\tilde u,\tilde v,\tilde w\}$ are quasi-steady and depend only on the ratio of wavenumbers 
$\tilde \kappa$ and not explicitly on the scaled spanwise wavenumber 
$\kappa _z$. Although the asymptotic solution is valid for 
$\kappa _z\gg 1$, the numerical calculations reveal the remarkable result that the algebraic growth of the quasi-steady asymptotic solution 
$\{\tilde u,\tilde v,\tilde w\}$ is indistinguishable from the full boundary-region solution even for 
$\kappa _z$ as low as 1. Figure 5 indeed shows that the trends of 
$\left |\tilde u\right |$ for different 
$\kappa _z\geq 1$ and the same 
$\tilde \kappa$ collapse onto one another when plotted as a function of 
$\tilde x$. It also means that the asymptotic solution describes the Klebanoff modes well even when the spanwise wavelength is comparable with the boundary-layer thickness, which is precisely the flow condition of interest in the experiments of MA01. Therefore, the asymptotic solution 
$\{\tilde u,\tilde v,\tilde w\}$ is utilized in the scaling analysis of § 4, where the collapse of the spectral distributions shown in figure 4 is obtained. Figure 5 also reveals that the initial growth of the disturbance is linear when 
$\kappa _z \geq 1$, that is, 
$|\tilde u|= G(\tilde \kappa )\left |\kappa _z\right | \sqrt {\bar {x}}$. The inset of figure 5 shows the slope 
$G(\tilde \kappa )$. The decay of 
$\tilde u$ as 
$\tilde \kappa \rightarrow \infty$, and therefore of 
$G(\tilde \kappa )$, is predicted by the asymptotic analysis because, in the limits 
$\kappa _z\gg 1$ and 
$\tilde \kappa \gg 1$, the solution can be written as 
$\breve u(\breve x,\eta )={\tilde \kappa }^2 \tilde u=\kappa _y^2 \bar {u}={O}(1)$, where 
$\breve x={\tilde \kappa }^2 \tilde x=\kappa _y^2 \bar {x}$.
Figure 5. Growth and decay of the scaled streamwise velocity component of the Klebanoff modes 
$|\tilde u|=\kappa _z^2\,|\bar {u}|$ at 
$\eta =1.46$ as a function of the scaled streamwise coordinate 
$\tilde x=\left |\kappa _z\right | \sqrt {\bar {x}}$ for different 
$\tilde \kappa$ values. The velocity is computed by solving numerically the boundary-region equations, found in LWG99. The straight solid lines denote the linear growth. The inset shows the slope of the linear growth, 
$G(\tilde \kappa )$.
4. Scaling of the Klebanoff modes
4.1. Variance of the boundary-layer streamwise velocity
The boundary-layer perturbations and the free-stream modes are related as (Hunt Reference Hunt1973; Hunt & Carruthers Reference Hunt and Carruthers1990)
\begin{equation} \hat{u}_i'(x,y,k_x,k_z) = \int_{-\infty}^\infty M_{ij}(x,y;k_x,k_y,k_z)\,\hat{u}_{\infty j}'(k_x,k_y,k_z) \,\mathrm{d}k_y, \end{equation}
where 
$M_{ij}$ is a tensor acting as a transfer function between the free-stream flow and the boundary-layer flow. The interest is in the correlation of the boundary-layer velocity components, delayed in time and 
$z$ (Batchelor Reference Batchelor1953),
\begin{equation} R_{ij}(x,y,r_z,\tau) = \epsilon^2 \left\langle {u^\prime_i(x,y,z+r_z,t+\tau)\,u^\prime_j(x,y,z,t)} \right\rangle_{zt}, \end{equation}which can be expressed as (refer to pp. 638–640 in Hunt Reference Hunt1973)
\begin{align} R_{ij}(x,y,r_z,\tau) &= \epsilon^2 \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \sum_{l=1}^3 \sum_{m=1}^3 M_{il}^{\dagger} M_{jm}\,\varPhi_{\infty lm}(\boldsymbol{k}) \nonumber\\ &\quad \times\exp[{\mathrm{i}(k_z r_z - k_x \tau)}]\,\mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z, \end{align}
where 
$\varPhi _{\infty lm}$ is the spectral tensor of the turbulence upstream of the flat plate, and the symbol 
${\dagger}$ indicates the complex conjugate. The focus is on the spectral properties of the mean-square streamwise velocity fluctuations, i.e. 
$i=j=1$, 
$r_z=\tau =0$ (LWG99),
\begin{equation} \epsilon^2 \langle {u^{\prime 2}} \rangle_{zt} = R_{11}(x,y,0,0) = \epsilon^2 \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \sum_{l=1}^3 \sum_{m=1}^3 M_{1l}^{\dagger} M_{1m}\,\varPhi_{\infty lm} (\boldsymbol{k})\,\mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
The relevant components of the transfer-function tensor 
$M_{ij}$ are
\begin{equation} M_{11} = \bar{u}^{(0)},\quad M_{12} = \frac{\mathrm{i} k_x}{\sqrt{k_x^2+k_z^2}}\,\bar{u}^{(0)} - \frac{k_z^2}{k_x \sqrt{k_x^2+k_z^2}}\,\bar{u}, \quad M_{13} = \frac{\mathrm{i} k_z}{k_x}\,\bar{u}, \end{equation}
where 
$\bar {u}^{(0)}$ is the next-order term of the expansion of 
$\bar {u}_0$ in (3.2) with respect to 
$k_x\ll 1$ (LWG99). By substituting (4.5a–c) into (4.4) and collecting the dominant terms 
${O}(k_x^{-2})$, the integrand in (4.4) becomes
\begin{equation} \sum_{l=1}^3 \sum_{m=1}^3 M_{1l}^{\dagger} M_{1m}\,\varPhi_{\infty lm}(\boldsymbol{k}) = \frac{k_z^2 \left|\bar{u}\right|^2}{k_x^2}\left(\frac{k_z^2}{\sqrt{k_x^2+k_z^2}}\,\varPhi_{\infty 22} + \varPhi_{\infty 33}\right)+ {O}(k_x^{{-}1}). \end{equation}As suggested by LWG99 on p. 187, an axial-symmetric turbulence model that describes free-stream turbulence is (Batchelor Reference Batchelor1953; Chandrasekhar Reference Chandrasekhar1950)
\begin{equation} \varPhi_{\infty ij} = \frac{ k_\perp^2 \delta_{ij}^\perp{-} k_{{\perp} i} k_{{\perp} j} }{k_\perp^2} \left( \varPhi_t - \frac{2 k_x^2}{k_\perp^2}\,\varPhi_x \right) + \frac{\varPhi_x}{k_\perp^2} \left( k_x^2 \delta_{ij}^\perp{-} k_x k_{{\perp} i} \delta_{i1} + k_\perp^2 \delta_{i1} \delta_{j1} \right), \end{equation}
where 
$k_{\perp i} = k_i - \delta _{i1} k_x$, 
$\delta _{i1}$ is the Kronecker delta, 
$\delta _{ij}^\perp = \delta _{ij} - \delta _{i1} \delta _{j1}$ is the cross-stream Kronecker delta, and 
$k_\perp =\sqrt {k_y^2+k_z^2}$. The functions 
$\varPhi _x=\varPhi _x(k_x,k_\perp )$ and 
$\varPhi _t=\varPhi _t(k_x,k_\perp )$ are the longitudinal and transverse spectra. In the limit 
$k_x\rightarrow 0$,
\begin{equation} \varPhi_{\infty 22} = \frac{k_z^2}{k_x^2+k_z^2}\,\varPhi_t,\quad \varPhi_{\infty 33} = \frac{k_x^2}{k_x^2+k_z^2}\,\varPhi_t. \end{equation}Substitution of (4.8a,b) into (4.6) and then into (4.4) leads to the variance of boundary-layer streamwise velocity,
\begin{equation} \langle {u^{\prime 2}} \rangle_{zt}(x,y) = \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \left(\frac{k_z}{k_x}\right)^2 \left|\bar{u}\right|^2(x,y)\, \varPhi_t(k_x,k_\perp)\,\mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
As discussed by LWG99, these results demonstrate that, at leading order, the growth and development of the Klebanoff modes is dictated by the transverse spectral function 
$\varPhi _t$ obtained by correlations of the velocity components perpendicular to the streamwise direction (refer to LWG99 on p. 188), and not by the longitudinal spectral function 
$\varPhi _x$, which is typically the object of experimental investigations of freely decaying grid-generated turbulence.
4.2. Free-stream turbulence spectrum
The axial-symmetric transverse turbulence spectrum 
$\varPhi _t(k_x,k_\perp )$ in § 4.1 is assumed to pertain to homogeneous turbulence and it is therefore independent of the streamwise direction (Hunt Reference Hunt1973; Hunt & Carruthers Reference Hunt and Carruthers1990). However, in a more general non-homogeneous case, the turbulence spectrum also depends on the position vector, 
$\varPhi _t(\boldsymbol {x},k_x,k_\perp )$, as for example discussed in Townsend (Reference Townsend1980). To the best of our knowledge, no detailed measurements of 
$\varPhi _t$ have been made, so our objective is to suggest a functional form for 
$\varPhi _t$ that is a satisfactory model for our problem.
Our choice of spectrum takes inspiration from the theory of temporally decaying turbulence discussed in Townsend (Reference Townsend1980) on p. 61. The results in the streamwise decaying case can be assumed to be qualitatively analogous to the temporally decaying case if the streamwise direction is considered in lieu of time for flows where the turbulence intensity is much smaller than the free-stream mean velocity, i.e. when Taylor's hypothesis is valid, as explained in Townsend (Reference Townsend1980) on p. 65. In the idealized limit of vanishingly small amplitude of free-stream turbulence generated by a grid swept through a still fluid, Batchelor (Reference Batchelor1953), on p. 93, shows that the time dependency is due solely to the viscous dissipation, and the temporal decay is exponential. However, Batchelor (Reference Batchelor1953) warns that this behaviour would occur only after a long time, and it would not apply to a real turbulent flow generated by a grid in a wind tunnel. The exponential decay would thus not pertain to locations relatively close to the turbulence-generating grid, which are certainly of interest in the study of the MA01 experimental results. Furthermore, if the turbulence spectrum 
$\varPhi _t$ were assumed to be independent of the streamwise direction, as in § 4.1, the streamwise evolution of the free-stream disturbance would affect the variance 
$\left \langle {u^{\prime 2}} \right \rangle _{zt}$ in the boundary layer only indirectly through the decaying free-stream wall-normal and spanwise velocity components because 
$\left |\bar {u}\right |$, the leading-order component in (4.9), vanishes as 
$y\rightarrow \infty$ (refer to (5.11) and (5.20)–(5.22) in LWG99). Neglecting the streamwise dependency of the free-stream spectrum would mean that the free-stream decay would be purely exponential because it is dictated by a linearized dynamics. Including the streamwise dependence in 
$\varPhi _t$ is therefore deemed to be more realistic, and it also serves the purpose of modelling mild effects of nonlinearity. Similar modelling of mild nonlinearity in a free-stream spectrum pertaining to realistic grid-generating turbulence has been proposed by LWG99 in their § 7.2.
Townsend (Reference Townsend1980) on p. 61 shows that the spectral function for decaying turbulence has the form
\begin{equation} E(k,t) = \langle {u'(t)^2} \rangle L(t)\,\mathcal{F}(k L(t)), \end{equation}
where 
$L(t)$ is an integral scale representing the free-stream isotropic turbulence, 
$\left \langle {\cdot } \right \rangle$ indicates spatial averaging, and 
$k$ is the wavenumber. The spectral function (4.10) is found by appropriate scaling of experimental data (Stewart & Townsend Reference Stewart and Townsend1996), as also discussed in Hinze (Reference Hinze1975) on p. 263. By substitution of (4.10) into the equation governing the rate of change of the turbulence spectrum, Townsend (Reference Townsend1980) finds
\begin{equation} \frac{\mathrm{d} \left\langle {u'(t)^2} \right\rangle}{\mathrm{d} t} \propto \frac{ \left\langle {u'(t)^2} \right\rangle^{3/2}}{L(t)}, \quad \frac{\mathrm{d} L(t)}{\mathrm{d} t} \propto \langle {u'(t)^2} \rangle^{1/2}, \end{equation}
as further explained in Batchelor (Reference Batchelor1953) on p. 103. The temporal decay of 
$\left \langle {u'(t)^2} \right \rangle$ that satisfies (4.11a,b) is
\begin{equation} \langle {u'(t)^2} \rangle \propto t^{-\gamma}, \end{equation}
which is consistent with numerous experimental data, for which 
$1.15<\gamma <1.45$ (refer to p. 160 of Pope Reference Pope2000), and with theoretical studies, which suggest 
$\gamma =1$ (refer to Tennekes & Lumley Reference Tennekes and Lumley1972) or 
$\gamma =3/2$ (refer to Davidson (Reference Davidson2004) on p. 407, where the Saffman spectrum is discussed). The decay constant 
$\gamma$ can then be assumed to be
\begin{equation} 1 \leq \gamma \leq 3/2. \end{equation}
The integral spatial scale 
$L$ is predicted to grow as
\begin{equation} L \propto t^\zeta, \quad 1/4 \leq \zeta \leq 1/2. \end{equation}Substitution of (4.12) and (4.14) into (4.10), and use of (4.13), lead to a simplified form of the spectrum
\begin{equation} E(k,t) \propto \frac{\mathcal{F}\left(k t^{d/2}\right)}{t^{c/2}}, \end{equation}for which the inequalities
\begin{equation} 1 \leq c \leq 5/2 \quad \mbox{and} \quad 1/2 \leq d \leq 1 \end{equation}
apply. Also, 
$c=3\gamma -2$ and 
$d=2-\gamma$, from which
\begin{equation} c=4-3d. \end{equation} The time-decaying isotropic spectrum (4.15) can now be used to obtain a spectrum that pertains to the grid-generated turbulence of interest in our problem. As the spectrum has to account for the streamwise decay of turbulence, the temporal dependence in (4.15) is converted to the streamwise dependence. The axial symmetry of the turbulence has to be modelled by including the effect of the cross-flow wavenumber 
$k_\perp$ because, as explained by Batchelor (Reference Batchelor1953), purely isotropic turbulence is extremely hard to obtain in the laboratory. Our axial-symmetric transverse spectrum therefore reads
\begin{equation} \varPhi_t(x,k_x,k_\perp;R_\lambda) = \frac{1}{k_\perp^b (k_x R_\lambda)^2 \delta^c}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{k_\perp^n}\right), \end{equation}
where, in lieu of 
$t$ in (4.15), we have introduced the streamwise coordinate 
$x$ and expressed this dependence through the boundary-layer thickness 
$\delta$ because 
$\delta \propto \sqrt {x}$. The dependence of the spectrum on 
$k_\perp$ is introduced inside and outside the function 
$\mathcal {F}$ to allow maximum generality. The spatial dependence of the spectrum (4.18) is mild compared with the long streamwise length scale of the Klebanoff modes because (4.18) is expressed as a function of 
$\delta =\delta ^*/\varLambda _z^*$, where 
$\delta ^*$ and 
$\varLambda _z^*$ are comparable. Consistently with the theoretical framework of § 3, the low-frequency assumption is adopted as the boundary layer acts as a low-frequency-pass filter. It is thus reasonable to consider a free-stream spectrum such as (4.18), dominated by low-frequency disturbances (
$k_x\ll 1$ with 
$k_x R_\lambda ={O}(1)$ or smaller).
The parameters 
$n,b,c,d$ in (4.18) are found by asymptotic analysis and by fitting the experimental data. The parameters 
$c$ and 
$d$ play analogous roles in (4.15) and (4.18).
4.3. Scaling of boundary-layer streamwise velocity spectra
By substituting the spectrum (4.18) into (4.9), the variance of the boundary-layer streamwise velocity becomes
\begin{equation} \langle {u^{\prime 2}} \rangle_{zt}= \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \left(\frac{k_z}{k_x}\right)^2\frac{\left|\bar{u}\right|^2}{k_\perp^b (k_x R_\lambda)^2 \delta^c}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{k_\perp^n}\right) \mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation} Expression (4.19) is used with (2.1) and (2.2a,b) to explain the scaling of the experimental results, shown in figures 2(c,d). The four parameters 
$n,b,c,d$ in (4.19) are found by using the following four conditions.
(i) In figure 2(c), the spectrum
$\hat E_\alpha$ depends only on the scaled streamwise wavenumber $k_x \delta _d$.(ii) In figure 2(d), the spectrum
$\hat E_\beta$ depends only on the spanwise wavenumber $k_z$ and is independent of the streamwise location.(iii) In figure 4(a), the best fitting of the experimental data leads to the power-law dependency (2.3) for
$\hat E_\alpha (k_x \delta _d)$.(iv) In figure 4(b), the best fitting of the experimental data leads to the power-law dependency (2.4) for
$\hat E_\beta (k_z)$.
4.3.1. Spectrum versus spanwise wavenumber
Motivated by the scaling of the spectrum 
$E_\beta$ by 
$Re_x$, given in the second expression in (2.2a,b), the variance (4.19) is rescaled by 
$Re_x$ as
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x}= \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \frac{k_z^2 \left|\bar{u}\right|^2}{C_e k_\perp^b (k_x R_\lambda)^3 \bar{x} \delta^c}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{k_\perp^n}\right) \mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
The streamwise velocity 
$\left |\bar {u}\right |$ is changed to 
$\left |\bar {u}\right |^2= \left |\tilde {u}\right |^2 k_x^2 R_\lambda ^2/k_z^4$, and the streamwise coordinate is eliminated by using 
$\bar {x}=\delta ^2 k_x R_\lambda /2$, to obtain
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x}= \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \frac{2 \left|\tilde{u}\right|^2}{C_e k_z^2 k_\perp^b (k_x R_\lambda)^2 \delta^{c+2}}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{k_\perp^n}\right) \mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
The asymptotic solution for 
$\kappa _z \gg 1$, i.e. 
$\left |\tilde {u}\right |^2 = (k_z \delta )^2 \left |G(\tilde \kappa )\right |^2/2$, shown in figure 5 and discussed at the end of § 3.2, is substituted into (4.21) to arrive at
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x}= \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \frac{\left|G(\tilde \kappa)\right|^2}{C_e k_\perp^b (k_x R_\lambda)^2 \delta^c}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{k_\perp^n}\right) \mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
The wavenumbers 
$k_\perp$ and 
$k_y$ are eliminated by using 
$k_\perp = |k_z|\,(1+\tilde \kappa ^2)^{1/2} = |k_z|\, K(\tilde \kappa )$ and 
$k_y = k_z \tilde \kappa$, and the integration limits are changed to 
$[0,\infty )$:
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x}= \int_0^\infty \int_0^\infty \int_0^\infty \frac{2^3 \left|G(\tilde \kappa)\right|^2}{C_e k_z^{b-1}(k_x R_\lambda)^2 K(\tilde \kappa)^b \delta^c}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{(k_z K(\tilde \kappa))^n}\right) \mathrm{d}k_x \,\mathrm{d}\tilde \kappa \,\mathrm{d}k_z. \end{equation}By using the rescaled (2.1),
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x} = \frac{C_\beta}{\epsilon^2} \int_0^\infty \hat E_\beta (k_z) \,\mathrm{d}k_z, \end{equation}we find
\begin{equation} \hat E_\beta (k_z)=\frac{2^3 \epsilon^2}{C_e C_\beta R_\lambda^2 k_z^{b-1} \delta^c} \int_0^\infty\frac{\left|G(\tilde \kappa)\right|^2}{K(\tilde \kappa)^b} \underbrace{ \int_0^\infty \frac{1}{k_x^2}\,\mathcal{F}\left(\frac{k_x R_\lambda \delta^d}{(k_z\,K(\tilde \kappa))^n}\right) \mathrm{d}k_x }_{I_\beta}\,\mathrm{d}\tilde \kappa. \end{equation}
By defining the integration variable 
$\sigma = k_x R_\lambda \delta ^d/[k_z\,K(\tilde \kappa )]^n$, the integral 
$I_\beta$ in (4.25) becomes
\begin{equation} I_\beta =\frac{R_\lambda \delta^d}{[k_z\,K(\tilde \kappa)]^n } \int_0^\infty \frac{\mathcal{F}(\sigma)}{\sigma^2}\,\mathrm{d}\sigma. \end{equation}Upon substitution of (4.26) into (4.25), we obtain
\begin{equation} \hat E_\beta (k_z)=\frac{2^3 \epsilon^2 \delta^{d-c}}{C_e C_\beta R_\lambda k_z^{b-1+n}} \int_0^\infty \frac{\left|G(\tilde \kappa)\right|^2}{[K(\tilde \kappa)]^{n+b}}\,\mathrm{d}\tilde \kappa \int_0^\infty \frac{\mathcal{F}(\sigma)}{\sigma^2}\,\mathrm{d}\sigma. \end{equation}
The key point here is that, as the function 
$\hat E_\beta$ must not depend on the streamwise direction, the dependence on 
$\delta$ must be eliminated. It follows that 
$c=d$.
The spectrum (4.27) becomes
\begin{equation} \hat E_\beta (k_z)= \frac{B_\beta G_\beta \varSigma_\beta}{{k_z}^{\tilde \beta}}, \end{equation}
where 
$\tilde \beta = b -1 + n$, and
\begin{equation} B_\beta = \frac{2^3 \epsilon^2}{R_\lambda C_e C_\beta},\quad G_\beta = \int_0^\infty \frac{\left|G(\tilde \kappa)\right|^2}{\left(1+\tilde \kappa^2\right)^{(n+b)/2}}\,\mathrm{d}\tilde \kappa, \quad \varSigma_\beta = \int_0^\infty \frac{\mathcal{F}(\sigma)}{\sigma^2}\,\mathrm{d}\sigma. \end{equation}
The algebraic decay emerging in (4.28) matches the behaviour of the experimental data in figure 4(b). At small 
$k_z$, the theoretical framework does not predict the trend of the data in figure 4(b), which is almost independent of 
$k_z$. At small 
$k_z$, the spanwise wavelength is larger than the boundary-layer thickness, the spanwise viscous effects are negligible, and the flow is ruled by the boundary-layer equations, as discussed in § 3.1. Our analysis instead hinges on the asymptotic solution of the boundary-region equations for which the spanwise wavelength and the boundary-layer thickness are comparable, i.e. the wall-normal and spanwise diffusion effects are both important (
$\kappa _z={O}(1)$ or larger). The same reasoning applies to the dash-dotted lines in figure 2(d), which do not collapse onto one another as they correspond to streamwise locations close to the leading edge, where spanwise-diffusion effects are negligible.
4.3.2. Spectrum versus streamwise wavenumber
Motivated by the scaling of the spectrum 
$E_\alpha$ by 
$Re_x^{3/2}$, given in the first expression in (2.2a,b), the variance (4.19) is rescaled by 
$Re_x^{3/2}$. By using 
$c=d$, found in § 4.3.1, we find
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x^{3/2}}= \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \frac{k_z^2 \left|\bar{u}\right|^2}{C_e k_x^{5/2} k_\perp^b R_\lambda^{7/2} \delta^c \bar{x}^{3/2}}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^c}{k_\perp^n}\right) \mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
The streamwise velocity 
$\left |\bar {u}\right |$ is changed to 
$\left |\bar {u}\right |^2= \left |\tilde {u}\right |^2 k_x^2 R_\lambda ^2/k_z^4$, and the streamwise coordinate is eliminated by using 
$\bar {x}=\delta ^2 k_x R_\lambda /2$, to find
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x^{3/2}} = \int_{-\infty}^\infty \int_{-\infty}^\infty \int_{-\infty}^\infty \frac{2^{3/2}\left|\tilde u\right|^2}{C_e k_z^2 k_x^2 k_\perp^b R_\lambda^3 \delta^{c+3}}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^c}{k_\perp^n}\right) \mathrm{d}k_x \,\mathrm{d}k_y \,\mathrm{d}k_z. \end{equation}
The changes of variable 
$k_z = k_\perp \sin \theta$, 
$k_y = k_\perp \cos \theta$, 
$k_\perp =k_o/\delta$ are used in (4.31) to find
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x^{3/2}}= \int_{-\infty}^\infty \int_0^\infty \int_0^{2 {\rm \pi}} \frac{2^{3/2}\left|\tilde u\right|^2}{C_e \delta^{c+3-b} k_o^{b+1} (\sin \theta)^2 k_x^2 R_\lambda^3}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^{n+c}}{k_o^n}\right) \mathrm{d}\theta \,\mathrm{d}k_o \,\mathrm{d}k_x. \end{equation}We substitute the asymptotic result
\begin{equation} \left|\tilde u\right|^2=\frac{k_z^2 \delta^2}{2} \left|G(\tilde \kappa)\right|^2 =\frac{k_o^2 (\sin \theta)^2}{2} \left|G(\cot \theta)\right|^2 \end{equation}into (4.32) to obtain
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x^{3/2}}= \int_{-\infty}^\infty \int_0^\infty \int_0^{2 {\rm \pi}} \frac{\sqrt{2} \left|G(\cot \theta)\right|^2 }{C_e \delta^{c+4-b} k_o^{b-1} k_x^2 R_\lambda^3}\, \mathcal{F}\left(\frac{k_x R_\lambda \delta^{n+c}}{k_o^n}\right) \mathrm{d}\theta \,\mathrm{d}k_o \,\mathrm{d}k_x. \end{equation}By using the rescaled (2.1),
\begin{equation} \frac{\left\langle {u^{\prime 2}} \right\rangle_{zt}}{C_e\,Re_x^{3/2}} = \frac{C_\alpha}{\epsilon^2} \int_0^\infty \frac{\hat E_\alpha (k_x \delta)}{\delta}\,\mathrm{d}(k_x \delta), \end{equation}
changing the limits of the integration along 
$k_x$ to 
$[0,\infty )$, and equating (4.34) and (4.35), we find
\begin{equation} \hat E_\alpha (k_x \delta)=\frac{2^{3/2} \epsilon^2 G_\alpha}{C_e C_\alpha R_\lambda^3} \int_0^\infty \mathcal{F}\left(\frac{k_x R_\lambda \delta^{n+c}}{k_o^n}\right) \frac{\mathrm{d}k_o}{k_x^2 \delta^{c+3-b} k_o^{b-1}}, \end{equation}where
\begin{equation} G_\alpha = \int_0^{2 {\rm \pi}} \left|G(\cot \theta)\right|^2 \mathrm{d}\theta. \end{equation}
We define the integration variable 
$\omega = k_x R_\lambda \delta ^{n+c}/k_o^n$ (
$n>0$) in (4.36) to obtain
\begin{equation} \hat E_\alpha (k_x \delta)=\frac{2^{3/2} \epsilon^2 G_\alpha}{n C_e C_\alpha R_\lambda^{3+(b-2)/n}} \int_0^\infty \frac{\mathcal{F}( \omega )\,\mathrm{d}\omega}{\omega^{1+(2-b)/n} \delta^{c+1+c(b-2)/n} k_x^{2+(b-2)/n}}. \end{equation}
By defining 
$\tilde \alpha = 2+(b-2)/n$ and 
$\tilde d = c + 1 + c(b-2)/n$, the spectrum becomes
\begin{equation} \hat E_\alpha (k_x \delta)=\frac{2^{3/2} \epsilon^2 G_\alpha}{n C_e C_\alpha R_\lambda^{3+(b-2)/n} \left(\delta^{\tilde d/\tilde \alpha} k_x\right)^{\tilde \alpha}} \int_0^\infty \frac{\mathcal{F}( \omega )\,\mathrm{d}\omega}{\omega^{1+(2-b)/n}}. \end{equation}
For the spectrum 
$\hat E_\alpha$ to depend only on 
$k_x \delta$, we set 
$\tilde \alpha = \tilde d$. It follows that 
$c=1$. The values 
$c=1$ and 
$d=1$ respect the inequalities (4.16a,b) and the relation (4.17) obtained in § 4.2 from Townsend's spectrum. The decay constant becomes 
$\gamma =1$, which also falls within the inequality range predicted by Townsend's theory and is consistent with theoretical and experimental studies (Tennekes & Lumley Reference Tennekes and Lumley1972; Fransson et al. Reference Fransson, Matsubara and Alfredsson2005).
By using the displacement thickness 
$\delta _d$ instead of 
$\delta$, as in the MA01 experiments, the spectrum (4.39) becomes
\begin{equation} \hat E_\alpha (k_x \delta_d)= \frac{A_\alpha G_\alpha \varOmega_\alpha}{{(k_x \delta_d)}^{\tilde \alpha}}, \end{equation}where
\begin{equation} A_\alpha = \frac{2^{3/2} \epsilon^2}{n C_e C_\alpha R_\lambda^{3+(b-2)/n}} \left(\frac{1.72}{\sqrt{2}}\right)^{\tilde \alpha},\quad \varOmega_\alpha = \int_0^\infty \frac{\mathcal{F}\left( \omega \right)}{\omega^{1+(2-b)/n}}\,\mathrm{d}\omega. \end{equation}
For 
$k_x \delta _d<0.04$, the experimental data shown in figure 4(a) decay algebraically at a smaller rate than at larger 
$k_x \delta _d$. For fixed 
$k_x$ and small 
$\delta _d$, the spanwise wavelength is larger than 
$\delta _d$, the spanwise diffusivity is negligible, and the flow is described by the boundary-layer equations. It is then expected that the spectrum behaves differently when the boundary-region equations, used in our theoretical framework, instead describe the flow.
4.4. Parameters of the transverse spectrum $\varPhi _t$
We use the exponents 
$\tilde \alpha =2.82$ and 
$\tilde \beta =1.55$ in (2.3), found from the best-fitting analysis in § 2.4, to solve the algebraic expressions 
$\tilde \beta =b-1+n$, found in § 4.3.1, and 
$\tilde \alpha =2+(b-2)/n$, found in § 4.3.2. The four coefficients of the transverse spectrum 
$\varPhi _t$ are
\begin{equation} c = d = 1, \quad n = \frac{\tilde \beta - 1}{\tilde \alpha - 1} = 0.302, \quad b = \frac{\tilde \alpha \tilde \beta + \tilde \alpha - 2 \tilde \beta}{\tilde \alpha - 1} = 2.248. \end{equation}
Table 1 presents the numerical values related to the energy spectra 
$\hat E_\alpha$ and 
$\hat E_\beta$.
Table 1. Numerical values of quantities related to the energy spectra 
$\hat E_\alpha$ and 
$\hat E_\beta$.
4.5. The spectral function $\mathcal {F}$
A spectral function 
$\mathcal {F}$ that satisfies the two integrals 
$\varSigma _\beta$ and 
$\varOmega _\alpha$, given in table 1, is now chosen. Inspired by Ishihara et al. (Reference Ishihara, Kaneda, Yokokawa, Itakura and Uno2005) and Sagaut & Cambon (Reference Sagaut and Cambon2008), we select 
$\mathcal {F}(\xi )= A_f \xi ^{a_1}\exp (- a_2 \xi ^{a_3})$, where the coefficients satisfy
\begin{equation} \varSigma_\beta =
\int_0^\infty
\frac{\mathcal{F}(\sigma)}{\sigma^2}\,\mathrm{d}\sigma =
A_f \int_0^\infty \sigma^{a_1-2} \exp ({-}a_2 \sigma^{a_3}) \mathrm{d}\sigma = \frac{A_f\,\varGamma\left(
\dfrac{a_1 - 1}{a_3}\right)}{a_3 a_2^{{(a_1-1)}/{a_3}}} =
0.488, \end{equation}
\begin{align} \varOmega_\alpha &=
\int_0^\infty
\omega^{\bar{\omega}}\,\mathcal{F}( \omega )\,\mathrm{d}\omega = A_f
\int_0^\infty \omega^{a_1+\bar{\omega}} \exp ({-}a_2
\omega^{a_3} ) \mathrm{d}\omega \nonumber\\
&= \frac{A_f\,\varGamma\left(
\dfrac{a_1 + \bar{\omega} +1}{a_3}\right)}{a_3 a_2^{{(a_1 +
\bar{\omega} + 1)}/{a_3}}} = 125.22,
\end{align}
with 
$\varGamma$ the Gamma function, and 
$\bar {\omega }=(b-2)/n-1=0.179$. We can find multiple combinations of 
$A_f$, 
$a_1$, 
$a_2$ and 
$a_3$ that satisfy 
$\varSigma _\beta$ and 
$\varOmega _\alpha$. Figure 6 shows an example of the spectral function 
$\mathcal {F}(\xi )$.
Figure 6. Spectral function 
$\mathcal {F}(\xi )$ for 
$A_f=0.03$, 
$a_1=3$, 
$a_2=0.3$ and 
$a_3=0.9$.
5. Conclusions and outlook
In this paper, we have continued our effort to obtain theoretical and numerical results that explain the experimental findings reported by Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001), one of the most important studies on the impact of free-stream turbulence on the growth and evolution of velocity perturbations in a flat-plate transitional boundary layer. In Ricco et al. (Reference Ricco, Luo and Wu2011), our theoretical framework and calculations reproduced the main features reported by Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001) on the initiation of nonlinear effects within the boundary layer, such as the enhancement of the wall-shear stress with respect to the laminar value, the growth of disturbances in the outer part of the boundary layer, and the motion of the peak fluctuations towards the wall. In the present paper, we have instead focused on the collapse of the energy spectral profiles, obtained by Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001) when appropriate rescaling was adopted.
The spectral theory of homogeneous temporal-decaying turbulence developed by Townsend (Reference Townsend1980) has been utilized to obtain a model spectrum for the streamwise-decaying axial-symmetric free-stream turbulence generated by Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001) by use of a grid located in the upstream section of their wind tunnel. Quasi-steady asymptotic solutions of the unsteady boundary-region equations, found by Leib et al. (Reference Leib, Wundrow and Goldstein1999), have been used in the analysis of the experimental results of Matsubara & Alfredsson (Reference Matsubara and Alfredsson2001). The quasi-steady approximation was justified by the established finding that the boundary layer acts as a low-frequency-pass filter on the free-stream fluctuations, i.e. low-frequency disturbances are amplified in the boundary layer, while high-frequency disturbances are less prone to reach the core of the boundary layer.
Further work should be directed at measurements of the cross-stream velocity components in the free stream to arrive at a functional form for the transverse spectrum, which is responsible for the generation of the low-frequency Klebanoff modes inside the boundary layer (Leib et al. Reference Leib, Wundrow and Goldstein1999). To the best of our knowledge, no experimental data of the free-stream transverse spectrum exist. These data would allow for a better understanding of the response of the boundary layer to the free-stream flow.
As our formulation considers quasi-steady components of the Klebanoff modes, more accurate models that would allow for comparison at any wavenumber and frequency should include perturbations at any value of the scaled wavenumber 
$\kappa _z$. The boundary-layer equations, valid near the leading edge where spanwise diffusion is negligible, should be solved for the cases with 
$\kappa _z\ll 1$. An evident complication is that the receptivity would then be dictated by the full free-stream spectrum (4.7), which is a combination of the streamwise and transverse spectra, and not only by the leading-order transverse spectrum (4.18). It also follows that velocity components of higher order (with respect to the frequency), such as those appearing in (4.5a–c) and the order-one components studied by Wu & Dong (Reference Wu and Dong2016), would have to be taken into account. These improvements could lead to better agreement between the theoretical results and the experimental data at small 
$k_x \delta _d$ and small 
$k_z$ in figure 4.
In our analysis, only a mild effect of free-stream nonlinearity has been included by modelling the streamwise dependency of the free-stream spectrum, along similar lines to the nonlinear model in § 7.2 of Leib et al. (Reference Leib, Wundrow and Goldstein1999). If the streamwise dependency of the free-stream spectrum had not been accounted for, the free-stream decay would have been exponential because dictated by a linearized dynamics and it would not have been representative of realistic turbulence generated by a grid in a wind tunnel (Batchelor Reference Batchelor1953). Lifting the assumption of low-amplitude disturbances would lead to a better understanding of the boundary-layer response to the free-stream perturbation flow during the nonlinear stages of transition, which may involve secondary instability and the formation of turbulent spots. An interesting line of research would be the quantitative comparison between such nonlinear receptivity results and experimental data during transition, such as those obtained, for example, by Verdoya et al. (Reference Verdoya, Dellacasagrande, Barsi, Lengani and Simoni2022).
Acknowledgements
I would like to express my deepest gratitude to Dr M.E. Goldstein. I give him full acknowledgement for initiating this project, for proposing to use the quasi-steady asymptotic solution of the boundary-region equations, and for carrying out the initial analysis on the scaling of 
$\hat E_\alpha$ and 
$\hat E_\beta$. I thank him for kindly suggesting that I should work on this problem. I am also grateful to Dr E. Marensi for her precious comments and suggestions. I also wish to acknowledge the support of the US Air Force through the AFOSR grant FA8655-21-1-7008 (International Program Office – Dr D. Smith) and EPSRC through grant no. EP/T01167X/1.
Declaration of interests
The author reports no conflict of interest.
Appendix A. Numerical procedures
The boundary-region equations, given by (5.2)–(5.5) on p. 180 in LWG99 and complemented by the free-stream and initial boundary conditions given by (5.28)–(5.31) on p. 183 and (5.25)–(5.27) on p. 182 in LWG99, are solved numerically. As the equations are parabolic along the streamwise direction, a streamwise marching scheme is employed. As shown in figure 7, a second-order implicit finite-difference scheme, central in 
$\eta$ and backward in 
$\bar {x}$, is adopted, where the derivatives of a velocity component are expressed as
\begin{align} \frac{\partial q}{\partial \eta} = \frac{q_{j+1} - q_{j-1}}{2\,{\rm \Delta} \eta},\quad \frac{\partial^2 q}{\partial \eta^2} = \frac{q_{j+1} - 2 q_{j} + q_{j-1}}{({\rm \Delta} \eta)^2},\quad \frac{\partial q}{\partial \bar{x}} = \frac{\tfrac{3}{2} q_{i,j} - 2 q_{i-1,j} + \tfrac{1}{2} q_{i-2,j}}{{\rm \Delta} \bar{x}}. \end{align}
If the pressure is computed on the same grid as the velocity components, a pressure decoupling phenomenon occurs. Therefore, the pressure is computed on a grid staggered in 
$\eta$:
\begin{equation} p = \frac{p_{j+1} + p_{j}}{2},\quad \frac{\partial p}{\partial \eta} =\frac{p_{j+1} - p_{j}}{{\rm \Delta} \eta}. \end{equation}
The pressure at the wall does not have to be specified and is calculated a posteriori by solving the 
$z$-momentum equation at 
$\eta = 0$. Due to the linearity of the equations, the system is in the form 
$\boldsymbol {A} \boldsymbol {x} = \boldsymbol {b}$. For a grid with 
$N$ points along 
$\eta$, 
$\boldsymbol {A}$ is an 
$(N-2) \times (N-2)$ block-tridiagonal matrix where each block is a 
$4 \times 4$ matrix associated with the four unknowns 
$\{\bar {u}, \bar {v}, \bar {w}, \bar {p}\}$. Therefore, the wall-normal index 
$j$ of the vectors and matrix runs from 1 to 
$N-2$. The numerical procedure used to solve the linear system is found in Cebeci (Reference Cebeci2002) on pp. 260–264.
Figure 7. Sketch of the regular grid (black circles) and staggered grid (grey circles) used for the numerical scheme, adapted from Viaro & Ricco (Reference Viaro and Ricco2019). BC stands for ‘boundary conditions’.
