The generation of Tollmien–Schlichting waves by free-stream turbulence in transonic flow

This paper studies the generation of Tollmien–Schlichting waves by free-stream turbulence in transonic flow over a half-infinite flat plate with a roughness element using an asymptotic approach. It is assumed that the Reynolds number (denoted Re) is large, and that the free-stream turbulence is uniform so it can be modelled as vorticity waves. Close to the plate, a Blasius boundary layer forms at a thickness of O(Re-1/2)\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$O(\mathrm{{Re}}^{-{1}/{2}})$$\end{document}, and a vorticity deformation layer is also present with thickness O(Re-1/4)\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$O(\mathrm{{Re}}^{-{1}/{4}})$$\end{document}. The report shows that there is no mechanism by which the vorticity waves can penetrate from the vorticity deformation layer into the classical boundary layer; therefore, a transitional layer is introduced between them in order to prevent a discontinuity in vorticity. The flow in the interaction region in the vicinity of the roughness element is then analysed using the triple-deck model for transonic flow. A novel asymptotic expansion is used to analyse the upper deck, which enables a viscous–inviscid interaction problem to be derived. In order to make analytical progress, the height of the roughness element is assumed to be small, and from this, we find an explicit formula for the receptivity coefficient of the Tollmien–Schlichting wave far downstream of the roughness.

Schematic of flow regime roughness elements found on aircraft wings. Solving the resulting linearised equations will enable an expression to be found for the amplitude of the Tollmien-Schlichting waves generated far downstream of the roughness element, at which point the objective of the current paper will have been met.

Problem formulation
Consider an ideal gas flow over a flat plate aligned with the free stream velocity. We assume that in the flow, there exist small amplitude vorticity waves (uniform free-stream turbulence). In addition, we shall have a roughness element on the flat plate a distant L from the leading edge of said plate. Placing the origin of our coordinate axes at the leading edge of the plate, it is natural to use (x,ŷ) as our coordinates. A schematic of the flow regime can be found in Fig. 1. We shall assume that the flow is two dimensional and that the body force is negligible. We define (û,v) as the velocity vector,p the pressure,ρ the density,ĥ the enthalpy, andμ the viscosity coefficient. We denote the free-stream velocity as U ∞ , and the free-stream pressure, density, and viscosity coefficient as p ∞ , ρ ∞ , and μ ∞ , respectively. In order to non-dimensionalise the compressible Navier-Stokes equations, we introduce the following variables: (2.1) Substituting these new variables back into our dimensional Navier-Stokes equations gives us the following nondimensional set of equations: Here, we have defined the Reynolds number (Re) and the free-stream Mach number (M ∞ ) as follows: respectively, where a ∞ = √ (γ p ∞ )/ρ ∞ is the speed of sound in the undisturbed flow upstream of the plate. In the analysis that follows, we will take the limit Re → ∞ and keep the Mach number as a constant which is close to 1.

Flow upstream of the roughness element
3.1 Free-stream flow (region 1) First, we shall consider the region far away from the plate (and upstream of the roughness). We look for solutions to the Navier-Stokes equations of the form where ε is a small parameter. Substituting these asymptotic expansions into the Navier-Stokes equations, we obtain It is clear that we can immediately eliminate h 1 from (3.2c) using (3.2e) giving We can look for wavelike solutions to this equation and (3.2a, 3.2b, 3.2d) of the form u 1 (ξ, y), v 1 (ξ, y), ρ 1 (ξ, y), Consequently, we write u 1 = U e i(αξ +βy) , v 1 = V e i(αξ +βy) , ρ 1 = Re i(αξ +βy) , p 1 = Pe i(αξ +βy) , (3.4) where α and β are the wave numbers in the x and y directions, respectively, and U , V , P, and R are constants. Substituting these representations into our equations results in a eigenvalue problem for c. In order for us to obtain non-trivial solutions, we require either (i) U = V = P = 0, R arbitrary or (ii) P = R = 0. Vorticity waves correspond to the case where R = P = 0, so this is the case we shall be concerned with. We will consider the case of one harmonic, so we define vorticity as follows: where ω 0 the amplitude of the wave. Substituting R = P = 0 into our governing equations, we see that the only equation not automatically satisfied is Combining this with our definition of vorticity, and without loss of generality setting ω 0 = 1, we obtain the following solutions for the velocity components: (3.7b) Of course due to Eq. (3.2e), we have that h 1 ≡ 0.

