NASA Technical Reports Server (NTRS) 20020078395: Transonic Flow Computations Using Nonlinear Potential Methods

60 Pages·4.8 MB·English

by  NASA Technical Reports Server (NTRS)

g
PROGRESS IN
PERGAMON Progress in Aerospace Sciences 36 (2000) 1-61
www.elsevier.com/locate/paerosci
Transonic flow computations using nonlinear
potential methods
Terry L. Holst *'1
NASA Ame.v Research Center, Mq_ett Fieid, Mail Stop T27B-/, CA 94035-1000, USA
Contents
3
1, Introduction ........................................................
4
2. Nonlinear potential governing equations ........................................
4
2.1. General .........................................................
5
2.2. Full or exact potential equation ...........................................
7
2.3. Boundary conditions and circulation ........................................
8
2.4. Nonconservative form .................................................
8
2.5. Transonic small disturbance (TSD) potential equation .............................
9
2.6. Small-disturbance boundary conditions ......................................
10
2.7. Mathematical properties ...............................................
tl
2.8. Transformation techniques ..............................................
13
2.9. Shock wave capture criteria .............................................
la
2.10. Conservative versus nonconservative forms ...................................
16
2.11. Full potential equation nonuniqueness ......................................
17
3. Algorithms and applications ................................................ 17
3.1 Earl',, TSD potential equation algorithms and applications ...........................
20
3.2. Early full potential equation algorithms and applications ...........................
21
3.3. Conservative full potential equation algorithms ..................................
21
3.3.1. Finite-volume schemes ............................................
23
3.3.2. Artificial density schemes ........................................... 24
3.3.3. Flux upwind schemes ............................................. 24
3.3.4. Entropy and vorticity corrections ......................................
25
3.3.5. Freestream consistency conditions ......................................
26
3.4. Approximate factorization iteration schemes ................................... 28
3.5. Convergence characteristics of SLOR, ADI and AF2 ..............................
3O
3.6. Multigrid iteration schemes .............................................
32
3.7. Other iteration schemes ................................................
32
3.8. Space marching schemes ............................................... 35
3.9. Time accurate schemes ................................................ 39
3.10. Design methods .................................................... 39
3.10.1. Indirect methods ...............................................
39
3.10.2. Inverse methods ................................................
41
3.10.3. Numerical optimization design (gradient methods) ...........................
43
3.10.4. Numerical optimization design (genetic algorithms) ...........................
43
311. Methods developed for complex geometry applications ............................
44
3.11.1. Zonal grid methods .............................................
45
3.t 1.2. Unstructured grid methods .........................................
48
3.11.3. Unstructured Cartesian grid methods ...................................
49
4. Concluding remarks ....................................................
Research Scientist, Computational Physics and Simulation Branch.
*Tel.: [phone]a-6032.
E-mail address' [email protected],nasa.gov (T.L. Hoist).
0376-0421.'00'$-see front matter Published by Elsevier Science Ltd. All rights reserved.
PlI: S03 76-042 11991000 10-X
TL
Hoist ,'Progrexs in Aeroxpace Sciences 36 (2000) I- 61
Nomenclature
Acronyms .',J,)', 2 Cartesian coordinates (physical domain)
AF approximate factorization X2, ':¢_t, "" physical domain transformation metrics
AF2 approximate factorization (scheme 2) intro- Y. angle of attack (°)or iteration scheme accel-
duced by Ballhaus and Steger [1] eration parameter
ADI alternating direction implicit central first-difference operators in the
CAD computer-aided design _,rt.... directions, respectively
CFD computational fluid dynamics 3:,3...... backward first-difference operators in the
CPU central processing unit _, i/.... directions, respectively,
FEM finite element method
forward first-difference operators in the
GA genetic algorithm _, q,... directions, respectively
GM gradient method As change in entropy across a shock wave
MFI OPS nil!lion floating-point operations per sec- 4, full or exact velocity potential
ond (p small-disturbance velocity potential
NSP number of supersonic points 7 ratio of specific heats = 1.4 (for air)
PDE partial differential equation P fluid density,
SIP strongly implicit procedure introduced by F circulation
Stone [2] f2 arbitrary closed control volume
SLOR successive line overrelaxation 0 angular coordinate in polar coordinate sys-
SOR successive overrelaxation tem (rad)
Notation spatial coordinates in computational do-
TSD transonic small disturbance main
a speed of sound computational domain transformation
metrics
A_-A_, transformation metrics
"r time in computational domain
,5" wing chord 172
CL lift coefficient Laplacian operator
CDp'Dr pmraetsesruirael coorefsfiucbiestnatntial derivative ?V gpraardtiiaelm deroipveartaivtoer operator
Lj, k unit vectors along x, y, z
J determinant of independent variable trans- Subscripts
formation Jacobian i,j, k grid indices corresponding to {, _,,_ com-
L lift putational coordinates
M_ freestream Mach number x, y,z indicates differentiation with respect Io
n unit outward normal vector
X, y, 2
P fluid pressure _,q, indicates differentiation with respect t_
stream and stream-normal coordinates
S boundary surface surrounding _ freestream condition
t physical time coordinate
It, F, W velocity components along .\, y, z
q velocity' vector magnitude Superscripts
q velocity vector n iteration index {nth iteration)
t physical time or time-like coordinate * critical condition (condition at M = I)
5. Recommendations for future work ............................................ 5t
Acknowledgements. 52
References ............................................................ 52
Abstract
This presentation describes the state of transonic flow simulation using nonlinear potential methods for external-
aerodynamic applications. The presentation begins with a review ofthe various potential equation forms (with emphasis
on the full potential equation) and includes a discussion of pertinent mathematical characteristics and all derivation
assumptions. Impact of the derivation assumptions on simulation accuracy, especially with respect to shock wave
capture, isdiscussed. Key characteristics ofall numerical algorithm types used for solving nonlinear potential equations,
including steady, unsteady, space marching, and design methods, are described. Both spatial discretization and iteration
scheme characteristics are examined. Numerical results for various aerodynamic applications are included throughout
the presentation to highlight key discussion points. The presentation ends with concluding remarks and recommenda-
tions for future work. Overall, nonlinear potential solvers are efficient, highly developed and routinely used in the
aerodynamic design environment for cruise conditions. Published by Elsevier Science Ltd. All rights reserved.
T.L.Holst ,'Progres_ inAerospace Sciences 36(2000) l-6I
1. Introduction tion is used so much in aerodynamic shape design for
transonic cruise conditions.
In the past three decades the field of computational The main emphasis of this review is to describe numer-
ical solution techniques for soh,ing transonic flow prob-
aerodynamics has evolved from a curious art barely
capable of predicting the inviscid flow over simple two- lems governed by the full potential equation. Because
dimensional shapes, such as airfoils, to the relatively' algorithms for solving the transonic small disturbance
mature current capability, which is capable of predicting (TSD) potential equation are very similar in nature, this
high-Reynolds-number viscous flows about complex topic is covered as well, but in less detail. In a general
three-dimensional shapes, including nearly complete air- sense, this presentation deals with relaxation schemes
suitable for the numerical solution of elliptic partial dif-
craft. The development of this scientific/engineering field
ferential equations. Of course, transonic flow is not pure-
has been paced by many milestones in a variety of areas,
including computational hardware, system software, ly' elliptic in nature, but consists of hyperbolic regions
computing environments, numerical algorithms, con> embedded in an otherwise elliptic domain. However, the
most successful numerical methods of solution for trans-
puter graphics, geometric modeling, flow solver algo-
rithms, etc. indeed, many volumes could be written in onic flow applications, at least for potential formulations,
have evolved from classical relaxation schemes asso-
describing the development of computational aerody-
namics. The purpose of this review is to describe a small, ciated with elliptic equations. Thus, most of the algo-
and yet very important aspect of computational aerody- rithms presented herein will have an elliptic-equation,
namics, that portion associated with nonlinear potential relaxation-algorithm flavor.
formulations. This area is important because numerical The transonic flow regime provides the most efficient
simulations based on nonlinear potential equations pro- aircraft cruise performance; hence, most large commer-
cial aircraft cruise in this regime. However, transonic flow
vide quantitative answers to aerodynamic questions in
a small amount of wall clock time. For aircraft design, fields tend to be sensitive to small perturbations in flow
reducing wall clock time is vitally important because it conditions or to slight changes in geometrical character-
means lower development costs and rapid product avail- istics, either of which can cause large variations in the
ability, which both contribute to larger market share. flow field. Large performance penalties can result be-
Potential equation numerical simulations are com- cause of relatively small perturbations away from desired
design conditions. Computational techniques, therefore,
putationally efficient because they involve the solution of
a simple scalar equation. The more complete formula- have enjoyed an increasing role in helping the aerody-
tions of computational fluid dynamics (CFD), the namics engineer find optimal design conditions, as well
Navier-Stokes and Euler. equations, consist of five as to evaluate design sensitivity. For more information
coupled equations. In addition, numerical iteration on how CFD methods, in general, are being used in (and
schemes for solving the potential equations typically con- have benefited) the aircraft design environment, the inter-
ested reader is referred to Rubbert [4].
verge in fewer iterations than iteration schemes for the
Euler or Navier-Stokes equations. Thus, potential sol- Transonic flow fields contain a variety of interesting
vers are typically an order of magnitude (or more) faster and unique characteristics. Typical airfoil and swept-
than Euler equation solvers on comparable grids E3]. wing flow fields are shown in Figs. 1 and 2. The outer
freestream flow is typically subsonic and elliptic in na-
The price for this extraordinary speed is limitation of
application. All potential formulations are inherently' ir- ture. Regions of supersonic flow usually exist on the
rotational and isentropic. These assumptions are gener- upper airfoil or wing surface and are generally' termin-
ally consistent with subsonic, transonic and supersonic ated by a weak "transonic" shock wave. For the case of
flows at or near cruise conditions providing all shock a swept-wing flow field, the shock wave may actually,
waves are weak. If strong shock waves exist in the flow consist of a system of shocks, as shown in Fig. 2. The
field, i.e., shocks with an upstream, normal-shock Mach first shock is swept and therefore has a supersonic
downstream Mach number. The aft shock is approxim-
number component at or above about 1.3, then the full
potential solution will be in error: the stronger the shock ately, normal to the local flow and therefore has a
subsonic downstream Mach number. Signals tend to
wave, the larger the error. A major ameliorating charac-
teristic of this situation is that for cruise conditions (asso- propagate very rapidly downstream in transonic flow
ciated with the transonic flow regime), the existence of fields where the propagation speed is u + a (local flow
strong shock waves is a very, undesirable characteristic. If speed plus speed of sound) and very slowly, upstream
acandidate configuration has a strong shock wave, anu- where the propagation speed isu - a. For a downstream
merical result does not have to be very accurate to disturbance to propagate upstream it must move around
eliminate it from further consideration in the design the supersonic zone, further increasing the difference
process. Ideally', as the configuration is refined, the shock between the upstream and downstream propagation
strength is reduced and the full potential equation accu- speeds. This situation tends to make transonic
numerical solution techniques, which depend on
racy is improved. This is why the full potential formula-
EL.
HoL;t 'Progress inAerospace Sciences 36 (2000) 1-61
Viscous effects are also extremely important in trans-
onic fows. This complex subject invoh, es four major
\
\ effects: (I) shock/boundary layer interaction effects,
\
(2) the decambering and thickness effects caused by the
addition of a simple displacement thickness, (3) trailing-
edge effects, and (4) near-wake effects. Although a dis-
cussion of viscous correction procedures is not within the
scope of this review, an ample number of references are
presented for those potential solvers that have viscous
correction procedures included.
This review begins with a discussion of the various
j-"
nonlinear potential formulations that have been utilized
'\ in the field of computational transonic aerodynamics
over the past two or three decades (chapter 2). Formula-
J
tion assumptions and limitations, nonconservative ver-
Fig. I. Mach number contours about an airfoil showing atypi- sus conservative forms, shock capturing capabilities and
cal two-dimensional transonic, inviscid flow field compuled us- nonuniqueness issues are discussed in detail in this chap-
ing a full potential algorithm.
ter. Next, in chapter 3, the presentation continues with
a review of past and present research activities involving
algorithm development and aerodynamic applications
with primary emphasis on the full potential formulation.
This chapter includes the milestone achievements that
have shaped the current state of the art in transonic
potential methods. Solution methods reviewed include
classical relaxation algorithms, time-accurate schemes,
supersonic space marching schemes and design methods.
Numerical result examples are included throughout
chapter 3 to highlight important discussion points. The
presentation ends with concluding remarks (chapter 4)
and recommendations for future work (chapter 5j.
Many additional review papers on this and other re-
lated topics are available. Afew of these include Hall [5]
and South [6] where a historical development of the
potential formulation in computational aerodynamics is
presented; Hotst et al. [7] Kordulla [8] and Nixon and
Fig. 2. Mach number contours on the upper surface ofa swept
wing showing a typical three-dimensional transonic, inviscid Kerlick [9] where a variety of transonic potential flow
flow field computed using a full potential algorithm. simulation surveys are presented; and the collected pa-
pers in Nixon [10] Caughey and Hafez [11] Zierep and
Oertel [12] Habashi [13] and Henne [14] where
a wealth of information about the more general topics of
physical-time-dependent algorithms, very stow (relative computational aerodynamics and transonic aerodynam-
to similar algorithms for subsonic or supersonic flow ics are presented. Finally, additional basic information
problemsJ. Such problems are said to be "stiff," and about numerical solution algorithms for nonlinear po-
require larger amounts of computer time. tential formulations is available in Hirsch [15] Anderson
Another characteristic of transonic flow is that it is et al. [16] and Pai and Luo [17"1.
governed by equations that are inherently nonlinear.
Linearization of these equations will remove the vital
flow field physics, which is responsible for the prediction 2. Nonlinear potential governing equations
of shock waves. In contrast, inviscid subsonic flow can be
linearized with good accuracy. The result is Laplace's 2.1. General
equation lot a relative thereof), which can be solved using
a direct method, i.e., a method without iteration. The There are several different potential equation formula-
inherent nonlinear behavior of transonic flow problems tions used in aerodynamic simulations. Although this
means that a direct solution is impossible. Thus, one presentation deals primarily with the full or exact velo-
basic feature associated with all transonic-flow numerical city potential formulation, it is of interest to review all
schemes is that they must be iterative. potential formulations to establish differences and
TL.
HMst / Progress inAerospace Sciences 36 (2000) 1-61
similarities. All potential flow formulations are based on where t is the physical time coordinate, p is the fluid
the ability to define a velocity potential, which requires density, f2 is an arbitrary closed control volume, S is the
an irrotational flow assumption. Thus, for a velocity boundary surface surrounding f2, and n is the unit out-
ward normal vector to the surface S. Eq. (3) states that
vector field defined by q, the requirement for a velocity
the time rate of change of the mass in an arbitrary fixed
potential to exist is
volume (first term) isbalanced by the net outflow of mass
Vx q = O. (1) leaving the same volume (second term). In order for Eq.
(3) to represent a closed-form description of a flow field,
If this condition holds in all locations of the flow field of
an algebraic expression for the density' in terms of the
interest, then a full or exact velocity potential function velocity potential must be utilized. Several forms of this
_pexists and is defined by expression will be discussed shortly.
Eq. (3) can be expressed in differential form by trans-
V4)= q. (2)
forming the surface integral term into a volume integral
Velocity components can be expressed in terms of partial using Gauss' Divergence Theorem:
derivatives of the velocity potential function. For Car- ffro
tesian coordinates this is given as n.p V4)dS = I7.p V4)dQ.
VO = q = ui + vj + wk = _b_i+ ¢,,j + 4):k
Using the fact that the control volume f2 is fixed with
or component by component as respect to time, the differentiation and integration asso-
ciated with the first term of Eq. (3) can be interchanged.
4)_,=_,, 4),=_', 4):=_,,,
Combining the two volume integrals into the same term,
where i, j and k are the standard unit vectors in the x, Eq. (3) becomes
y and z directions, respectively; and u, v and w are the fIIo( )
Cartesian velocity vector components, also along the x,
--_ + V.p V4) &Q=0.
y and z directions, respectively. In the above expressions,
the quantity 4)_(for example) is used to indicate apartial
Since the control volume is arbitrary, the integrand in
derivative of 4) with respect to the spatial coordinate x.
the above equation vanishes everywhere, which results in
The velocity potential function has a spatial variation
the desired differential form of the unsteady full potential
which is independent of path. It can be defined for incom-
equation.
pressible or compressible flows that are either steady or
unsteady. It is restricted, however, because of the irrota-
_p
tional-flow assumption, to flows without viscous effects, -- + V'pV4) =0. 141
&
i.e., potential flows are inherently inviscid in nature. Of
course, viscous boundary layer corrections can be in- The above integral and differential forms of the full
cluded quite easily (at least for attached flows) by solving potential equation still need an additional relation to
a potential formulation in conjunction with the bound- complete the formulation. In particular, a relation that
ary layer equations, but this is beyond the scope of this expresses the fluid density p as a function of the velocity
presentation and will not be discussed further. Other components 4):,, 4)_.and 4':, is required. A suitable ap-
types of rotational flow corrections can be included, e.g., proach for this derivation starts with the inviscid mo-
flows with circulation and/or vortices, but this requires
mentum equation given by (this is actually a form of the
additional (often empirical) modeling. More on circula- Euler momentum equation)
tion modeling will be presented in the section on bound-
ary conditions. Next, the discussion turns to the various
gq __VP = O,
governing equation forms that utilize a velocity potential _ + q. vq -7
function.
where p is the fluid pressure. The second term in the
2.2. Full oi exact potential equation above equation can be reduced to a convenient form
using the Lagrange acceleration formula [18] given by
The most general form of the full potential equation is
derived from the mass continuity equation using the
Dt -& +V -qx(Vxqj,
definition of the velocity potential given above l-Eq. (2)].
This equation, written in integral form, is given by
where qin the second term on the right-hand side repres-
!-fo,,-;
ents the magnitude of the velocity vector. The D/Dr
=- dr2 + n' PV4)dS = 0, (3) notation used on the left-hand side stands for the
EL. Hoist ProgressinAerosl)aceScience_36(2000)1-61
material or substantial derivative defined by speed of sound definition. Thus, the final form of the
integrated unsteady Bernoulli equation isgiven by
Dq dq
D--t= ,at + q"P'q
a_4_ - + =c(o. (7)
8t ,-1
For irrotational fluids the following equation is easily
obtained
The ultimate goal of this derivation is obtaining a rela-
tionship between the fluid density and derivatives of the
velocity potential, i.e., the fluid's velocity components,
thus allowing closure ofthe fullpotential equations given
Using this relation the momentum equation becomes above by Eqs. (3) or (4). This is accomplished using
Eq. (7), the speed of sound definition, and the isentropic
density-pressure relation [Eq. (6)], yielding the finaI de-
_'-7+ r" + _' sired relation
Substituting Vdofor q and using the fact that
p = 1+ .,__._1(M_,_- 2do,- cb_- do_- ¢_')
vi dp _ Vp
(8)
for a barotropic fluid: yields In this equation the density p and the velocity compo-
nents ¢_,, Cy and 4): are nondimensionalized by the
v =-+ + =0. freestream values of the density p_ and the freestream
, ,! p value of the speed of sound a,, respectively. The Car-
The integration of this expression along an arbitrary line tesian coordinates x, y and z and the time t are non-
in the flow domain yields dimensionalized by a characteristic length, usually the
airfoil or wing chord c and the quantity a:,/c, respect-
ively.
?-?-qt5+ _ + q-72 = Cit), (5) Several other forms of the full potential equation with
different types of nondimensionalization have been used
where C(0 is a constant of integration that m general is
for numerical computations. For example, the steady
a function of time but not space. Eq. (5)is the unsteady
flow version of the differential form of the full potential
Bernoulli equation. It and various related forms of
equation in which the density is nondimensionalized by
the Bernoulli equation are used in many areas of fluid
the stagnation density p_,_ and all velocity components
dynamics.
are nondimensionalized by the critical speed of sound
The integral in the unsteady Bernoulli equation must
a* isgiven by
be evaluated before this equation can be used further.
This is accomplished using (pdo_)_+ (pdo,.)),+ (pdo-)_-= 0. (9)
P
-- = const, (6)
p7
where ;' isthe ratio of specific heats (equal to 1.40for air)
Alllength scales are stillnondimensionalized by the same
and "const" is a constant that can be evaluated when
characteristic length. Use ofthe above nondimensionaliz-
anondimensionalization ischosen. Eq. (6)isthe standard ation creates the following useful conditions. At stagna-
densit)-pressure relationship foran isentropic flow. With tion points
this relation the integral in Eq. (5)can easily be evaluated
yielding p=l, do.,-= doy= do:= 0,
and at sonic lines
CdP_ ;-' P a:
. ;,-lp 7-1"
where a is the fluid speed of sound. The last equality in
the above equation is obtained using the perfect gas
=0.633938145... (7= 1.40).
In addition, either of the two nondimensionatizations
2A barotropic fluid isone in which the density can be ex-
given above can be used to evaluate the constant in the
pressed solely as a function of pressure, i.e.. p =O(P).For
example,afluidundergoing anisentropic process isabarotropic isentropic density-pressure relation [Eq. (6)] or the con-
fluid. stant in the steady Bernoulli equation, which is Eq. (7)
EL. Hoist .'Progress in Aerospace Sciences 36 L2000) 1-61
with the time terms removed. For example, using the stated, is given by
second nondimensionalization from above (p_,,_ and a*),
x2+y:+z2--, _, q_--+qS_,
the following results are obtained:
Isentropic density-pressure relation where 4'_ is the freestream distribution of the velocity
potential, usually uniform flow. The latter two boundary
p _7+1 (11) conditions, symmetry planes and geometric surfaces, are
rOY O_, both treated in the same manner, i.e., with a flow tan-
gency assumption given by
Stead), Bernoulli's equation:
q-n=0,
q2 a2 t 7 + 1
+ (12)
2 7-I 2",:-1 where n is a unit vector normal to the geometry of
interest. More on flow tangency boundary conditions is
Another widely used potential equation form is obtained
presented in Section 2.8 where transformation techniques
from Eq. (9) by making the additional assumption of for the full potential equation-are discussed.
incompressible flow. This yields the familiar Laplace's For aerodynamic applications to be useful the numer-
equation ical formulation must be able to predict aerodynamic
loads, e.g., lift. The Kutta-Joukowski theorem says
v% =4,_.,+ 4,,,,.+ _:: = 0. (13
L = rO<q_F,
Unlike the nonlinear full potential equation given by
Eqs. (3), (4) or (9) Laplace's equation is linear. Although where L is the lift, p:_ is the freestream fluid density, q:_ is
this further limits the flow field physics that can be the freestream fluid velocity magnitude, and F is the
simulated, e.g., shocks cannot be captured with Laplace's circulation. The circulation around (for example) an air-
equation, the linear nature lends itself to the powerful foil is mathematically defined as
method of superposition. A numerical approach that
utilizes superposition to solve Laplace's equation is typi-
call)' called a panel method. A key panel method charac-
F = qDq' dl,
teristic is that only the geometric surface of interest need
t
be discretized. In contrast, all numerical methods used to
solve nonlinear governing equations require a field or where I is arts, closed path surrounding the airfoil for
volume discretization, a feature that greatly increases the which the velocity vector field is defined. Using Stokes'
numerical algorithm complication. In a panel method, Theorem it can be seen that circulation is inherently tied
flow tangency along the aerodynamic surface is obtained to vorticity, that is
by solving for source, vortex or doublet distribution
strengths for the surface's discretized elements. This op-
eration requires and is computationally paced by the F=_qdl= ffsVXq'ndS'.
inversion of a large matrix whose rank is equal to the
number of boundary elements. Thus, a key aspect of ans,
where S is the surface constructed such that its boundary
panel rnethod implementation is the judicious selection
is Iand nis the unit outward normal vector to S. Thus, it
(both in number and placement) of an appropriate sur-
face element discretization. Panel method details are can be seen from the above equation that an irrotational
beyond the scope of this presentation and will not be velocity vector field, such as that predicted by the futl
discussed further. The interested reader is referred to potential equation, is not capable of supporting circula-
Anderson et al. 1-16] for basic information on panel tion, i.e., there is no lift.
This situation can be corrected by adding a linear
methods; Smith 1-191 for a historical presentation of
potential vortex solution to the nonlinear potential that
panel method development; and Hess 1,,201Hoeijmakers
surrounds the airfoil. This is accomplished by modifying
1,21] and Roggero and Larguier [22] for information on
current panel method applications. the freestream boundary condition as follows
4_ob= 4',. + _1_,',
2.3. Boundary conditions and circulation
where qSobis the new outer boundary condition, 0_ is the
To complete the full potential governing equation spe- usual uniform-flow velocity potential solution, and qS,is
cification, boundary conditions are required along all the newly added potential vortex solution gixen by
boundaries. Specifically, these boundaries fall into three
F
categories: freestream, symmetry planes, and geometric ¢b.... 0.
surfaces. The freestream boundary condition, simply 27:
"I_L.Hoist/Progress in Aerospace Sciences 36 (2000) 1-61
In the above equation, F/2n is a constant representing This is the steady', nonconservative full potential equa-
the vortex strength and 0 is the usual angular coordinate tion written in three-dimensional Cartesian coordinates.
associated with a traditional polar coordinate system The velocity components 05x, 4_, and 42=have been re-
centered inside the airfoil. The vorticity associated with placed by u, v and w, respectively.
the resulting potential solution iszero everywhere except The unsteady version of the nonconservative full po-
at the center of the vortex where it is infinite. The circula- tential equation can be derived in a similar fashion and is
tion is a constant (equal to F) for all integration given by
curves/that include the airfoil and zero for all integration
cur\es that do not include the airfoil. Because of the (a2-- u2)05xx + (a2 -- F2)05,,r-.b(a2 -- n,2)05z:
periodic nature of O,the 42vfunction is double-valued at -- 2uv42_y -- 2uw42= -- 2vw¢_.:
0 = 0,taking on values of 0and 2n. In other words, along
some "cut" in tile velocity potential solution, usually = 42,, + 2u42._,+ 2v42y, + 2w42=,. (15)
emanating from the airfoil trailing edge to downstream
More on the characteristics of the full potential equation,
infinity along the airfoil wake, the velocity potential
both conservative and nonconservative forms, is pre-
"jumps" from its 0 = 0 value to its 0 = 2rc value. The sented in Section 2.10.
magnitude of this velocity potential jump (easily derived
by looking at the definition of 05,.)is equal to F. In effect,
2.5. D'ansonic small disturbance (TSD) potential equation
forcing the 05, vortex strength to be equal to the airfoil
trailing edge velocity potential jump, is like a Kutta
Another potential equation formulation used for com-
condition that forces the airfoil-surface upper and lower
puting transonic flows about aircraft is the transonic
pressures to match.
small-disturbance (TSDI potential equation. Many of the
numerical algorithm breakthroughs realized in solving
2.4. Noncon.vercative _rm
the full potential equation were first developed using the
simpler TSD potential formulation. The TSD potential
Eqs. (4) and (9) are two forms of the full potential
equation is derived from the full potential equation by
equation which are commonly used in numerical applica-
first defining a small-disturbance velocity potential O
tions, especially when shock waves are expected to be
captured in the solution. These versions of the full poten- V_o= q - q_,
tial equation are written in the so-called conservative
where q is the usual local velocity" vector and q_. is the
form, which is characterized by having all variables inside
freestream velocity vector, which is assumed to be alig-
the outer-most differentiation. The steady full potential
ned with the x direction. It is defined by
equation has also been solved in nonconservative form for
subsonic and transonic applications. For subsonic ap- q_ = u_i.
plications the nonconservative and conservative solutions,
With the above definitions the small-disturbance velo-
assuming equivalentl3 small levels of numerical error, are
city components of Vcp are given by
virtually identical. For transonic applications involving
captured shock waves the nonconservative and conserva-
(Px = U -- U_ (Oy _ F, Cpz = W.
tive solutions are different. At shock waves the nonconser-
vative approach produces an error in the form of amass Derivation of the TSD potential equation begins by
substituting these small-disturbance velocity compo-
source that causes an error in shock position and strength.
The consequence of this error is discussed in detail in nents into the full potential equation [either Eq. (14) for
Section 2.10. The nonconservative form of the full poten- steady flows or Eq. (15) for unsteady flows]. Then after
tial equation is presented here because of its historical neglecting small terms according to the small-distur-
importance and because it provides a useful framework bance assumptions given by
for analyzing the conservative full potential equation.
The nonconservative form of the full potential equa-
tion is derived from the conservative form by using the
chain rule to expand derivatives. Expressions for the
density derivatives are obtained from the density expres-
the three-dimensional unsteady TSD potential equation
sion [Eq. (10)], the speed of sound definition, and the
becomes
isentropic density-pressure relation [Eq. (11)]. Substitu-
tion of these derivatives into the expanded full potential
I1- M_ -- M2 (;,+ 1}cPx-],p.,., -: (p,., + (p::
equation yields
L [Iz
(a2 - u-' cb_ + la: - v:142_.,+ (a: - w2to:= 1
= _-(o, + 2u_cPxt). (16)
- 2uvqS_.r - 2uw05._: - 2vw42>- = 0. {14) a_
EL Hoist /Progress inAerospace Sciences 36 (200(111-61
One term containing a small-disturbance quantity' (the bance term. the unsteady small-disturbance equation
first term) survives the small-disturbance analysis. This is given by
because for transonic flow
1
(l - M2 )_p,.,.+ _p.,. +_p:_ = c-_7(__p,, + 2u_ _p._,) (20)
(O-
1 _ Mk' _...(M?2 + I)-22-,
U_ is obtained. This equation is linear and valid for either
subsonic or supersonic flow, but not transonic flow. As in
Eq. (16) is valid for subsonic, supersonic and transonic
the TSD potential equation case. a steady version of this
flows that satisfy the original full potential equation
equation is obtained by neglecting all time terms, which
assumptions (inviscid, isentropic and irrotational flow)
and that are a "small disturbance" away from freestream. yields
An essential ingredient of Eq. (16), and all transonic (1- M_)_p,._+ <, + _p=:=0. (21)
governing equations, is that they are nonlinear, i.e., the
O:,(P_,._term is nonlinear. This term is required to predict This is the famous Prandtl-Glauert equation and can be
shock waves, which are inherently a nonlinear phe- used to describe steady, small-disturbance potential flow.
nomena. Linearization of a transonic flow governing Lastly, if the flow is assumed to be incompressible
equation removes the essential mathematics required to (i.e., a_ --, _, M_ _ 0) the small-disturbance potential
predict shock waves. equation, valid for either steady or unsteady flow, be-
Another variation of Eq. (16) that has been used in comes
many' applications (see, e.g., [1,2311) is the low-frequency
q)x_ + _P_.y+ _P:: = 0, (22l
TSD equation. Ifthe frequency' of oscillation of the prob-
lem under consideration is small enough, then the which is another version of Laplace's equation. This
(p, term can be neglected yielding version is based on the small-disturbance potential func-
tion and differs from Eq. (13) in the boundary conditions
2tl_
1-_v,"; - ,¥I5(: + .,,j_(Ojox7,, + r&-r+ O= = _ (&,,. that are applied, both in the freestream and at the air-
foil/wing surface.
(17)
2.6. Small-disturbance boundam. _conditions
A convenient nondimensional measure of a problem's
frequency is obtained using the reduced frequency' de-
The freestream boundary' condition consistent with
fined by k = o)c/u_, where ca is the problem's physical
small-disturbance theory is that all disturbances must
oscillation frequency measured in cycles per second.
vanish in the freestream, i.e.,
When the reduced frequency of a specific unsteady prob-
lem is less than about 0.2 (and all other TSD assumptions ¢p._.:,o= ¢p_._ = _0:,_ = O.
are valid), then the low-frequency TSD equation is usu-
ally considered to be valid. The flow tangency boundary condition for a typical
Of course, the steady TSD equation is obtained by "thin" airfoil or wing used in conjunction with any of the
neglecting all time derivatives in Eq. (16), which yields small disturbance formulations presented in the last sec-
tion is generally derived as follows: The standard tan-
gency condition can be implemented at the airfoil or wing
1-:'v/_- M_:(;'+ 1) q_xT_p.,.,. + cp,._.+ _.-_-=0. (18 surface using
/'/_ _]
Eqs. (16)-(18t are in nonconservative form, i.e., not all ¢p ± cp_ dg ±(x, y)
variable coefficients are inside the outer differentiation.
To transform these equations into conservative form is
where x is aligned with the freestream direction, y is in
easily accomplished and yields [e.g., for Eq. (18)]
the span direction and - is in the vertical direction. The
subscript "ws" indicates that the boundary condition is
(1-._4_)_p_-M_e,, + L_2T,.j._+ <.,.+ _P.-=_o. (19) applied at the wing surface. The functions _l-(x, y) and
g-(x, y) define the upper and lower wing surfaces, respec-
gqs. (16)-(19) represent classical forms of the TSD equa- tively. The small disturbance version of this boundary
tion. but several other forms exist, including forms de- condition is obtained by making two simplifications.
rived to better predict transonic flows on swept wings. First, q0_is neglected relative to u_. in the middle term
For more information on other TSD equation forms see denominator. Second, the flow tangency boundary con-
van der Vooren et al. [24] or Slooff [25]. dition is applied at the airfoil slit, i.e.. at z = 0, instead of
Additional useful equations can be obtained from the airfoil surface. This latter approximation greatly sim-
Eq. (16!. For example, by' neglecting the last small-distur- plifies the volume grid generation process for TSD
10
TL. Holst/ Progress in Aerospace Sciences 36(2000) 1-61
potential applications because geometrical surfaces do istics are real and distinct; i.e., if the discriminant of Eq.
not have to be fitted with grids. This is a key reason that (25) is greater than zero (B2- 4AC > 0), then the equa-
this formulation was so widely used for three-dimen- tion is hyperbolic; parabolic ifthe characteristics are rea;
sional transonic flow applications in the early years of and coincidental (B2- 4AC = 0); and elliptic ifthe char-
CFD development. The final small-disturbance flow-tan- acteristics are complex and distinct (B" - 4AC < 0).
gency boundary condition with these simplifications be- By using the discriminant test described above, it can
comes be shown that the TSD equation given by Eq. (18) (two-
dimensional version) is hyperbolic when
dy_-!x, y)
(a-(x, y. 0+-)= u_--
dx
_>
u_ (;,+ l)M_'
This expression approximates the required flow tangency
boundary' condition at the airfoil surface to an accuracy and elliptic when
consistent with small-disturbance theory. A key flaw in
small-disturbance theory is displayed in this boundary
<:
condition at (for example) an airfoil leading edge where u.. (v+ 1)MU
the slope of the surface becomes infinite, and accurate
boundary condition implementation is impossible. This In other words, the sign of the first term coefficient
difficulty is a symptom of the breakdown in small-distur- determines the equation type. If the coefficient is positive,
bance theory at stagnation points. The streamwise the local flow is subsonic; ifit is negative, the local flow is
velocity component perturbation becomes large and is supersonic. The nonlinearity of the first term is essential
actually equal to the freestream velocity at the stagnation for describing the mixed character of transonic flow and
point. This fundamental limitation in the TSD potential is the mechanism by which shock waves are formed.
equation approach is the primary reason that its use has The characteristic directions associated with the TSD
declined in recent years. potential equation are given by
dv = + 1-M_-M_b+ll
2.7. Mathematical properties
1,2
The primary motivation for studying the nature of
Notice that these characteristic slopes are symmetric
partial differential equations (PDEs) in the present con-
about the x-axis regardless of the local velocity vector
text is to gain insight into the physics they describe and orientation; i.e., the characteristics are not a function of
to develop guidelines for the implementation of numer-
the y component of the velocity %.. A sketch of the
ical solution procedures. Different equation types gener-
steady, two-dimensional TSD potential equation charac-
ally require different solution algorithms. With this
teristics for a typical supersonic point is presented in
purpose in mind, consider the following general quasi-
Fig. 3.This situation, which isin dramatic contrast to the
linear, second-order PDE:
full potential or Euler formulations, where the character-
istics are symmetric about the local stream direction, has
Au_.,. 4- Bu,._ + Cu;_ = F. (23)
certain simpli_'ing implications regarding spatial discret-
where u is an arbitrary dependent variable and A, B, ization approximations for the TSD potential equation.
C and F are (at most) functions ofx, y, u, ux and u_..This In particular, it is much easier to construct a spatial
equation can be studied and classified by considering discretization scheme for the TSD potential equation
the corresponding characteristic equation given by (for with a numerical domain of dependence that is compat-
a derivation of the characteristic equation and additional ible with the PDE's mathematical domain of dependence.
discussion on this topic see Ames [26] and Mitchell 1-27]) As seen in Fig. 3, as v grows relative to u, the x-axis
symmetry condition for the characteristics becomes more
"dy'.t-".-, / ., =o. and more nonphysical. This is a direct result of the
(24)
small-disturbance assumption that requires v to remain
small in order to keep the TSD potential equation valid
Using the quadratic formula, the two characteristic direc- Classification of the stead 3, full potential equation
tions associated with Eq. (23) are given by given by Eqs. (9) and (10) is difficult because this set of
equations is not in the standard form given by Eq. (23).
(dr" I + B ___- \ B"_- 4AC (25) However, the nonconservative full potential equation
\_.x I , 2.4 [Eq. (14)] is ideally suited for this purpose. This equation
written in two dimensions is given by
The nature of these characteristics determines the equa-
tion classification. Eq. I23) is hyperbolic ifthe character- (a 2 -- uZ)_bxx - 2uv4)._, + (az - v:)4)_,_.= O. (26)