1. Introduction

Despite the intensive development of computer technologies and numerical methods for the Navier-Stokes and Reynolds equations, problems of the three-dimensional boundary layer are of significant interest in the fluid dynamics. So far these problems have been little studied as a result of objective difficulties related with the large dimensionality and complexity of equations. Therefore, analytic results in this field can play an important role in the depth understanding of fluid dynamics phenomena and their study. In this part, some modern results in the three-dimensional boundary layer theory are discussed.

The small perturbation theory for inviscid flows is well developed and widely applied to estimate aerodynamic characteristics of real flight apparatus. Also it has been attempted to develop such theory for the boundary layer [1]. However, the zero approximation (“flat plate” approximation, zero cross-flow approximation) only given a rational contribution and were used in calculations. Equations for perturbations were complex. They required a numerical solution that was not much simpler than the full equation system. Father investigations of three-dimensional effects in the boundary layer theory became possible only after developments of computers with the enough power, numerical methods, and turbulence models [2].

Another approach was developed on the base of the rational perturbation theory including the first-order approximation [3, 4, 5, 6, 7, 8, 9, 10] for some class of flows, such as flows over aircraft wings and fuselages at small angle of attack, which have high importance as for the theory and the practice. In this case, zero-order approximation functions do not depend on the cross coordinate. Equations of the first-order approximation reduce to a two-dimensional system by introducing a new variable. The cross coordinate is included to this system as a parameter. This property of the self-similarity simplifies the solution procedures allowing to apply two-dimensional numerical methods and to reduce computing resources.

The singularity in the solution of 2D steady boundary layer equation is well known as the separation. Singularities arising in solutions of unsteady or 3D laminar boundary layer equations are not related directly with the flow separation and are slightly studied due to difficulties of analytical investigations of complex equations and uncertainty of numerical result treatments. However, this task is of interest for the mathematical physics and for numerical modeling of aerodynamic applications. For the first time, a singularity was found in the solution of 2D unsteady BL equations for the flow around the flat plate impulsively set into motion [12]. The singularity of the similar type was discovered on the side edge of a quarter flat plate in a uniform freestream [13] and at a collision of two jets [14]. In Ref. [15], necessary conditions were formulated for a singularity formation in self-similar solutions of the unsteady model and 3D incompressible laminar boundary layers on a flat surface with pressure gradients. Sufficient conditions and singularity types were not studied, and real flow conditions were not considered. Singularities of numerical solutions (the nonuniqueness or the absence of a solution) were found for the laminar boundary layer in the leeward symmetry plane on a round cone at incidence [16, 17, 18]. Similar results were obtained inside the computation region of the 3D turbulent boundary layer on the swept wing [19]. The singular behavior of boundary layer characteristics (the skin friction tends to the infinity in the symmetry plan) was found for the boundary layer on the small span delta wing [8, 10]. The explanation of these phenomena was found on the base of analytical solutions of laminar boundary layer equations on conical surfaces [10, 21, 22, 23, 24]. The asymptotic flow structure on the base of Navier-Stocks equations in the singularity vicinity is constructed.

The problem of the flow separation control using plasma actuators on the base of the electrical discharge is assumed as a perspective aerodynamic instrument [26, 27, 28]. It is considered as a one method for the control of the separated flow asymmetry near the nose part of aircrafts. The problem was complicated by the absence of an adequate model for the boundary layer-discharge interaction and a criterion for flow asymmetry arising. The use as a criterion numerical results and experimental data is restricted as a result of the high sensitivity of the asymmetry origin to different parameters [29]. Solution of these problems was obtained with the development of new models [30, 31, 32, 33, 34].

2. Small perturbation theory for three-dimensional boundary layer

As follows from the cross-flow impulse equation in biorthogonal coordinates [2], the necessary conditions for a small cross velocity (|w| << 1) are the relations

1 λ H 2 ∂ p ∂ z ∼ cos θ ∼ k 1 ∼ ε λ << 1 . E1

The small parameter ε characterizes the gradient of the pressure p t x z with respect to transverse nondimensional coordinate z; t and x here are dimensionless time and longitudinal coordinate, λ is body span, H 2 is metric coefficient, k 1 is longitudinal coordinate line curvature, and θ is the angle between coordinate lines on the body surface. Using these conditions flow parameters in the 3D boundary layer are presented by asymptotic expansions:

v w = v w 0 t s + ε v w 1 t s z , h w = h w 0 t s + ε h w 1 t s z , w = ε λ w 1 t s n z , V = u v h ρ μ κ = V 0 t s n + ε V 10 t s n z + ε 1 V 11 t s n z . E2

Here, s t x z is a dimensionless length of the coordinate line z = const measured from the critical point x c t z ; n is normal coordinate transformed with Dorodnitsyn transformation; v w t x z and h w t x z are blow (suction) velocity and the surface temperature; u and v , h , ρ , κ and μ are dimensionless longitudinal and normal velocities, enthalpy, density, thermal conductivity, and viscosity. The parameter ε 1 << 1 is not known a priori, it describes own flow perturbations inside the boundary layer. The following is found from the analysis of equations: for thin wings ε 1 = ε / λ 2 , for slightly asymmetric bodies ε characterizes the asymmetry and ε 1 = α ∗ / λ , where λ is relative body thickness [4, 5, 6, 7, 8, 9, 10].

To calculate boundary layer characteristics, the equation system for the composite solution incorporated in all terms of asymptotic expansion (2) was derived:

ρ ∂ u ∂ t + u − βw λ H 2 ∂ u ∂ s + v ∂ u ∂ n + ∂ p ∂ s = ∂ ∂ n μ ∂ u ∂ n , ρ ∂ h ∂ t + u − βw λ H 2 ∂ h ∂ s + v ∂ h ∂ n = ∂ ∂ n κ Pr ∂ h ∂ n + γ − 1 M 2 ∂ p ∂ t + u − βw λ H 2 ∂ p ∂ s + μ ∂ u ∂ n 2 , ρ ∂ w ∂ t + u ∂ w ∂ s + v ∂ w ∂ n + k 1 u 2 − k 2 uw = ∂ ∂ n μ ∂ w ∂ n − 1 λ H 2 ∂ p ∂ z − εβ ∂ p 1 ∂ s + cos θ ∂ p ∂ s , ρ ∂ q ∂ t + u ∂ q ∂ s + v ∂ q ∂ n + ∂ k 1 ∂ z u 2 − k 2 uq = ∂ ∂ n μ ∂ q ∂ n − − ∂ ∂ z 1 λ H 2 ∂ p ∂ z − εβ ∂ p 1 ∂ s − cos θ ∂ p ∂ s , q t s n z = ∂ w ∂ z , ∂ ρ ∂ t + ∂ ρu ∂ s + ∂ ρv ∂ n + 1 λ H 2 ρq − β ∂ ρw ∂ s − k 2 ρu = 0 , n = 0 : u = w = 0 , v = v w t x z , h = h w t x z ∂ h ∂ n = 0 , n → ∞ : u = u e t x z , w = w e t x z , h = h e t x z . E3

Eq. (3) is not true in the vicinity of the wing leading edge, where the pressure perturbation has the singularity. Using the asymptotic theory, singular regions near blunted and sharp leading edges were analyzed. It was found that the boundary layer in these regions is described by equations for the boundary layer on the sweep parabola or wedge. On a body the boundary layer begins in the critical point.

The system (2) was applied to the solution of different problems for wings and bodies [4, 5, 6, 7, 8, 9, 10]. To illustrate the developed approach in Figures 1 and 2, calculations of displacement thicknesses (Figure 1) and skin frictions (Figure 2) on the wind tunnel model of the US Air Force fighter TF-8A supercritical wing at Mach numbers M = 0.99 and 0.5 are presented. Solid lines correspond to solutions of Eq. (3) for the wing model ( Re = 2.246 ⋅ 10 6 ); dotted lines on Figure 4 are results for full scale wing ( Re = 2.58 ⋅ 10 7 ); symbols present solutions of full 3D boundary layer equations [11].

Figure 1.
Displacement thickness distributions on the model of supercritical wing; M = 0.99 , Re = 2.246 ⋅ 10 6 , and α ∗ = 3.12 ° .

Figure 2.
Skin friction distributions on the model of supercritical wing; M = 0.5 , α ∗ = 12.09 ° , Re = 2.246 ⋅ 10 6 (solid lines), and Re = 2.58 ⋅ 10 7 (dotted lines).

These figures demonstrate that the asymptotic solution very well reproduce numerical results as for the skin friction and for displacement thicknesses in the large parameter diapason.

3. Singularities in solutions of three-dimensional boundary layer equations

The laminar boundary layer problem on a thin round cone with the half apex angle δ c << 1 at the angle of attack α* depends on the parameter k = 4 α ∗ / 3 δ c only. Firstly, analytical results about singularities were obtained for outer BL part for a such cone. It is understood from previous works [15, 16, 17, 18, 20], the singularity can arise when two subcharacteristic (streamlines) families collided—this is a necessary condition. Such situation arises usually in the leeward symmetry (runoff) plane over a body of revolution at an angle of attack. Unusual properties in numerical solutions of self-similar equations in this plane for a round slender cone in supersonic freestreams were studied in many works due to the practical interest of the heat exchange on flying vehicles head parts [16, 17, 18, 20]. In this case, one parameter defines the flow. Two solutions were found in the windward symmetry (attachment) plane and at small angles of attack ( k ≤ k c ) in the leeward symmetry plane. In this plane, no solutions were obtained at moderate angles of attack ( k c ≤ k ≤ 2 / 3 ) and many solutions at larger incidences up to BL separation ( 2 / 3 ≤ k < 1 ). Full BL equation solutions with initial conditions in the windward symmetry plane fixed the violation of symmetry conditions in the runoff plane, a velocity jump through this plane in the angle of attack diapason, when the self-similar solution has been absent [10, 21]. The task for the cone was solved numerically on the base of parabolized Navier-Stokes equations, without the streamwise viscous diffusion [20]. However the problem is retained since the flow structure and reasons of unusual BL properties have not been explained.

Analytical solutions of full equations for the outer BL part on the slender round cone with initial conditions in the windward symmetry plane showed the singularity presence in the leeward symmetry plane of the logarithmic type at k = 1 / 3 and of a power type at k > 1 / 3 [10, 21]. It had been shown numerical solutions provided incorrect results near the singularity due to the accuracy loss. Similar but more complex results were obtained for arbitrary cones; they allow defining the sufficient conditions of the singularity arising [10, 22]. The asymptotic flow structure at large Reynolds number near the singularity on the base of Navier-Stokes equations was constructed, and analytical solutions in different asymptotic regions were obtained, which were matched with BL solutions. The analysis of the viscous-inviscid interaction region, in particular, revealed that the singularity can arise not only in self-similar but in full 3D BL equations [10, 22]. The theory showed that the singularity appearance relates with eigensolutions of the BL equations appearing near the runoff plane; it also explained numerical modeling results on the base of parabolized Navier-Stokes equations.

In the outer BL part, the theory gives the critical angle of attack for the singularity appearance k c = 1 / 3 . However calculations showed that this parameter is a function on numbers of Mach M and Prandtl Pr and the wall temperature h w , k c = k c M ∞ Pr h w [10]. This indicates that a singularity can arise in the near-wall region. The series decomposition of the near-wall solution in the runoff plane showed the presence of a parameter α , the linear combination of skin friction components, and the sign change of which leads to the change of the physical flow topology near this plane [24]. The analysis of BL equations in the near-wall region showed that α = 0 corresponds to the critical value k c , and it was confirmed by all published numerical calculations [16, 17, 18, 20]. In the runoff plane, the new power type singularity in solutions of full BL equations revealed that it is related with the eigensolutions appearing near this plane. Calculation results for BL on delta wing confirm the singularity presence.