Vorticity deformation layer (region 2)
The work of [10] showed that introducing free-stream turbulence to flow over a flat plate causes a new layer of flow to be formed at y ∼ O(Re −1/4 ). With that in mind, we begin conducting analysis of this region by introducing the scaled variables: We look for perturbed solutions to the Euler equations of the form: (3.9) Eliminating pressure we find that Using this in conjunction with the definition of vorticity, we can see that vorticity does not change but is merely transported downstream with the mean free-stream velocity. As the vorticity upstream of the plate is given by (3.5), we can write where ξ = x − t. This equation suggests that u 2 and v 2 are functions of ξ and y only. Looking for solutions to the Euler equations of this form, we find Taking into account that there are no pressure perturbations in the free stream, we see that p 2 = 0. Consequently, the energy equation becomes In the previous section, we found that upstream of the leading edge of the plate h 1 = 0. By the above equation, we must, therefore, have h 2 ≡ 0 in the entire vorticity deformation layer, from which we can deduce that ρ 2 ≡ 0. The continuity equation is then reduced to (3.14) Using this, we can introduce a stream function ψ 2 which is found to be where C and D are constants of integration. Therefore, we have (3.16b) Clearly, we require C = 0 otherwise u 2 will grow exponentially with y. We can calculate D using the impermeability condition at the plate from which we can deduce that (3.17)

Steady boundary layer flow (region 3)
When the y coordinate becomes of order (Re −1/2 ), a viscous boundary layer forms close to the plate. As there are no leading order pressure perturbations in the vorticity deformation layer, free-stream turbulence does not penetrate into the boundary layer. Therefore, the flow in the boundary layer is steady. Assuming that that x remains an O(1) quantity, and that Y = Re 1/2 y is an O(1) quantity as Re → ∞, we look for solutions to the compressible Navier-Stokes equations of the form: (3.18) Substituting these representations into the governing equations and letting Re → ∞, we obtain (to leading order) The boundary conditions of the system are as follows. At the leading edge of the plate, we have (3.20) On the surface of the plate, we have the no-slip and impermeability conditions, which are expressed as follows: The matching condition at the outer edge of the boundary is We will also choose a thermally insulated wall, which has a constant temperature. We write this condition as follows: Self similar solutions to this system of the variable η = Y/ √ x exist. We seek solutions of the form (3.24) In the limit η → ∞, we have the solutions: where A, B, and C are constants.
In the limit η → 0, we consider the flow in the neighbourhood of the roughness element and Taylor expands our solutions in this limit around x = 1. We find to leading order (3.26) where λ is the usual Blasius constant.

Transition layer (region 4)
Outside of the boundary layer, we can see that there are no pressure perturbations to leading order. As a consequence, there is no mechanism by which our vorticity waves can penetrate into the boundary layer. However, we cannot have a discontinuity in vorticity at the outer edge of the boundary layer. This leads us to introduce a transition region between the classical boundary layer and the vorticity deformation layer, where the vorticity waves decay. We shall seek the solution to the Navier-Stokes equations in this region in the form: where we assume that the perturbations of pressure are small enough to be neglected. The perturbations of u and v are assumed to be of the same order as those in the vorticity deformation layer, and the order of σ is to be found. The coordinates (x, Y ) are the same as those in the boundary layer, and ξ is as in the vorticity deformation layer. We define Y through where we have that as Re → ∞, (Re) → ∞, δ(Re) → 0. Rearranging for Y , we have . (3.29) Substituting the expansions into the Navier-Stokes equations, we find that If this is constant, then the periodic solution in ξ would be a superposition of exponential functions to the power of Y . This would make it impossible to satisfy both boundary conditions for vorticity at each edge of the transition layer, so in order to prevent this, we set In what follows we will take our leading order expansions U 4 , V 4 , ρ 4 , and μ 4 from the outer edge of the Blasius boundary layer studied in the previous section. To avoid degeneration of the first term in the energy equation, we will set This reduces our equations to (3.33d) By the principle of least degeneration, we want to keep as many of the terms in these equations as possible. To this end, we set which is consistent with δ → 0 as Re → ∞, as well setting σ = ε. Assuming that the solutions are periodic in ξ and matching with the perturbations in both the Blasius boundary layer and vorticity deformation layer, we find where J 0 is the zeroth-order Bessel function of the first kind.

