A k-ω Turbulence Model for Quasi-Three-Dimensional Turbomachinery Flows
Rodrick V. Chima* NASA Lewis Research Center Cleveland, Ohio 44135

Abstract
A two-equation k-ω turbulence model has been developed and applied to a quasi-three-dimensional viscous analysis code for blade-to-blade flows in turbomachinery. The code includes the effects of rotation, radius change, and variable stream sheet thickness. The flow equations are given and the explicit Runge-Kutta solution scheme is described. The k-ω model equations are also given and the upwind implicit approximate-factorization solution scheme is described. Three cases were calculated: transitional flow over a flat plate, a transonic compressor rotor, and a transonic turbine vane with heat transfer. Results were compared to theory, experimental data, and to results using the Baldwin-Lomax turbulence model. The two models compared reasonably well with the data and surprisingly well with each other. Although the k-ω model behaves well numerically and simulates effects of transition, freestream turbulence, and wall roughness, it was not decisively better than the Baldwin-Lomax model for the cases considered here.

The Baldwin-Lomax model is popular because it is easy to implement (at least in 2-D) and works fairly well for predicting overall turbomachinery performance. However, the model has both numerical and physical problems. Numerical problems include awkward implementation in 3-D, difficulty in finding the length scale [2], and slow convergence if the length scale jumps between grid points. Physical problems include a crude transition model and the neglect of freestream turbulence, surface roughness, and mass injection effects which are often important in turbines. These effects are sometimes added to the Baldwin-Lomax model using techniques developed for boundary layer codes [3, 4]. Physical problems also include poor prediction of separation [5], which is important in compressors, and underprediction of wake spreading [2]. A few researchers have used other turbulence models for turbomachinery problems. Choi et. al. have used the qω model [6], Hah (who participated in the blind test case) used a k-ε model [7], and Kunz and Lakshminarayana used an algebraic Reynolds stress k-ε model [8]. Unfortunately none of these researchers have used a BaldwinLomax model in the same code for comparison. Ameri and Arnone have compared the q-ω, k-ε, and BaldwinLomax models for turbine heat transfer problems [9, 10]. Two papers have compared the k-ω and BaldwinLomax models for turbomachinery problems. Bassi, et. al. examined a film-cooled turbine cascade [11], and Liu et. al. examined a low pressure turbine cascade [12]. Both papers compared the computed results primarily with experimental pressure distributions. In the present work the k-ω model developed by Wilcox [13] was incorporated in the author’s quasi-threedimensional (quasi-3-D) turbomachinery analysis code [14]. The code includes the effects of rotation, radius change, and stream surface thickness variation, and also includes the Baldwin-Lomax turbulence model. The k-ω model was chosen for several reasons. First, the effects of freestream turbulence, surface roughness, and mass injection are easily included in the model [13]. Second, transition can be calculated using the low-Reynolds-number version of the model [15]. Third, Menter has shown that the k-ω model does well for flows with adverse pressure gradients [5 ,16]. Finally, the k-ω model should behave well numerically since it avoids the use of the distance to the wall and complicated damping functions.

Introduction
A large percentage of computational fluid dynamics (CFD) analysis codes for turbomachinery use the Baldwin-Lomax turbulence model [1]. This was evident in the results of the blind test case for turbomachinery codes sponsored by ASME/IGTI at the 39th International Gas Turbine Conference held in The Hague in June of 1994. The results have not yet been published. Of the 12 participants, nine used the Baldwin-Lomax turbulence model, one used an algebraic mixing length model, and two used k-ε models. One of the objectives of that test case was to investigate the effects of turbulence models. However, because of differences in grids, large variations between the computed solutions, and lack of experimental measurements in the boundary layers, it was not possible to draw any conclusions regarding the effect of turbulence models.
*Aerospace Engineer, Associate Fellow AIAA Copyright 1995 by AIAA, Inc. No copyright is asserted in the United States under Title 17, U.S. Code. The U.S. Government has a royalty-free license to exercise all rights under the copyright claimed herein for Governmental purposes. All other rights are reserved by the copyright owner.

1

q = J – 1 ρ, ρu, ρvr, e

T T

θ m r(m)

h(m)

= J – 1 0, K 2, 0,

0

(3)
T

S = J – 1 0, S 2, S 3, S 4 ρU ρuU + ξ m p

E = J

–1

( ρvU + ξ θ p )r ( e + p )U + ξ θ rΩp ρV ρuV + η m p ( ρvV + η θ p )r ( e + p )V + η θ rΩp

(4)

z
Figure 1. Quasi-three-dimensional stream surface for a compressor rotor.

This paper describes the quasi-3-D flow equations and the explicit Runge-Kutta scheme used to solve them. The paper also gives the k-ω equations written in a quasi-3-D form, gives the boundary conditions, and describes the implicit upwind ADI scheme used to solve the turbulence model equations. The model was tested on three cases and compared to the Baldwin-Lomax model and to experimental data. The cases included a flat plate boundary layer with transition, a transonic compressor rotor, and a transonic turbine vane.

