See discussions, stats, and author profiles for this publication at: https://www.researchgate.

net/publication/260310402

Improving the CFD Predictions of Airfoils in Stall

Conference Paper · January 2005

DOI: 10.2514/6.2005-1227

CITATIONS READS
15 721

4 authors, including:

Marilyn J. Smith
Georgia Institute of Technology
199 PUBLICATIONS   1,483 CITATIONS   

SEE PROFILE

Some of the authors of this publication are also working on these related projects:

Nonlinear Transients and Gusts View project

Dynamic stall in helicopter rotor environment View project

All content following this page was uploaded by Marilyn J. Smith on 10 May 2016.

The user has requested enhancement of the downloaded file.

 Improving the CFD Predictions
of Airfoils in Stall

Andrew Shelton∗, Jennifer Abras†, Robert Jurenko‡, and Marilyn J. Smith§

Georgia Institute of Technology, Atlanta, GA 30332-050

The inconsistent ability of Reynolds-Averaged Navier-Stokes (RANS) simulations to

predict the aerodynamic characteristics and flow field physics of airfoils and wings at and
near stalled conditions is well known. At stall, various turbulence models can give widely
varying results, while at lower angles of attack the results are essentially invariant. The
Spalart-Allmaras turbulence model, in conjunction with Michel’s Criterion to determine
transition, provides the most consistent correlations with experimental data. However, the
addition of a transition model does not guarantee good post-stall comparison with
experimental data. Some preliminary modifications to the Spalart-Allmaras turbulence
model appear to improve post-stall correlation with experimental data and warrant further
investigation.

Nomenclature
cb1 = Spalart-Allmaras turbulence model constant, 0.135
cb2 = Spalart-Allmaras turbulence model constant, 0.622
Cd = section drag coefficient
Cl = section lift coefficient
Cm = section moment coefficient
Cp = pressure coefficient
ct1 = Spalart-Allmaras turbulence model trip constant, 1.0
ct2 = Spalart-Allmaras turbulence model trip constant, 2.0
ct3 = Spalart-Allmaras turbulence model constant, 1.2
ct4 = Spalart-Allmaras turbulence model constant, 0.5
cv1 = Spalart-Allmaras turbulence model constant, 7.1
cw1 = Spalart-Allmaras turbulence model constant, cb1 κ 2 + (1 + cb 2 ) σ
cw2 = Spalart-Allmaras turbulence model constant, 0.3
cw3 = Spalart-Allmaras turbulence model constant, 2
d = distance to the closest wall
f = function
~
g,r, S = intermediate variables
H32 = shape factor
r̂ = surface roughness
Rex = Reynolds number based on local surface length
Reθ = momentum thickness-based Reynolds number
ue = boundary layer edge velocity
γ

= intermittency factor
κ = Kármán constant, 0.41
λθ = pressure gradient parameter
ν = kinematic viscosity

∗
Graduate Research Assistant, School of Aerospace Engineering, M/S 0150, AIAA Member.
†
Graduate Research Assistant, School of Aerospace Engineering, M/S 0150, AIAA Student Member.
‡
Undergraduate Research Assistant, School of Aerospace Engineering, M/S 0150, AIAA Student Member.
§
Associate Professor, School of Aerospace Engineering, M/S 0150, AIAA Associate Fellow.

1
American Institute of Aeronautics and Astronautics
νt = eddy viscosity
σ = turbulent Prandtl number, 2/3
χ = ν~ / ν

subscripts
m = momentum
tr = transition

I. Introduction