3.1 Self-similar boundary layer on a cone

The 3D laminar boundary layer on a conical surface in the orthogonal coordinate system xyφ (Figure 3) is described by following self-similar equations and boundary conditions [10, 22]:

u yy = Awu φ + vu y + A 1 w u − w , w yy = Aww φ + vw y + w 2 3 u + Kw − h 2 3 + K , h yy = Awh φ + vh y − M e u y 2 + 3 2 A 1 w y 2 , ρh = 1 , y = ε 3 ρ e u e 2 x μ e ∫ 0 y ∗ ρd y ∗ l , Re = ε − 2 = ρ ∞ u ∞ l μ ∞ , f y = u , g y = w , v = − f − K − 1 2 A ln ρ e μ e / u e φ g − Ag φ , y = 0 : u = v = w = 0 , h = h w h y = 0 ; y = ∞ : u = w = h = 1 . E4

Figure 3.
Solutions of boundary layer equations (dotted lines) and parabolized Navier-Stokes equations (solid lines).

Equation coefficients are defined by expressions

M e φ = γ − 1 M ∞ 2 u e 2 h e , K φ = 2 w eφ 3 Ru e , A φ = 2 w e 3 Ru e , A 1 φ = 2 3 w e u e 2 E5

In these equations, to reduce formulas, Pr = 1 and the linear dependence of the viscosity on the temperature ( ρμ = 1 ) are assumed. Indexes y and φ denote derivatives with respect to the corresponding variables; x is the distance from the body nose along the generator referenced to the body length l ; y is the Dorodnitsyn variable; y ∗ is normal to the body surface; φ is the transversal coordinate, and it can be the polar angle for a round cone; f y φ and g y φ are longitudinal and transverse stream functions; v y φ is transformed normal velocity; and R φ is the metric coefficient. The density ρ , the enthalpy h , the viscosity μ , the longitudinal u , and transversal w velocities are referenced to the values at the outer boundary indexed by e , which are normalized to their freestream values indexed by ∞ ; they are functions of φ only. The transversal velocity on the outer boundary layer edge w e = 0 in the initial value plane (the attachment plane) φ = 0 , in which K 0 > 0 , and in the runoff plane φ = φ 1 , in which K φ 1 = − k < 0 , and two boundary layer parts that came from different sides of the attachment plane collided. For the round cone, φ 1 = π .

Eq. (4) is simplified for slender bodies since in this case, u e = ρ e = μ e = 1 , A 1 << 1 . Neglecting proportional to A 1 terms in (4), we obtain the Crocco integral for the enthalpy and momentum equations in the form

h = h w + h r u − 1 2 M e u 2 , h r = 1 − h w + 1 2 M e , M e = γ − 1 M 2 , v = − f + Kg + Ag φ , u yy = Awu φ + vu y , w yy = Aww φ + vw y + w 2 3 u + Kw − h 2 3 + K . E6

For the slender round cone with the apex half angle δ c << 1 at the angle of attack α ∗ , simple expressions for outer functions are

w e = 2 α ∗ sin φ , K φ = k cos φ , A φ = k sin φ , k = 4 α ∗ 3 δ c E7

3.2 Singularities in the outer boundary layer region

In the outer boundary layer region, y ≫ 1 , ( y is the Dorodnitsyn variable normal to the wall), flow functions are represented as [17, 36]

u = 1 + U η φ , w = 1 + W η φ , η = y − δ / a φ , h = 1 + H = 1 − 1 2 M e + h w − 1 U − 1 2 M e U 2 E8

Here δ φ is the displacement thickness defined by the equation of F. Moore [6], the function a φ is found from the local self-similarity condition, and U << 1 and W << 1 are velocity perturbations with respect to boundary conditions, which in the first-order approximation satisfy to equations [10, 21]

U ηη + η U η − aAU φ = 0 , W ηη + η W η − 2 3 a 3 2 AW φ + 1 + 3 K W = 2 3 ap φ U . E9

These equations have solutions:

U η φ = C 1 erfc η / 2 W η φ = − b φ U , W 1 η φ = − b φ U + B 1 k V η φ E10

Constants C 1 and B 1 are calculated from matching condition with a numerical solution inside the boundary layer. These solutions satisfy to initial conditions in the attachment plane and must tend to zero at η → ∞ . The function V η φ is the solution of the homogeneous equation for the cross-velocity perturbation, when the right-hand side equals to zero; it is expressed by Veber-Hermite functions [21]. The coefficient B 1 ∼ 1 / K 0 has the singularity at K 0 → 0 . For the round cone this limit corresponds to zero angle of attack; in this case, the analytical expression for W 1 η φ shows the presence of the power type singularity in the leeward plane φ = φ 1 [10, 21]. The first solution W η φ is regular in this limit, and its behavior is defined by functions a φ and b φ , which satisfy to equations [10, 21, 22]

w e b φ + 2 1 + M w eφ b = 2 pMw eφ , p φ = 1 + 1 + 3 2 K 1 2 M e + h w − 1 , w e a φ + 2 N + 1 w eφ a = 2 Nw eφ , N φ = 3 M φ = K − 1 . E11

Solutions of these equations with initial conditions in the attachment plane are represented in integral forms in the general case and have analytical expressions for the round cone [10, 21]. Their properties near the leeward plane, at ζ = φ 1 − φ << 1 , are represented by expressions

w e = 3 2 kRζ , k = − K φ 1 , R = R φ 1 , p 1 = p φ 1 , n = 3 m = − K − 1 φ 1 , m ≠ 1 : b = mp 1 m − 1 − b m ζ 2 m − 1 , m = 1 : b = − 2 p 1 ln ζ + b 1 , n ≠ 1 : a = n n − 1 + a n ζ 2 n − 1 , n = 1 : a = − 2 ln ζ + a 1 E12