F = J

–1

(5)

1 e = ρ C v T + -- ( u 2 + v 2 ) 2 is the total energy per unit volume, 1 p = ( γ – 1 ) e – -- ( u 2 + v 2 ) 2 is the pressure, rm hm 1 dr 1 dh ----- = -, and ----- = -r dm hdm r h are derivatives of the streamtube geometry, and

(6)

(7)

Quasi-3-D Navier-Stokes Equations
The Navier-Stokes equations have been developed in an ( m, θ ) coordinate system as shown in figure 1. Here m is the arc length along the surface, dm = dz + dr
2 2 2

(8)

(1)

and the θ-coordinate is fixed to the blade row and rotates with angular velocity Θ. The radius r and the thickness h of the stream surface are assumed to be known functions of m. The equations have been mapped to a body-fitted coordinate system, simplified using the thin layer approximation, and nondimensionalized by arbitrary reference quantities ρ 0 , c 0 , and µ 0 . The Reynolds number Re and Prandtl number Pr are defined in terms of these reference quantities. The final equations are given in [14] and are summarized below. ∂ t q + ∂ ξ E + ∂ η ( F – Re S ) = K where
–1

rm hm –1 –1 K 2 = ( ρv 2 + p – Re σ 22 ) ----- + ( p – Re σ 33 ) ----- (9) r h The viscous fluxes are given by S 2 = η m σ 11 + η θ σ 12 S 3 = ( η m σ 12 + η θ σ 22 )r µ 2 - 2 S 4 = --------------------- ( η m + η θ )∂ η a 2 + uS 2 + vS 3 ( γ – 1 )Pr a 2 = γ p ⁄ ρ is the speed of sound squared. Using Stokes’ 2 hypothesis, λ = – -- µ , the shear stress terms are given by 3

(10)

(2)

2

σ 11 = 2µ∂ m u + λ∇ ⋅ V 2µ σ 22 = ------ ( ∂ θ v + ur m ) + λ∇ ⋅ V r hm σ 33 = 2µu ----- + λ∇ ⋅ V h rm 1 - σ 12 = µ  ∂ m v – v ----- + -- ∂ θ u   r r r m hm 1 2 λ∇ ⋅ V = – -- µ ∂ m u + u  ----- + -----  + -- ∂ θ v r 3 h r U and V are relative contravariant velocities U = ξm u + ξθ w V = ηm u + ηθ w

A spatially-varying time step and implicit residual smoothing are used to enhance convergence. Details of the solution scheme used here are given in (18) and are described briefly below. (11) The discrete equations are solved using the explicit multistage Runge-Kutta scheme developed by Jameson, Schmidt, and Turkel [19]. A four-stage scheme is used. For efficiency, physical and artificial dissipation terms are computed only at the first stage. The Baldwin-Lomax model is updated every five time steps. The k-ω model is usually updated every two time steps with twice the ∆t of the flow solution. The spatially-varying time step is calculated as the harmonic mean of inviscid and viscous components in each grid directional. Artificial dissipation consisting of blended second and fourth differences is added to prevent point decoupling and to enhance stability. Eigenvalue scaling, as introduced by Martinelli and Jameson [20] but modified by Kunz and Lakshminarayana [8], is used to weight the artificial dissipation in each direction. The scaling is based on a blend of the one-dimensional time step limits at each point. The artificial dissipation is also reduced linearly by grid index near the wall and wake centerline to minimize the effects on the boundary layer. The explicit four-stage Runge-Kutta scheme has a Courant stability limit of about 2.8. Implicit residual smoothing introduced by Jameson and Baker in [21] can be used to increase the time step, and hence the convergence rate, by a factor of two to three. On high aspect ratio grids the stability limit is dominated by the grid spacing in the finest direction, and it is sufficient to use implicit smoothing in that direction only. The stability analysis given in [21] is used to calculate the smoothing parameter required at each point, then the same Eigenvalue scaling used for the artificial dissipation is used to reduce or eliminate the smoothing parameter in grid directions where it is not needed. The use of Eigenvalue scaling for both the artificial dissipation and implicit smoothing greatly increases the robustness of the numerical scheme.

(12)

where w is the relative tangential velocity, w = v – rΘ . The θ-metrics are scaled by 1 ⁄ r and the Jacobian is scaled by rh . The metric terms are found using central differences and θ –mη ⁄ r ξm ξθ = J η ηm ηθ –θξ mξ ⁄ r J = [ rh ( m ξ θ η – m η θ ξ ) ] – 1 The effective viscosity is µ = µL + µT (14)

(13)