Flow in the interaction region
In order to analyse the flow around the roughness element, we shall use the same triple-deck structure as [8,9]. We require a time scale of t = O(Re −2/9 ) in order for the resonance conditions with Tollmien-Schlichting waves to be satisfied, and we shall further assume that the Mach number is of the form: where Q ∞ is the Kármán-Guderley parameter.

Lower deck (region 5)
The horizontal extent of the roughness element is assumed to be of the order x = O(Re −1/3 ), and the vertical extent of the lower deck is y = O(Re −11/18 ). We shall, therefore, assume that the shape of the roughness is of the form: It is natural to introduce re-scaled coordinates (x * ,y 5 ) and time t * as follows: We seek solutions to the Navier-Stokes equations in this region in the form: where we assume that close to the roughness element enthalpy, viscosity, and density are constants, and ε << Re −1/9 to prevent the expansion becoming disordered. The first terms in these expansions correspond to the steady flow solution, whereas those that are O(ε) are due to time-dependent perturbations. Substituting these expansions into the non-dimensional Navier-Stokes equations yields to O(1): Our boundary conditions for these equations come from our no-slip conditions on the roughness element, and matching u with flow upstream of the roughness as x * → −∞. These are Clearly, we can see that the pressure p 5 does not change across the lower deck with respect to the y-direction. Using the continuity equation we find that in the limit y 5 → ∞ where the function A * is still to be found. We can perform a similar analysis on the O(ε) terms in our expansions. The boundary conditions for this new problem arẽ Substituting into the Navier-Stokes equations, we obtain the O(ε) equations: In this case, we find that as y 5 → ∞ In order to understand, the physical meanings of A * andÃ * consider the streamline slope angle θ at the outer edge of region 5. Taking the limit as y 5 → ∞, we have Due to this relationship between θ , A * , andÃ * we call A * the O(1) displacement function andÃ * the O(ε) displacement function.

Middle deck (region 6)
The middle deck region is a continuation of the boundary layer upstream of the roughness. As a consequence, it has vertical extent y = O(Re −1/2 ) which leads us to re-scale our y-coordinate as follows: Relating this to the solution for the lower deck by considering the matching conditions for u and v as y 6 → 0, we can see that in this limit the solution in the middle deck is Consequently, we shall seek solutions in this region of the form: Substituting these expansions into the Navier-Stokes equations and using elimination to find u and v, we find the streamline slope angle to be where K (x * ),K (x * , t * ) are unknown functions. As θ has no dependence on y 6 we can match this with streamline slope angle from the lower deck yielding

