Asymptotic approximation of the eigenvalues and the eigenfunctions for the Orr-Sommerfeld equation on infinite intervals

Asymptotic eigenvalues and eigenfunctions for the Orr-Sommerfeld equation in two and three dimensional incompressible flows on an infinite domain and on a semi-infinite domain are obtained. Two configurations are considered, one in which a short-wave limit approximation is used, and another in which a long-wave limit approximation is used. In the short-wave limit, WKB methods are utilized to estimate the eigenvalues, and the eigenfunctions are approximated in terms Green's functions. The procedure consists of transforming the Orr-Sommerfeld equation into a system of two second order ordinary differential equations for which eigenvalues and eigenfunctions can be approximated. The approximated eigenvalues can, for instance, be used as a starting point in predicting transitions in boundary layers with computer simulations (computational fluid dynamics). In the long-wave limit approximation, solutions are expressed in terms of generalized hypergeometric functions.


Introduction
Stability and transition in shear flows and in boundary layers are mechanisms which need to be understood importantly for applications in mechanical and aerospace Engineering and in the atmospheric sciences [10,17]. Both theoretical and experimental studies have been carried out by different researchers in order to improve our knowledge on the properties of mechanisms that govern transitions and instabilities in fluid flows and in boundary layers, and so, different types of instabilities that include linear inviscid, viscous and nonlinear instabilities, and transitions in fluid flows and in boundary layers were investigated and many more [7,9,15,17] and references therein.
With that goal in mind, analytical and numerical methods have been used to solve the Orr-Sommerfeld equation that governs the mechanisms of linear stabilities of fluid flows on finite domains (channels) , semi-infinite domains (e.g. boundary layers) and infinite domains (e.g. wakes) [7,9,17]. The task becomes even more complicated when nonlinearities and turbulences are taken into consideration. Due to the evolution of computers, this may be accomplished by performing computer simulations. However, obtaining initial input that allow computer simulations to converge to the correct solutions remains an important challenge.
For instance, Gregory et al. [8] performed numerical simulations and analyzed the linear inviscid stability of the three-dimensional boundary layer and applied it to the rotating disk flow but their results were biased. Brown [3] extended the work by Gregory et al. [8] and included the viscous effects and applied the Orr-Sommerfeld equation to the rotating disc and sweptback wing. He used the temporal instability theory but his results did not match the observed values. Cebeci and Stewartson [5] used the spatial stability theory and solved the Orr-Sommerfeld equation on rotating disc profiles and their results were almost the same as Brown [3]. Cooper and Carpenter [6] investigated the stability of the rotating-disc boundary-layer flow over a compliant wall, and analyzed the the so called Type I and II instabilities. Reid [16] derived an exact solution to the Orr-Sommerfeld equation for the plane Couette flow. Walker et al. [18] investigated a physically-based computational technique intended to estimate an initial guess for complex values of the wavenumber of a disturbance leading to the solution of OrrSommerfeld equation in a shear flow on a semi-infinite domain. However, the problem of stability and transition in shear flows and in boundary layers still is a challenging problem [15,17].
Here, a different procedure which is rather analytical than numerical as in Walker et al. [18] is proposed. Two approaches are considered, the shortwave limit approximation and the long-wave limit approximation. As it will be shortly shown (section 3), in the short-wave limit approximation, the wavelengths in the spanewise and streamwise directions are relatively short while in the long-wave limit approximation, on the other hand, the wavelengths in the spanewise and streamwise directions are relatively long.
In the short-wave approximation, the Orr-Sommerfeld equation is written as a system of two second order ordinary differential equations, and thus the eigenvalues are approximated using the WKB method and the eigenfunctions are approximated in terms of Green's functions. Their corresponding outer solutions are investigated since their properties (e.g. behavior at infinity) are easier to analyze than those of the solutions obtained by means of Green's functions. For instance, instead of solving the Orr-Sommerfeld equation in the limit of infinite Reynolds numbers and using the fact that the mean flow velocityŪ (y) = 0 in the boundary layer near the wall as described in Schmid and Henningson [17], where critical layers and singularities are introduced in the equations while the Orr-Sommerfeld equation does not have any, we rather assume thatŪ (y) = y, 0 < y < ∞, approximate the eigenvalues using WKB methods and derive the approximate solutions in terms of Green's functions.
On one hand, the advantage of this procedure is that it works for any value of the Reynolds number, while on another hand, it can always be used to approximate eigenvalues and eigenfunctions of the Orr-Sommerfeld equation for any type of the mean flow velocity profile in three dimensions (3D) that can be approximated in two dimensions (2D) using the Squires' Theorem (transformation) asŪ(y) = ay 2 + by + c, where a, b and c are arbitrary constants. The approximated eigenvalues may be used as a starting point in predicting transitions in shear flow (e.g. boundary layers) in two or three dimensions with computer simulations (CFD).
In the long-wave limit, the Orr-Sommerfeld is reduced to a second order ordinary differential equation, and the solutions are written in terms of the generalized hypergeometric function. The obtained results can be usefull in geophysical fluid dynamics (GFD) where the size of the fluid disturbances can be of the same order as the radius of the earth or greater, for example, in the atmospheric boundary layer [10].