where the molecular (laminar) viscosity µ L is evaluated using a power law function of the temperature, and the turbulent viscosity µ T is evaluated using either the BaldwinLomax model [1] or Wilcox’s k-ω model [13, 15]. Minor modifications to the coefficients and blending functions used in the Baldwin-Lomax model are described in [2]. Boundary Conditions At the inlet the total pressure, total temperature, and tangential velocity component are specified and the upstream-running Riemann invariant based on the axial velocity is extrapolated from the interior. At the exit, three of the four conserved variables are extrapolated and the θaveraged pressure is specified using the method described by Giles in ref [17]. Periodic boundaries between the blades are solved like interior points using a dummy grid line outside the domain. Multistage Runge-Kutta Scheme The flow equations are discretized using finite differences and solved using an explicit Runge-Kutta scheme.

k-ω Turbulence Model
The k-ω turbulence model was first postulated by Kolmogorov in 1942 and later independently by Saffman in 1970 (see Wilcox’s book [13] for references.) It has been under development by Wilcox for many years and is described in detail in [13]. The model solves two turbulence transport equations for the turbulent kinetic energy k and the specific dissipation rate ω. The model has a basic formulation for fully turbulent flows that satisfies the law of the wall without knowledge of the distance to the wall or complicated near-wall damping terms. There is also a low-Reynolds-number formulation used for modeling

3

transition [15]. Boundary conditions can be specified to simulate mass injection or surface roughness. Most of Wilcox’s development of the model used boundary layer codes, but recently Menter has shown several applications to Navier-Stokes codes [16]. Menter found that the model exhibited strong dependence on freestream values of ω and proposed a somewhat ad hoc fix. In this work many of Menter’s suggestions for numerical implementation of the model have been used, but his fix for the problem of freestream dependence has not. The quasi-3-D form of the k-equation has been derived by writing the m- and θ-momentum equations in non-conservative form, multiplying each by its fluctuating velocity component, and Favre averaging. The usual turbulence modeling approximations are made, i. e., the Boussinesq model is used for the Reynolds stress terms, pressure work, diffusion, and dilatation are all neglected, and turbulent dissipation is taken to be proportional to k × ω . The production term is written in terms of the vorticity magnitude using Menter’s suggestion [16]. Source terms that arise from the quasi-3-D equations are neglected. The ω-equation is derived from the k-equation by dimensional considerations. Wilcox’s constants are used without modification. The final form of the model equations is as follows: J 1 ∂ t q + U ∂ ξ q + V ∂ η q – Re– 1 -- G = -- ( P – D ) ρ ρ where q = k, ω
T

* D --- = β ωk ρ βω 2

(21)

The baseline k-ω model has five coefficients: β = 3 ⁄ 40 , β * = 9 ⁄ 100 , σ = 1 ⁄ 2 , σ *= 1 ⁄ 2 , α = 5 ⁄ 9 , and the trivial constant α * = 1 . The low-Reynolds-number model replaces three of the constants with the following bilinear functions of the turbulence Reynolds number Re T : β * = ( 9 ⁄ 100 )F β , α = ( 5 ⁄ 9 ) ( F α ⁄ F µ ) , and α * = F µ , where 5 ⁄ 18 + ( Re T ⁄ R β ) 4 F β = ---------------------------------------------1 + ( Re T ⁄ R β ) 4 α 0 + Re T ⁄ R ω F α = --------------------------------1 + Re T ⁄ R ω
* α 0 + Re T ⁄ R k F µ = -------------------------------1 + Re T ⁄ R k

(22)

ρk Re T = ---------µL ω
* with α 0 = 1 ⁄ 10 , α 0 = β ⁄ 3 = 1 ⁄ 40 ,

(15)

and R β = 8 , R ω = 27 ⁄ 10 , R k = 6 . (16) Boundary Conditions (17) At the inlet the turbulence intensity Tu and turbulent viscosity µ T are specified. Then k and ω are found from 3 2 k = -- Tu 2 U in 2 ρk ω = α * ----µT

ρk µ T = α * ----ω

The molecular plus turbulent diffusion terms G are written using the thin-layer approximation giving
2 2 ( η m + η θ ) ( µ + σ * µ T )∂ η k G = -----------------------J ( µ + σµ T )∂ η ω

(23)

(18)

Menter’s form of the production terms is used [16]. µT 2 1 -Ω Re– 1 ----- 2 – -- k  ∂ m u + -- ∂ v ρ 3  r θ  α α*Ω2 2 1 – -- ω  ∂ m u + -- ∂ v 3  r θ 

where α * = 1 for the baseline model or α * = F µ for the low-Reynolds-number model. Substituting equation (22) (19) for α * into equation (17) for µ T gives a quadratic for Re T . The solution is µT µT 1 * * * Re T = -- – α 0 R k – -----  + (α 0 R k ) 2 + 2R k ----- (2 – α 0 ) 2  µL µ L (24) (20) and ω may be found from

P -- = ρ

where rm 1 Ω = ∂ m v – -- ∂ θ u + v ----r r is the vorticity. The destruction terms are given by

4

ρk ω = ---------------µ L Re T