Here a n and b m are known coefficients [10, 21]. These formulas are true for non-slender bodies also [10, 22, 23].

These results show the presence in the outer BL part of two singularity types in the leeward plane related with properties of functions a φ and b φ . For k < 1 the function U η ζ exists at ζ = 0 but reaches this limit irregularly; its behavior is studied analytically in details for the slender round cone [10, 21]. For k ≥ 1 the function U η ζ is singular at ζ → 0 since a ζ → ∞ and the BL thickness tend to infinity as a ζ : the logarithmic singularity type takes place at k = 1 , and it is of the power type at k > 1 . At k ≥ 1 the flow separation is observed in experimental and numerical studies; this phenomenon changes not only the outer part but also the inner boundary layer structure. It should be noted that such behavior of velocity viscous perturbations near the BL outer part at the separation development is a new property in the comparison with the 2D flows.

The function W η ζ has irregular but finite limit in the leeward plane for ζ → 0 at k < 1 / 3 . This limit is singular at k ≥ 1 / 3 : the singularity has the logarithmic or power type, if k = 1 / 3 or k > 1 / 3 . At 1 / 3 ≤ k < 1 the singularity is related with the behavior of cross-flow velocity only. This singularity leads to the longitudinal vortex component strengthening in the outer part of the viscous region. The singularity takes place, if the pressure gradient is negative ( k ≤ 2 / 3 ) or positive ( k > 2 / 3 ). It is formed by BL proper solutions, which have homogeneous conditions on both boundaries and arise near the runoff plane. The critical value k c = 1 / 3 for the outer BL part is undependable on the wall temperature and Mach and Prandtl numbers; however the considered singularities define the real flow structure near the leeward plane at k ≥ 1 / 3 [17, 36, 37].

3.3 Asymptotic flow structure near the singularity

Due to the irregularity of solutions already at k ≥ 1 / 6 ( m ≤ 2 ), the vortex boundary region near the runoff plane is formed with transverse dimension ζ ∼ ε 1 2 − m ; at m ∼ 1 this value is of the order of the BL thickness ε . In this region, the transverse diffusion is the effect of the first order, and to describe it we introduce the following variables:

ε 1 = 3 2 Re ρ e φ 1 u e φ 1 / μ e φ 1 − 1 2 . z = kx Rζ / ε 1 , u = u y z , h = h y z , w = w y z E13

Using these variables from Navier-Stokes equations at ζ ∼ ε 1 << 1 for this region, we derive self-similar equations, which in its outer part, at y ≫ 1 , reduce to the form

U yy + kU zz + 1 − k yU y + kzU z = 0 , W yy + kW zz + 1 − k yW y + 2 z + kz W z + 2 k m − 1 W + 2 3 p 1 U = 0 E14

For k < 1 these equations have the solution corresponding to the regular at k → 0 solution of BL equations:

U y z = C 1 erfc y 1 − k / 2 erf z / 2 , W = − B z C 1 erfc y 1 − k / 2 , B zz + 2 z + z B z − 2 m − 1 B = − 2 mp 1 F z , F z = erf z / 2 E15

The function B z is expressed by Kummer’s function Φ a b x [10, 22, 23]:

B = mp 1 B 0 z + B m Φ 1 − m 3 2 − 1 2 z 2 , B m = b m R kx / ε 1 2 1 − m . E16

B 0 z is a particular solution of the inhomogeneous equation; the coefficient B m is determined from matching condition.

In Figure 4, comparisons of solutions of BL (dotted lines) and Navier-Stokes (solid lines) equations for m = 1 / 2 (curves 1 and 2) and m = 1 (curves 3 and 4) are presented. It is seen that regular solutions of Navier-Stokes equations are converged quickly to singular solutions of BL equations.

Figure 4.
The general flow scheme and the coordinate system for a cone.

Another effect generated by the singularity at k ≥ 1 / 3 due to the BL growth at ζ → 0 is the viscous-inviscid interaction. This effect is important in the region, where the inviscid and induced cross velocities have same orders; this condition defines the transverse dimension of the region Δ φ and the velocity scale:

Δ φ ∼ ε x − 1 4 , w e ∼ kRu e ε x − 1 4 . E17

In this region, the flow has the two-layer structure. Assuming the potential flow in the outer inviscid region, the solution here is presented by the improper integral from the displacement thickness δ x s . In the boundary layer, the flow is described by full 3D equations:

s = Rζ ε , w e = 3 2 u e ε W e x s W e x s = − ks 1 + r , r = 4 m π ∂ ∂ x ∫ 0 ∞ δ x t dt s 2 − t 2 v = f + Kg + Ag s + 2 3 xf x , h = h w + h r u − 1 2 M e φ 1 u 2 u yy = W e wu s + vu y + 2 3 xuu x w yy = W e ww s + vw y + w 2 3 u + W es w − h 2 3 + W es + 2 3 xuw x E18

For these equations boundary conditions have the form (1). A solution of these equations will be matched with the boundary layer solution at s → ∞ . Initial conditions are needed at some streamwise location x = x 0 , which can be obtained from a solution of Navier-Stokes equations near the body nose; this feature does the problem more complicated. Obtained equations allow a self-similar solution for hypersonic flows at some additional assumptions.

The solution in the outer boundary layer part, at y ≫ 1 , is described by formulas

t = y / d x s , u = 1 + U x t s , w = 1 − c x s U U = C 1 erfc t / 2 , p 0 = 3 2 1 2 M 0 + h w − 1 1 + r sd s − 2 mxd x − 2 n − 1 − r s d = − 2 n 1 + r sc s − 2 mxc x − 2 m − 1 − r s c = − 2 m p 1 − qp 0 E19

Along characteristics ξ x s = const , which are streamlines of the inviscid flow, the equations for functions d = d ξ s and c = c ξ s are integrated. At s → 0 these functions are represented in the form