T he efficient prediction of aircraft wing and helicopter rotor performance, vibratory loads, and aeroelastic
properties is a major concern to the aerospace community. The prediction of these characteristics is only as
accurate as the weakest component of an overall analysis that comprises aerodynamics, structural mechanics, and
dynamics. The majority of these simulations rely heavily on the use of Reynolds-Averaged Navier-Stokes (RANS)
Computational Fluid Dynamics (CFD) for performance analysis, as well as more efficient simplified versions of the
aerodynamics equations for many of the aeroelastic predictions. These simplified aerodynamics modules often
utilize look-up tables to provide two-dimensional aerodynamic characteristics, which are then corrected by a number
of theoretical and empirical factors for sweep, unsteady aerodynamics, finite wing tip effects, etc. Typically these
tables are comprised of a combination of wind tunnel data, empirical data and numerical analyses. The potential to
rely more heavily on CFD for these tables, as well as on RANS simulations for the complex three-dimensional
simulations can be realized with the advent of faster computers and more sophisticated physical models.
The issue of what is the best turbulence model and how accurate are CFD predictions are not new questions.
There exists a number of papers in the literature where researchers have applied a number of RANS models to both
static and unsteady configurations with experimental data, for example, Srinivasan, Ekaterinaris and McCroskey1,
Davidson2, Ekaterinaris3, Szydlowski and Costes4, and Rumsey and Gatski5. Many RANS codes do not have
transition models, relying on fully turbulent or user-defined transition locations, though recent research in the field
has shown a greater interest in transition modeling. Here, for example, He et al6, Kim, Zaman and Panda7, Lee and
Chen8, Dorsey9, and Xu and Sankar10 have studied current transition on a number of different applications.
This research began as an effort to evaluate the accuracy of airload predictions currently utilized in fixed-wing
and rotorcraft applications11. From this study, it was observed that for the linear aerodynamic regime of an airfoil,
the choice of turbulence model was not critical, as long as the grid was sufficient to resolve the suction peak. When
non-linear effects became important – as in transonic flows and the stall region – the choice of the turbulence model,
along with the fidelity of the grid, resulted in widely varying results. As a side note, while the lift and drag
predictions remained within the experimental error limits, the combination of lower lift and higher drag resulted in
poor L/D predictions. From these results, the next step has been an evaluation of current RANS CFD methods to
predict the two-dimensional characteristics of two additional airfoils with varying levels of experimental data. The
ability of several different CFD codes to predict the lift, drag and pitching moment of these airfoils has been
examined and extensive correlations of these simulations with the results of wind tunnel tests have been performed.
Next, an evaluation of the impact of transition models on the different turbulence models, particularly in the stall
regime, has been completed. Finally, some modifications to the Spalart-Allmaras turbulence model were evaluated.
These simulations occurred over a range of Mach numbers and Reynolds numbers to make this information useful
for full-scale aircraft and rotorcraft applications, as well as the emerging smaller Unmanned Aerial Vehicle (UAV)
applications.

II. Test Cases

A. Symmetric Airfoil at Moderate Reynolds Number

The NACA0015 has been widely tested with respect to applications in rotorcraft. The airfoil has been
extensively tested for both static and dynamic stall by Piziali12. Both tripped and untripped data is available for a
Mach number of approximately 0.3 and a Reynolds number of about 2 million. The data includes lift, drag and
moment data, as well as pressure coefficient data at select angles of attack. Piziali’s static data was obtained by
moving the airfoil at very slow changes in pitch, hence there are both “up” and “down” data points.

2
American Institute of Aeronautics and Astronautics
B. Cambered Airfoil at Moderate to High Reynolds Number
The SC1095 is the primary airfoil used in constructing the UH60A Blackhawk helicopter main rotor. It has been
tested in multiple wind tunnels, the data of which has been collected and evaluated by Bousman13. This data
includes composite performance characteristics, although most individual test data are not available. The Mach
numbers range from 0.2 to 1.0 and Reynolds numbers of 1 to 10 million. The data is primarily untripped.

C. Symmetric Airfoil at Low Mach Numbers and Reynolds Numbers

The NACA0012 has been tested by Maresca14 for Reynolds numbers ranging from 10,000 to 100,000 and free
stream velocities of 5 – 10 m/s. These data are well-suited to evaluate applications for small unmanned vehicles.
Data available include forces and moments, as well as boundary layer velocity profiles and pressure coefficient
distributions. These data are from untripped configurations.

III. Description of the RANS Methodologies

A. CFL3D
CFL3D, developed at NASA-Langley, solves the RANS equations on structured grid topologies for both two-
and three-dimensional geometries15. Spatial discretization is an upwind-biased scheme for the advection terms and a
central scheme for the diffusion terms. Temporal discretization is implicit with subiterations. Numerous turbulence
models are available, including Spalart-Allmaras, Wilcox k-ω, and Menter SST models, and an option for user-
defined transition is included.

B. DSS2
DSS2 is an in-house Georgia Tech RANS code for two-dimensional applications, developed by Professor
Lakshmi Sankar and his students. It has been utilized for a variety of applications, including dynamic stall of
rotorcraft and wind turbines10. It utilizes a fourth-order spatial central difference scheme coupled with a second-
order time integration. It also includes a user-defined transition for use with the Baldwin-Lomax and Spalart-
Allmaras turbulence models.

IV. Transition Modeling

In combination with the turbulence models available for the methodologies in Section III, many transition onset
models are available. For this study, simple models employed by a wide range of RANS codes were evaluated.
These include Michel’s criterion16, Eppler’s criterion17, and Abu-Ghannam and Shaw’s criterion18.
One of the more widely used methods is Michel’s Criterion. This method correlates the transition Reynolds
number to the momentum thickness Reynolds number, and is of the form given by

 22400  0.46