(25)

Using first-order backward Euler time differencing and linearizing the right hand side about the previous time step gives   J 1  I + ∆t U ∂ ξ + V ∂ η – Re – 1 -- ∂ η G′ – -- ( P′ – D′ ) ∆q ρ ρ   1 J = – ∆t U ∂ ξ q + V ∂ η q – Re – 1 -- ∂ η G – -- ( P – D ) ρ ρ (31)

A turbulent length scale can be defined using (17) as ρk α* k µ T = α * ----- = ρ k ------------- = ρ k l  ω  ω (26)

The effects of varying the inlet values of ω or l is discussed with the results. On solid walls k = 0, and ω is set using Wilcox’s roughness model.
2 µτ ∂u ω = ----- S R = S R ν ∂ y wall

where primes indicate Jacobians. The advective terms are approximated using first-order upwind differences. U ∂ ξ q ≈ U ( q i, j – q i – 1, j ) – U ( q i, j – q i + 1, j )
+ -

(27)

1 U ± = -- ( U ± U ) 2 The diffusive terms (18) can be written as

(32)

where   -----  +  kR SR =   100  -------k  R   50 2  + k R < 25    + k R ≥ 25   

(28)

2 2 gk ∂η k ( η m + η θ ) ( µ + σ * µ T )∂ η k G = -----------------------= gη ∂η ω J ( µ + σµ T )∂ η ω

(33)

+ and k R is the equivalent sand grain roughness height in + turbulent wall units. For all results shown here k R was set

Since σ * = σ = 1 ⁄ 2 , g k = g ω = g , and the Jacobian G v ′ is simply ∂ ∆k ∂G ∆q = g η ∂q ∂ η ∆ω

to 5, giving a hydraulically smooth surface. To avoid numerical difficulties near leading edges where µ τ becomes large, an upper limit was imposed on ω using a boundary condition suggested by Menter [16]. 10 6ν ω max = ----- × -----------Re β∆y 2 where ∆y is the grid spacing at the wall. k and ω were extrapolated at the exit and treated as periodic across trailing edge wake cut lines and between blade rows. ADI Solution Scheme An alternating direction implicit (ADI) scheme was used to solve the k-ω equations. An implicit scheme was chosen so that the equations could be updated less often than the flow equations without stability problems, and also because the k-ω equations are dominated by complex source terms which are evaluated only once in an ADI scheme, but potentially at every stage of a Runge-Kutta scheme. Equation (15) may be written as J 1 ∂ t q = – U ∂ ξ q + V ∂ η q – Re – 1 -- ∂ η G – -- ( P – D ) ρ ρ (30) (29) G′∆q =

(34)

The diffusive terms are approximated by second-order central differences. 1 ∂ η g∂ η k ≈ -- [ ( g i, j + 1 + g i, j ) ( k i, j + 1 – k i, j ) 2 + ( g i, j – 1 + g i, j )( k i, j – 1 – k i, j )] The difference approximations (32) and (35) are diagonally dominant and have zero row sum, which according to Baldwin and Barth [22], makes the implicit operator an M-type matrix with a non-negative inverse. Menter’s linearization is used for the source terms [16]. The production terms are treated explicitly, i.e., P′ = 0 . The destruction terms are linearized using
* * D′ ----- = β ω β k ρ 0 2βω

(35)

(36)

The term in the upper right corner β *k is the only coupling between the k and ω equations. Here the term has been neglected and the equations are solved uncoupled from each other and from the flow equations. Equation (31) is solved using an approximate factorization.

5

laminar 1/5 power law

1000 x C f

turbulent 1/7 power law

u+

k-ω model Baldwin-Lomax model 10 0 Rex 10 1

Spalding formula k-ω model Baldwin-Lomax model 10 2 y+ 10 3 10 4

Figure 2.

Theoretical and computed skin friction for a flat plate boundary layer. D′ I + ∆t  U ∂ ξ + -----  ρ J I + ∆t  V ∂ η – Re – 1 -- ∂ η G v ′ ∆q   ρ

Figure 3.

Theoretical and computed velocity profiles for a flat plate boundary layer.

Results
Flat Plate Boundary Layer

J 1 = – ∆t U ∂ ξ q + V ∂ η q – Re – 1 -- G v – -- ( P – D ) ρ ρ

(37)

The destruction terms D′ are included with the streamwise ξ operator. Treating the destruction terms implicitly improves the diagonal dominance but gives the implicit operator a nonzero row sum. The Baldwin-Lomax and k-ω models have both been coded fairly efficiently. A turbulent flow solution updating the Baldwin-Lomax every time step model takes 1.53 times as long as a laminar solution. Past experience has shown that it is sufficient to update the Baldwin-Lomax model every five time steps, reducing the CPU time to about 1.1 times that of a laminar solution. A solution updating the k-ω equations takes 1.6 times as long as a laminar solution. Some cases have been successfully computed updating the k-ω equations every five time steps (with five times the ∆t of the flow solution), making the kω model nearly as fast as the Baldwin-Lomax model. Other cases failed to converge unless the k-ω equations were updated every other time step, and that strategy has been used for all results reported here. The net result is that a flow solution with the k-ω model takes about 1.18 times as long as a solution with the Baldwin-Lomax model.