The three-dimensional linear stability model
We consider an incompressible parallel flow in three dimensions with velocity components w(x, y, z, t) =w(y) + w ′ (x, y, z, t), and pressure p(x, y, z, t) =p(y) + p ′ (x, y, z, t), where the terms with a bar represent the mean flow quantities, while the terms with a prime represent small perturbation quantities or waves. We then substitute then (1), (2), (3) and (4) in the Navier-Stokes equation and neglect nonlinear terms (products of perturbation quantities). This yields the linearized Navier-Stokes equations, and where ν is the fluid kinematic viscosity. The full model describing the linear stability of the three-dimensional incompressible fluid flow comprises equations (5), (6) and (7) complemented by the continuity equation [17] ∂u ′ ∂x Applying the divergence operator to (5), (6) and (7) gives We now eliminate the pressure term in (6) by applying the Laplacian operator and combining with (9). We then obtain We make all variables in (10) non-dimensional with respect to a typical reference speed V , a typical length scale H in the y-direction and a typical length scale L in both the streamwise direction (x-direction) and spanewise direction (z-direction). Thus, we find that equation (10) does not change but the Laplacian operator takes the form where is the spatial aspect ratio, while on the other hand, the non-dimensional kinematic viscosity is given by where ν * is the dimensional kinematic viscosity and R is the Reynolds number. And the Reynolds number R can now be written in terms of the aspect ratio r 2 as 3. The Orr-Sommerfeld equation: from three dimensions to two dimensions We discuss the three dimensional model and reduce it to a two dimensional one. In stability theory of fluid flow, this procedure is known as Squire's theorem (or transformation) [17]. We consider that the perturbations are small wavelike perturbations propagating in xz-plane with amplitude φ(y), and assume these perturbations are plane waves. Thus where α and β are wavenumbers in streamwise and spanwise directions, and ω is the transient frequency of the waves. Substituting (15) in (10) gives where the subscript y stands for differentiation with respect to y, and φ, for instance, satisfies the boundary conditions φ(0) = φ(∞) = 0 and φ y (0) = φ y (∞) = 0. (17) in the boundary layer. Or φ satisfies the boundary conditions This boundary condition, for example, shall be important in investigating the stability of the two dimensional wake. We importantly observe that (16) is the famous Orr-Sommerfeld equation but slightly modified by the aspect ratio r 2 . Settingα = rα,β = rβ,R = rR andω = rω gives which is the Orr-Sommerfeld equation.
Moreover settingŪ =ū cos θ +v sin θ andŪ yy =ū yy cos θ +v yy sin θ in (19) gives the two dimensional Orr-Sommerfeld equation We consider two configurations as mentioned before. In one configuration, we assume that the spacial aspect ratio is large r ≫ 1. This configuration is representative for waves with relatively short wavelengths in the xz-plane since the wavelength in the streamwise direction λ x = 2π/α = 2π/(rα) and that in the spanewise λ z = 2π/β = 2π/(rβ) become shorter and shorter as the aspect ratio becomes larger and larger (r → ∞). And so, this is a short-wave limit approximation. In that case, it is possible to write the two-dimensional Orr-Sommerfeld equation (20) in terms of a system of two second ordinary differential equations. We will shortly see (section 4) that WKB methods can be used to approximate eigenvalues.
In the other configuration, the aspect ratio is small, 0 < r ≪ 1 (r → 0 + ). And so λ x and λ z become larger as the aspect ratio becomes small (r → 0 + ). Therefore, this configuration is representative for waves with long wavelength in the xz-plane. Hence, this a long-wave limit approximation, and we will illustrate it (section 5) with two examples in which analytical solutions can be obtained.