Re ,tr = 1.1741 +  (1)
 Re  Rex,tr
θ

 x,
tr 

This method relies only on the boundary layer thickness to calculate transition onset. It is based on data gathered
using a flat plate with a small pressure gradient. Comparison with experimental data shows that Michel’s Criterion
has good agreement for small angles of attack, but for higher angles of attack the accuracy of this criterion starts to
degrade19.
Eppler’s method not only relates the momentum thickness, but also relates the energy to the momentum
shape factor at transition using the relation

u θ 
log e  > 18.4 H 32 − 21.74 − 0.34rˆ (2)
 υ 

The Abu-Ghannam and Shaw relation provides an expression to correlate free stream turbulence and pressure
gradients with the momentum thickness Reynolds number at the start of transition. This relation can be used even if
there is no pressure gradient or free stream turbulence. The expression is

3
American Institute of Aeronautics and Astronautics
  ) 
λ

F(
Re ,tr = 163 + exp F ( 
θ

)−
λ

(3)
 6.91 
θ

where from experimental data for adverse pressure gradients,

2
F (λ ) = 6.91 + 12.75λ + 63.64λ λ <0 (4)
θ θ θ θ

or from experimental data for favorable pressure gradients,

2
F (λθ ) = 6.91 + 2.48λθ + 12.27λθ λθ > 0 (5)

The Chen and Thyson modification20 is provided as an intermittency model in order to avoid abrupt transition. It
is governed by

 u3 
 e Re − 1.34 
( )
3
γ = 1 − exp − x,tr  (6)
 213 log Re x,tr − 4.7323 ν 2 
 

V. Results

A. Grid Studies
Before beginning the computational evaluation, an extensive grid study was performed to determine an optimal
grid with which to proceed. The NACA0015 configuration was used for this grid study. Several different grid
parameters were varied for simulations in the post-stall regime. There are several recommendations for grids that
are to be utilized for simulations in the stall/post-stall region.
The suction peak is very sensitive to the grid, in particular the normal spacing next to the surface. Normal
surface spacing at the leading edge should be greatly refined, compared with pre-stall simulation requirements. For
example, a simple refinement from a surface normal of 1.0 X 10-6, which provided excellent linear range correlation,
to 1.0 X 10-7 resulted in a pressure coefficient drop at the suction peak of +1.0 in the direction towards experimental
data. At 1.0 X 10-7, additional resolution of the normal spacing did not significantly change the peak further.
Trailing edge vortex shedding or separation can also be greatly enhanced with streamwise grid refinement near
the trailing edge and in the near wake region. For 2nd order spatially-accurate schemes, an additional 50 to 60 points
in the last 10% of the airfoil and first 10% of the wake resulted in the capture of wake physics – in particular the
shed wake vorticity. It was interesting to note that the Splarat-Allmaras model wake results showed very little
difference between the RANS and DES formulations, but the Menter SST DES model did show an improved wake
correlation. Both simulations used an infinite wing assumption.
The leading edge streamwise grid was also refined where leading edge separation was expected, however, the
separation and subsequent leading edge vortex shedding were not captured. This was true for both conventional
turbulence modeling (Spalart-Allmaras, k-w, Menter SST) and two DES models (Spalart-Allmaras DES, and Menter
SST DES). This implies that transition modeling will be needed for airfoils that experience leading edge separation.
While these grid variations were studied independently, their combinatory effect was also verified. Assuming
that the grids include minimized skewness, the aspect ratio of the grid cells at the leading edge also plays an
important role. Grids with aspect ratios of greater than 200 had poorer correlation compared with their lower aspect
ratio counterparts.
Using these criteria – it was found that a grid of size 545 X 201 was within 1 – 3 % of the force and moment
values of a grid double that size. This was primarily because the normal spacing characteristics of the grids at the
leading edge were very similar. While all of the leading edge suction peaks were over predicted, the force and
moment predictions were within 5% of each other.

4
American Institute of Aeronautics and Astronautics
 2.0
Experimental Clmax B. The Impact of Different Turbulence Models
Range (Ref. 13)
1.5
Using two different codes, the differences of fully
turbulent simulations was evaluated for the three
1.0
airfoils. For clarity, only partial results are shown so
that the figures remain readable. Differences between
Lift Coefficient

0.5
the code performances will be noted. Consider first the
SC1095 airfoil. Different fully turbulent simulations
0.0
were performed with a number of different grids and
-10 -5 0 5 10 15 20 codes in Ref. 11, and simulations with the same grid
-0.5 and code show similar trends, as seen in Fig. 1.
Experimental data from a number of independent wind
CFL3D Spalart Allmaras
-1.0 CFL3D Wilcox K-omega tunnel tests for this airfoil, as documented in Ref. 13 are
CFL3D Menter's SST shown, where available. As in Ref. 11, the lift curve
-1.5 slopes and zero lift angles of attack for the different
Angle of Attack, deg turbulence models are almost identical and fit well
a) Lift Coefficient within experimental limits. The fully turbulent runs
0.35
1.5
0.30