A turbulent flat plate boundary layer test case was used to verify the code. The flow had a Mach number of 0.3, a Reynolds number based on plate length of 10 7 , and a freestream turbulence of one percent. The grid had 127 points streamwise with 12 points upstream of the leading edge, and was stretched exponentially to resolve the transition region. There were 41 points above the plate with a wall spacing of y + = 2 to 4 . Figure 2 compares the computed skin friction using both the Baldwin-Lomax and k-ω models to a 1/5 power law correlation for laminar flow and a 1/7 power law correlation for turbulent flow. The k-ω results agree closely with the turbulent correlation and the Baldwin-Lomax results are slightly low. The Baldwin-Lomax model was run with the transition model given in their original paper [1], that is, the turbulent viscosity was calculated for each grid line but set to zero for the whole line if the maximum was less than cmutm = 14 times the laminar viscosity. The BaldwinLomax model transitions at Re x = 3 ×10 , while the k-ω model transitions at Re x = 5 ×10 . Both results are within the experimentally observed range. The k-ω model transition point was found to be strongly dependent on the freestream turbulence level, but fairly independent of the freestream value of ω. However, the transition point calculate for a turbine vane shown later was found to be very sensitive to the freestream value of ω. The BaldwinLomax transition point was strongly dependent on the
5 5

6

0.

60

0.

80

0.

80
0. 90

1.

00
0. 80

0.

80

0. 90

0.

1.

40
0. 70

0.

70

0.

80
0 1.4

0.

80

0.

00

90

0. 90

1.
1.4

40

computed Baldwin-Lomax / k-ω

measured

Figure 4. Computed and measured contours of relative Mach number for a transonic compressor

rotor.

Baldwin-Lomax model

k-ω model

Figure 5.

Turbulent viscosity contours around a transonic compressor rotor.

7

0. 80

1.

0. 80

0. 80

60

0. 80

0. 70

0.

70

0. 80

90 0.

00 1.

40 1.
1.4 0

0 1.4

0 1.4

1.50

0 1.3
1.50
1. 4 0
0 1.4
1.4

30 1.

90 0.

00 1.

0

1.50

0

0 1.3
1.50
1.4 0

1. 30 1. 30

1.4 0

1.0

1.2

0.8

1.0

M rel

0.6

0.4 experimental data k-ω model Baldwin-Lomax model 0.2 0 40 80 120 percent pitch 160 200

M rel

0.8

0.6 experimental data k-ω model Baldwin-Lomax model 0.4 0 40 80 120 percent pitch 160 200

Figure 6.

Computed and measured near wake Mach number profiles for a transonic compressor rotor.

Figure 7.

Computed and measured far wake Mach number profiles for a transonic compressor rotor.

value of the cutoff parameter cmutm, suggesting that this parameter could be calibrated to simulate freestream turbulence effects. Figure 3 compares computed velocity profiles located at the end of the plate to Spalding’s composite sublayer/ law-of-the-wall profile. Both results agree closely, but the Baldwin-Lomax profile has a slightly stronger wake component. Transonic Compressor Rotor The transonic compressor rotor described by Suder et al. [3] was used to test the quasi-3-D effects in the model. The rotor was tested experimentally at NASA Lewis Research Center using both laser anemometry and conventional aero probes. A section of the rotor at 70 percent span was analyzed. The radius was specified as a line 70 percent of the way between the hub and shroud, and the stream surface thickness was specified as the local distance between the hub and shroud, normalized to one at the inlet. A meridional view of the streamtube is shown in figure 1. A C-type grid was used with 319 points around the blade and 45 points away from the blade. The grid spacing gave y + < 2 over most of the blade. The calculations were run 2000 iterations, which took about 3.25 minutes for the k-ω model on the Cray C-90 called eagle at NASA Ames Research Center. Figure 4 compares contours of relative Mach number computed with the k-ω model (the Baldwin-Lomax model gives identical contours) to contours measured experimentally using laser anemometry. The inlet Mach number is