Short-wave limit approximation (r ≫ 1)
According to (14), the Reynolds number R is proportional to r 2 . Therefore, the last term in the coefficient for φ in (20) which is proportional to r 2 is negligible compared to the other terms which are proportional to r 4 . And so, it can be dropped in the short-wave limit configuration where r ≫ 1. We then have We note that the assumption r ≫ 1 is quite important (see section 4.1) since it helps us to write (21) in a form allowing us to make use of WKB methods to approximate solutions for (20). We further observe that (20) will become (21) ifŪ (y) is a linear function of y, e.g. Couette flow. Now let us consider the differential operator its expansion is Applying this differential operator to φ gives φ yyyy + (P + Q)φ yy + 2P y φ y + (P yy + QP )φ = 0.
Comparing this with (21) yields and This gives and where according to (14), χ = R/r 2 = LV /ν * . Hence, setting in (24) implies that Ψ has to satisfy We use the following theorem to establish the boundary conditions for (33). Proof.

(a)
We first observe that the solution φ for (32) is φ(y) = φ h (y) + φ p (y) = Ae rky + Be −rky + φ p (y), where A and B are constants and φ h and φ p are the homogeneous solution and particular solution respectively. Applying the method of undetermined coefficients to (32), the particular solution has to take the form φ p (y) = f (y)Ψ(y) where f (y) is a function that has to be chosen in order to make φ h and φ p independent. But the form of P (y) and that of Q(y) given by (29) and (30) respectively, indicate that φ h and φ p will always be independent for all y. And so φ p = εΨ(y) where ε is some constant. Moreover A must vanish in order φ(y) to be finite as y → ∞. And so φ(y) = Be −rky + εΨ(y). Therefore, A, B and ε are constants. And so B must vanish in order φ(y) to be finite as y → −∞. This gives φ(y) = Ae rky +εΨ(y). Therefore, φ(−∞) = lim y→−∞ (Be rky +εΨ(y)) = ε lim y→−∞ Ψ(y) = 0. Hence, Proof. Corollary 2's proof follows from Theorem 1. To prove Corollary 2, we set a = 0, A = B and let δ 1 → 0 in Theorem 1 so that we obtain Ψ(−∞) = Ψ(∞) = 0.

WKB approximation for eigenvalues on a semi-infinite domain in the
short wave-limit approximation Using Theorem 1, in the boundary layer, we solve and with boundary conditions and where δ 1 is a very small positive constant (δ 1 → 0 + ), and as before, Following Theorem 1, δ 2 may be set to zero. In that case, we can restrict the domain to R + = (0, ∞), and thus use the WKB method described in Appendix B to approximate the eigenvalues k. We write Q as, where λ(k, ω) = (iχω−k 2 )/iχk. Now, (35) can be rewritten as the Schrodinger equation [2], where E = iχkλ(k, ω), V (y) = iχkŪ (y) and ǫ 2 = 1/r 2 → 0 is a small constant since r 2 ≫ 1. It is shown in the Appendix B that using the WKB method, the eigenvalues do satisfy where the limits of integration y = 0 and y = b are the turning points of V − E. On substituting E = iχkλ(k, ω) and V (y) = iχkŪ (y) in (40) gives We consider two configurations and approximate the eigenvalues using WKB method [2]. In the first configuration (configuration 1), the background mean flow is given byŪ (y) = by + c, 0 < y < ∞, where b and c are constants, while in the second configuration (configuration 2), the background mean flow is given byŪ(y) = ax 2 + by + c, 0 < y < ∞, where a, b and c are constants. Hence, If the flow is steady (e.g. laminar boundary layer), ω = 0 and the eigenvalues k n satisfy which can be solved using basic numerical methods. In the special case wherē U(y) = y (b = 1 and c = 0), (44) can be explicitly solved to obtain If for example θ = π/6, then the streamwise wavenumber is while the spanewise wavenumber is