0.25
Drag Coefficient

1.0
0.20
Lift Coefficient

0.15 0.5

0.10
0.0
0.05
-10 -5 0 5 10 15
0.00 -0.5 CFL3D Spalart Allmaras
-10 -5 0 5 10 15 20 CFL3D Wilcox K-omega
Angle of Attack, deg
CFL3D Menter's SST
b) Drag Coefficient -1.0
0.04 Angle of Attack, deg

0.02
Figure 2. Comparison of different turbulence
model simulations for the SC1095 airfoil at a Mach
0.00 number of 0.7.
Moment Coefficient

-10 -5 0 5 10 15 20
-0.02 show a Clmax that is higher than and typically occurs
later (at higher angles of attack) than experimental data.
-0.04 The break in the drag and moment coefficients is
-0.06
likewise delayed. Similar scatter in the CFD results
occur at Mach 0.7 (Fig. 2), though no reliable
-0.08 experimental data was available at that Mach number.
-0.10
Similar results are seen for the NACA0015 airfoil in
Fig. 3, when compared with Piziali’s experimental
-0.12 Angle of Attack, deg data12. The maximum lift coefficient is over predicted
c) Moment Coefficient by all of the turbulence models, and the angle of attack
Figure 1. Comparison of different turbulence at which the maximum lift coefficient occurs is
model simulations using CFL3D for the SC1095 typically about 2 to 3 degrees above the experimental
airfoil at a Mach number of 0.4. data. Similarly, the drag and moment coefficient breaks
are also seen to occur at higher angles of attack than
expected.
Finally, the low Reynolds number NACA0012 is examined. It should be recognized that at these low Reynolds
number, the fully turbulent assumptions are least appropriate for all of these test cases. From Tables 1 and 2, it is
seen that at 9 degrees angle of attack – just before stall – the computational lift correlates very well with the
experimental lift, but for post-stall conditions at 14 degrees angle of attack, the computational lift is over predicted.

5
American Institute of Aeronautics and Astronautics
 2 Piziali (Up) Table 1. Comparison of Turbulence Models:
Piziali (Down) NACA0012 Airfoil at 9 Degrees Angle of Attack
1.8 CFL3D Spalart Allmaras
1.6 CFL3D Menter's SST Cl Cd Cm
CFL3D Wilcox K-omega
1.4 Experiment 0.95 0.0786 -0.0247
Lift Coefficient

DSS2 Spalart Allmaras

1.2 DSS2 Baldwin Lomax Laminar 0.93 0.0305 0.0140
1 Spalart-Allmaras 0.94 0.0312 0.0114
0.8 k-ω 0.94 0.0324 0.0096
0.6 Menter SST 0.93 0.0330 0.0114
0.4
SA Doubled Grid 0.90 0.0822 -0.0022
0.2
0
0 5 10 15 20
Angle of Attack, Deg Table 2. Comparison of Turbulence Models:
a) Lift Coefficient NACA0012 Airfoil at 14 Degrees Angle of Attack
Cl Cd Cm
0.25
Experiment 1.10 0.2949 -0.159
0.2 Laminar 1.40 0.0828 0.011
Drag Coefficient

Spalart-Allmaras 1.38 0.0736 0.014

0.15
k-ω 1.36 0.0720 0.015
0.1 Menter SST 1.38 0.0771 0.016

0.05 Comparatively, using the same grid with a fourth-

order scheme, the drag is predicted much better,
0 indicating that this test case, a finer grid is necessary
0 5 10 15 20 for second-order schemes.
Angle of Attack, Deg
b) Drag Coefficient
0.06 C. Evaluation of Transition Models
Given the over predictions of the fully-turbulent
0.04 simulations, it is apparent that transition may play a
Moment Coefficient

0.02 critical role in determining the airfoils’ stall