about 1.4. The flow passes through a weak upstream-running wave system, then through a strong normal passage shock, and leaves the rotor at about Mach 0.8. The computed Mach number behind the normal shock is somewhat lower than the measured values. This is partly due to the assumed stream surface but may also be partly due to the inability of either turbulence model to capture the shockboundary layer interaction correctly. One motivation for using a turbulent transport model was that algebraic models frequently fail to find the correct length scale and thus give nonsmooth values of µ T [2]. Occasionally nonconvergence or instability problems can be traced to poor numerical behavior of the turbulence model. Figure 5 shows contours of µ T computed using the Baldwin-Lomax model (left) and the k-ω model (right.) Figure 4 also shows that the computed results seriously underpredict the wake spreading. This is shown quantitatively in figures 6 and 7. Figure 6 compares computed and measured near-wake profiles about 0.28 chords downstream of the trailing edge while figure 7 compares far-wake profiles about 2 chords downstream. The two turbulence models give surprisingly similar results, especially considering the erratic behavior of the BaldwinLomax model seen in figure 5. The computed wakes are both narrower and deeper than the measured wakes. Wilcox has shown in [13] that his model gives the best prediction for planar wake spreading as ω in → 0 . This limit corresponds to µ T → ∞ , which seems unreasonable for an inlet value. Varying ω in by five orders of magnitude had very little effect on the

8

max.

.2

.3 Baldwin-Lomax .8 .9 1. residual rms k-ω

.2

.3 .9

Baldwin-Lomax k-ω .8 .9 1. Iterations .9 .8 .9 .9 .8 .9

Figure 9.

Residual histories for the VKI turbine vane computations.

from a Mach number of about 0.15 at the inlet to about 0.9 at the exit. Surface heat transfer was converged to plotting accuracy in 3000 iterations in fully turbulent regions. Laminar parts of the flow took longer to converge, so all calculations were run 5000 iterations. Solution times were about eight minutes per case on the Cray C-90 computer. A typical residual history is shown in figure 9. The k-ω calculations converged monotonically, but with the BaldwinLomax model the maximum residual oscillated near the round trailing edge. Computed distributions of isentropic surface Mach number are compared to experimental data for 2,is = 0.875 and 1.02 in figure 10. The Baldwin-Lomax and k-ω models give identical results and are not shown separately. The subsonic results agree very well with the experimental data. The transonic results slightly underpredict the Mach number on the rear (uncovered) part of the suction surface. All subsequent results are for the subsonic case. Computed wake profiles located 43 percent of axial chord downstream of the trailing edge are compared to the experimental data (digitized manually from [23]) in figure 11. Again the computed wakes are narrower and deeper than the measured wakes, but here the k-ω results are slightly better than the Baldwin-Lomax results. Figures 12 - 15 show the effects of various parameters on surface heat transfer coefficient H [ W/(m K ) ] . In each figure the abscissa is the arc length S [mm] along the vane surface. Figure 12 compares computations using the baseline k-ω model, the low-Reynolds-number k-ω model, and the Baldwin Lomax model to experimental data for
2

Figure 8.

Computed Mach number contours for the VKI turbine vane.

computed wake spreading. The freestream turbulence was set to three percent, and doubling it had little effect on the wake spreading. Transonic Turbine Cascade A transonic turbine vane tested by Arts, et al. [23] was computed as a third test case. The vane was tested experimentally in the Isentropic Light Piston Compression Tube Facility at the von Karman Institute. The facility has independent control over the exit Reynolds number Re 2,is , the exit Mach number, M 2,is , and the inlet turbulence intensity Tu. Surface pressures were measured with static taps, and wake total pressure profiles were measured with a high-speed traversing probe. The vanes were initially at 300 K and the freestream temperature was 415 K. Unsteady blade surface temperatures were measured during a run using platinum thin film gauges, then converted to heat fluxes using a one-dimensional semi-infinite-body model. For the computations a C-type grid was used with 383 points around the vane and 49 points away from the vane. The grid spacing gave y + < 1.5 over most of the vane. Blade shapes and computed Mach contours for 2,is = 0.9 are shown in figure 8. The flow accelerates

9

experimental data Baldwin-Lomax model k-ω model

M 2,is 1.02 0.875 0 20 40 60

expt. comp.

80

100

Figure 10. Computed and measured distributions of

isentropic Mach number for the VKI turbine vane.

1-P /P 02 01

M isen

Figure 11. Computed and measured total pressure profiles

0.43 chords behind the VKI turbine vane. µ T,in /µ ref 7.25e-2 7.25e-4 2.50e-4 6.25e-5

baseline k-ω low Re k-ω Baldwin-Lomax expt., Tu = 4% H ( W / m2 K) H ( W / m2 K)

l in (mm) 1.0e 1 1.0e-1 3.5e-2 9.0e-3

comp.

pressure surface

S (mm)

suction surface

pressure surface

S (mm)

suction surface

Figure 12. Surface heat transfer coefficient predicted by

three turbulence models.

Figure 13. Effects of inlet turbulent length scale or

viscosity on predicted heat transfer coefficient. comp. Re2,is expt. 2. M 1. M .5 M H ( W / m2 K) comp.

Tu (%) expt. 6 4 1 H ( W / m2 K)

pressure surface

S (mm)

suction surface

pressure surface

S (mm)

suction surface

Figure 14. Effects of inlet turbulence intensity on heat

Figure 15.

transfer coefficient.

coefficient

Effects of Reynolds number on heat transfer

10