Green's function approximation for eigenfunctions
In this section, we approximate the eigenfunctions φ n by means of Green's functions. We first solve (35) for Ψ and then use Green's functions to solve the inhomogeneous equation (34). In that case, solutions to (34) will be given by where G(y, ξ) are Green's functions. Green's functions associated with the boundary value problem defined by (34) and (36) are given by (A.8) (see Appendix A).

Configuration 1 (r ≫ 1):Ū (y) = by + c, 0 < y < ∞
We consider that the background flow is linear,Ū = by + c, where b and c are constants as before. We first solve, for Ψ, subject to boundary conditions as before. We make the change of variable η = (ir 2 χkb) 1/3 (y −λ), wherẽ λ = −ir 2 χk[c − λ(k, ω)]/b. Rearranging terms gives the Airy equation [1] Ψ where D 1 and D 2 are constants [1]. A solution satisfying the boundary conditions (57) is thus given by which has eigen-solutions where H n are Hermite polynomials of n [1]. The boundary condition (63) shall be satisfied if and only if n is odd. In that case, solutions to (62)-(63) take the form where H 2m+1 are Hermite polynomials of 2m+1. To obain the eigen-solutions φ m (y), we apply (55) as as in configuration 1.

Outer solution approximation for the eigenfunctions
The analysis of the eigenfunctions obtained by means of Green's functions does not seem to be an easy task due to the complexities involved in the computations of the integrals. In order to understand the behaviors of these solutions we can instead look at the outer solutions. Outer solutions are valid when ǫ = 1/r → 0 is a small parameter and give insight into the behavior of the solution for large argument y ≫ 1 (Bender and Orszag [2]). In the short-wave limit approximation, the homogeneous solution is proportional to e −rky = e −ky/ǫ . Therefore, we expect the homogeneous solution to quickly vanish with y as ǫ → 0. In that case, the outer solution would be accurate even for small values of y of order O(ǫ/k), ǫ → 0. Hence, outer solutions shall be valid on R + as long as ǫ → 0.
For small enough ǫ, we obtain from (32) that P φ ∼ Ψ. This gives This means that the solution φ of the Orr-Sommerfeld equation (16) is driven by Ψ whenever ǫ → 0 or the aspect ratio is large (r ≫ 1) since the homogeneous solution which is proportional to e −rky = e ky/ǫ rapidly vanishes with y. Thus, the constant ε in Theorem 1 should be ε = r 2 /k 2 . We look at two special cases of configurations 1 and 2 (see section 4.2), for which the outer solutions and the eigenvalues can explicitly be obtained.
4.3.1. Case 1:Ū (y) = y, 0 < y < ∞ This case corresponds to the configuration 1 in section 4.2 with the constants b = 1 and c = 0. We then haveλ = λ in (60). We also consider that ω = 0 (laminar flow) and let the phase velocity angle orientation θ = 0 so that k = α the streamwise wavenumber. Using (60) and (67) gives the outer solution where λ n (see section 4.1) is given by with, see equation (45), Some plots of the outer eigen-solution φ n as a function of y are shown on Figure 1. It is seen that the amplitude of φ n increases as the Reynolds number R decreases.
(72) Some results are shown in Figure 2. It is seen that the amplitude of φ n increases as the Reynolds number R decreases as in case 1 whereŪ(y) = y.

Stability of the two-dimensional wake in the short-wave limit
We write our solutions in terms of the hypergeometric function 2 F 1 . For reference, we shall first define the generalized hypergeometric function. Functions of this type are also used in section 5.