characteristics. Three of the more popular algebraic
0 models have been evaluated: Michel, Eppler, and Abu-
0 5 10 15 20
-0.02 Ghannam and Shaw. The impact of the Chen-Thyson
intermittency modification is also evaluated in
-0.04 conjunction with the Michel Criterion. Again, all
-0.06 three airfoils are evaluated to provide consistency
across geometry, Reynolds number, and Mach number
-0.08
Angle of Attack, Deg
combinations.
Figure 4 illustrates a typical set of force and
c) Moment Coefficient moment predictions using the aforementioned
Figure 3. Comparison of different turbulence transition models. Compared to the data from fully
model simulations for the NACA0015 airfoil at a turbulent simulations, it is seen that the only transition
Mach number of 0.291. model that significantly improves the simulation is the
Michel Criterion model. It has moved the maximum
lift coefficient prediction to the extreme corner of the experimental data. Similarly, the breaks for the drag and
moment data occur much earlier, more in keeping with the experimental results.
The NACA0015 results show similar trends, as seen in Fig. 5. The use of Michel’s Criterion improves the
predictions, though it does not correct all of the problems. Comparison is given with both the tripped and untripped
experimental data, which show some differences in the post-stall region.

6
American Institute of Aeronautics and Astronautics
 2 1.6 Piziali Tripped (Up)
Experimental Clmax Range (Ref. 13) Piziali Tripped (Down)
1.5 1.4 Piziali Untripped (Up)
Piziali Untripped (Dn)
1.2 Fully Turbulent
1

Lift Coefficient
User-Defined @ 0.05c
Lift Coefficient

1.0 Michel's Criterion

0.5 Eppler
0.8
0
0.6
-10 -5 0 5 10 15 20
-0.5 Fully Turbulent 0.4
Michel
-1 Michel with Chen Thyson 0.2
Eppler
-1.5 AbuGhannam and Shaw
0.0
Angle of Attack, deg 0 5 10 15 20
Angle of Attack, Deg
a) Lift Coefficient
a)Lift Coefficient
0.4
0.25
0.35
0.3 0.20
Drag Coefficient

Drag Coefficient
0.25
0.15
0.2
0.15 0.10
0.1
0.05 0.05

0
0.00
-10 -5 0 5 10 15 20
Angle of Attack, deg 0 5 10 15 20
Angle of Attack, Deg
b) Drag Coefficient
b) Drag Coefficient
0.06
0.06
0.04
Moment Coefficient

0.04
0.02
0.02
Moment Coefficient

0
-10 -5 -0.02 0 5 10 15 20 0
-0.04 0 5 10 15 20
-0.02
-0.06
-0.04
-0.08
-0.06
-0.1
-0.12 -0.08
-0.14 -0.1
Angle of Attack, deg Angle of Attack, Deg

c) Moment Coefficient c) Moment Coefficient

Figure 4. Comparison of different transition model Figure 5. Comparison of different transition model
simulations for the SC1095 airfoil at a Mach simulations for the NACA0015 airfoil at a Mach
number of 0.4. number of 0.291.

Comparisons of the pre- and post-stall angles of attack for the NACA0012 (Tables 3 and 4) show a slightly
different result than the previous two airfoil results. Pre-stall, there is no change in the lift predictions. For the 14
degree post-stall case, the implementation of Michel’s Criterion does not help to improve the correlations. Using the
4th order option within the DSS2 code, however, indicates that the higher order scheme does a much better job at
predicting the lift, getting within 2% of the value obtained by integrating the pressures. As with the fully turbulent
case, a fine (doubled) grid for the 2nd order simulation changes the drag error from an under prediction of 61% to an
over prediction of 4%. The moment also moves in the appropriate direction. These results again signify the need
for refined grids or higher-order methods for this airfoil test case.

7
American Institute of Aeronautics and Astronautics
 Table 3. Comparison of Transition Models: Table 4. Comparison of Transition Models:
NACA0012 Airfoil at 9 Degrees Angle of Attack NACA0012 Airfoil at 14 Degrees Angle of Attack
Cl Cd Cm Cl Cd Cm
Experiment 0.95 0.0786 -0.0247 Experiment 1.10 0.2948 -0.1567
Laminar 0.93 0.0305 0.0140 Laminar 1.40 0.0828 0.0112
Fully Turbulent 0.93 0.0312 0.0114 Fully Turbulent 1.38 0.0736 0.0138
Michel’s Criterion 0.94 0.0312 0.0114 Michel’sCriterion 1.38 0.0739 0.0139
Michel’s Criterion 1.38 0.07412 0.0140
Eppler 0.94 0.0312 0.0113 with Chen-Thyson
Abu-Ghannam and 0.94 0.0310 0.0112 Eppler 1.38 0.0736 0.0138
Shaw Abu-Ghannam and 1.38 0.0738 0.0138
Michel’s Criterion 0.89 0.0819 0.0003 Shaw
Doubled Grid 4th order – 1.08 0.1075 0.0089
Michel’s Criterion