Re 2, is = 1 ×10

6

and Tu = 4 percent. Triangles show the

experimental data. The baseline k-ω solution is fully turbulent on the suction surface giving high values of H. The pressure surface has a highly favorable pressure gradient and acts laminar over part of the chord, transitioning near S = -20. The low-Reynolds-number k-ω solution remains laminar on the pressure surface and transitions near the measured transition point on the suction surface; however, the transition point was forced by choice of µ T,in , as discussed later. The Baldwin-Lomax solution agrees closely with the low-Reynolds-number k-ω solution. As discussed with figure 14, the laminar parts of the flow have augmented heat transfer due to freestream turbulence that none of the models predict. Figure 13 shows the effects of µ T,in or the corresponding length scale l in on H for the same flow conditions. The inlet turbulent viscosity was varied by about 3 orders of magnitude to produce solutions that ranged from fully laminar to almost fully turbulent. Only a small range of values of µ T,in gave transition near the measured location. Corresponding turbulence length scales are shown on the figure. The length scale that gives the best transition location, l in = 3.5 ×10 , is about 1/1600 times the pitch or about 7 times the grid spacing at the wall. All subsequent calculations were run with this length scale. Although this strong dependence on µ T,in is disconcerting since appropriate values are not known in advance, Wilcox points out that transition in real flows is not simply a function of Tu, but is also frequency dependent [15]. He also suggests that two coefficients in the low-Reynolds-number model could be adjusted to better match other flows, although this has not been attempted here. The effects of freestream turbulence intensity Tu are shown in figure 14. The effect on suction surface transition location is modeled reasonably well by the low-Reynolds-number k-ω model. The experimental data shows a strong augmentation of heat transfer at the leading edge and on the pressure surface as Tu is increased. Although the low-Reynolds-number k-ω model depends directly on Tu, in laminar regions the model gives values of turbulent viscosity that are much too small to affect the heat transfer. Boyle has added algebraic correlations to the BaldwinLomax model to simulate freestream turbulence effects [4]. These correlations could be added to the k-ω code, but they do involve the distance from the wall. The effects of Reynolds number e 2,is are shown in figure 15 for Tu = 4 percent. The data shows a large increase in heat transfer with e 2,is . The k-ω results show qualitative agreement, but the magnitude of the heat transfer is underpredicted due to the failure of the model to pre–2

dict freestream turbulence effects. Effects on suction surface transition location are overpredicted, and predicted transition is too abrupt. The data shows transition on the pressure surface at the highest Reynolds number that is not predicted, although a small change in µ T,in does cause pressure surface transition, as shown in figure 13.

Concluding Remarks
Wilcox’s k-ω turbulence model has been added to a quasi-3-D Navier-Stokes analysis code for turbomachinery. The code includes the effects of rotation, radius change, and stream sheet convergence, and also included the Baldwin-Lomax turbulence model. The quasi-3-D flow equations and boundary conditions were described. An explicit multistage Runge-Kutta scheme with spatiallyvarying time step and implicit residual smoothing was used to solve the flow equations. The quasi-3-D form of the k-ω model equations and boundary conditions were also described. An upwind implicit ADI scheme was used to update the turbulence model equations uncoupled from the flow equations. The numerical scheme was quite robust, but about 18 percent slower than the BaldwinLomax model. Calculations were made for three test cases: a flat plate boundary layer with transition, a transonic compressor rotor with significant quasi-3-D effects, and a transonic turbine vane. The flat plate calculations agreed very well with theory for both turbulence models. Transition predictions were reasonable for both models and suggested that the Baldwin-Lomax transition model could be calibrated to simulate free stream turbulence effects. The compressor rotor calculations showed very close agreement between the two turbulence models, but both models failed to capture the measured wake spreading. The turbine calculations showed very good agreement with measured surface pressures for both turbulence models. Predicted wake profiles were thinner and deeper than measured profiles, although the k-ω model gave marginally better results. The Baldwin-Lomax model, the baseline k-ω model, and the low-Reynolds-number k-ω model were compared for heat transfer calculations. The Baldwin-Lomax model did reasonably well considering the simple transition model used. The low-Reynolds-number k-ω model showed a high sensitivity to inlet values of ω, expressed as an inlet turbulent viscosity or length scale which are not generally known. The k-ω model was able to capture the effects of inlet turbulence intensity on transition but not on augmentation of heat transfer in laminar regions. It may be possible to model this effect with a simple algebraic model. Effects of Reynolds number were predicted qualitatively. The k-ω model exhibited some attractive numerical properties, but for the cases considered here, predictions were not decisively better than those made with the Baldwin-Lomax model. Other test cases may identify areas

11

where the k-ω model is significantly better than the Baldwin-Lomax model. Future work will extend the k-ω model to three dimensions where algebraic models are poorly defined and difficult to implement.