c = Cs L + m p 1 + p 0 r s m − 1 − r , L ξ s = m − 1 − r s 1 + r ; d = Ds I + n n − 1 − r s , I ξ s = n − 1 − r s 1 + r , C = b m ε m − 1 , D = a n ε n − 1 . E20

Coefficients C and D are obtained by matching d ξ s and c ξ s at s → ∞ in relation to a ζ and b ζ at ζ → 0 [17, 36]. The logarithmic singularity appears in these functions at I = 0 or L = 0 . At L ξ 0 < 0 or I ξ 0 < 0 , the singularity is of the power type.

Following from presented results, in contrast with the 2D separation, the viscous-inviscid interaction does not eliminate the singularity in 3D boundary layer; this effect moves only the critical value of k c .

3.4 Singularities in the boundary layer near-wall region

The singularity in the outer BL part gives the critical value k c = 1 / 3 , although calculations show k c = k c M ∞ Pr h w . This indicates on the possibility of singularity arising in the near-wall region. To study this possibility, at the first, we study the solution behavior of Eq. (6) at y << 1 in the runoff plane φ = φ 1 where the solution is presented in the form

τ 0 = du 0 0 dy , θ 0 = dw 0 0 dy , u 0 = τ 0 y + U 0 y , w 0 = θ 0 y + W 0 y , v 0 = − α y 2 − F 0 + kG 0 , α = 1 2 τ 0 − k θ 0 E21

Second terms of these decompositions can be presented by series

U 0 = F 0 y = ∑ i = 0 α i y i + 4 i + 4 ! , F 0 = ∑ i = 0 α i y i + 5 i + 5 ! , W 0 = G 0 y = ∑ i = 0 β i y i + 2 i + 2 ! , G 0 = ∑ i = 0 β i y i + 3 i + 3 ! . E22

First three coefficients of these series are defined by relations

α 0 = − 2 τ 0 α , α 1 = k τ 0 β 0 , α 2 = k τ 0 β 1 β 0 = − ph w , β 1 = − p τ 0 h r , β 2 = 1 3 τ 0 − 3 k θ 0 θ 0 + 2 pM e τ 0 2 E23

Using these decompositions we can study qualitatively a dependence of the flow structure near the runoff plane from parameters by analyzing the subcharacteristic behavior. The transformed normal to the body surface v and transverse w velocities at ζ << 1 and y << 1 in the first-order approximation are represented in the form

v = v 0 = − α y 2 − 1 6 k β 0 y 3 = − 1 6 k β 0 y 2 y + y c y c = − 6 α k β 0 = 6 α kph w , w = − w 0 = − k θ 0 ζy E24

In the plane ζ = 0 , the cross-flow velocity w = 0 due to symmetry conditions. Here two critical points, in which v = 0 , can be. The first point locates on the cone surface y = 0 , and the second one y = − y c appears in the physical space at α < 0 , if p > 0 ( k < 2 / 3 ), that corresponds to small angles of attack for the round cone and at α > 0 , if p < 0 . Commonly, the critical value of the cross-flow velocity gradient k c ≤ 1 / 3 corresponds to the negative cross-flow pressure gradient p > 0 , the transverse skin friction in this region θ 0 > 0 .

Using these expressions, the equation for the subcharacteristics is obtained in the form

y c dy y y + y c = β dζ ζ , β = α k θ 0 ; α ≠ 0 : y = y c y 0 s β y c + y 0 1 − s β , s = ζ ζ 0 , α = 0 : y = y 0 1 − y 0 d ln s , d = ph w 6 θ 0 E25

Here y 0 and z 0 define the initial point in the cross-plane.

The subcharacteristic behavior is shown in Figure 5a and b for p > 0 . At α > 0 velocities v < 0 and w < 0 ; the only critical point node is in the coordinate origin, and subcharacteristics go to it from the region ζ ≠ 0 (Figure 5a). At α = 0 y c = 0 and the point ζ = y = 0 is double critical point of the type saddle node: the saddle is in the lower half-plane, i.e., out of the physical space. The node is in the upper half-plane, and the subcharacteristic pattern retains the same as at α > 0 . At α < 0 the node drifts in the point ζ = 0 , y = y c > 0 , and the coordinate origin becomes by the saddle point (Figure 5b). In this case, at y > y c the normal velocity v < 0 and at 0 < y < y c v > 0 ; v = 0 on the line y = y c .

Figure 5.
Subcharacteristics in the cross-plane at α ≥ 0 (a) and α < 0 (b); p > 0 .

This analysis shows that at the parameter α sign change, the physical flow structure varies qualitatively, and the value α = 0 is a criterion of the new flow property appearance. It should be noted that in solutions of Navier-Stokes equations for similar problems near the coordinate origin z = y = 0 in the leeward symmetry plane, the streamwise-oriented vortex arises, and the flow is not described by the BL theory since the viscous diffusion inside the vortex is distributed along the radius from its axis, but not along the normal to the body surface. On the base of this qualitative analysis, it is supposed that the critical value k c h w M is defined by the relation

2 α k c = τ 0 k c − k c θ 0 k c = 0 E26

To support this hypothesis, equations for functions U 0 y and W 0 y are analyzed by substituting near-wall decompositions to Eq. (6). Considering functions U 0 y and W 0 y as perturbations, we can linearize resulting equations and obtain in the first-order approximation:

U 0 yy + α y 2 U 0 y + τ 0 F 0 − kG 0 = − ατ 0 y 2 , W 0 yy + α y 2 W 0 y − 2 3 τ 0 − 3 θ 0 yW 0 + θ 0 F 0 − kG 0 = β 0 + β 1 y + 1 2 β 2 y 2 + 2 3 θ 0 y − p h r − 2 M e τ 0 y U 0 E27

At y → 0 U 0 y and W 0 y are expressed by above series, and in order to match them with the solution of full Eq. (6) in the main BL part, it is required that these functions will grow at y → ∞ not faster than a power function. To study their solution behavior at y → ∞ and α ≠ 0 , we introduce the new variable:

ξ = − α y 3 / 3 , y = − 3 ξ / α 1 3 . E28

At the limit ξ → ∞ , previous equations are reduced in the first-order approximation to the form

ξ ∂ 2 U 0 ∂ ξ 2 + 2 3 − ξ ∂ U 0 ∂ ξ = − τ 0 3 3 ξ α 1 3 , c = 2 τ 0 − 3 k θ 0 9 α ξ ∂ 2 W 0 ∂ ξ 2 + 2 3 − ξ ∂ W 0 ∂ ξ + cW 0 = − β 1 3 α + β 2 6 α 3 ξ α 1 3 − 2 9 α θ 0 + 3 M e τ 0 U 0 E29

Solutions of these equations can be represented as

U 0 = A 00 ∫ 0 y e − 1 3 α s 3 ds + τ 0 3 ξ α 1 3 , W 0 = B 00 Φ − c 2 3 ξ + B 01 ξ 1 3 Φ 1 3 − c 4 3 ξ − 3 β 1 2 τ 0 − 3 k θ 0 − 3 β 2 τ 0 − 9 k θ 0 3 ξ α 1 3 − θ 0 + 3 M e τ 0 τ 0 − 3 k θ 0 U 00 E30

First terms of these expressions are solutions of homogeneous equations, with zero right-hand sides; A 00 , B 00 , and B 01 are constants; Φ a b x is Kummer’s degenerate hypergeometric function, which has asymptotes at ξ → ∞ :

α > 0 , ξ < 0 : Φ ∼ − ξ c ; α < 0 , ξ > 0 : Φ ∼ e ξ ξ c − 2 3 E31

Solutions grow exponentially at α < 0 and p > 0 ; they cannot be matched with the solution in the main BL part. Therefore, at these conditions a solution of BL equations cannot exist. This conclusion and also the criterion (26) for the boundary of the existing leeward symmetry plane solution are confirmed by numerical calculations for the slender round cone at an angle of attack [25, 26, 27, 28, 29, 30, 31, 32, 37], a part of which is presented in Figure 6. In this figure, symbols correspond to calculations of limit values α k c for the solution existing at different boundary conditions in the diapason of Mach numbers from 2 to ∞ at the Prandtl number 1 for different surface temperatures. At k < 1 / 3 data are grouped near the value α = 0 in accordance with the criterion (26). The data scatter is, apparently, due to the decrease of the calculation accuracy at the approach to the critical value k c and also with errors of data copying from papers. At k > 1 / 3 , all calculations are finished with α > 0 , since the solution existing in this region is determined by singularities in the outer BL part, but not in the near-wall region.

Figure 6.
The boundary of the solution existing in the leeward symmetry plane of the slender round cone at the angle of attack and Pr = 1 in the dependence of the critical value k c : ▲, [28]; ■, [29]; and ○, [37].

Then we consider the solution behavior of full BL equations in the near-wall region beside the runoff plane at ζ << 1 . 3D BL equations have the parabolic type, and their solution before the runoff plane knows nothing about the solution in this plane; however, in order for the first solution to move smoothly into the last one at α > 0 , the first will be locally self-similar. Due to this condition, the streamwise τ ζ and cross-flow θ ζ friction stresses and the self-similar variable η at ζ << 1 will be defined by expressions

τ ζ = τ 0 a ζ , θ z = θ 0 a ζ , η = y a ζ E32

The function a z at α ≥ 0 will satisfy to the condition a 0 = 1 . In this case, flow functions in the boundary layer near the wall can be represented in the form

f η ζ = a ζ τ 0 η 2 2 + F η ζ , u η ζ = f η = τ 0 η + U η ζ g η ζ = a ζ θ 0 η 2 2 + G η ζ , w η ζ = g η = θ 0 η + W η ζ v = a α − 1 2 θ 0 kζ a ζ a η 2 + F − kG 1 + kζ a ζ a − kζ G ζ − k ζη ζ W E33

Substituting these expressions to Eq. (6) and linearizing the result with respect to disturbances, we obtain the first-order approximation for the flow in the near-wall region beside the runoff plane:

U ηη + αη 2 U η + a 2 k θ 0 ζη U ζ + τ 0 F − kG 1 + a ζ a − kζ G ζ = − ατ 0 η 2 W ηη + αη 2 W η + a 2 k θ 0 ζη W ζ + θ 0 F − kG 1 + ζ a ζ a − kζ G ζ − 3 α c η W = − αθ 0 η 2 + a 2 β 0 + β 1 η + 1 2 β 3 η 2 + 2 3 θ 0 + 3 pM e τ 0 η − ph r U E34

Here β 3 = 2 3 τ 0 θ 0 − k θ 0 2 + pM e τ 0 2 . Due to local self-similarity at α ≥ 0 , we define the function a ζ as

α a 2 − 1 2 k θ 0 ζ aa ζ = α , a 2 = 1 + C ζ q , q = 4 α k θ 0 E35

The constant C is found from a comparison with numerical calculations. It follows from this relation at α ≥ 0 and q < 2 the solution of Eq. (6) in the near-wall region at ζ << 1 can find in the form of the series:

F η ζ = F 0 η + ζ q F q η + … , U η ζ = U 0 η + ζ q U q η + … , G η ζ = G 0 η + ζ q G q η + … , W η ζ = W 0 η + ζ q W q η + … , E36

The first term of this expansion is the solution for the runoff plane but depends on the self-similar variable. Second terms define the proper solution of BL Eq. (6) at ζ << 1 , which at η → ∞ , has the form [37]

ξ = − αη 3 3 , U q ξ = A q 0 Φ 4 3 2 3 ξ + A q 1 ξ 1 3 Φ 5 3 4 3 ξ W q ξ = B q 0 Φ 4 3 − c 2 3 ξ + B q 1 ξ 1 3 Φ 5 3 − c 4 3 ξ + 9 β 1 2 τ 0 + 3 β 2 11 3 τ 0 − k θ 0 3 ξ α 1 3 + θ 0 + 3 M e τ 0 2 τ 0 U 0 − θ 0 + 3 M e τ 0 τ 0 − 3 k θ 0 U q . E37