D. Turbulence Model Modifications

With the addition of transition models, the prediction of the basic forces and moments for the three airfoils were
improved, but they did not completely address all of the problems. The common error for all of the simulations
appears to be an over prediction of the stall, including both the magnitude of the lift and the angle of attack where
the break first occurs. Figure 6 illustrates the turbulent kinetic energy contained within the boundary layer at the
leading edge for a high angle of attack of 16 degrees. Fully-turbulent and computed transition via Michel’s
Criterion simulations are compared with a completely laminar computation. It can be seen that the energy contained
in the boundary layers are significantly higher at the leading edge, so that while the laminar boundary layer
simulation has separation, the other two simulations
remain fully attached. Comparisons between
computational and experimental velocity profiles for
Fully Turbulent the NACA0012 (Fig. 7) indicate similar disparities in
the boundary layer energy and character.
Thus, it appears that, while the turbulence models
do well in predictions at low angles of attack, they
Tripped tend to break down at higher angles of attack. This is
not surprising as the turbulence models were primarily
developed and correlated using flat plate and other
data with less dramatic pressure gradients. Further
Michel Criterion consideration of the turbulence models show that a
number of constants are typically utilized within each
Figure 6. Comparison of turbulent kinetic energy model, and these constants, usually derived from
near leading edge of NACA0015 airfoil at 16 degrees experimental data, are based on lower pressure
angle of attack gradients. It is hypothesized that at higher angles of
attack, the constants utilized within the turbulence
models may be causing the simulations to generate boundary layers with too much energy in them.
To investigate this, the Spalart-Allmaras turbulence model was chosen for trial exploration. This turbulence
model was chosen because of all of the turbulence models tested, it tended to provide forces and moments closest to
the experimental data, particularly in conjunction with the Michel Criterion. Consider the Spalart-Allmaras
turbulence model formulation, as per Ref.21 where the eddy viscosity is defined as

χ3
ν t = ν~f v1 , f v1 = (7)
χ + cv31
3

and the transport equation is

8
American Institute of Aeronautics and Astronautics
 Dν~
[ ]  ν~ 
2

= c b1 [1 − f t 2 ]S ν~ + ∇ • ((ν + ν~ )∇ν~ ) + c b 2 (∇ν~ ) − c w1 f w − b21 f t 2    + f t1 ∆U 2
~ 1 2 c
(8)
Dt σ  κ  d 

where

1/ 6
χ  1+ c 6 
~ ν~ f w = g  6 w63  ν~
S =S+
κ 2d 2
f v2
,
f v2 =1−
1 + χf v1
,  g + c w3  ,
r= ~ 2 2
(
S κ d , g = r + c w2 r − r
6
) (9)

Because it appears that too much turbulent energy is being created, the first two terms of the equation associated
with the near-wall turbulence for finite Reynolds numbers were evaluated. These terms affect the eddy viscosity
where there are steep gradients in the flow, i.e., close to the wall. The value of cv1 regulates the intercept of the
viscous sublayer with the log law of the wall. A lower value will cause a lower intercept. The impact of changing
this variable is evaluated. As seen in Figs. 8 and 9, the reduction of the cv1 term results in an overall shift of the
forces and moments towards the experimental values. Notice that the higher-order scheme and fully turbulent
results show more improvement than the lower-order and computed transition results.
A second parameter,σ, in the viscous
0.05 diffusion term was also examined for its
Experiment
O(4), Spalart Allmaras impact on the boundary layer
Distance Normal to Airfoil

O(2), Laminar predictions. This parameter was chosen

0.04 O(2), Spalart-Allmaras
O(2), Wilcox k-omega because, per Ref. 21, “we pay little
0.03
O(2), Menter's SST attention to a factor σ multiplying it [the
diffusion]” because “It [the diffusion
term] vanishes in the ideal near-wall
0.02 solution”. This is no longer true for
large pressure gradient flows, and this
0.01 term may aid to identify the important
budget terms in these stalled regions. It
0 is observed in Fig. 8 that this parameter
also plays a role in the modifying the
-0.01 0 0.01 0.02 0.03 0.04
Streamwise (u) Velocity
boundary layer for the stalled
conditions, though it appears to have
Figure 7. Comparison of streamwise (u) velocity profile at x/c=0.3 less of an effect than the cv1 term. The
for the NACA0012 airfoil at 14 degrees. combination of the best solutions from
both terms does not result additional
improvements, however, indicating that further study is needed to determine optimum constant values.
The impact of the cv1 term on the flow field is readily apparent in Figs. 10 and 11, which illustrate the pressure
distribution and the upper surface streamlines, respectively. The streamlines indicate that a lower cv1 term results in
the trailing edge separation moving forward between 0.15c and 0.20c in the stalled region. The pressure distribution
is particularly promising – the suction peak reduces significantly with the reduction in the cv1 term and the pressure
distribution matches the experimental data much better over the upper forward half of the airfoil. Comparisons of
the results of this grid with the modified constant are comparable to results obtained by the authors and Szydlowski
and Costes4 with significantly refined grids.