Upper deck (region 7)
In order to analyse the upper deck, we are required to re-scale y again. As we have that the vertical extent of the upper deck is y = O(Re −5/18 ), we shall introduce a new variable y 7 defined as follows: We then look for solutions to the Navier-Stokes equations in this region of the form: where u α is a time-dependent perturbation due to vorticity which is known to be of the form e iαt * , and u st , v st , p st , ρ st , and h st are the solutions in the steady flow. Derivatives of u α will be cancelled out by next order expansions of the perturbation terms, which have been omitted in the asymptotic expansion of u. Substituting these expansions into the Navier-Stokes equations, we obtain to O(1): Using (4.19c) in (4.19d) and comparing this to (4.19a), we can see that this system is degenerate. It is, therefore, necessary to retain the next order terms, from which we obtain the equation We, therefore, have that the governing equation for p st is given by This is accompanied by the boundary condition that the perturbations are expected to decay far away from the roughness We now turn to the unsteady terms in our expansions (4.18). Our O(ε) equations are This system is again degenerate, so it is necessary to consider the next order equations which are Making use of Eqs. (4.19) and (4.23), this problem can be reduced to the following equation for the unsteady pressure perturbation in the upper tier and of course, we expect a similar decay of these perturbations as in the O(1) case. In order to close our viscous-inviscid interaction problem for p * , we require an additional boundary condition at y 7 = 0 for both steady and unsteady flow. To do this, we consider our streamline slope angle θ for region 6. So far, we have not considered the continuation of our transition layer into the triple-deck model. However, using a similar argument to [12], we can see that this is not necessary. The displacement thickness of the boundary layer has only been affected by the viscous lower deck, where the velocity was of small magnitude in comparison with other regions. We expect the velocity to be continuous across the transition layer of the triple deck, so it would be natural for it to be of the same order as the velocity in regions 6 and 7. Therefore, we can reasonably assume that the streamline slope angle will not be affected by this layer, as it is unaffected by region 6 which has the same magnitude of velocity and is much thicker. Consequently, the streamline slope angle at the bottom of region 7 is the same as it is in region 6 so we have From this, we can derive the boundary condition on p st and p * at y 7 = 0 which are respectively. For steady flow, we, therefore, have For unsteady perturbations, we have

Linear receptivity
In order to make analytic progress with our systems of equations (both steady and non-steady), we assume that the height of the roughness is small in comparison to the lower deck. We, therefore, write the equation of the roughness to be Using this assumption, we can linearise our equations for the lower deck. In this region, we, therefore, write the expansions for our velocities terms and pressure as follows: where σ is a small parameter. We shall also perform similar expansions for the displacement function A, and upper deck pressure p of the form: A = δ A * (x * ) + · · · + δσÃ * (x * , t * ) + · · · , (5.3a) p = δp st (x * , y 7 ) + · · · + δσ p * (x * , t * ) + · · · . (5.3b) We now substitute these representations into our equations for the lower deck. First, we shall consider the problem for steady flow.

Steady problem
Our linearised equations for steady flow to O(δ) read We also have the boundary conditions: We also need to relate P 5 to p st . As the pressure distribution does not change across the middle deck, we can relate P 5 back to the pressure in the lower deck by the equation: This system of equations will be solved in the usual way of transforming our variables into Fourier space. We define the Fourier transforms as follows: Defining a new variable z = (ik) 1 3 y 5 , we find the solution for subsonic flow (Q ∞ < 0) to be where p st is the Fourier transform of p st , and f is the Fourier transform of the roughness shape. To find the solution for supersonic flow (Q ∞ > 0), we have to take into account the causality condition, which yields the solutions At this point, we note that the branch cut of (ik) 1 3 in the complex k-plane is chosen in the usual way to be the positive imaginary semi-axis.

Unsteady problem
In order to solve our unsteady problem, we consider the terms of O(σ δ) in our linearised unsteady lower deck equations which give The boundary conditions for this problem arẽ In order to make analytical progress in finding solutions to these equations, we shall assume that our solutions are periodic in time, with wavenumber α corresponding to the wavenumber of u α (t). This means we can write them in the form: where the functions on the right-hand side have no time dependence, and c·c denotes the complex conjugate. Substituting these expressions into (5.10) and (5.11) yields Performing the same substitution in the boundary value problem (4.29) derived for the pressure in the upper deck p * gives In similar fashion to how we approached the steady problem, we begin by Fourier transforming our variables in x * , so that the above problem for p u becomes where p u and A u are the Fourier transforms of p u and A u , respectively, and p st was found in our analysis of the steady problem. The solution to this problem of course depends on the nature of flow (i.e. whether it is subsonic or supersonic) due to the fact that the forcing on the right-hand side of (5.15a) depends on the sign of Q ∞ . In addition, as we need to use an inverse Fourier transform to return to physical space, we require a solution for all k. Therefore, for each sign of Q ∞ , we have to split our solution into 3 different regions. For Q ∞ < 0, we have We begin by considering the case of subsonic flow.