Here A q 0 , A q 1 , B q 0 , and B q 1 are constants. These relations show that the proper solution in near-wall BL region near the runoff plane is nonzero. It is irregular at α ≥ 0 and it is singular at α < 0 . The logarithmic singularity is not in this case, and the solution of BL equations exists at the critical value k c in contrast to the outer region.

In the work of [15], at the analysis of perturbations in the boundary layer related with the angle of attack, it was found that they lead to infinite disturbances in the symmetry plane, although equations have no visible singularities contained. In this case, the first-order approximation is described by the Blasius solution for the delta flat plate. In Figure 7, dimensionless longitudinal and transverse skin friction distributions f 1 ʹʹ z and g 1 ʹʹ z , induced by the second order BL approximation (Figure 7a) and the angle of attack (Figure 7b) are presented in dependence on transverse coordinate z = 1 − Z / X , where X and Z are Cartesian streamwise and transverse coordinates. By approaching the symmetry plane ( z = 1 ), skin friction perturbations infinitely grow. Detailed investigation of equations for these functions showed that in these cases singularities take place as in the near-wall and outer BL parts. In the outer part, the singularity corresponds to values of the parameter m = 3/4 and 7/8 in relation to cases а and b, respectively. The longitudinal velocity perturbation singularity is related only with the near-wall singularity.

Figure 7.
Skin friction distributions on the small aspect ratio delta wing at М = 2 related with (a) second BL approximation and (b) angle of attack.

Near-wall singularities generate the flow structure including three asymptotic sublayers describing the viscous-inviscid interaction similar as near the 2D separation point. However, the viscous-inviscid interaction is not enough to remove the singularity of the obtained type. Near the wall sublayer close to the symmetry plane the fourth region is formed, in which the flow is described by the parabolized Navier-Stocks equations similar to the above case of the outer singularity.

4. Studies of the symmetric flow instability over thin bodies and the control possibility on the base of the interaction model of 3D boundary layer with the electrical discharge

The electric discharge is considered as one of effective methods for control of the flow asymmetry over bodies [23, 24, 25, 26, 27]. However, to select optimal control parameters, it needs to have a reasonable criterion for the asymmetry origin and a possibility for fast estimation of the control effect. For the second problem, the model of the boundary layer and discharge interaction is proposed. The scheme of this model is shown in Figure 8 [28, 29, 30, 31, 37].

Figure 8.
A scheme of discharge interaction with the boundary layer.

It is assumed the plasma discharge effect can be modeled by the heat source in the boundary layer. The effect of gas ionization is neglected since the ionization coefficient is of the order of 10⁻⁵. This source in the energy equation is presented by formulas:

Q = Q ∗ xl h ∞ u ∞ = Q 0 y 2 exp − y − y c φ 2 σ , y c = 2 y 0 φ − φ 1 φ 2 − φ E38

Here Q ∗ is a dimensional source intensity; Q 0 is a maximum of dimensionless heat-release intensity; σ characterizes the discharge width; y c φ is a centerline of the discharge that is approximated by the parabola; y 0 is a maximum distance from the discharge centerline to the wall; and the angles φ 1 and φ 2 determine the electrode locations.

Calculations of the turbulent boundary layer characteristics were conducted using the method [10] for a slender cone of half-apex angle δ c = 5 ∘ at the angle of attack α = α ∗ / δ c = 3.15 . Other parameters are: l = 1 m, T ∞ = 288 K , u ∞ = 10 m / s , σ = 1 , and y 0 = 1 ; the center between electrodes is located at ϕ 0 = 0.5 ϕ 1 + ϕ 2 = 1.714 rad 98.25 ° , ϕ 1 = ϕ 0 − 3 Δ ϕ , and ϕ 2 = ϕ 0 + 3 Δ ϕ , where Δ ϕ = 0.0314159 is the integration step of the finite-difference approximation.

In Figure 9, the dimensionless enthalpy (Figure 9a) and circumferential velocity (Figure 9b) profiles across the boundary layer are shown as functions of η for Q 0 = 200 and for different polar angles φ . These profiles are similar to the source heat intensity distribution across the boundary layer. The temperature reaches the maximum value near the rear electrode, φ = φ 2 = 1.809 . Behind the heat source region, the temperature maximum decreases and moves toward the upper boundary-layer edge due to the heat diffusion. The station φ = 1.87 is located just after the separation point

Figure 9.
Profiles of temperature (а) и and circumferential velocity (b) across the boundary layer near the heat-release region.

Figure 10a demonstrates the plasma discharge effect on the separation point. As the heat source intensity increases from 0 to 400, the separation angle, φ s , decreases from 133° to about 105°. It is seen that the plasma heating is more effective in the range Q 0 < 100 , where the slope d φ s / dQ 0 is relatively large.

Figure 10.
Discharge effect on the separation angle (a) and flow state (b).

Figure 10b illustrates feasibility of the vortex structure control using a local boundary-layer heating on the base of the developed criterion of symmetric flow stability (solid line). Due to the heat release, the flow configuration changes from the initial asymmetric state ( φ s ≈ 133 ° , symbol 1) to the symmetric state with θ s ≈ 120 ° (symbol 2). This requires a nondimensional heat source intensity Q 0 ≈ 30 that corresponds to the total power which is approximately equal to 480 W. This example indicates that the method is feasible for practical applications of the global flow structure control.

The method of the global flow stability was developed [27, 28, 29, 30, 31] using the asymptotic approach for the flow over slender cones, the separated inviscid flow model [34] and the stability theory of autonomous dynamical systems [35]. Comparison of the calculated criteria for different elliptic slender cones with experimental data for laminar and turbulent boundary layers sowed its efficiency.

