Linear stability of high-speed mixing layers

Applied Numerical North-Holland

Mathematics

7 (1991) 93-127

93

Linear stability of high-speed mixing layers * Michele NASA

G. Macaraeg

and Craig

Langley Research Center, Hampton,

L. Streett VA 23665-5225,

USA

Received August 1988 Revised December 1988

Abstract Macaraeg, M.G. and C.L. Streett, 7 (1991) 93-127.

Linear stability

of high-speed

mixing layers, Applied

Numerical

Mathematics

A newly developed spectral compressible linear stability code (SPECLS) is presented for analysis of shear flow stability and applied to high-speed boundary layers and free shear flows. The formulation utilizes the first application of a staggered mesh for a compressible flow analysis by a spectral technique and a multi-domain spectral discretization (MDSPD) option to resolve highly irregular structures. An order of magnitude less number of points is needed for equivalent accuracy in disturbance growth rates compared to calculations by a finite difference formulation, and a factor of three fewer points is required by the MDSPD relative to a single-domain spectral discretization (SPSPD). The mean flow in this study consists of two parallel gases, one of which is quiescent. A sudden appearance of multiple supersonic modes is found to occur near Mach 3. Their existence has not been observed in past studies. This occurrence is attributed to the boundary conditions which impose zero perturbations (reflecting boundary conditions) of all disturbances in the far field. The condition is imposed sufficiently far that subsonic disturbances are unaffected (i.e., results match the case of an unbounded free shear layer). However, at approximately Mach 2 there are radiating solutions (supersonic disturbances) both above and below the shear region. The imposed “wall” far-field condition causes a multiplicity of higher supersonic modes at Mach numbers exceeding two. The issue of relevance in studying the stability of a free shear flow is the impact of transition on fuel-air mixing efficiency in scramjet combustors. The bounded free shear layer in this study is the proper analogue of this situation. In particular the observation of a pronounced drop in mixing efficiency as Mach number is increased may be related to these multiple modes which have very low growth rates relative to the subsonic disturbances which predominate at low Mach numbers. At M, = 3 it is found that increasing the temperature of the quiescent gas relative to the injected gas in the mixing layer inhibits the development of these supersonic modes for low streamwise wave numbers. Verification of the method will be given for the case of boundary layers using an existing finite difference compressible stability code. and by comparison with analytical results obtained for free shear flows.

1. Introduction There is an enduring interest in the linear stability of shear flows. It can be attributed to the fact that currently this theory (used in conjunction with some semi-empirical procedure such as an N-factor method) is virtually the only means of predicting the location of the transition or indeed of determining whether a given laminar flow will become turbulent or not. At the present * Presented

at SAE Aerotech

0168-9274/91/$03.50

‘88, Anaheim,

0 1991 - Elsevier

CA, October

Science Publishers

3-7, 1988. B.V. (North-Holland)

94

M. G. Macaraeg, C.L. Streett / Linear stability

time, there is no prospect of a unified theory of transition even in low speed flows (where some of the underlying mechanisms are relatively well-known), let alone in hypersonic flows. As a short-term goal, it is imperative then to obtain linear stability results. The implicit assumption is that supersonic-hypersonic transition has its origin in linear instability (not that Morkovin’s bypasses [ll], i.e., cases where a flow transitions without any initial linear instabilities, are inapplicable to the hypersonic regime). The disturbances in this regime are small and exponentially growing, and the disturbance with the highest growth rate is assumed to cause a laminar flow to transition. Analysis of these disturbances can be handled by solving the linearized Navier-Stokes equations about a specified mean flow. Discretizing this equation set results in a complex eighth-order algebraic eigenvalue problem. The solution gives the structure of the disturbance (complex eigenfunction) and its associated phase speed and growth rate (complex eigenvalue). A spectral collocation technique on a pressure staggered mesh has been incorporated into a spectral compressible linear stability code (SPECLS) to solve the above system. The newly developed pressure staggered mesh is shown to perform very well for spectral solutions of compressible flow problems; and the multi-domain spectral discretization in SPECLS [8] allows for efficient resolution of flows with very irregular structures or widely disparate scales [6]. The focus of this study is the stability of a viscous, compressible free shear flow or mixing layer. It is well known that the initial stages of shear flow instabilities are driven by linear mechanisms. The growth and propagation of these disturbances in these early stages must be understood before transition can be accurately modelled or predicted. An issue of relevance in studying the stability of free shear flows is the impact of transition on fuel-air mixing efficiency in scramjet combustors. The eventual goal is to be able to manipulate the downstream evolution of free shear flows to increase mixing efficiency. Since it is well known that turbulent mixing is several orders of magnitude stronger than laminar, the ability to induce transition is a very desirable goal. It is of concern that mixing efficiency is decreased fourfold in the Mach number range of one to four [12], since this drop encompasses the Mach number range of interest. The physical mechanism responsible for this trend is unknown. The linear stability analysis to be presented has provided a clue regarding the underlying mechanism possibly responsible for this decrease. The supersonic disturbances which predominate at high Mach numbers are shown to have very small growth rates relative to subsonic disturbances at the low streamwise wave numbers (CX) studied. These disturbances have a nontrivial oscillatory structure which must be resolved for accurate studies. In the past, supersonic disturbances have not been studied in any depth [4] and have never been adequately resolved numerically. The MDSPD incorporated in SPECLS makes an analysis of these disturbances possible. The high accuracy in the SPECLS formulation has allowed for the first detailed look at the oscillatory structures of these supersonic modes which dominate the initial stages of instability for high speed mixing layers [ll]. The mean flow in this study consists of two parallel gases, one of which is quiescent. A sudden appearance of multiple supersonic modes is found to occur near Mach 3. Their existence has not been observed in past studies. This occurrence is attributed to the boundary conditions which impose zero perturbations (reflecting boundary conditions) of all disturbances in the far field. ’ The condition is imposed sufficiently far (approximately 1500 momentum thicknesses in the far ’ Dr. Leslie Mack, private

communications.

M. G. Macaraeg, C. L. Streett / Linear stahilit,v

95

field, on either side of the shear) that subsonic disturbances are unaffected (i.e., results match the case of an unbounded free shear layer). However, at approximately Mach 2 there are radiating solutions (supersonic disturbances) both above and below the shear region. The imposed “wall” far-field condition causes a multiplicity of higher supersonic modes at Mach numbers exceeding two. As mentioned earlier, the issue of relevance in studying the stability of a free shear flow is the impact of transition of fuel-air mixing efficiency in scramjet combustors. The bounded free shear layer in this study is the proper analogue of this situation. The observation of a pronounced drop in mixing efficiency as the Mach number is increased may be related to these multiple modes which have very low growth rates at low (Yrelative to the subsonic disturbances which predominate at low Mach numbers. It is found that by increasing the temperature of the quiescent gas relative to the injected gas at M, = 3 and these low values of a, the development of multiple supersonic modes is inhibited. A very different trend is found for high values of IX, which is addressed in [7]. Verification of the method will be given for the case of boundary layers using an existing finite difference compressible stability code [lo] and by comparison with analytical results obtained for free shear flows.

2. Basic equations for parallel shear flows 2.1. Mean j7ow The mean flows utilized in these studies are obtained from a spectral similarity solution for compressible shear flows. The basic spectral collocation technique is given in [15]. For the free shear flow analysis two gases are considered, with one of these gases taken as quiescent. The effect of relative temperature is studied by varying the parameter &= TJT,, where T2 refers to the free-stream temperature of the quiescent gas, and T, the free-stream temperature of the injected gas. A representative mean flow streamwise velocity profile is given in Fig. 1. 2.2. Compressible

stability equations

The basic equations governing the flow of a viscous compressible fluid are the Navier-Stokes equations. These equations in terms of a three-dimensional Cartesian coordinate system for a

1.0 .8 “0

1; .2 0

Fig. 1. Streamwise

velocity profile of mean flow for free shear layer. Mm = 3, pT = 1. T, = 520 o R.

96

M. G. Macaraeg,

heat-conducting

perfect

aii* J++u

gas in dimensional

_* au;

C. L. Streett

/ Linear stability

form are:

1 a?:

(1)

*=“gqY J ax,

ap* a __ -at* + -(p*u,“) ax:

= 0,

(2) t _-

P*C,* i

J* -+,,.jh,

a(u,*p*) aji* at* + ax;

7

(3) (4)

where (5) (6) and G: =(u*, u*, w*), x,* =(x*, y*, z *), and i, j = 1, 2, 3, according to the summation convention. Asterisks denote dimensional quantities and overbars denote time-dependent quantities. In these equations K* is the coefficient of thermal conductivity; R *, the gas constant; c:, the specific heat at constant pressure (assumed constant); and A*, the second coefficient of viscosity. In this study all velocities are scaled by U,, the main stream velocity, and all lengths are scaled by 6 *, a modified displacement thickness of the velocity profile in the x direction. The Reynolds number and Mach number are given by Re = U,S */v, .

(7)

M = U,//yR*T, >

(8)

where v, and T, are the kinematic viscosity and mean temperature in the free stream and y is the ratio of specific heats. The Prandtl number, u, is assumed to be 0.70 and ,ii* is evaluated by Sutherland’s law. If we assume that the base flow is a locally parallel boundary layer or free shear layer, then u(x,

y, 2, t) = U,(y)

+ ii(y)

u(x,

y, 2, t) = u”(y) ei(ax+/3r--wf),

00)