11. Basi, F., Rebay, S., and Savini, M., “A NavierStokes Solver with Different Turbulence Models Applied to Film-Cooled Turbine Cascades,” in Heat Transfer and Cooling in Gas Turbines, AGARDCP-527, Oct. 1992, pp. 41.1 to 41.16. 12. Liu, F., and Zheng, X., “Staggered Finite Volume Scheme for Solving Cascade Flow with a k-ω Turbulence Model,” AIAA Journal, Vol. 32, No. 8, Aug. 1994, pp. 1589-1597. 13. Wilcox, D., C., Turbulence Modeling for CFD, DCW Industries, Inc., La Canada, CA, 1994. 14. Chima, Rodrick V., “Explicit Multigrid Algorithm for Quasi-Three-Dimensional Viscous Flows in Turbomachinery,” AIAA Journal of Propulsion and Power, Vol. 3, No. 5, Sept.-Oct. 1987, pp. 397-405. 15. Wilcox, D., C., ‘‘Simulation of Transition with a Two-Equation Turbulence Model,’’ AIAA Journal, Vol. 32, No. 2, Feb. 1994, pp. 247-255. 16. Menter, F., R., ‘‘Improved Two-Equation k-ω Turbulence Model for Aerodynamic Flows,’’ NASA TM-103975, Oct. 1992. 17. Giles, M. B., “Nonreflecting Boundary Conditions for Euler Equation Calculations,” AIAA Journal, Vol. 28, No. 12, Dec. 1990, pp. 2050-2058. 18. Tweedt, D. L., and Chima, R. V., ‘‘Rapid Numerical Simulation of Viscous Axisymmetric Flow Fields,’’ submitted for presentation at AIAA 34th Aerospace Sciences Meeting, Jan. 15-18, 1996. 19. Jameson, A., Schmidt, W., and Turkel, E., “Numerical Solutions of the Euler Equations by Finite Volume Methods Using Runge-Kutta Time-Stepping Schemes,” AIAA Paper 81-1259, June 1981. 20. Martinelli, L., and Jameson, A., “Validation of a Multigrid Method for the Reynolds Averaged Equations,” AIAA Paper 88-0414, Jan. 1988. 21. Jameson, A., and Baker, T. J., “Solution of the Euler Equations for Complex Configurations,” AIAA Paper 83-1929, July 1983. 22. Baldwin, B. S., and Barth, T. J., ‘‘A One-Equation Turbulence Transport Model for High Reynolds Number Wall-Bounded Flows,’’ NASA TM 102847, Aug. 1990. 23. Arts, T., Lambert de Rouvroit, M., and Rutherford, A., W., ‘‘Aero-Thermal Investigation of a Highly Loaded Transonic Linear Turbine Guide Vane Cascade,’’ von Karman Institute Technical Note 174, Sep. 1990.

References
1. Baldwin, B. S., and Lomax, H., “Thin-Layer Approximation and Algebraic Model for Separated Turbulent Flows,” AIAA Paper 78-257, Jan. 1978. 2. Chima, R. V., Giel, P. W., and Boyle, R. J., “An Algebraic Turbulence Model for Three-Dimensional Viscous Flows,” NASA TM-105931, Jan. 1993. 3. Suder, K. L., Chima, R. V., Strazisar, A., J., and Roberts, W., B., ‘‘The Effect of Adding Thickness and Roughness to a Transonic Axial Compressor Rotor,’’ ASME Paper 94-GT-339, June 1994. 4. Boyle, R. J., ‘‘Navier-Stokes Analysis of Turbine Blade Heat Transfer,’’ ASME Paper 90-GT-42, June 1990 (also NASA TM 102496). 5. Menter, F., R., ‘‘Performance of Popular Turbulence Models for Attached and Separated Adverse Pressure Gradient Flows,’’ AIAA Journal, Vol. 30, No. 8, Aug. 1992, pp. 2066-2071. 6. Choi, D., and Knight, C. J., “Computations of 3D Viscous Flows in Rotating Turbomachinery Blades,” AIAA Paper 89-0323, Jan. 1989. 7. Hah, C., and Wennerstrom, A. J., “Three-Dimensional Flowfields Inside a Transonic Compressor With Swept Blades,” ASME Paper 90-GT-359, June 1990. 8. Kunz, R. F., and Lakshminarayana, B., “Computation of Supersonic and Low Subsonic Cascade Flows Using an Explicit Navier-Stokes Technique and the k-ε Turbulence Model,” in Computational Fluid Dynamics Symposium on Aeropropulsion, NASA CP-10045, Apr. 1990. 9. Ameri, A., and Arnone, A., “Navier-Stokes Turbine Heat Transfer Predictions Using Two-Equation Turbulence Closure,” AIAA Paper 92-3067, July 1992 (also NASA TM-105817). 10. Ameri, A., and Arnone, A., “Prediction of Turbine Blade Passage Heat Transfer Using a Zero and a Two-Equation Turbulence Model,” ASME Paper 94-GT-122, June 1994.

12