Long-wave limit approximation (r → 0 + ) on a semi-infinite domain
In this section, we consider the semi-infinite domain 0 < y < ∞. As seen in section 4.3, the analysis of the solutions obtained by means of Green's functions is not an easy task, here we use a different procedure. However, in the long-wave limit approximation, WKB methods cannot be applied, but (20) can be reduced to a form that allows us to readily obtain solutions in terms of hypergeometric functions whose properties are known.
with boundary conditions in the boundary layer.
Here, we consider two velocity mean profiles, the linear velocity mean profileŪ (y) = by + c and the quadratic velocity mean profileŪ(y) = δy 2 + by + c and solve the boundary value problem (91) -(92). For the quadratic mean flow profile, we assume that δ is a small constant. This implies that U yy = δ is also small, and consequently the third term involvingŪ yy = δ may be dropped.

Configuration 4 (r → 0 + ):Ū(y) = δy 2 + by + c, δ = a small constant
We consider the mean flow profile given by U(y) = δy 2 + by + c, where δ is a small constant. We also consider that U yy = δ ∼ O(r 2 ) to make sure the third term in (83) is negligible compared to the other terms and can therefore be dropped. This gives φ yyyy − [2r 2 k 2 + ir 2 χk(δy 2 + by + c − ω/k)]φ yy = 0. (90) gives Hermite equation as in section 4.2. We then obtain, in terms of the Reynolds number, where as before H n are Hermite polynomials of order n.
We note that we have to chose Hermite polynomials with even order in order the solution to satisfy the boundary condition φ y (0) = φ y (∞) = 0, and then integrate (91) twice to approximate the eigen-solutions for (83)-(84). The dispersion relation can be approximated as λ(k m , ω m ) = 2m+(1/4), m = 0, 1, 2, · · · . In that case, in the laminar boundary layer, the eigenvalues can be approximated by while we have (93) Example 1. If m = 0, for example, H 0 (y) = 1. Then using Proposition 1 in Nijimbere [12] gives where the constant of integration is set to zero in order to satisfy the boundary condition φ y (0) = φ y (∞) = 0. Integrating (94), see equation (C.3) in Appendix C, gives where the constant of integration is set to zero in order to satisfy the boundary condition φ(0) = φ(∞) = 0.
Making non-dimensional all variables and parameters, a spacial aspect ratio was introduced. This mainly allowed us to consider two configurations, the short-wave limit approximation which is obtained by letting the aspect ratio taking large values, and the long-wave limit approximation in which the aspect ratio takes small values in the Orr-Sommerfeld equation. This also allowed us to utilize analytical and asymptotic methods to obtain asymptotic solutions of the Orr-Sommerfeld equation and their corresponding eigenvalues. Most importantly, the procedure used in the present paper works regardless of the value of the Reynolds number.
In the short-wave limit approximation, the dispersion relation, the asymptotic eigenvalues and their corresponding asymptotic eigenfunctions were derived for configurations where the velocity mean flow profiles can, using Squire's transformation, be represented either as linear function or as a quadratic function. The eigenvalues were approximated using WKB method. Asymptotic eigenvalues and asymptotic eigenfunctions were also derived for the two-dimensional wake. The eigenfunctions were written in terms of Green's functions, and their corresponding outer approximate solutions were obtained as well. The results showed that the amplitude of the wave become larger as the Reynolds number becomes small. This is in agreement with the fact that small viscosity induces viscous instabilities. In the long-wave limit approximation, solutions were derived in terms of hypergeometric functions whose properties are known.
Due to the evolution of computer technology, Computational Fluid Dynamics (CFD) should help more in predicting transition from laminar flows to turbulent flows in three-dimensional shear flows. But starting simulations within a good range of eigenvalues allowing convergence of simulations to correct solutions remains a challenging problem. For this reason, we have approximated eigenvalues in the present paper which may, for instance, be used as a starting point in CFD simulations.

Appendix A. Green's functions
We want to solve, using Green's functions, the equation In that case, Green's functions are some functions G(y, ξ) that solve the equation where δ is the delta Dirac function, and the functions G satisfy the boundary conditions Green's functions are given by   in terms of the aspect ratio.

Appendix C. Some useful integrals
Some useful integrals involving the generalized hypergeometric function (73) are evaluated here. The method used is similar to that in Nijimbere [12].