w(x,

y, 2, t) = W,(y)

(11)

p(x,

y, 2, t) = P,(y)

+p(y)

ei(ax+pr-wt),

(12)

T(x,

y, z, t) = T,(y)

+ i(y)

ei(ax+pr-wr),

(13)

+ G(y)

ei(ax+pz-wr),

e’(ax+pr-wt),

(9)

where U,, W,, and To represent the steady unperturbed parallel shear layer (mean flow), and quantities with tildas denote complex disturbance amplitudes. V,(y) is assumed zero since the flow is parallel, and P,(y) is assumed constant across the shear layer. (Y and fi are x and z disturbance wave numbers, respectively, and w is the complex frequency. Equations (9)-(13) are substituted into equations (l)-(6), the mean flow terms subtracted out, and the terms which are

M. G. Macaraeg

C. L. Streett

/ Linear stahilit?

91

quadratic in the disturbance amplitudes are neglected. The resulting system is the linearized compressible Navier-Stokes equations for the disturbance quantities given as follows [lo]: D2(c& + /GIG) +

afi+PG)

+i

T

- 1

I(a2 + ,&‘)Da

(15) D~+i(c&+/3,i)-~iii,M2(~Uo+~Wo-o)~-~(aLI,+PWo-w)~=O. 0

0 (16)

D2?+ 2( y - 1)M2u (&+Pw,‘) ff2 + p2 +2(y-1)M2a

(WA%,‘)

(“2+p2J

D(aii+/36)+---

2T,’ dP0 /_L~ dTo Di

Db+W)