Subsonic flow
In the case of subsonic flow (Q ∞ < 0), the solution to the system (5.15) is and K 1 is defined through Relating P u to p u using the condition (5.6), we have that We now return to our linearised lower deck equations. After performing Fourier transforms and eliminating v u (the Fourier transform of v u ) in a similar fashion to the analysis performed for the steady problem, we obtain the boundary value problem for the linearised lower deck variables in Fourier space We now make the transformations At this point, we note that we are required to introduce a branch cut of (ik) 1 3 , which we will again place along the positive imaginary semi-axis in the k-plane. Our boundary problem for u u is, therefore, The general form of the solution to the ODE above is where q 1 = ρ w λμ w , and Ai and Bi are the Airy functions. In order to satisfy (5.22c), we require D = 0; otherwise, the derivative of u u will be infinite while u u is tending to a finite value which is a contradiction. Differentiating this solution once and matching it with (5.22b), we have In addition, we can integrate (5.23) and considering (5.22d), we obtain In the limit ζ → ∞, we can deduce from (5.22c) that Eliminating C using (5.24), we have that A u can be expressed in the following form: where we have defined the parameters: Substituting this into (5.19), we find It is now necessary to return to physical space. Performing an inverse Fourier transform for P u yields In order to evaluate this integral, we are going to deform the integration contour into the complex plane. In preparation for this, we shall evaluate the poles of the integrand. These occur when we have This is a dispersion relation for k. In the limit of k → 0, if we assume that ζ 0 remains finite (i.e. α ∼ k 2/3 ), then this equation becomes Similar to Eq. (5.31), this equation has infinite roots. It is known from properties of the Airy function that these come in complex conjugate pairs and all lie in the negative real half of the ζ 0 -plane. It can be seen that the first root of (5.31) behaves differently in the ζ 0 -plane. As k → ∞, instead of tending to a finite point which is a root of (5.32), it instead tends to infinity. This root also behaves differently in the k-plane. As α goes from zero to infinity, all other roots remain in the second quadrant of the k-plane for all α. The first root, however, starts in the second quadrant and crosses the negative real semi-axis at a critical point k * < 0 which has a corresponding critical frequency α * . For all α > α * , this root remains in the third quadrant. This is the root which corresponds to the wavenumber of the Tollmien-Schlichting wave. Trajectories of the first five roots in both the k and ζ 0 planes can be found in Fig.  2, where the root corresponding to the Tollmien-Schlichting wave is in bold.
The process used to calculate the behaviour of the roots as α increased was as follows. First, Eq. (5.30) is written in terms of ζ 0 only, using (5.21). Then, Eq. (5.31) is used to find the position of each root at α = 0. For each subsequent value of α, the previous value of ζ 0 was used as an initial guess, and the dispersion relation was solved using Newton Iterations.
Returning to (5.29), we deform our contour of integration in the k-plane as shown in Fig. 3. We have had to split our contours into the regions of k < 0 and k > 0 due to the branch cut we introduced along the positive imaginary semi-axis. We, therefore, split our integral up in the following way:  , Here, it is important to note that the receptivity coefficient is a function of the Kármán-Guderley parameter Q ∞ and frequency α only, and does not depend on the roughness shape. Numerical calculations of the receptivity coefficient can be found in Fig. 4.

Supersonic flow
For supersonic flow, the analysis proceeds in much the same way. The solution to the problem (5.15) in this case is where is now defined to be and K 2 is Performing an identical analysis to the subsonic case, we find that for Q ∞ > 0  where in this case, the receptivity coefficient is with the coefficients having the subscript defined through While similar to the subsonic case, we note that the receptivity coefficient for supersonic flow is not identical to (5.36). The results of the numerical calculations for receptivity coefficient in supersonic flow can be found in Fig.  6.