5. Investigations of abnormal features of the heat transfer and the laminar-turbulent transition for hypersonic flows around flat delta wing with blunted leading edges

Although found in the experimental zones of abnormal high heat fluxes on the windward flat surface of the half cone with blunted nose and delta wings with blunted leading edges, the phenomenon of the early laminar-turbulent transition [38, 39, 40, 41, 42, 43, 44, 45, 46] cannot be explained in frameworks of the boundary layer theory and on the base of solutions of parabolized Navier-Stocks equations. Only detailed flow simulations using full Navier-Stocks equations allowed to find reasons of such anomalies [46, 47, 48].

Figure 8 shows the comparison of calculated (the upper part) and experimental (the lower part) heat flux distributions on the delta wing with the leading edge sweep angle χ = 75°, the bluntness radius of cylindrical edges and the spherical nose R = 8 at the angles of attack α = 0°, M = 6, unit Reynolds numbers Re₁ = 1.1556 × 10⁶ m⁻¹ [47, 48]. Similar patterns were obtained in numerical simulations for different Reynolds numbers and Mach numbers up to 10.5 [46]. At moderate Mach numbers, a flow on such simple surface outside the nose and leading edge regions is described very well by the flat plate approximation and has no anomalies.

At hypersonic speeds, high heat flux regions, which is present in Figure 11, are observed in the middle wing span and near the symmetry plane. It is seen that the experimental middle high heat flux streak is finished by the turbulent wedge. Calculations were conducted only for the laminar flow.

Figure 11.
Comparison of numerical (the upper part) and experimental (the lower part) specific heat flux distribution on the wing surface.

To understand the reason for the heat flux anomaly, the cross-flow pattern helps (Figure 9). Three longitudinal vortexes are in this flow. The largest vortex is in the inviscid region above shock (the dark layer) and boundary (the light layer) layers. Vortex near the symmetry plane and in the middle of the span occupies both layers. Its mutual location depends on the blunt radius, Mach, and Reynolds numbers [43, 46]. For the considered case, the middle vortex is above the high heat flux region that is shown below the cross-flow pattern (Figure 12).

Figure 12.
The cross-flow structure above the wing in the section Х = 0.1 m.

The analysis shows that high heat flux streaks are formed by the convective transfer of heat gas from the shock layer to the wing surface by the gas rotation inside the vortex. In the considered case, the middle vortex is formed before the symmetry plane vortex near the nose in the narrowing flow region between the head shock and the leading edge due to the cross-flow acceleration near the leading edge and the induced pressure gradient related with the domed flow structure near the symmetry plane.

In considered conditions, the middle vortex also is the reason for the laminar-turbulent transition. Formed along the vortex center, streamwise velocity profiles have inflection points that lead to the Rayleigh instability development. Transverse velocity profiles along this line have the S-shaped form that leads to the cross-flow instability. Both these processes result to the more early transition than Tollmien-Schlichting wave evolution.

6. Conclusions

In this work, the short review of researches on the study of BL equation singularities, which are formed when two streamline families are collided, is presented. This phenomenon can arise only in unsteady and 3D problems and has no analogue in 2D flows. A typical example of such problem is the flow around a slender cone in the vicinity of the runoff plane. In this case, solutions are found in the analytical form that allows to analyze explicitly the singularity character.

The analysis of solutions for the outer flow part revealed two singularity types. One type is in streamwise and cross-velocity viscous perturbations; it arises at values of relative cross pressure gradient k ≥ 1 and leads to the exponential disturbance growth as the runoff plane is approached. At k = 1 the singularity is logarithmic and at k > 1 it is power; its appearance is correlated with the BL separation appearance. Another singularity type at smaller values of k ≥ 1 / 3 in the first-order approximation leads to the infinite growth of transverse velocity perturbations only and is not related directly with the flow separation; at k = 1 / 3 the singularity is logarithmic, and at k > 1 / 3 it is power. These BL singularities correspond to some asymptotic flow structure at Re ≫ 1 . This structure includes the boundary region with the dimension of the order of the BL thickness, in which the viscous transverse diffusion effect smoothes the singularity. The comparison of obtained parabolized Navier-Stokes equation solutions describing the flow in the boundary region with BL equations solutions confirms this conclusion. Second region induced by the viscous-inviscid interaction effect has the transverse dimension of the order of square root from the BL thickness and the two-layer structure. For the potential flow in the outer inviscid subregion, the integral solution representation is found on the base of the slender wing theory. The inner subregion is described by full 3D BL equations, the solution of which is obtained for the outer viscous subregion part. It was shown that the viscous-inviscid interaction does not eliminate the singularity but drifts it in the parametric space. To eliminate the irregularity, the boundary region is needed.

To find the dependence of the critical parameter of the singularity appearance k c on Mach and Prandtl numbers and the wall temperature BL equations, solutions are studied in the near-wall region beside the runoff plane. Equation subcharacteristic (streamlines) analysis showed the presence of one parameter α , the sing of which defines the qualitative change of the streamline topology and, consequently, the physical flow structure. It is shown and confirmed by comparison with all available calculations that the boundary of the solution which exists in the runoff plane corresponds to the criterion α k c = 0 . The analysis of BL equation solutions near the runoff plane revealed the presence at α ≥ 0 of irregular and at α < 0 singular proper solutions. This is confirmed by numerical calculations of the flow around slender delta wing with the small aspect ratio. Singularities in the near-wall region generate the some flow structure in its vicinity, the study of which is out of this paper framework. Presented results do not depend on outer boundary conditions and are true for the full freestream velocity diapason including hypersonic flows.

Presented research allows concluding that the flow in symmetry planes, for example, on wings, has the complex structure, which is needed to take into account the numerical modeling in order to eliminate the accuracy loss. Regular flow function decompositions commonly used at solutions of BL equations are not applied near this plane, and it cannot be considered as a boundary condition plane due to a possible solution disappearance.

Rajeswari R.