(y-

l)M2((YU0 + pw, - +

(17)

(18)

98

M.G. Macaraeg, C. L. Streett

where primes and D indicate

differentiation

/ Linear stability

with respect

to Y. The equation

of state is (19)

where p. is the mean conditions are: .6 CXctp6 a+-pii

flow density

=0

at

and

y=o

b the complex

boundary

iy + - m

bounded

density

disturbance.

The

boundary

layer, free shear flow.

(20)

t”

3. Solution technique Equations (14)-(18) and (20) constitute an eigenvalue problem for the complex frequency w, once the disturbance wave numbers (Y and p are specified. Discretization of these equations in the Y direction forms a generalized matrix eigenvalue problem, suitable for computer solution. Equations (14)-(18) are essentially an eighth-order system; thus the eight boundary conditions (see (20)) are sufficient for solution, and no boundary condition is applied to the disturbance pressure. Whatever discretization scheme is used must respect this arrangement, if a truly accurate solution is to be expected. Since stability analyses are inherently much more sensitive to inaccuracies relative to a large-scale aerodynamic calculation, an artificial boundary condition could potentially impact stability results. The discretization scheme used here is a spectral collocation technique, using Chebyshev polynomials as basis functions. The nodes of the variables 6, q+= aii + p@, q-= at? - /Iii, and i are located at the Gauss-Lobatto points (the extrema of the last Chebyshev polynomial retained in the expansion [12]); the energy and momentum equations are collocated at these points. Thus discrete boundary conditions may be imposed for these variables at both endpoints of the domain. The pressure nodes are located at the Gauss points (the zeroes of the first neglected polynomial) of a Chebyshev series one order less than that used for the other variables; the continuity equation is collocated at these points. Since no Gauss points fall on the boundary, we are free of any requirement of providing an artificial numerical boundary condition for pressure, and we have the proper balance of number of equations and unknowns. The far-field boundary of the discretized domain is placed at a finite distance, typically 20-100 6*-units from the wall or shear layer centerline. Extensive sensitivity studies were performed to determine the effect of this finite domain truncation. Stretching is employed in the discretization to improve resolution near the wall or centerline. In the case of the boundary layer, either of two stretching forms are used: y, =Y,,,(c*

- l)c'Yc/~c2

-YV

(21)

or Yp = Yln&YA1

+ C, - Y, ) 3

(22)

M.G. Macaraeg,

C. L. Streett / Linear stability

99

where yr, is the coordinate in the physical space [0, y,,,], _yc is the computational coordinate, y, E [0, 11, and C,, C,, and C, are adjustable constants. In (21) C, controls the amount of stretching and C, the rate. The smaller the quantity (C, - 1) the stronger the stretching. In (22) C, controls both the rate and amount of stretching. The smaller C, the stronger the stretching. In this work, C, is either 4 or 6, or C, ranged from 1.2 to 2.0; C, between 0.01 and 0.03 is used. For the shear layer, equation (21) is used as the stretching, yielding a physical space of [-y,,,, y,,,] from ~~y,E[-l, 11. Using standard spectral collocation discretization formulas [14,15], matrix differentiation operators are formed for both the Gauss-Lobatto (6, q+, q-, r”) and the Gauss (p) grids, incorporating the selected stretching function. Mean flow quantities from the spectral boundary layer code of [15] are spectrally interpolated onto the new mesh, and derivatives of these quantities obtained using the differentiation operators. The generalized matrix eigenvalue problem which results from this discretization is of the form:

(23) for the momentum

and energy equations,

and

for the continuity equation, where A, B, C, D, and E are matrix coefficients derived from equations (14)-(18) L denotes a spectral differentiation operator, the unknown vector $I = pressure p_ Subscripts GL and G (v”, q+, q-, t”)T, and th e vector P contains the disturbance denote location at or operation on Gauss-Lobatto and Gauss point grids, respectively; ZzL and ZgL are spectral interpolation matrices, from Gauss to Gauss-Lobatto points and vice versa. The unknown vectors $ and P are collected into a single vector, and the matrix equations (23) and (24) are assembled into a large generalized matrix eigenvalue problem for input into a standard library routine for solution. A complex modified QZ algorithm [16] is used to obtain the eigenvalues of the system directly; this is referred to as a global search. The most unstable eigenvalue is then selected and used as an input to an inverse Rayleigh method to purify the eigenvalue of the effects of round-off error and to obtain the solution eigenvectors. In all cases, the global and local (Rayleigh-iterated) eigenvalues agreed to better than eight decimal places. The codes used in this study were not optimized, relying on LINPACK and EISPACK [16] routines. A global search with 91 points requires about 120 seconds on a Cray 2, whereas a local search with 199 points takes 17 seconds. 3. I. Spectral multi-domain

technique

The above discretization scheme is utilized in the spectral multi-domain technique developed by Macaraeg and Streett [5]. This technique was formulated to handle both advectionand diffusion-dominated flow problems. Extremely large differences in discretization across an interface, through domain size, number of points and stretchings, have been shown to not disrupt exponential-order accuracy [S]. These advantages are crucial for solving problems with widely disparate scales, which is the case for flows undergoing transition or involving chemical reaction [41-

100

M.G. Macaraeg,

C.L. Streett / Linear stability

A simple one-dimensional two-region example will serve to illustrate the present method for interfacing two collocation-discretized regions. Given the second-order, potentially nonlinear boundary value problem: [WJ)l,-~u,,=W), U(-1)

XE 1-L

11, (25)

U(l) = 6,

=a,

we wish to place

an interface at the point x = m, and have independent collocation discretizaE [m, 11. Even though the point x = m is an interior tions in the regions x(l) E [ - 1, m] and xC2) point to the problem domain, simply applying a collocation statement there, utilizing a combination of the discretizations on either side, will not work; the resulting algebraic system is singular. This is because the spectral second-derivative operator has two zero eigenvalues; thus the patching together of two spectrally discretized domains yields potentially four zero eigenvalues in the overall algebraic system. Two of these eigenvalues are accounted for by imposition of boundary conditions, and one by continuity of the solution at the interface, leaving one zero eigenvalue in the system. To alleviate this difficulty, a global statement of flux balance is used. Rewriting (25) as: (26)

[G(u)] x = S(u)> where the flux is G(U) then integrating

= F(U) (26) from

- vUX,

(27)

- 1 to 1 results in

G(U)I.=,-G(I/)l.=~,+[G]l.=,=IIIS(U)dx. If the jump in flux at the interface,

G(U) I

x=

-I+

(28)

[G], is zero, then (28) may be written:

/_“:S(U) dx = G(U)

1x=,

-

/‘,S( U) m

dx.

The statement of global flux balance across the two regions, along with the assumption that the solution is continuous, provides the condition necessary to close the equation set which results from spectral discretization of (25) in two regions. Note that the left-hand side of (29) involves side involves the region the discretization in the region x(i) E [ - 1, m], while the right-hand x(*) E [m, 11. Since spectral collocation discretization strongly couples all points in their respective regions, equation (29) couples all points in both discretizations. Note also that no statement is made concerning whether or not (25) is advectionor diffusion-dominated. Equation (25) is considered. a scalar equation here, although the above is extendable to a system. A thorough treatment of the spectral multi-domain technique for compressible flow stability can be found in [8].

4. Verification 4.1. Boundary

layer

For verification, calculations are performed for the stability analysis of compressible dimensional similarity boundary-layer profiles. A spectral mean flow is utilized. The

twobasic

M. G. Macaraeg, Table

101

C. L. Streett / Lmear stability

1

Calculation

of

temporal

eigenvalues

for

different

grids

using

SPECLS

(M,

= 0.00001,

T, = 520”R.

Re = 2200,

(Y= 0.2, p = 0) w

Number of points 35

(5.97289555Ox1O-2,

4.163678340x10-‘)

39

(5.982453898

x 10-2,

4.051750872

x 1O-3)

41

(5.984158831

x lo-*,

4.016799466

x lo-‘)

45

(5.983655920xW2,

65

(5.983575125

x 10m2, 4.023290767

x 10m3)

95

(5.983575084

x 10

x 10 - s)

200

(5.983575084x

4.022189152x10-‘) 2, 4.023291203

10p2,

4.023291203

x 10m3)

technique is given in [15]. The cases to be presented are for M, = 0.00001, 4.5, and 10, assuming adiabatic flow. Comparison with second-mode calculations by Mack [9] will also be presented. Initially, a low Mach number (A4, = 0.00001, T, = 520 o R, Re = 10,000) essentially incompressible case is analyzed. A resolution study using SPECLS is shown in Table 1 for the eigenvalues versus number of grid points. A similar study is given for COSAL, a finite difference compressible stability code [lo], in Table 2. Accuracy of the growth rate to three significant digits is obtained with 45 points in the spectral code; COSAL requires approximately 500 points for equivalent accuracy. It is of interest to note that converged values of growth rate differ in the third decimal place between the two codes. An independent calculation by Mack agreed to four decimal places with SPECLS. Corresponding eigenfunctions for the amplitude of B are given in Figs. 2(a)-(b) for both codes. The normal coordinate, y, is in units of 6*. The profiles are virtually identical. A higher Mach number boundary-layer profile perturbed by a three-dimensional disturbance is given next. Conditions are J4, = 4.5, T, = 520 o R, Re = 10,000, (Y= 0.6, p = 1.0392. A resolution study for the eigenvalue computations is given in Tables 3 and 4 for SPECLS and COSAL, respectively. Accuracy to three significant digits for values of growth rate are obtained with the spectral code using 81 points; in Table 4 COSAL is seen to require approximately 800 points for equivalent accuracy. Eigenfunctions for D are again shown in Figs. 3(a)-(b) from each

Table

2

Calculation

of temporal

eigenvalues

for different

grids using COSAL

(M,

= 0.00001,

P = 0) w

Number of points 100

(5.983647668

X 10-3)

(5.981935522

x 10-2, x 10-2,

3.96449684

200

4.004225579

x 10p3)

300

(5.981623402

x 10p2,

4.01151043

X 10p3)

500

(5.981464629x

10-2,

4.015224481

X10-‘)

x 10-2,

4.016810997

X 10m3)

4.017305399

X 10p3)

1025

(5.981396790

1200

(5.981375400x

(with Richardson

extrapolation)

lo-*,

c = 520 o R, Re = 2200, (Y= 0.2,

102

M. G. Macaraeg, C.L. Streett / Linear stability Mm = .OOOl , adiabatic

a

M

= 4.5 , adiabatic

-.4 U-

Fig. 2. ii velocity fluctuation amplitude function for (Y= 0.2, /3 = 0, M, = 0.00001, Re = 2200, 7” = 520 o R, adiabatic. Calculation by (a) SPECLS, (b) COSAL.

Fig. 3. ii velocity fluctuation amplitude function for a = 0.6, fi = 1.0392, h4, = 4.5. Re = 10,000. T, = 520 o R, adiabatic. Calculation by (a) SPECLS, (b) COSAL.

code. A multi-domain spectral discretization (MDSPD [5]) in SPECLS with two domains obtains equivalent accuracy with the single-domain spectral discretization (SDSPD) with one-third of the number of points as illustrated in Table 5. The savings is significant considering that the number of operations in SPECLS goes as the cube of the number of points. Insensitivity of these solutions to the locations of the far-field boundary is next illustrated with the MDSPD. The compressible case discussed above utilized an outer extent of 30 6 *-units. This

Table 3 Calculation of temporal ,B = 1.0392)

eigenvalues

for different

grids using SPECLS

Number of points

w

45 51 65 81 95 120 151 200

(0.49578803, (0.49591572, (0.496217679, (0.496093364, (0.496106105, (0.496104826, (0.496104845, (0.496104846,

(Mm = 4.5, T, = 520 o R, Re = 10,000, (Y= 0.6,

4.36603960 3.749838454x 3.727396074 3.761833505 3.763494422 3.764796852 3.764820857 3.764820881

x 10 - 3, 10-3) x lo-‘) x 10-3) x 10-3) x 10-3) x 10-3) x lo-‘)

M. G. Macaraeg Table 4 Calculation of temporal j3 = 1.0392)

eigenvalues

for different

C. L. Streett / Linear stahlit)

grids using COSAL

Number of points

w

211 513 813 1025 1500 (with Richardson

(0.496035465. (0.495990048, (0.495984488, (0.495983122, (0.495981891,

25/17 31/25 37/25

= 4.5, T, = 520 o R, Re = 10,000, (Y= 0.6,

3.720663427 x lo-“) 3.688194736 x lo- ‘) 3.684224749 x lo- 3, 3.683247436x10-s) 3.682367822 x 10 -s)

extrapolation)

Table 5 Calculation of temporal (Y= 0.6, ,B = 1.0392) Two domains

(M,

103

eigenvalues

for different

(Total number

grids using MDSPD

of points)

on temporal

(M = 4.5, 7;, = 520 o R, Re = 10,000,

w (0.496102327, (0.496104982, (0.496104984,

(42) (56) (62)

Table 6 Effect of far-field boundary location Re = 10,000, (Y= 0.6, p = 1.0392)

in SPECLS

eigenvalues

using MDSPD

Two domains

Far-field boundary (S *-units)

w

37/21 37/31 37/51

15 30 60

(0.496104813, (0.496104828, (0.496104848,

3.764737236 x 10ml) 3.764889142x10-s) 3.764890672 x lo-‘)

in SPECLS

(Mm = 4.5, T, = 520 o R,

3.764828351~10~~) 3.764822175 x 10 - ‘) 3.764817702 x 10p3)

same case was run halving and doubling this far-field boundary distance. Results are given in Table 6. A factor of four change in outer extent yields basically no change in o: The phase speed is maintained to seven significant digits, and the growth rate is unchanged for five significant digits. First- and second-mode eigenfunctions for &I, = 10, T,= 90 o R, Re = 100,000 are compared with compressible inviscid modes of neutral subsonic solutions [9]. Plotted in Figs. 4(a)-(b) are the pressure amplitude functions of the first mode from SPECLS and [9], respectively. The value of cy is 1.844. Figures 5(a)-(b) display similar plots for the second mode, with (Y= 4.9877. Nondimensionalizations and properties differ for these comparison cases; of relevance are the distinctive shapes of these jj modes, i.e., that the number of zeroes in jj is one less than the mode number [9]. The spectral code utilized 99 points for this case. A further comparison by Mack and SPECLS for equivalent flow properties and viscous flow at IV, = 4.5, and T,= 300 o R, gave agreement to three significant figures in growth rate, and five significant figures in phase speed.

hf. G. Macaraeg,

104 M

C. L. Streett / Linear stability

= 10 , adiabatic M

a

= 10 , adiabatic

,lL

,.!,1.

1,

1.2

.6 Y

Y

-1 .o

b 1.0

Fig 4. Pressure fluctuation amplitude function for the first mode of neutral subsonic solution at (Y= 1.844, p = 0, A4, = 10, T, = 90”R, adiabatic. Calculation (a) by SPECLS (Re = lOO,OOO),(b) as in [9] (inviscid).

4.2. Shear

b -2.0 /;; Fig. 5. Pressure fluctuation amplitude function of the second mode of neutral subsonic solution at a = 4.9877, /3 = 0, A4, = 10, 7;, = 90” R, adiabatic. Calculation (a) by SPECLS (Re = lOO,OOO),(b) as in [9] (inviscid).

layer

Further vertification of SPECLS is presented for parallel shear flows. The normal spatial scale (y) used for plotting disturbance eigenfunctions is in units of the modified “displacement thickness” 6 *. 6 * ranges from 14.498 for incompressible flow to 15.505 for a Mach 3 shear layer. The corresponding nondimensionalized momentum thickness is 0.8757 for incompressible and 0.8188 for the Mach 3 shear layer. In other words, multiplying the normal coordinate by roughly a factor of 18 will give the far-field extent in terms of momentum thickness. Conversion of nondimensionalized spatial wave numbers similarly follows. In the case of the tanh mean flow, S * is assumed to be one. Results are compared with those obtained by Blumen [l], Gropengiesser [4] and Ragaab and Wu [13]. The velocity profile whose stability is investigated for comparison with Blumen is given by U, = tanh(y).

(30)

The basic thermodynamic state is assumed constant [l]. Table 7 displays maximum growth rates for a range of Mach numbers as determined by Blumen [l]. The corresponding growth rates using SPECLS (Re = 10,000) are also displayed and were obtained using at most 81 collocation points. To determine the maximum growth rates using SPECLS a plot of growth rate versus (Yis displayed in Fig. 6 for A4, = 0.5 and A4, = 0.9. The peaks occur at the growth rates given in Table 7 for A4, = 0.5 and A4, = 0.9. Note that the lower Mach number flow is unstable to a larger range of wave numbers.

M.G. Macaraeg, C. L. Streett / Linear stahilit,v Table 7 Maximum

values of wI for Mach numbers

105

> 0

Mach number

(Y

0, HI (tanh( Y ))

o, (SPECLS)

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

0.445 0.433 0.426 0.417 0.409 0.397 0.370 0.326 0.279 0.208 0.000

0.190 0.187 0.181 0.171 0.158 0.141 0.122 0.101 0.078 0.055 0.000

0.18954 0.18736 0.18112 0.17105 0.15760 0.14101 0.12180 0.10012 0.07760 0.05450 0.00000

(tanh( y))

It is interesting to note the effect of Reynolds number on the maximum growth rate. Figure 7 displays a plot of (Y versus w for M, = 0.9 and Re = 100. The maximum growth rate is seen to increase with increasing Reynolds number indicating that this shear layer profile has an inviscid instability.

.08 Oi .06 ,,-’

_-.--_

‘\

‘\\

/’

‘1\ ‘1 \\ \

.02

0

I

0

.l

I

.2

I .3

I

\\ \

\\ \\

I \

.4

I

I

.5

I

I

.6

I

I .7

I

I, .8

1 .9

1.0

a Fig. 6. Growth

rate versus wave number,

(Y,for the hyperbolic

tangent

velocity profile for M, = 0.5 and

= 0.9.

M. G. Macaraeg,

106

C. L. Streett

/ Linear stability

.lO -

.08 -

-

Re=lOO

-----

Re = 10,000

.06 -

_C---__

Oi

Fig. 7. Growth

rate versus wave number,

(Y,for the hyperbolic

tangent

profile

for Re = 10,000, Re = 100, M, = 0.9.

The real and imaginary parts of the eigenfunctions ii, 6, and j! for maximum growth rate are given in Figs. 8(a)-(c), for M, = 0.9, Re = 10,000. Each eigenfunction plot is normalized to one. Note that the real part of ii displays a small structure about the origin. It is found that this feature disappears as the neutral curve is approached as illustrated in Figs. 9(a)-(f). The manifestation of this structure is reflected in plots of magnitude and phase for ii as given in Figs. 10 and 11. Note the pronounced peak in the magnitude which decreases radically as the structure of ii about the inflection point disappears (Fig. 10). The plots of Fig. 11 exhibit a change of phase when the structure is prominent ((Y = 0.1) as opposed to the smooth phase plot

Fig. 8. Complex

disturbance

amplitudes

versus

y for the hyperbolic Mm = 0.9, Re = 10,000.

tangent

profiles

at the maximum

growth

rate:

M.G. Macaraeg C. L. Streett / Linear stability

0

-25

ia

107

f

~

25 -25

0

Y

25

Y

Fig. 9. E velocity fluctuation amplitude functions for increasing values of (Y for the hyperbolic tangent profile; M, = 0.9, Re =lO,OOO; (a) (Y= 0.1, w = 0.0428, (b) (Y= 0.208, w = 0.0544. Cc) (Y= 0.35. w = 0.0313, (d) (Y= 0.40, w = 0.0317, (e) (Y= 0.42, w = 0.00605, (f) (Y= 0.43, b =0.00216

associated with CY= 0.4, when the structure has disappeared. appearance of this structure as exhibited in Fig. 12, which increasing (Y(Re = 100).

1.0

k a

-1.0 1.0

k

b

-1.0

10

L &

lid 0

II4 0 2.0

UR

II611 -1 .o k -20

0

Y Fig. 10. Plots of // 111)for increasing

L!L-

2.0

1.0

C

Viscosity appears to widen is a progression of ii plots

20

F c

0 -20

I

0

20

Y

a. M, = 0.9. Re = 10,000; (a) a = 0.1, (b) a = 0.35, (c) a = 0.4.

the for

M.G. Macaraeg

108

CL.

Streett

/ Linear stability

cl=O.l

-----

.4

a = 0.4

I

$h

I I I I I

-.4

I -50

-40

-30

-20

-10

I

I

I

I

I

0

IO

20

30

40

I

50

Y

Fig. 11. Phase plots of ii corresponding

to a = 0.1and a = 0.4 cases given in Fig. 9.

More recent studies have utilized similarity solutions for the mixing layer mean flow [4,13] and carried out a spatial stability analysis. In these studies the Prandtl number is assumed to be one, and p = T. A comparison with these studies is given in Fig. 13, in terms of spatial growth rate (CT*) versus frequency (w), for Mm = 0 and fi = 0. It is interesting to note that the Gaster transformation used to obtain the spatial eigenvalues from the temporal results compare well with the full spatial analysis performed in [4,13]. For comparison purposes the spatial nondimensionalization for Fig. 13 is the momentum thickness utilized in [4,13].

Re = 100

c

rJz+ -25

d

0

25

I+, -25

Y Fig. 12. ii eigenfunction

plots

for increasing a for the hyperbolic tangent profile; a = 0.01, (b) a = 0.1, (c) a = 0.2, (d) (Y= 0.35.

0

25

Y M, = 0.9, Re = 100, /3 = 0; (a)

M. G. Macaraeg, C. L. Streett / Linear stahilitJ ___

p,= .6. SPECLS ~~

__.

109

&= 1.0, SPECLS &= 2 0, SPECLS a &= .6. 15, 131

.I 4 r

Q &= 1.0, [5. 131 ,121

Fig. 13. Comparison

4.3. Sensitivity

/

,_A -1. A\

7 p,= 2 0. [5, 131

(SPECLS with Gaster transformation) with Gropengiesser stability analysis). Spatial growth rate versus frequency.

[4]. and Ragaab M, = 0, B = 0.

and Wu [13] (spatial

analysis

When interpreting the results of a linear stability analysis one must realize the sensitivity of these results to the accuracy of the mean flow with respect to the numerical scheme as well as properties in the flow. To elucidate these points we choose to scrutinize a case by examining the sensitivity of calculated growth rates to accuracy of the mean flow and assumed thermodynamic properties. Consider a boundary-layer flow with M, = 4.5, T, = 300 “R, h = - $E(, a = 0.6 and ,l3= 1.0392. This flow defines Case 1 in Table 8. (A lower T, than in Table 3 is chosen to emphasize the extreme sensitivity of results to thermodynamic properties even at these relatively low temperatures.) A less than 3% change in the Prandtl number for this flow causes a 20% change in the calculated growth rate as illustrated in Table :3. Even more dramatic is the effect of changing the free-stream temperature (T,) by 42%, which causes a 48% change in growth rate. In other words, a small deviation in the velocity or temperature profile has significant effect on disturbance growth rates. The accuracy of assumed mean flows is a primary concern in studies of stability and transition. The preceding observation emphasizes the importance of this concern. A lesser effect, but nonetheless important for code validation purposes, is the effect of bulk viscosity illustrated by comparing Cases 4 and 1. The effects of these properties are of course intensified at higher temperatures. Further points of interest are the effects of the numerical accuracy of the mean flow, as well as the interpolation scheme employed to transfer mean flow values onto SPECLS grid. These latter effects, although minimal, will cumulatively cause a difference in the second significant figure for the growth rate of Case 1.

110

M.G. Macaraeg

C.L. Streett / Linear stability

Table 8 Sensitivity analysis: comparing Case 1 with varying boundary-layer mean flow properties (A4, = 4.5, Re = 10,000, (Y= 0.6, p = 1.0392; SPECLS discretization: SDSPD, 81 points); bold quantities denote a property changed relative to Case 1 Case

Mean i-low

T,(OR)

Pr

h

Interpolation of mean flow onto SPECLS grid

spectral

300

0.72

- :p

spectral

0%

300

0.72

- 4~

cubic spline

0.56%

300

0.72

- ip

linear

0.57%

UOI spectral

300

0.72

$/.L

spectral

0.13%

u51 spectral

300

0.7

- +j~ spectral

20%

520

0.72

- $p

spectral

30%

520

0.7

- $r~ spectral

48%

w

(0.497254134,

7.260461085

x

10m3)

% difference with Case 1

1151 (0.497308042, 7.219787293x 10e3)

2nd-order finite-difference

1101 (0.497308964,

(0.49725909, (0.496759786,

7.218772535

7.269563811

X

X

10m3)

10m3)

5.883195138 X 10m3)

Znd-order finite-difference

1151 (0.496601687, 5.149929437X 10m3)

spectral

1151 (0.496093364, 3.761833505X 10m3)

spectral

1151

It is of relevance to note the effect the Prandtl number has on the generalized inflection points (defined by (C&‘/T’,) = 0) for shear flows, since these flows are characterized by inflectional instabilities. For demonstration a comparison is given in terms of a plot of (C&‘/T,) versus y in Fig. 14 corresponding to two different shear layer mean flows. The solid line represents a mean flow with u = 0.7, Sutherland’s law viscosity, and & = 0.6, representative of a typical profile in the present report. Both flows have a free-stream Mach number of 5. The dashed curve is for u = 1, p = T, and & = 0.6, representative of a profile from [4]. Note that the number and occurrences of inflection points in these two cases are vastly different, which has important implications for studies involving overly simplified constitutive relations. The simplified flow

Fig. 14. Plot of the function

(U,,/r,)’

I

I

I

I

2

4

6

8

I

10

versus y for two shear layer mean flows (Mm = 5, & = 0.6).

M. G. Macaraeg, C. L. Streett

/ Linear stability

111

( CI= 1, p = r) predicts three generalized inflection points, whereas the more realistic flow predicts only one. Two inflection points are found not to occur for as high a Mach number as 10 assuming the latter constitutive relations.

5. Stability of a compressible

shear layer

The spectral stability code is used to analyze the stability of a compressible shear layer. This mean flow is a spectral collocation similarity solution for free shear flows. Issues of relevance in this study are understanding the impact of transition on fuel-air mixing efficiency in scramjet combustors. It is well known that turbulent mixing is many orders of magnitude faster than laminar. Ideally, one would like to be able to manipulate the downstream evolution of shear layers to enhance mixing. The initial stages of shear flow instabilities are driven by linear mechanisms. Understanding the growth and propagation of the disturbance in these early stages will allow not only a better understanding of the onset of transition but will allow initiation of the transition process in a numerical model so that the physics might be systematically studied. The study investigates a range of Mach numbers, gas temperatures and disturbance wave numbers. The mean flow in all cases involves a jet being injected into a quiescent gas. In addition the experimental observation that mixing efficiency is decreased fourfold in the range of Mach 1 to 4 [2,12] is implied in some new results from linear theory using SPECLS. 5. I. Effect of Mach number and ,l&- on stability For the cases under study, it is observed that lower Mach numbers correspond to a more unstable mean flow relative to higher Mach number cases. This point is illustrated in Fig. 15 which displays growth rate versus (Y for a range of Mach numbers. The wave propagation angle ( 8) is 60 O. This angle is defined as 8 =

cos-l

Ca’

+;y

(31)

.

.40 -

-MM,=1 -__--Mao=2 ---Ma=3

0

.5

1.0

1.5

2.0

2.5

3.0

3.5

I 4.0

Q

Fig. 15. Growth

rate versus wave number

a for a range of Mach numbers;

0 = 60 O, Pr = 1. T, = 520 o R

M. G. Macaraeg,

112

.I5 -

C.L. Streett

/ Linear stability

~

Cold into hot gas (pT = 5)

-----

Coldgas(pT=l)

---

Hot into cold gas (pT = 0.2)

.I0 w .05 -

J 0

.5

1.0

1.5

2.0

2.5

3.0

Q

Fig. 16. Growth

rate versus wave number

a for various temperatures;

M, = 2, 0 = 60 o

This angle is chosen since it is known [6] that disturbances which are propagated at an angle of 50 o to 60 o relative to the direction of flow experience the greatest amplification. As can be seen in Fig. 15, the lower the Mach number the wider the band of wave numbers which can cause the flow to go unstable. The stability characteristics of the shear flow are also quite sensitive to the temperature difference between the injected and the quiescent gas. Figure 16 illustrates this sensitivity. The plot displays growth rate versus (Y for a wave propagation angle of 60’ and Mach 3 flow. The three temperature cases given in the plot are cold injection (500”R) into a hot quiescent gas (2500’R), & = 5; cold injection into a cold gas (500”R), &- = 1; and hot injection (2500” R) into a cold quiescent gas (500 o R), & = 0.2. The condition of relevance to the scramjet is & = 5, since in this context a cold fuel is injected into a much hotter air stream [3]. Mixing efficiency is greatly enhanced by transition, so the greater instability of the shear layer at this temperature is a favorable scenario for increased air-fuel mixing. Cold injection into a cold quiescent gas is less unstable, and hot injection into a cold gas is the more stable of these cases. It seems that if a fuel temperature increases relative to an air stream its stability will be enhanced for the parameter range studied. This observation has important implications. For example, cooling an aircraft by running fuel under the skin will hinder mixing efficiency since heating the injected fuel will decrease the range of wave numbers which can induce an instability necessary for the flow to transition. 5.2. Effect of Mach number, &. and wave number on disturbance

structure

5.2.1. 30 modes It is not surprising that distinct differences in the shapes of temperature disturbance eigenfunctions occur as the gas temperatures are varied. To illustrate this point four cases are displayed in Fig. 17 corresponding to disturbance eigenfunctions of mean flows which differ only in the temperature of the injected gas with respect to the quiescent gas. The Mach number of the differences injected gas is 2, and 6’= 60 O, with (Y= 0.0862 and /? = 0.1493. The qualitative essentially reflect the effects of viscosity which is a consequence of the different temperature distributions for each case.

113

M. G. Macaraeg, C. L. Streett / Linear stability

1.0

-1.0

-10

‘.‘.‘*‘*‘.’ 0 Y

I d. F.I.I.ISI10 -10

1,

, , , , , , , , ,

0 Y

10

Fig. 17. Effect of temperature on i disturbance eigenfunctions; Mm = 2, 0 = 60 O, a = 0.0862; (a) cold injection (500 o R) into a hot gas (2500 o R), & = 5; (b) cold injection into a cold gas (520 o R), & = 1; (c) hot injection into hot gas (2500 o R), ,& = 1; (d) hot injection (2500 o R) into cold gas (500 o R), &- = 0.2.

Y

Fig. 18. Effect of increasing

Y

(Yon i? and d eigenfunctions; 0 = 60 O, M, = 3, & = 5, TI = 520 o R; (a) (Y= 0.25; (b) a = 1.0; (c) (Y= 1.25.

114

M.G. Macaraeg,

C. L. Streett / Linear stability

i;

b

-10

0

10

-10

Y

Fig. 19. Effect

of increasing

a on ,G and i eigenfunctions; f3 = 60 O, M-=3, (Y= 0.8; (c) cx= 1.25.

0

10

Y

&=l,

T,=520”R;

(a) ol=O.4;

(b)

It is further observed that the disturbance eigenfunctions significantly tighten in structure as (Y increases which puts a greater demand on the resolution required for a SDSPD. A progression of eigenfunction plots for increasing (Y(8 = 60 o ) is given in Figs. 18 (pr = 5, Tr = 500 o R) and 19 (& = 1, Ti = 520 o R) illustrating this observation. Difficulties in resolution similarly occur for higher Mach number disturbances. The SDSPD in SPECLS had difficulty resolving eigenfunctions beyond 44, = 3.75 for &- = 5 and as early as M, = 3 for & = 0.2. Restrictions on the allowable stretching for single-domain spectral methods contributes to this difficulty. The MDSPD in SPECLS is especially useful in these cases. Examples of its usage for cases requiring very severe stretchings will be given later. To illustrate an effect of increasing Mach number on the three temperature cases of Fig. 16, plots of pressure disturbances are given in Figs. 20-21. The wave propagation angle is 60” and corresponds to comparable locations on each case’s respective w versus (Y plot, so that growth rates are approximately 10e2. Note that for & = 0.2 (Fig. 20 (left)) an additional lobe develops at M, = 3. Beyond M, = 3 the single-domain code has difficulties resolving disturbance eigenfunctions. Contrast the left-hand side of Fig. 20 with the right-hand side, PT = 5. The jj eigenfunctions maintain a single lobe and are one-signed. Cases are resolved up to M, = 3.75 (note that the spatial scale is for the interval [ - 10, lo], though the actual extent is [ - 100, loo]). Similarly, PT = 1 (Fig. 21) is resolvable to Mm = 3.75 although the real part of jj is different from the preceding case.

M. G. Macaraeg, C. L. Streett / Linear stability Hot

Cold

into Cold

115 into Hot

i

a -1.0 1.0

i;

b -1.0 1.0

"p C -1.0 -100

0 Y

100

-10 \:

Fig. 20. p disturbance

eigenfunctions for increasing Mach number. Left-hand side: & = 0.2, T, = 2500 o R, 13= 60 O; (a) M, = 1.0, (b) A4, = 2.0, (c) Mm = 3.0. Right-hand side: /& = 5, T, = 520 o R, 0 = 60 O; (a) M, = 1.0, (b) Mm = 2.5, (c) Mm

Y

Fig. 21. jj disturbance

= 3.75.

Y

eigenfunctions for increasing Mach number; & = 1, @= 60”, TI = 520° R; (a) M, =l.O, I& = 2.0, (c) M, = 2.5, (d) IV, = 3.0, (e) IV& = 3.5, (f) M, = 3.75.

(b)

116

M.G. Macaraeg,

C. L. Streett

0

100

/ Linear stability

I

-100

Y Fig. 22. Effect

5.2.2.

20

of increasing

CYon

-10

0

10

Y

j3 and i disturbance eigenfunctions; M, = 1, /3 = 0, Re = 10,000, T,=2500°R; (a) 0=1,(b) &=5,(c) a=6.

& = 0.2,

modes

The above studies involve three-dimensional disturbances with a wave propagation angle of 60”. Preliminary results indicate that two-dimensional disturbances (/? = 0) exhibit similar trends; however, difficulties in resolving eigenfunctions with a SDSPD occur at lower Mach numbers. Again, at higher (Y, disturbance eigenfunctions tighten radically beyond the limit at which a single-domain spectral method can resolve, and cases become harder to resolve as & decreases (cooling). This observation is demonstrated in Fig. 22 which plots a sequence of p and i eigenfunction for increasing (Y, p = 0. The mean flow corresponds to M, = 1, Re = 10,000 with &- = 0.2. Note the extremely tight structures in p and t” as (Yincreases. Again note that while t” is plotted on the interval [ - 10, lo], the actual spatial extent is [ - 100, 1001. Numerical oscillations in p for (Y= 5 and 6 are quite pronounced, indicating breakdown of resolution. A MDSPD is ideal for flows with radical scale differences [4]. To illustrate this point d and 6 eigenfunctions obtained from a multi-domain discretization are displayed in Fig. 23 for the Re = 10,000, (Y= 5 case discussed previously. Adjacent to the MDSPD results are the SDSPD solutions for p and u”.The multi-domain solution remains oscillation-free. The number of points in each discretization is roughly 100; however, the multi-domain utilizes three domains, with 41 points in the center domain between - 0.5 and 0.5, 25 points in the left domain between - 100, and -0.5, and 37 points in the right domain between 0.5 and 100. The plot of the pressure disturbance for the entire spatial extent [ - 100, 1001 as obtained from the multi-domain solution is given in Fig. 24, to illustrate the fineness of the structure which is resolved.

M.G. Mucaraeg,

-25

117

C.L. Streett / Linear stability

0

25

0

-25

25

Y

Y

Fig. 23. u” and jj eigenfunctions obtained with MDSPD as opposed to SDSPD in SPECLS for the case described Fig. 22(b), (Y= 5; single-domain: (a) 6, (b) p; multi-domain: (c) 6, (d) jj.

in

The preceding case ((Y = 5) is found to be resolvable by a SDSPD but only after quite a bit of trial and error stretching of the mesh for a variety of resolutions. The important point is that the MDSPD is quite robust and gives accurate eigenvalues over a wide range of stretching parameters. This observation is illustrated in Table 9 which lists a range of stretching parameters and corresponding phase speeds and growth rates for this case. Both the SDSPD and the MDSPD for this illustration utilized 99 points with the same outer extent. Note that for the SPSPD cases changing the stretching parameter by approximately 25% causes about a 10% change in phase speed and over a 50% change of value in growth rate. Contrast this sensitivity with the MDSPD. The stretching parameter is allowed to vary 100% in the center domain (three domains are utilized with interfaces at f 1). The phase speed has changed by less than 1% and the growth rate by only 6%. This robustness is extremely important since one usually has no idea

1.0

r

_--.____. p,

I

PR

1

6

II

-1 .o -100 Fig. 24. p disturbance eigenfunctions case of Fig. 23, a = 5. Discretization:

I

1

I

I

I

I

I

I,

I

0 Y

I.

I.

I

.,

100

obtained by MDSPD in SPECLS displayed on entire computational domain for center domain on [ -0.5, 51, 41 points; left domain on [ - 100, -0.51, 21 points; right domain on [0.5, 1001, 35 points.

M.G. Macaraeg

118 Table 9 Effect of stretching

parameter:

SDSPD

CL. Streett

versus MDSPD

/ Linear stability

in SPECLS;

Mm = 1.0, Re = 10,000, & = 0.2, a = 5, fi = 0

Single-domain 99 points on [ - 100, 1001

Multi-domain 17/65/17 (99 points)

Stretching parameter

0

Stretching parameter

w

1.1 1.2 1.3 1.4

(1.47179, 0.272236) (1.49883, 0.159461) (1.49806, 0.158226) (1.57814, 0.266635)

1.6 2.0 2.2 4.0

(1.497268, (1.496913, (1.497146, (1.496359,

on [ - lOO/-

l/1/100]

0.167511) 0.159498) 0.158781) 0.156986)

of the value of the phase speed or growth rate. In addition, one needs to determine whether the disturbance mode is spurious-which is measured by its persistence over a wide range of resolutions. 5.3. Supersonic

disturbances

As mentioned earlier the disturbance eigenfunctions become increasingly complex as the Mach number is increased. These higher Mach number cases are unresolvable by a SDSPD. To illustrate, a series of J and u” eigenfunctions calculated by the SDSPD stability code are displayed in Figs. 25 and 26, respectively, for a disturbance wave angle of 60 O. Note that at M, = 3.5 the injected gas side of the disturbance begins to take on an oscillatory nature; at M, = 3.75 these oscillations are more pronounced. It is well known that for supersonic disturbances the eigenfunction structure will be oscillatory [4]. (A supersonic disturbance occurs

1.0

M,=2.0

j ..:. *- :

a

b

1.0

b

M,=3.0

: d

6

-100

0

+_._

r. I. I. I. I. I. I, 1.1.1, -100 0

100

Y

Y Fig. 25. p disturbance

eigenfunctions

I

100

for increasing

Mach

number

Tl = 2500 o R; (a) M, = 2.0, (b) Mm = 3.0,

(c)

and a SDSPD IV,

=

3.5,

(d)

in SPECLS; M,

=

3.75.

0 = 60 o, & = 0.2,

M. G. Macaraeg,

b

CL.

Streett

/ Linear stability

;.~c_~, -100

0

100

d

~~~~~

-100

Y Fig. 26. 6 disturbance

119

0

100

Y

eigenfunctions for increasing Mach number and a SDSPD in SPECLS; T, = 520 “R; (a) Mm = 2.0, (b) Mm = 3.0, (c) M, = 3.5, (d) Mm = 3.75.

0 = 60 O, & = I,

when the wave velocity of the disturbance relative to the local flow, in the direction of wave propagation, has a magnitude greater than the speed of sound.) The MDSPD is able to capture the structure of these supersonic modes with relative ease. Figure 27 displays two unstable supersonic modes ( CI= 0.30416, 8 = 33.5 o ) which are associated with the instability of a Mach 4,

Fig. 27. ji disturbance eigenfunctions for two unstable supersonic modes corresponding 0 = 60 O, (Y= 0.30416. MDSPD in SPECLS with three domains, interface

to M, = 4, Re = 10,000, at + 1.

120

M.G. Macaraeg,

-10

0

CL. Streett

/ Linear stability

10

-10

Y Fig. 28. i disturbance eigenfunctions

for the two unstable

supersonic domain

[- 100, 1001and expanded

0 Y

modes of Fig. 27 plotted [ - 10,101.

10 on the full domain

& = 1, free shear flow. This pair of competing unstable modes are linearly independent solutions of the complex eigenvalue problem. The multiplicity of unstable modes occurs for supersonic shear layers, and not for subsonic cases. The MDSPD involves three domains: 105 points on the oscillatory side, 41 points in the inner domain and 25 points in the outer domain where the profile is smooth. Interface locations are +l. Note further the level of complexity of the 2 eigenfunctions of this case. Figure 28 is a plot of t” for both modes on a full scale + 100 and greatly expanded scale + 10. The center structure is an added complexity to the structure which also requires adequate resolution, and further illustrates the necessity of a flexible discretization scheme like a MDSPD. 5.4. Higher supersonic modes Associated with these supersonic modes are higher supersonic modes which have a very tight oscillatory structure and are much lower in growth rate than subsonic disturbances. These higher modes are easiest to identify from plots of growth rate versus wave number. In Fig. 29 an additional lobe at approximately (Y= 2.6 is seen for this growth rate versus wave number plot for M, = 2, &- = 1 and p = 0. A plot of phase velocity for these disturbances depict these higher supersonic modes with initially high phase speed followed by a decrease in speed to one which is commensurate with the primary mode (Fig. 30). Figure 29 also depicts a family of slightly unstable modes, between (Y= 0.8 and (Y= 1.6, which also have phase speeds different from the primary mode (Fig. 30). The MDSPD in SPECLS allows for the resolution of these tightly oscillatory modes. Their development structurally is given in Fig. 31, which depicts pressure disturbance eigenfunctions for the flow in Fig. 29 as a function of increasing wave number.

M. G. Macaraeg,

C. L. Streett / Linear stabihty

*

M,

=

2,

2D

121

rrodes

.16

Fig. 29. Appearance

of a higher supersonic

Past studies of shear flow linear analysis reveals a multiplicity of &=land/3=O.Th e associated The plot is not complete, but the

mode:

M, = 2, pT = 1. I”, = 520 o FL P = 0.

stability have found at most one higher mode [6]. The present higher supersonic modes occurring at approximately A4, = 3, growth rate curves for these higher modes are given in Fig. 32. complex trend of multiple mode development is evident. The

*

M,

=

2,

2D modes

.6

phase

v

%. .4

tt ot~r~l~I~l,l,I,r,l,I,r 0

.4

.8

1.2

1.6

2.0

2.4

Fig. 30. Phase speed plot for modes depicted

2.8

in Fig. 29.

3.2

3.6

4.0

122

M. G. Macaraeg,

C.L. Streett

/ Linear stability

a = .9

-.6

-1.0

1

III

c

I

I

-40

-60

-20

I

I

20

0

I

1

I,

I

40

60

.6

-.6 1

-1.ot

+I -60

1

-40I

/

-20I

I

0I

I

20I

40I

,

1

I

60

-60

1

I

,

-40

of a higher supersonic

,

I 0

I

I 20

,

I 40

/

1, 60

Y

Y Fig. 31. Formation

I -20

mode as a function

*

of increasing

M,

=

3,

LYfor the disturbances

2D

described

modes

.36c

a Fig. 32. Appearance

of multiple

higher supersonic

modes:

Mm = 3, & = 1, T, = 520 o R, /3 = 0.

by Fig. 29.

M. G. Macaraeg,

C. L. Streett / Linear stability

*

IIt,

1

I

3

.25

.50

M,

=

3,

2D

123

modes

I .75

1.30

1.25

1.53

Fig. 33. Phase speed plot for modes depicted

1

1.75

I

I

2.00

2.25

2.53

in Fig. 32.

accompanying phase speed plot is given in Fig. 33. These modes are not found for unbounded free shear flows. There are basically three types of supersonic modes found in this analysis: modes which are supersonic with respect to the injected gas (Fig. 34(a)), supersonic with respect to the quiescent gas (Fig. 34(b)), and supersonic with respect to both gases (Fig. 34(c)). To emphasize the complex structure associated with these modes, Figs. 35(a)-(b) depict the temperature disturbance for the doubly supersonic mode, on both a full and expanded spatial scale. These modes are found to have lower growth rates than the one-sided supersonic disturbances, which is expected since in the former case the mode is radiating energy in both directions. Recall the experimental observation which finds a strong decrease in mixing efficiency for increasing Mach number. In this regime supersonic and higher supersonic modes begin to predominate with their associated very low growth rates. Ideally one would like to render these flows as unstable as possible so that transition is encouraged which of course enhances mixing. For the flow in Fig. 32, this destablization is accomplished by increasing PT, which is known to destabilize the free shear flow at this Mach number and for this range of (Y.By increasing pr by a factor of five, multiple supersonic mode formation is suppressed and the growth rate is increased by a factor of three (Fig. 36). At high values of (Y the trend is quite different [7]: increasing &- stabilizes the mixing layer.

6. Conclusions The first spectral collocation linear stability code for compressible flow (SPECLS) is utilized for high-speed flows. The technique employs the first pressure staggered grid for a spectral discretization applied to compressible flows. This accurate discretization can be employed in

124

M.G. Macaraeg,

-1.3

;a__

-100

1.0

C.L. Streett / Linear stability

-50

I

.6L

lb,-6

-1.3 -100

l

Y 18

I

-60

-40

-20

3

20

40

60

-60

-40

-73

0

70

40-

60

,

I

I 80

100

80

100

Li____J -100

-SC

Fig. 34. Pressure disturbance eigenfunctions for higher supersonic modes described by Fig. 32: (a) supersonic respect to the injected gas, a = 2.06, (b) supersonic with respect to the quiescent gas, (Y= 0.925, (c) supersonic respect to both gases, a = 0.925.

with with

M. G. Macaraeg,

C. L. Streett / Linear stability

125

Y

Fig. 35. Temperature

disturbance

eigenfunction

for mode of Fig. 34(c): (a) full spatial scale: [ - 100, 1001, (b) expanded spatial scale: [ - 10, lo].

nonlinear simulations. Verification cases for high-speed boundary layers indicate an order of magnitude reduction in the number of points required to obtain equivalent accuracy in growth rates with a finite difference formulation. In addition, a multi-domain spectral discretization (MDSPD) in SPECLS is found to require a factor of three less points, which is significant since the operation count of the spectral formulation goes as the cube of the number of points. The highly irregular structure of the supersonic modes are easily resolved by a MDSPD which is shown to be highly robust over a wide range of stretching parameters and resolutions. The stability analysis of a compressible shear flow is presented. The study indicates that stability of mixing layers is enhanced by increasing Re, increasing Mach number and higher temperatures of the injected gas. The appearance of a multitude of higher supersonic modes is found to occur at about A4, = 3, /IT = 1 and p = 0, and is attributed to the far-field reflecting

126

M.G. Macaraeg, CL. Streett

/ Linear stability

.06

.04

a Fig. 36. Inhibition

of higher supersonic

mode formation;

M, = 3, T, = 520 o R, p = 0.

boundary conditions (bounded free shear flow). The effect of a finite outer extent is found to have a significant impact on supersonic disturbances, despite the fact that the far field extends roughly 1500 momentum thicknesses in distance. At M, = 3, formation of these multiple supersonic modes is suppressed and growth rate greatly increased by increasing & for these low values of streamwise wave numbers. Controlling the development of supersonic mode formation may be a key in inducing transition at high Mach numbers which will greatly enhance mixing efficiency.

Acknowledgement The authors would like to thank M.Y. Hussaini, Institute for Computer Applications in Science and Engineering, for helpful discussion during the course of this work. The authors are also deeply appreciative for the assistance of L.M. Mack, Jet Propulsion Laboratory, during the verification stages of the spectral code, and for providing insight into the results obtained. Further acknowledgment is given to B. Sekar for performing the computations for Fig. 13.

References [l] W. Blumen, Shear layer instability of an inviscid compressible fluid, .I. Fluid Mech. 40 (1970) 769-781. [2] G.L. Brown and A. Roshko, On density effects and large structure in turbulent mixing layers, J. Fluid Mech. 64 (1974) 775-816. [3] J.P. Drummond, Numerical simulation of a supersonic, chemically reacting mixing layer, NASA TM 4055 (1988). [4] H. Gropengiesser, Study on the stability of boundary layers and compressible fluids, NASA TTF-12, 786 (1969). [5] M.G. Macaraeg and C.L. Streett, Improvements in spectral collocation discretization through a multiple domain technique, Appl. Numer. Math. 2 (1986) 95-108.

M. G. Macaraeg,

C. L. Streett / Linear stability

121

[6] M.G. Macaraeg and C.L. Streett, A spectral multi-domain technique for viscous, compressible, reacting flows, Presented at 5th International Conference on Numerical Methods in Laminar and Turbulent Flow, Montreal, Que. (1987); Internat. J. Numer. Methods Fluids (1988). [7] M.G. Macaraeg and C.L. Streett, New instability modes for bounded, free shear flows, Phys. Fluids A (submitted). [8] M.G. Macaraeg, C.L. Streett and M.Y. Hussaini, A spectral multi-domain technique applied to high-speed chemically reacting flows, in: Proceedings Second International Symposium on Domain Decomposition Methods, Los Angeles, CA (1988). [9] L.M. Mack, Boundary layer stability theory, NASA Cr-131591, Jet Propulsion Lab. (1969). [lo] M.R. Malik and S.A. Orszag, Linear stability analysis of three-dimensional compressible boundary layers, J. Sci. Comput. 2 (1) (1987) 77-97. [ll] M.V. Morkovin, Transition at hypersonic speeds, NASA Contract No. NASI-18107 (1987). (121 D. Papamoschou and A. Roshko, Observation of supersonic free shear layers, AIAA Paper 86-0162 (1986). [13] S.A. Ragaab and J.L. Wu, Instabilities in the free shear layer formed by two supersonic streams, AIAA Paper 88-0038, presented at the AIAA 26th Aerospace Sciences Meeting, Reno, NV (1988). [14] C.L. Streett, Spectral methods and their implementation to solution of aerodynamic and fluid mechanic problems, Presented at the Sixth International Symposium on Finite Element in Flow Problems, Antibes, France (1986). [15] C.L. Streett, T.A. Zang and M.Y. Hussaini. Spectral methods for solution of the boundary-layers equations, AIAA Paper 84-0170, presented at the 22nd AIAA Aerospace Sciences Meeting, Reno, NV (1984). [16] J.H. Wilkinson, The AIgebraic Eigenualue Problem (Clarendon Press, Oxford, 1965).