VI. Conclusions
The prediction of forces and moments over airfoils in stall and post-stall conditions has been studied with
different turbulence models, transition models, and computational methods. Refined grids are required to capture
some of more complex flow field features found at higher angles of attack, including vortex shedding and
separation. The Spalart-Allmaras turbulence model performed the best consistently across the airfoils, in particular
with the algebraic transition model, Michel’s Criterion. It was seen that the boundary layer contained too much
energy at these angles of attack, causing the simulations to miss the maximum lift coefficient magnitude and
location. Preliminary investigations of the Spalart-Allmaras turbulence model are encouraging, as the modification

9
American Institute of Aeronautics and Astronautics
 Sigma
0 0.5 1 1.5 2 1.6
1.6 1.4
1.4 1.2

Lift Coefficient
Lift Coefficient

1.2 1
1 0.8
0.8 4th-order
2nd-order 0.6
0.6 2nd-order, sigma=1/2
4th order, sigma=1/2 0.4 4th-order, Michel's Criterion
0.4 4th-order, varying sigma Piziali (Ref. 13), Untripped Upper
Piziali (Ref. 13), Tripped Upper 0.2
0.2 Piziali (Ref. 13), Tripped Lower
Piziali (Ref. 13), Untripped Lower
0 0
0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8
Cv1 Cv1

a) Lift Coefficient a) Lift Coefficient

Sigma 0.14
0 0.5 1 1.5 2
0.14 0.12

Drag Coefficient
0.12 0.1
Drag Coefficient

0.1 0.08
0.08
0.06
0.06
0.04
0.04
0.02 0.02

0 0
0 1 2 3 4 5 6 7 8 0 1 2 3 4 5 6 7 8
Cv1
Cv1
b) Drag Coefficient b) Drag Coefficient
Sigma 0.06
0 0.5 1 1.5 2
0.04 0.04
Moment Coefficient

0.03 0.02
0.02 0
Moment Coefficient

0.01 -0.02
0 -0.04
-0.01 -0.06
-0.02 -0.08
-0.03 -0.1
-0.04 0 1 2 3 4 5 6 7 8
Cv1
0 1 2 3 4 Cv1 5 6 7 8
c) Moment Coefficient
c) Moment Coefficient
Figure 9. Variation in NACA0015 forces and
Figure 8. Variation in NACA0015 forces and
moments with Spalart-Allmaras turbulence model
moments with Spalart-Allmaras turbulence model
modifications for a fully turbulent simulation.
modifications for a fully turbulent simulation.

of the production and diffusion empirical constants showed that results comparable to finer grids can be obtained
numerically with the moderate grid size. These results will be significant if consistent improvements for moderate
three-dimensional grids can be verified.
Future work in this area is continuing, studying additional terms in the Spalart-Allmaras model and expanding
the correlations with more test cases and with large eddy simulation computational results. Additionally, a method
to determine how to best modify the constants in the turbulence model as angle of attack is changed is necessary.

10
American Institute of Aeronautics and Astronautics
 -10

-8 Experiment
Pressure Coefficient

Cv1=3.5
Sigma=1/3
-6 Original

-4

-2

2
0 0.2 0.4 0.6 0.8 1
x/c
Figure 10. Pressure coefficients for the NACA0015 airfoil at 16 degrees angle of attack.

Figure 11. Comparison of the flow separation areas for the upper surface between original Spalart-Allmaras
model and the modified Spalart-Allmaras model (cv1=3.5). Upper is orginal, lower is modified.

11
American Institute of Aeronautics and Astronautics
 Acknowledgments
This work is sponsored by the National Rotorcraft Technology Center (NTRC) at the Georgia Institute of
Technology. Dr. Yung Yu is the technical monitor of this center. Computational support for the NTRC was provided
through the DoD High Performance Computing Centers at ERDC and NAVO through an HPC grant from the US
Army, S/AAA Dr. Roger Strawn. The computer resources of the Department of Defense Major Shared Resource
Centers (MSRC) are gratefully acknowledged.