Summary
This paper studied the receptivity of the Blasius boundary layer to free-stream turbulence in transonic flow. As the receptivity analysis is local, it is, therefore, applicable to an arbitrary boundary layer, for example, a wing surface. We have seen that, unlike acoustic waves, vorticity waves do not produce pressure perturbations to leading order in the free stream or the vorticity deformation layer. Therefore, the vorticity waves cannot penetrate into the boundary layer, so the process by which Tollmien-Schlichting waves are generated differs greatly from the acoustic case. In the present case, Tollmien-Schlichting waves arise due to the interaction between the free-stream turbulence and a roughness element on the plate surface. In order for this interaction to take place, we required that the strict double resonance condition associated with Tollmien-Schlichting waves be followed, i.e. the frequency of the vorticity wave is close to the frequency of the Tollmien-Schlichting wave, and in the Fourier spectrum of the perturbations produced by the roughness element is a wave number matching the wave number of the Tollmien-Schlichting wave. Unlike in the case of the acoustic wave, the interaction of the vorticity wave with steady perturbations (produced by the roughness element) takes place in the upper tier of the triple-deck structure. After first non-dimensionalising the governing equations, solutions in the flow upstream of the roughness element were derived. Due to the lack of leading order pressure perturbations, it was necessary to introduce a transition layer at the outer edge of the boundary layer in which the velocity perturbations caused by the free-stream turbulence decayed. The analysis of the boundary layer was, therefore, simple, as there was no need to introduce a Stokes layer close to the plate as is necessary in the case of acoustic waves.
To analyse the flow around the roughness, we used the scalings corresponding to transonic viscous-invisicid interaction theory, which enabled equations governing the velocity and pressure distributions in the three regions of the triple-deck structure to be derived. We note at this point that while the analysis performed on the lower and middle decks is similar to previous research, the inclusion of free-stream turbulence meant that the asymptotic expansions of the upper deck were novel. The external perturbations present in this deck interacted with the perturbations caused by the roughness element to produce the Tollmien-Schlichting wave.
When finding the amplitude of the Tollmien-Schlichting wave, it was necessary to assume that the height of the roughness element was small compared to LRe − 11 18 in order to make analytic progress. This assumption allowed the equations for the lower deck of the triple-deck model to be linearised, so the resulting viscous-inviscid interaction problems could be solved analytically. The behaviour of Tollmien-Schlichting waves far downstream of the roughness element is of primary interest when modelling laminar-turbulent transition as noted in the Introduction. By first solving the steady problem, and then assuming the unsteady solutions were periodic in time before making use of Fourier transforms, elementary results from complex analysis, and Watson's Lemma, a formula for the amplitude of Tollmien-Schlichting waves in this limit was found, depending only on the frequency of the vorticity wave and the Kármán-Guderley parameter. It is interesting that unlike the acoustic case studied by [7], the receptivity coefficient did depend on the sign of the Kármán-Guderley parameter.
To conclude, in this paper, we have derived an explicit equation for the amplitude of Tollmien-Schlichting waves generated by free-stream turbulence in transonic flow far downstream of roughness. This first required solving the Navier-Stokes equations in the various regions upstream of the roughness, before introducing the triple-deck model of [8,9], the use of which allowed the equations governing the flow around the roughness element to be derived. The receptivity coefficient for both subsonic and supersonic flow were found and calculated numerically, and depended only on the frequency of the vorticity wave and the Kármán-Guderley parameter.
Of course, the generation of the Tollmien-Schlichting waves can be also studied by means of the direct numerical simulations of the Navier-Stokes equations. However, a clear advantage of the asymptotic theory presented here is that it provides an explicit dependence of the amplitude of the Tollmien-Schlichting wave on the shape of the wall roughness. We see that in order to delay the laminar-turbulent transition, one needs to minimise the Fourier transform of the wall roughness function at a particular value of the wave number that corresponds to the Tollmien-Schlichting wave. This result may be directly used to improve the wing design. For example, of particular interest is a small gap between the nasal retractable section of the wing with the main body of the wing. While it cannot be avoided without increasing the landing speed of an aircraft, the effect of this gap can be minimised by simply changing the gap width.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.