References
1
Srinivasan, G.R., Ekaterinaris, J.A., and McCroskey, W.J., “Evaluation Of Turbulence Models For Unsteady Flows Of An
Oscillating Airfoil,” Computers & Fluids, Vol. 24, No. 7, Sept. 1995, pg 833-861.
2
Davidson, L., and Rizzi, A., “Navier-Stokes Stall Predictions Using An Algebraic Reynolds-Stress Model,” Journal of
Spacecraft and Rockets, vol. 29, no. 6, Nov-Dec, 1992, pg. 794-800.
3
Ekaterinaris, J.A., and Menter, F.R., “Computation Of Oscillating Airfoil Flows With One- And Two-Equation Turbulence
Models,” AIAA Journal, vol. 32, no. 12, Dec. 1994, pg. 2359-2365.
4
Szydlowski, J., and Costes, M., “Simulation of Flow Around a NACA0015 Airfoil for Static and Dynamic Stall
Configurations using RANS and DES,” presented at the 4th Decennial Specialist’s Conference on Aeromechanics, San Francisco,
CA, January, 2004.
5
Rumsey, C.L., and Gatski, T.B., “Recent Turbulence Model Advances Applied To Multielement Airfoil Computations,”
Journal of Aircraft, Vol. 38, No. 5, September/October, 2001, pg. 904-910.
6
He, X., Senocak, I., Shyy, W., Gangadharan, S. N., and Thakur, S., “Evaluation of laminar-turbulent transition and
equilibrium near wall turbulence models, Part A: Applications”, Numerical Heat Transfer; , Vol. 37, No. 2, Feb. 2000, pg. 101-
112.
7
Kim, S.-W., Zaman, K.B.M.Q.; Panda, J., “Numerical Investigation Of Unsteady Transitional Flow Over Oscillating
Airfoil,” Journal of Fluids Engineering, Vol. 117, No. 1, March 1995, pg 10-16.
8
Lee, S. C., and Chen, C. R., “Navier-Stokes Computations Of Transition To Turbulent Flow Around Airfoils,” SAE
Technical Paper Series, 1990, 901808.
9
Dorsey, D. J., “A Comparative Study of Four Algebraic Transition Models,” SAE Technical Paper Series, 1994, 941142.
10
Xu, G., and Sankar, L. N., “Effects of Transition, Turbulence and Yaw on the Performance of Horizontal Axis Wind
Turbines,” AIAA 2000-0048, 2000.
11
Smith, M. J., Wong, T.-C., Potsdam, M., Baeder, J. and Phansee, S., " Evaluation of CFD to Determine Two-Dimensional
Airfoil Characteristics for Rotorcraft Applications." Presented at the AHS 60th Forum, Baltimore, MD, June 2004.
12
Piziali, R., A., "2-D and 3-D Oscillating Wing Aerodynamics for a Range of Angles of Attack Including Stall." NASA TM
4632, September 1994.
13
Bousman, William G., "Aerodynamic Characteristics of SC1095 and SC1094 R8 Airfoils," NASA TP-2003-212265, also
AFDD TR-04-003, December 2003.
14
Maresca, C., Berton, E, Favier, D., and Benyahia, A., “Study of 2D and 3D Boundary Layter of Moving Walls by
Embedded LDV Measurements,” Final Technical Report, European Research Office of the U.S. Army, March 2002.
15
Rumsey, C., Biedron, R., and Thomas, J., CFL3D: Its History and Some Recent Applications, NASA TM-112861, May
1997; presented at the Godunov's Method for Gas Dynamics Symposium, Ann Arbor, MI, May 1-2, 1997.
16
Michel, R., “Etude de la Transition et Calcul de le Trainée des Profils d’aile en Incompressible, ” ONERA Publ. 58, 1951.
17
Eppler, R., Airfoil Design and Data, New York, NY, Springer-Verlag, 1990, Pg. 562.
18
Abu-Ghannam, B. J., and R. Shaw, “Natural Transition of Boundary Layers-The Effects of Turbulence, Pressure Gradient,
and Flow History,” Journal of Mechanical Engineering Science, Vol. 22, No. 5, 1980, Pgs. 213-228.
19
Johansen, J., “Unsteady Airfoil Flows with Application to Aeroelastic Stability,” RisØ National Laboratory, RisØ-R-
1116(EN), Oct 1999.
20
Chen, K. K., and Thyson, N. A., “Extension of Emmons’ Spot Theory to Flow on Blunt Bodies,” AIAA Journal, Vol. 9,
1971, Pg. 821-825.
21 Spalart, P. R., and Allmaras, S. R., “A One-Equation Turbulence Model for Aerodynamic Flows,” La Recherche
Aerospatiale, No. 1, 1994, pg. 5 – 21.

12
American Institute of Aeronautics and Astronautics

View publication stats