- Email: [email protected]
Physics of forced unsteady flow for a NACA 0015 airfoil undergoing constant-rate pitch-up motion∗
Fluid Dynamics North-Holland
Research
10 (1992) 351-369
Physics of forced unsteady flow for a NACA 0015 airfoil undergoing constant-rate pitch-up motion * K.N. Ghia a, J. Yang a, G.A. Oswald
a and U. Ghia b
* Department of Aerospace Engineering and Engineering Mechanics, Computational Fluid Dynamics Research Laborator)), Unicersity of Cincinnati, Cincinnati, ON 45221, USA h Department of Mechanical, Industrial und Nuclear ~ng~neerjng~ ~omputatjo~lal Fluid Dynamics Laboratory, Unicersity of Cincinnati, Cincinnati, OH 45221, USA Received
6 May 1992
Abstract. The unsteady Navier-Stokes (NS) analysis of Gsswald, Ghia and Ghia in velocity-vorticity variables is modified to study the dynamic stall phenomenon for a N’ACA 0015 airfoil undergoing constant R,, pitch-up maneuvers at Reynolds number Re = 10000 and 45000. The use of third-order accurate biased upwind differenciny for the nonlinear convective terms in the vorticity transport equation removes the spurious oscillations observed in the earlier studies by the authors for these values of Re. The fully implicit and vectorized ADI-BGE method of the authors is used to solve the unsteady NS equations. Instantaneous inertial surface vorticity, which is an invariant of the choice of reference frame selected, is employed to determine the location of separation of the boundary-layer flow on the suction surface; also a separation bubble embedded within the boundary layer is observed for both cases somewhere between the leading edge and the quarter-chord point. Primary, _ secondarv. tertiatv . and auarternarv vortices have been observed before . the dynamic-stall vortex evolves and gathers its maximum strength.
1. Introduction
Specifically, the large amplitude rapid pitching motion associated with the initiation of a high angle-of-attack maneuver typically leads to the generation of a dynamic stall vortex whose evolution results in large transient lift, drag and moment that can, for short periods of time, produce loadings significantly larger than those expected during either steady or quasi-steady flight. Indeed, the successful completion of an abrupt, drastic maneuver can depend upon the ability of holding the dynamic stall vortices in place, at least for the duration of the maneuver, and subsequently bleeding the excess accumulated vorticity in a controlled manner into the wake. Abrupt shedding of large amounts of locally concentrated vorticity can so rapidly alter the lift distribution on a body that a tumbling loss-of-control incident can occur, as the associated rapid changes in moment distribution cannot be tolerated. Carr (1988) has comprehensively reviewed the literature on the dynamic stall phenomenon and has also articulated the effect of key parameters on this phenomenon. Helin (1989) has Correspondence to: K.N. Ghia, Department of Aerospace Engineering and Engineering Mechanics, Cincinnati, 759E Baldwin, MW70, Cincinnati, OH 45221-0070, USA. * This research is supported, in part, by AFOSR Grant No. 90-0249, and by the Ohio Supercomputer No. PESO80. 0169-5983/92/$05.75 All rights reserved
0 1992 - The Japan
Society
of Fluid Mechanics
/ Elsevier
Science
Publishers
University Center
B.V
of
Grant
352
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
also highlighted recent advances on the subject, while stressing the importance of unsteady aerodynamics for highly maneuverable and agile aircraft. In addition, he has raised the important issue of the effect of flow separation on the formation of the energetic dynamic stall vortex. These two reviews adequately point out some of the unresolved issues associated with the problem of dynamic stall. In light of the review article by Carr (19881, no attempt has been made in the present paper to review some of the noteworthy contributions made to this subject in the past. Ghia et al. (1992) recently presented some preliminary results of their unsteady NS analysis using the vorticity stream function formulation. A brief discussion of their analysis is given in the next two sections. Their results for low Reynolds number, Re = 1000, for a NACA 0015 airfoil undergoing’a constant 0, pitch-up motion, showed that the eruption of the fluid, in the form of a secondary vortex from the unsteady boundary layer, rotating in counterclockwise direction near the quarter-chord point (QCP) triggers the formation of the dynamic stall vortex. Subsequently, Osswald et al. (1990) presented results for Re = 10000, a test configuration for which Visbal and Shang (1989) have presented detailed numerical results. At the rapid pitch rate of 0, = 0.2 selected, the central-difference scheme (CDS) with a grid of (566 x 101) showed spurious oscillations in the vorticity field at 13> 32”. Past experience suggests that this type of oscillation can be eliminated using further clustering and/or increasing the grid size in critical areas. In addition, the recent success achieved by Blodgett et al. (1990) in the use of third-order biased upwind differencing scheme (UDS) to obtain results for unsteady separated flow past an elliptic cylinder at higher Re motivated the use of UDS in the present study. Hence, to eliminate the spurious oscillations in the vorticity field, the nonlinear convective terms are approximated using a biased third-order UDS. This then leads to the scope of the present study. The objectives of the present study are as follows: (1) To extend the earlier work of the present investigators to higher Re flows. (2) To understand and quantify the role of unsteady separation in the evolution of the dynamic stall vortex. (3) To determine the physical mechanism involved in moment stall, with the ultimate goal of revealing aerodynamic mechanisms to delay or mediate its effect.
2. Governing
equations
and boundary
conditions
Osswald et al. (1990) modified their earlier analyses (Ghia et al., 1985; Osswald et al., 1985) for fixed body flows to now permit arbitrary six degree-of-freedom maneuvers. The general 3D analysis of Osswald et al. (1990) is in terms of velocity and vorticity variables. This analysis is used here in its 2D form and is briefly presented next; for additional details, the reader is advised to refer to earlier work of the authors. An arbitrary maneuver can be completely defined by specifying the trajectory rn,r(t) of some point B fixed on the body (fig. 1) with respect to an inertial observer, together with the specification of the instantaneous angular velocity n,(t) of the body. Kinematically, the translational velocity of the origin B of the body-fixed frame is then Vn,r(t) = drn,,/dt, the acceleration of the translational acceleration is an,,(t) = d2rn,r/ dt2, while the angular body-fixed reference frame is aa = da,/dt. In the present analysis, these functions, which define a specific maneuver, are assumed to be explicitly prescribed known functions of time and will be given in a later section. With the apparent position vector r being given as
(1)
K.N. Ghia et al. / Forced unsteady flow on NACA 001.7 airfoil
353
Body Fixed Observer (x*,x*)
iii
e
. .....w ..’
- - *..
Body Fixed General
(5’,5*)
Inertial Obierver (x,y) Fig. 1. Inertial and body-fixed coordinate
systems for arbitrary maneuvering body. f&(t) = - R&G
= - t&(r)i,
a,(t)= a&M = -B&M.
the apparent velocity V is kinematically related to the inertial velocity V, as I/= V, - V,,,(t)
Xr.
-O,(t)
(2)
For the 2D simulations of interest in the present study, it is computationally more efficient to solve for a d~~z~rb~~cestream function rather than for inertial velocity directly, where the disturbance stream function @,P is defined as the deviation from uniform flow. Thus, the inertial stream function $r can be written as $1 =Y
+ @o(f)
+ ILP(S’,
f2,
q,
(3)
where y is the vertical inertial coordinate and I&, is the constant of integration which represents the instantaneous displacement of the zero inertial streamline at infinity owing to the inviscid generation of lift (inviscid circulation generated) by the body. The inertial velocity then becomes
where i*is the unit base vector of the inertial Cartesian reference frame and e, and e2 are the covariant base vectors of the generalized-coordinate (
g = (g,,g22-g,,g,,),
where
g,,=
axk axk
c : Y kc, a’$’ a<’
(5)
and (x1, x2> represent the body-fixed Cartesian coordinate system shown in fig. 1, whose unit base vectors are e^, and e^,. For the maneuvering airfoil shown in fig. 1, the unsteady Navier-Stokes equations in generalized orthogonal non-inertial body-fixed coordinates are as shown below where Re = &c/u, with U, being the free-stream reference velocity, c the chord of the airfoil being the reference length, and Y the kinematic viscosity.
354
K. N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil Elliptic
stream
by the following
function problem. The disturbance elliptic equation
stream
function
distribution
is governed
(6) subject
to the boundary $p = 0
conditions
at infinity,
(7)
-tf43w[wf along the body surface. and J/,1(t) This elliptic
+ (x212])
Here,
e,(t)
= 5,&>4
problem
-x’[
- cos 0,(t)]
*P = (x2[ G,,(t)
+
V,‘,,( t) - sin t9,( t)] (8)
+ h(f)
is the instantaneous
pitch angle of the airfoil,
see fig. 1,
(9)
VBZ,lWZ.
for the disturbance
stream
function
is coupled
with the following:
Temporally parabolic uorticity transport problem.
(10) subject
to the boundary w1 = 0
conditions
(11)
at infinity,
(12) along the body, together
with the constraint
(g,,/&)a#:‘/ac*= -([ cos +[sin
e,(t)
O,(t)
along the body. In eq. (lo), the contravariant &v’
= (gl,/G)P’( +[sin
O,(t)
dw = (622/G)-'( +[sin Finally, non-inertial
e,(t)
[ cos O,(t) - Vi,,(t)
[cos
e,(t)
- Vi,I(t)
-
VA,,< t) SX2R,( t)]&
-Vi,,(t)
-x’~,(t)]e^,)
components -l&,(t)
of the apparent
(13) velocity
are
+x2n,(t)]$, *e, +a*P/ag2,
-x’R,(t)]e,) -L&(t)
*e,
+x”l&(t)]e^,
-.x’f2B(t)]e^2)
*e2-a$Clp/a[*.
(15)
the pressure is computed by solving the Neumann-Poisson problem. In the reference frame, the linear momentum equation in terms of disturbance velocity,
K.N. Ghia et al. / Forced unsteady flow on NACA 001.5 airfoil
355
is given as V.(aVD/at+v.VvD+a,(t)
xYD+Re-’
VXW,)=
-V+Vp
(16)
with the disturbance velocity being given as VD = (e,/&)
YfDM2
- (eJ&)
(17)
a$D/%‘,
and subject to the Neumann boundary conditions given as Vp*e,=
-(WD/&fV~VVD+~,(t)
xVD+Re-’
VXw,)*e,,
(18)
VP-~,=
-(WD/at
XVD+Re-’
VXo,)*e,.
(19)
+‘I/.VVD+R,(t)
All terms in this problem are well behaved at infinity. It remains only to prescribe the maneuver and the formuiation
is complete.
3. Pitch-up maneuver and grid generation Visbal and Shang (1989) have investigated a rapid pitch-up maneuver for a NACA 0015 airfoil using a compressible Navier-Stokes code based on the solution algorithm of Beam and Warming (1978). Their maneuver is defined as %/i(t) e,(t)
= 0, [ 1 - exp( - (46,‘t,)t]
= L’,{t - (rJ4.6)
1.
Visbal and Shang (1989) investigated several cases, with 0, ranging from 0.2 to ranging from 0.5 to 1.0. Pitch-axis locations of leading edge, trailing edge, and point were investigated. In that study, one set of parameters used were L?, = 0.2, t, the airfoil pitched about the axis through the l/4-chord point. Specifically, the currently studied is such that rn,i(t)
=O,
n,(t)
= -a,[1
vu,,(t)
= 0,
au(t)
= -(4.6/t0)f10
a,,,(t)
= 0,
R, = 0.2,
(20) 0.6 and t, l/4-chord = 0.5 with maneuver
-exp(-(4.6/t,)t],
t, = 0.5.
exp(-(4.6/t,)t), (21)
The governing equations (6), (10) and (161, with the boundary conditions given in eqs. (7), (81, @l&(15) and (181, (191, are to be solved numerically using a suitable orthogonal grid. A generalized Schwarz-Christoffel grid generation procedure (Osswald et al., 1989) is empioyed to map a C-grid topology of infinite extent to a unit square computational domain. A portion of a typical (444 x 101) clustered orthogonal body-fixed grid distribution is shown in fig. 2 for a NACA 0015 airfoil with a rounded trailing edge. Next, brief comments are made on the numerical method used to solve the discretized problem.
4. Numerical method for the discrete problem In the earlier studies (Ghia et al., 1992; Osswald et al., 1990) by the authors, all spatial derivatives were approximated using central differences. In the present study, the convective derivatives in the vorticity transport equation (10) are appro~mated using a third-order biased UDS, with all other spatial derivatives still being approximated using CDS. Following
356
Y
0.50
-
0.40
-
0.30
-
0. to
-
0.10
-
0.00
-
-a
10
-0.20 -0.30
-
-0.40
-
-0.50
-
f
’ -0.50
0. 00
0.
so
1. 00
X
Fig. 2. Grid distribution
for a NACA
0015 airfoil with (444~
loll
points.
Blodgett et al. f1990), the vorticity transport equation is solved using a vectorized form of the alternating direction imphcit CADI) procedure of Douglas and Gunn ~19~4~. This procedure has been shown by Osswald et ai. (1987) to extend to three-Dimensions easily for the velocity-vort~c~ty formuiatio~ and, when cast into delta form, it is computatjona~ly shghtly more efficient. As described by Osswald et al. f1987), second-order accurate time extrapolations are employed to predict the boundary vorticity at the new time level. These predicted boundary vorticities are combined with the mid time-level stream functions to produce a spatially and temporally second-order accurate solution for the interior vorticity field at the new time level. The boundary condition for the vorticity equation is solved using a second-order accurate expression, since it leads to a dissipative derivative as the leading truncation error term and is therefore consistent with the biased third-order upwind differencing scheme used in this study. The elliptic stream function equation is solved by a direct block Gaussian eIimination (BGEI technique. This ADI-BGE fully implicit numerical simulation technique of the authors has been previousfy described by Osswaid et al. 09851 for a generafized coordinate inertial reference frame. Owing to the ~~~~~-~~~~~~~~~~ of the velocity-vorticity forrmdation, only minor modifications were required in the program to implement the arbitrary maneuvering problem. The Neumann-Poisson problem is solved using a direct solver for which the methodology used is described by Blodgett (1990). Next, the results obtained in this study are discussed.
5. Resuits
and discussion
A test configuration was seiected with Re = 10000, constant pitch-up rate 0, = 0.2 and the pitch axis located at l/4 chord. Osswald et al. (1990) had obtained results for this configura-
tion with a C-grid of (566 x 101) mesh points using CDS. Their detailed results were compared with those obtained by the present method wherein the convective terms are differenced using third-order biased UDS. The results are very similar for 0” < 6 5 17”, beyond which the earlier CDS scheme results show a slight time lag for 17” < 8 I 32”. Beyond 0 = 32”, no comparison is attempted since CDS had spurious oscillations. The medium size C-grid of (444 x 101) points place 204 points on the surface of the airfoil, a number consistent with that used by Visbal and Shang (1989) for this configuration and one that was considered adequate by them from their grid-independence study. The grid distribution for this C-type (444 X 101) mesh was arrived at by examining the disparate length scales of this separated flow problem given by the asymptotic theories. These include the boundary-layer scale of UfRe-i? for attached flow, and the streamwise and normal scales of 0(Re-3’“) and OfRe-s/s) respectively, near flow separation, the separated free shear layer scale of the order of the shear layer thickness, the massively separated zone of O(l) and the scaIe of OfRe-“‘ii) for the boundary-layer eruption near the unsteady separation point; this last scale has been discussed by Smith (1990). An attempt is made to honor all of these scales. Instead of using a further refined grid, results were obtained for this same configuration with a C-grid of (376 x 81) mesh points. The detailed vorticity field computed with this grid was compared with the corresponding results of the earlier grid and it was observed that, in the pitch-up motion for 0” zz B I 22”, the two results were very similar. For the range 22” < 8 _<40”, the fine-grid solution lagged as much as 3” by the time when B = 40” was reached. Finally, for the range 40” < @I 54”, it was felt that there were differences in the dynamics of the evohing vortices. To keep this paper of manageable size, these results are not shown here. At this juncture, it wouid have been appropriate to obtain results using a finer grid, but with fimited sup~~computer resources, it was decided to use the (444 x 101) grid with 44844 mesh points in the caiculations. It should be noted that, at one point away from the body in the normal direction, a third-order accurate forward scheme is used for the vorticity, rather than the CDS used by Blodgett et al. (1990). With the grid size of (444 X 1011, the numerically permissible value of At was 0.05. The minimum grid spacings in the streamwise and normal directions near the leading edge (LEl were 0.003 42 and 0.000 838, respectively. To resolve the temporal scales of this unsteady separated flow problem and to achieve consistency of the solution, a value of At = 0.001 was chosen based on the numerical experiments.
Detailed results are obtained for two cases, fi) with Re = 1~~00 and (ii) for Re = 45flOO. For both of these cases, pitch-up rate Q, = 0.2 was used and the pitch axis was located at the quarter-chord point. The overall fully implicit methodology used for generating the results is highly vectorizable and currently achieves a computational index of 7 ps per time step per mesh point, using a single processor on the CRAY Y-MP 8/864 at the Ohio Supercomputer Center. Vectorization has been achieved satisfactorily, with the code currently performing in the range of 130-180 Mflops. The resufts for the two cases studied are discussed next. 5.3. Case ii): Re = 10 000
Comprehensive results for this case are presented in figs. 3-6. The inertia1 vorticity shown in figs. ?a--3h depicts the formation and subsequent shedding of the dynamic stall vortex, whose evolution in time and space wilt now be discussed. Additional rest&s avaiIabIe in the form of animation for this pitch-up motion are also employed in describing this dynamic stall event. The NACA 0015 airfoil used is started impulsively from rest at 19= 0” and calculations
K. N. Ghta et al. / Forced unsteady flow on NACA 0015 airfoil
0 BU
II 8” C’
Fig.
m
. OM
3. Instantaneou<
OX
100
150
2.m
z.91
.
. 50
tm
. LldO
1.00
150
2.m
vorticity for NACA 0015 airfoil undergoing constant d pitch-up motion: H(deg) = 0.0 (a). 11.4 (h). 14.8 Cc). 22.X Cd). 27.4 CC), 334.3 (f). 36.6 (g), 54.9 (h).
L.SO
Re = 10000.
Fig. 4. Instantaneous
stream
function
for NACA
0015 airfoil
undergoing
constant
i pitch-up
B(deg) = 0.0 (a), 11.4 (b), 14.8 cc), 22.8 cd), 27.4 (e), 34.3 (f), 36.6 (g), 54.9 (h).
motion;
Re = 10000.
(a)
(0 \
l-7
(cl
(h)
Fig. 5. Instantaneous
C, distribution
for NACA
0015 airfoil undergoing
constant
H pitch-up
motion:
Re = IOOO~I.
361
K N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
THETA-4.502
THETA- 18.241
-1
(9
_‘~:___$L_-_ -2mo.0
-.zaoD
-
THETA-7.928
wmo-
-4um.o -I -0.10
(4
O.&l
THETA-21.678
THEM- 11.36 mao-
do
Pi0
0.b
o.ia da CHORD
Fig. 6. Wall vorticity
0.43
for NACA
o.jo
0.b
o.hl
1.60
--a
0015 airfoil undergoing
THETA-25.1 16
(h)
a10
O.&a do
constant
oh
ah
b pitch-up
O.&l oh CHORD
motion:
ah0
da
0.b
Re = 10000.
ah
1.60
362
K. N. Ghiu rt ul. / Forced unsteady ,flow on NACA 0015 airfoil
are carried out until an asymptotic unsteady state is reached. The unsteady shedding process is set up quickly and it is observed that, between 20 and 25 characteristic times, the flow field results do not show any change, there by establishing that an asymptotic state has been achieved. The pitch-up motion is then initiated; separation occurs near the trailing edge (TE) at 6’ = 1.25” and propagates upstream towards the QCP by /3 = 7.93”. An extremely thin separation bubble is simultaneously observed at 0 = 7.93” near QCP; this is barely visible in the inertial surface vorticity plot of fig. 6. As the pitch-up motion continues to 6’ = 12.52”, the LE shear layer develops a kink and the energetic clockwise spinning fluid sets up a pressure gradient, which results in formation of a counterclockwise spinning vortex within the boundary layer and its strength continues to intensify near 0 = 15.95”. At H = 17.1”, a tertiary vortex with clockwise spinning fluid is formed underneath the secondary vortex. As 0 increases further, all of the vortices gain strength and, around 0 = 22.80”, a quaternary vortex emerges near the surface directly underneath the tertiary and secondary vortices. Also, these vortices slowly displace towards QCP; this is also evident in fig. 5 for C,,. It was observed by Smith and Walker (1992) and Walker et al. (1987) that, as these vortices intensify to their maximum strength, strong updrafts begin to develop and, at 0 = 28.56”, the unsteady boundary layer erupts, leading to unsteady inviscid-viscid interactions and the boundary-layer fluid is ejected into the inviscid flow, with some of this fluid being entrained in the formation of the dynamic stall vortex. As 19increases further. The stall vortex slowly convects over the airfoil surface, as shown in fig. 4g; the complete evolution of the stall vortex and its passage on the suction surface is seen in figs. 4c-4h. The various inviscid-viscid interactions near TE are not discussed here. It should be noted that, at higher values of 0, the maximum positive vorticity occurs near TE. 5.4. Case (ii): Re = 4.5 000 This case is simulated to compare the results of the present study with the experimental data of Walker et al. (1985) and, as such, the configuration parameters Re = 45 000, n,, = 0.2, with pitch axis at QCP are used in this study. For this case also, the flow was initiated from impulsive start at 19= 0” and an asymptotic state was attained with the results not varying between characteristic time 20 and 25. A key observation from the results not shown here is that the vorticity field did show the shedding process at early time, but as the time increased the shedding process ceased. Also, at 0 = O”, the vorticity field near TE showed a symmetric needle-type structure with opposite-sign vortices embedded in the flow field. An effort is underway to critically examine these features. The pitch-up motion is then initiated starting with 6’ = 0” and the vorticity contours for this case are depicted in figs. 7a-7h. The discussion that follows is also partly inferred from the animation sequence available to the authors. At 0 = 6.78”, separation is seen on the suction surface, extending between 0.6 and 0.75 c (c = chord length) and at 0 = 7.93” this range extends further as well as a separation bubble first appears just downstream of the LE on the suction surface. As the pitch-up motion continues, the leading-edge separation bubble embedded in the boundary layer intensifies and the upper-surface boundary layer beyond 0.3 c is lifted off up to 0.85 c, as seen from surface vorticity predictions in fig. log. At 0 = 14.80”, a tertiary vortex with clockwise rotating fluid is formed underneath the secondary vortex at about 0.1 c from the LE. With further increase in the value of 0. the secondary bubble grows abruptly in height around 0 = 17.1” and, by t) = 19.39”, the height decreases as the secondary vortex is being squeezed by a second clockwise spinning vortex formed just upstream of it. As the pitch-up motion continues, a quarternary vortex emerges at 0 = 21.68” and the eruption of the secondary vortex is initiated. The process is not as simple as in the case of Re = 10000, due to the presence of a second clockwise spinning vortex. The eruption process appears to
Fig. 7. Instantaneous
vorticity for NACA 0015 airfoil undergoing constant 4 pitch-up motion: (a), 10.2 (b), 12.5 Cc), 20.5 td), 22.8 (e), 34.3 (f), 40.0 (g), 44.6 (h).
Re = 45000. e(deg) = 0.0
Fig. 8. Instantaneous
stream function for NACA 0015 airfoil undergoing constant d pitch-up motion; fKdeg) = 0 (a), 10.2 (b), 12.5 cc). 20.5 cd), 22.8 (e), 04.3 (f), 40.0 (g), 44.6 (h).
Re = 45000.
KN. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
365
(h)
Fig, 9. Instantaneous
C, distribution for NACA 0015 airfoil undergoing constant i pitch-up motion; Re = 45000. 8(deg)= 0.0 (a), 10.2 (b), 12.5 CC), 20.5 cd), 22.8 (e), 34.3 (f), 40.0 (g), 44.6 (h).
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
366
THETA-14.603
CHORO
THETA-18.241
META=
(4
Fig. 10. Wall vorticity
,676
THETA-11.36
for NACA
0015 airfoil undergoing
constant
0 pitch-up
motion;
Re = 45 000.
K.N. Gizia ei al. / Forced unsteady flow on NACA 001.7 airfoil
367
36X
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
have been completed by 8 = 26.27” when the LE shear layer is pinched off by the secondary vortex and subsequently roils into a dynamic stall vortex. The subsequent passage of the dynamic stall vortex over the suction surface is depicted well in the stream function contours in figs. 8f-8h. The Cp distribution for this case is shown in fig. 9, which also shows the movement of the complex vortex structure towards QCP before the secondary vortex erupts. Finally, for this case, the streaklines were compared with the experimental data of Walker et al. (198.5) and the agreement as seen in fig. 11 is quite satisfactory. Rumsey and Anderson (1988a, bl have used a thin-layer compressible Navier-Stokes formulation to give results for this same NACA 0015 airfoil. Their results also showed qualitative agreement with the same experimenta data of Walker et ai. (198.5). In addition, they had also obtained results for a compressible turbulent flow configuration using Baldwin-urns as well as Johnson-~ng models; the latter model compared better with the experimental data used.
6. Conclusions
The unsteady NS analysis of Osswald et al. (1990) is modified suitably to treat higher-Re, constant R,-pitch-up motion for a NACA 0015 airfoil. The third-order accurate biased UDS used appears to be satisfactory for obtaining results up to Re = 45 000. A consistent secondorder accurate wall vorticity boundary condition is used to provide an accurate vorticity field. The numerical method used is a fully vectorized ADI-BGE technique of the authors. From the numerical simulation results for Re = 10000 and 45 000, the structure of the boundary layer is depicted and the separation process is partly quantified. It appears that both cases exhibit a leading-edge separation bubble embedded in the bounda~ layer that eventually leads to a complex vortex structure with four vortices. When the secondary vortex structure erupts, it sets up strong unsteady inviscid-viscid interactions and the stall vortex is created from the LE shear layer fluid. It appears that, if the separation bubble is controlled through suction appropriately between LE and QCP, the dynamic stall may be delayed and/or prevented, depending on the parameters of the problem. In the immediate future, it is planned to carry out numerical experiments to control the dynamic staI1 vortex. Also, software has been developed to visualize streaklines as well as the structure of the boundary layer; this should aid in comparing the unsteady computed results with the time-dependent flow visualization data more readily. Currently, effort is also under way to extend this work to study dynamic stall on wings.
Acknowledgement
The authors are sincerely indebted to Thomas Rohling for creating the animation for this study, to Keith Blodgett for providing help in various phases of the work and to Paul Melde, Christopher No11 and Hugh Thornburg in the preparation of the manuscript using LaTeX.
References Beam, R.M. and R.F. Warming (197X) An implicit factored scheme for the compressible Navier-Stokes equations, AI“54 L 16, 393401. Blodgett, K.J. (1990) Unsteady separated flow past an eliiptic cylinder using the ho-dimensional incompressible Navier-Stokes equations, nit..% ‘2%&s, University of Cincinnati. BIodgett, K., K.N. Ghia, G.A. Osswald and U. Ghia (19901 Analysis of unsteady separated flow past an elliptic cylinder, BUN. Am. Phys. Sot. 35. 2244.
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
369
Carr, L.W. (1988) Progress in analysis and prediction of dynamic stall, J. Aircraft 25, 6-17. Douglas, J. and J.E. Gunn (1964) A general formulation of alternating direction methods, Part I - parabolic and hyperbolic problems, Num. Math. 6, 428-453. Ghia, K.N., G.A. Osswald and U. Ghia (1985) Analysis of two-dimensional incompressible flow past airfoils using unsteady Navier-Stokes equations, in: Numerical and Physical Aspects of Aerodynamic Flows, Vol. III, ed. T. Cebeci (Springer, Berlin) pp. 318-338. Ghia, K.N., U. Ghia and G.A. Osswald (1992) Characterization of dynamic stall phenomenon using two-dimensional unsteady Navier-Stokes equations, Proc. NASA /AFOSR / AR0 Workshop Physics of Forced Unsteady Separation, Moffett Field, California, April; to appear as NASA CP3144, ed. L. Carr, pp. 129-147. Helin, H.E. (1989) The relevance of unsteady aerodynamics for highly maneuverable and agile aircraft, Proc. 4th Symp. Numerical and Physical Aspects of Aerodynamic Flows, Long Beach, CA (January 1989). Osswald, G.A., K.N. Ghia and U. Ghia (1985) An implicit time-marching method for studying unsteady flow with massive separation, AIAA CP 854, pp. 25-37. Osswald, G.A., K.N. Ghia and U. Ghia (1987) A direct algorithm for solution of incompressible three-dimensional unsteady Navier-Stokes equations, AIAA CP 874, pp. 408-421. Osswald, G.A., K.N. Ghia and U. Ghia (1989) Analysis of potential and viscous flows past arbitrary two-dimensional bodies with sharp trailing edges, AZAA CP 895, pp. 668-677. Osswald, G.A., K.N. Ghia and U. Ghia (1990) Simulation of dynamic stall phenomenon using unsteady Navier-Stokes equations, Proc. IMACS Isf Int. Conf: Comp. Physics (June 1990) Boulder CO, pp. 252-261; Comput. Phys. Commurr. 65 (1991) 209-218. Rumsey, CL. and W.L. Anderson, (1988a) Some numerical and physical aspects of unsteady Navier-Stokes computations over airfoils using dynamic meshes, AIAA Paper No. 88-0329. Rumsey, C.L. and W.L. Anderson (1988b) Parametric study of grid size, time step and turbulence modelling on Navier-Stokes computations over airfoils, AGARD CP 437, Validation of Computational Fluid Dynamics Vol. I, Paper No. 5. Smith, C.R., J.D.A. Walker (1992) Some aspects of unsteady separation, Proc. NASA/AFOSR/ARO Workshop on Physics of Forced Unsteady Separation, Moffett Field, CA (April 1990), NASA CP3144, ed. L. Carr, pp. 203-208. Smith F.T. (1990) Finite time singularities in unsteady separation, Proc. Workshop Analytical Methods in Unsteady Separation, Columbus, OH (January 1990) pp. l-22. Visbal, M.R. and J.S. Shang (1989) Investigation of the flow structure around a rapidly pitching airfoil, AU.4 J. 27, 104441055; Proc. AFSOR /FJSRL /DFAM/ UC Workshop II on Unsteady Separated Flows, FJSRL-TR-88-084, pp. 91-108. Walker, J.M., H.E. Helin and J.H. Strickland (1985) An experimental investigation of an airfoil undergoing large amplitude pitching motions, AHA J. 23, 1141-t 142; AIAA Paper No. 85-0039. Walker, J.D.A., C.R. Smith, T.L. Doligalski and A.W. Cerra (1987) Impact of a vortex ring on a wall, J. Fluid Mech. 1x1. 90- 140.
Research
10 (1992) 351-369
Physics of forced unsteady flow for a NACA 0015 airfoil undergoing constant-rate pitch-up motion * K.N. Ghia a, J. Yang a, G.A. Oswald
a and U. Ghia b
* Department of Aerospace Engineering and Engineering Mechanics, Computational Fluid Dynamics Research Laborator)), Unicersity of Cincinnati, Cincinnati, ON 45221, USA h Department of Mechanical, Industrial und Nuclear ~ng~neerjng~ ~omputatjo~lal Fluid Dynamics Laboratory, Unicersity of Cincinnati, Cincinnati, OH 45221, USA Received
6 May 1992
Abstract. The unsteady Navier-Stokes (NS) analysis of Gsswald, Ghia and Ghia in velocity-vorticity variables is modified to study the dynamic stall phenomenon for a N’ACA 0015 airfoil undergoing constant R,, pitch-up maneuvers at Reynolds number Re = 10000 and 45000. The use of third-order accurate biased upwind differenciny for the nonlinear convective terms in the vorticity transport equation removes the spurious oscillations observed in the earlier studies by the authors for these values of Re. The fully implicit and vectorized ADI-BGE method of the authors is used to solve the unsteady NS equations. Instantaneous inertial surface vorticity, which is an invariant of the choice of reference frame selected, is employed to determine the location of separation of the boundary-layer flow on the suction surface; also a separation bubble embedded within the boundary layer is observed for both cases somewhere between the leading edge and the quarter-chord point. Primary, _ secondarv. tertiatv . and auarternarv vortices have been observed before . the dynamic-stall vortex evolves and gathers its maximum strength.
1. Introduction
Specifically, the large amplitude rapid pitching motion associated with the initiation of a high angle-of-attack maneuver typically leads to the generation of a dynamic stall vortex whose evolution results in large transient lift, drag and moment that can, for short periods of time, produce loadings significantly larger than those expected during either steady or quasi-steady flight. Indeed, the successful completion of an abrupt, drastic maneuver can depend upon the ability of holding the dynamic stall vortices in place, at least for the duration of the maneuver, and subsequently bleeding the excess accumulated vorticity in a controlled manner into the wake. Abrupt shedding of large amounts of locally concentrated vorticity can so rapidly alter the lift distribution on a body that a tumbling loss-of-control incident can occur, as the associated rapid changes in moment distribution cannot be tolerated. Carr (1988) has comprehensively reviewed the literature on the dynamic stall phenomenon and has also articulated the effect of key parameters on this phenomenon. Helin (1989) has Correspondence to: K.N. Ghia, Department of Aerospace Engineering and Engineering Mechanics, Cincinnati, 759E Baldwin, MW70, Cincinnati, OH 45221-0070, USA. * This research is supported, in part, by AFOSR Grant No. 90-0249, and by the Ohio Supercomputer No. PESO80. 0169-5983/92/$05.75 All rights reserved
0 1992 - The Japan
Society
of Fluid Mechanics
/ Elsevier
Science
Publishers
University Center
B.V
of
Grant
352
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
also highlighted recent advances on the subject, while stressing the importance of unsteady aerodynamics for highly maneuverable and agile aircraft. In addition, he has raised the important issue of the effect of flow separation on the formation of the energetic dynamic stall vortex. These two reviews adequately point out some of the unresolved issues associated with the problem of dynamic stall. In light of the review article by Carr (19881, no attempt has been made in the present paper to review some of the noteworthy contributions made to this subject in the past. Ghia et al. (1992) recently presented some preliminary results of their unsteady NS analysis using the vorticity stream function formulation. A brief discussion of their analysis is given in the next two sections. Their results for low Reynolds number, Re = 1000, for a NACA 0015 airfoil undergoing’a constant 0, pitch-up motion, showed that the eruption of the fluid, in the form of a secondary vortex from the unsteady boundary layer, rotating in counterclockwise direction near the quarter-chord point (QCP) triggers the formation of the dynamic stall vortex. Subsequently, Osswald et al. (1990) presented results for Re = 10000, a test configuration for which Visbal and Shang (1989) have presented detailed numerical results. At the rapid pitch rate of 0, = 0.2 selected, the central-difference scheme (CDS) with a grid of (566 x 101) showed spurious oscillations in the vorticity field at 13> 32”. Past experience suggests that this type of oscillation can be eliminated using further clustering and/or increasing the grid size in critical areas. In addition, the recent success achieved by Blodgett et al. (1990) in the use of third-order biased upwind differencing scheme (UDS) to obtain results for unsteady separated flow past an elliptic cylinder at higher Re motivated the use of UDS in the present study. Hence, to eliminate the spurious oscillations in the vorticity field, the nonlinear convective terms are approximated using a biased third-order UDS. This then leads to the scope of the present study. The objectives of the present study are as follows: (1) To extend the earlier work of the present investigators to higher Re flows. (2) To understand and quantify the role of unsteady separation in the evolution of the dynamic stall vortex. (3) To determine the physical mechanism involved in moment stall, with the ultimate goal of revealing aerodynamic mechanisms to delay or mediate its effect.
2. Governing
equations
and boundary
conditions
Osswald et al. (1990) modified their earlier analyses (Ghia et al., 1985; Osswald et al., 1985) for fixed body flows to now permit arbitrary six degree-of-freedom maneuvers. The general 3D analysis of Osswald et al. (1990) is in terms of velocity and vorticity variables. This analysis is used here in its 2D form and is briefly presented next; for additional details, the reader is advised to refer to earlier work of the authors. An arbitrary maneuver can be completely defined by specifying the trajectory rn,r(t) of some point B fixed on the body (fig. 1) with respect to an inertial observer, together with the specification of the instantaneous angular velocity n,(t) of the body. Kinematically, the translational velocity of the origin B of the body-fixed frame is then Vn,r(t) = drn,,/dt, the acceleration of the translational acceleration is an,,(t) = d2rn,r/ dt2, while the angular body-fixed reference frame is aa = da,/dt. In the present analysis, these functions, which define a specific maneuver, are assumed to be explicitly prescribed known functions of time and will be given in a later section. With the apparent position vector r being given as
(1)
K.N. Ghia et al. / Forced unsteady flow on NACA 001.7 airfoil
353
Body Fixed Observer (x*,x*)
iii
e
. .....w ..’
- - *..
Body Fixed General
(5’,5*)
Inertial Obierver (x,y) Fig. 1. Inertial and body-fixed coordinate
systems for arbitrary maneuvering body. f&(t) = - R&G
= - t&(r)i,
a,(t)= a&M = -B&M.
the apparent velocity V is kinematically related to the inertial velocity V, as I/= V, - V,,,(t)
Xr.
-O,(t)
(2)
For the 2D simulations of interest in the present study, it is computationally more efficient to solve for a d~~z~rb~~cestream function rather than for inertial velocity directly, where the disturbance stream function @,P is defined as the deviation from uniform flow. Thus, the inertial stream function $r can be written as $1 =Y
+ @o(f)
+ ILP(S’,
f2,
q,
(3)
where y is the vertical inertial coordinate and I&, is the constant of integration which represents the instantaneous displacement of the zero inertial streamline at infinity owing to the inviscid generation of lift (inviscid circulation generated) by the body. The inertial velocity then becomes
where i*is the unit base vector of the inertial Cartesian reference frame and e, and e2 are the covariant base vectors of the generalized-coordinate (
g = (g,,g22-g,,g,,),
where
g,,=
axk axk
c : Y kc, a’$’ a<’
(5)
and (x1, x2> represent the body-fixed Cartesian coordinate system shown in fig. 1, whose unit base vectors are e^, and e^,. For the maneuvering airfoil shown in fig. 1, the unsteady Navier-Stokes equations in generalized orthogonal non-inertial body-fixed coordinates are as shown below where Re = &c/u, with U, being the free-stream reference velocity, c the chord of the airfoil being the reference length, and Y the kinematic viscosity.
354
K. N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil Elliptic
stream
by the following
function problem. The disturbance elliptic equation
stream
function
distribution
is governed
(6) subject
to the boundary $p = 0
conditions
at infinity,
(7)
-tf43w[wf along the body surface. and J/,1(t) This elliptic
+ (x212])
Here,
e,(t)
= 5,&>4
problem
-x’[
- cos 0,(t)]
*P = (x2[ G,,(t)
+
V,‘,,( t) - sin t9,( t)] (8)
+ h(f)
is the instantaneous
pitch angle of the airfoil,
see fig. 1,
(9)
VBZ,lWZ.
for the disturbance
stream
function
is coupled
with the following:
Temporally parabolic uorticity transport problem.
(10) subject
to the boundary w1 = 0
conditions
(11)
at infinity,
(12) along the body, together
with the constraint
(g,,/&)a#:‘/ac*= -([ cos +[sin
e,(t)
O,(t)
along the body. In eq. (lo), the contravariant &v’
= (gl,/G)P’( +[sin
O,(t)
dw = (622/G)-'( +[sin Finally, non-inertial
e,(t)
[ cos O,(t) - Vi,,(t)
[cos
e,(t)
- Vi,I(t)
-
VA,,< t) SX2R,( t)]&
-Vi,,(t)
-x’~,(t)]e^,)
components -l&,(t)
of the apparent
(13) velocity
are
+x2n,(t)]$, *e, +a*P/ag2,
-x’R,(t)]e,) -L&(t)
*e,
+x”l&(t)]e^,
-.x’f2B(t)]e^2)
*e2-a$Clp/a[*.
(15)
the pressure is computed by solving the Neumann-Poisson problem. In the reference frame, the linear momentum equation in terms of disturbance velocity,
K.N. Ghia et al. / Forced unsteady flow on NACA 001.5 airfoil
355
is given as V.(aVD/at+v.VvD+a,(t)
xYD+Re-’
VXW,)=
-V+Vp
(16)
with the disturbance velocity being given as VD = (e,/&)
YfDM2
- (eJ&)
(17)
a$D/%‘,
and subject to the Neumann boundary conditions given as Vp*e,=
-(WD/&fV~VVD+~,(t)
xVD+Re-’
VXw,)*e,,
(18)
VP-~,=
-(WD/at
XVD+Re-’
VXo,)*e,.
(19)
+‘I/.VVD+R,(t)
All terms in this problem are well behaved at infinity. It remains only to prescribe the maneuver and the formuiation
is complete.
3. Pitch-up maneuver and grid generation Visbal and Shang (1989) have investigated a rapid pitch-up maneuver for a NACA 0015 airfoil using a compressible Navier-Stokes code based on the solution algorithm of Beam and Warming (1978). Their maneuver is defined as %/i(t) e,(t)
= 0, [ 1 - exp( - (46,‘t,)t]
= L’,{t - (rJ4.6)
1.
Visbal and Shang (1989) investigated several cases, with 0, ranging from 0.2 to ranging from 0.5 to 1.0. Pitch-axis locations of leading edge, trailing edge, and point were investigated. In that study, one set of parameters used were L?, = 0.2, t, the airfoil pitched about the axis through the l/4-chord point. Specifically, the currently studied is such that rn,i(t)
=O,
n,(t)
= -a,[1
vu,,(t)
= 0,
au(t)
= -(4.6/t0)f10
a,,,(t)
= 0,
R, = 0.2,
(20) 0.6 and t, l/4-chord = 0.5 with maneuver
-exp(-(4.6/t,)t],
t, = 0.5.
exp(-(4.6/t,)t), (21)
The governing equations (6), (10) and (161, with the boundary conditions given in eqs. (7), (81, @l&(15) and (181, (191, are to be solved numerically using a suitable orthogonal grid. A generalized Schwarz-Christoffel grid generation procedure (Osswald et al., 1989) is empioyed to map a C-grid topology of infinite extent to a unit square computational domain. A portion of a typical (444 x 101) clustered orthogonal body-fixed grid distribution is shown in fig. 2 for a NACA 0015 airfoil with a rounded trailing edge. Next, brief comments are made on the numerical method used to solve the discretized problem.
4. Numerical method for the discrete problem In the earlier studies (Ghia et al., 1992; Osswald et al., 1990) by the authors, all spatial derivatives were approximated using central differences. In the present study, the convective derivatives in the vorticity transport equation (10) are appro~mated using a third-order biased UDS, with all other spatial derivatives still being approximated using CDS. Following
356
Y
0.50
-
0.40
-
0.30
-
0. to
-
0.10
-
0.00
-
-a
10
-0.20 -0.30
-
-0.40
-
-0.50
-
f
’ -0.50
0. 00
0.
so
1. 00
X
Fig. 2. Grid distribution
for a NACA
0015 airfoil with (444~
loll
points.
Blodgett et al. f1990), the vorticity transport equation is solved using a vectorized form of the alternating direction imphcit CADI) procedure of Douglas and Gunn ~19~4~. This procedure has been shown by Osswald et ai. (1987) to extend to three-Dimensions easily for the velocity-vort~c~ty formuiatio~ and, when cast into delta form, it is computatjona~ly shghtly more efficient. As described by Osswald et al. f1987), second-order accurate time extrapolations are employed to predict the boundary vorticity at the new time level. These predicted boundary vorticities are combined with the mid time-level stream functions to produce a spatially and temporally second-order accurate solution for the interior vorticity field at the new time level. The boundary condition for the vorticity equation is solved using a second-order accurate expression, since it leads to a dissipative derivative as the leading truncation error term and is therefore consistent with the biased third-order upwind differencing scheme used in this study. The elliptic stream function equation is solved by a direct block Gaussian eIimination (BGEI technique. This ADI-BGE fully implicit numerical simulation technique of the authors has been previousfy described by Osswaid et al. 09851 for a generafized coordinate inertial reference frame. Owing to the ~~~~~-~~~~~~~~~~ of the velocity-vorticity forrmdation, only minor modifications were required in the program to implement the arbitrary maneuvering problem. The Neumann-Poisson problem is solved using a direct solver for which the methodology used is described by Blodgett (1990). Next, the results obtained in this study are discussed.
5. Resuits
and discussion
A test configuration was seiected with Re = 10000, constant pitch-up rate 0, = 0.2 and the pitch axis located at l/4 chord. Osswald et al. (1990) had obtained results for this configura-
tion with a C-grid of (566 x 101) mesh points using CDS. Their detailed results were compared with those obtained by the present method wherein the convective terms are differenced using third-order biased UDS. The results are very similar for 0” < 6 5 17”, beyond which the earlier CDS scheme results show a slight time lag for 17” < 8 I 32”. Beyond 0 = 32”, no comparison is attempted since CDS had spurious oscillations. The medium size C-grid of (444 x 101) points place 204 points on the surface of the airfoil, a number consistent with that used by Visbal and Shang (1989) for this configuration and one that was considered adequate by them from their grid-independence study. The grid distribution for this C-type (444 X 101) mesh was arrived at by examining the disparate length scales of this separated flow problem given by the asymptotic theories. These include the boundary-layer scale of UfRe-i? for attached flow, and the streamwise and normal scales of 0(Re-3’“) and OfRe-s/s) respectively, near flow separation, the separated free shear layer scale of the order of the shear layer thickness, the massively separated zone of O(l) and the scaIe of OfRe-“‘ii) for the boundary-layer eruption near the unsteady separation point; this last scale has been discussed by Smith (1990). An attempt is made to honor all of these scales. Instead of using a further refined grid, results were obtained for this same configuration with a C-grid of (376 x 81) mesh points. The detailed vorticity field computed with this grid was compared with the corresponding results of the earlier grid and it was observed that, in the pitch-up motion for 0” zz B I 22”, the two results were very similar. For the range 22” < 8 _<40”, the fine-grid solution lagged as much as 3” by the time when B = 40” was reached. Finally, for the range 40” < @I 54”, it was felt that there were differences in the dynamics of the evohing vortices. To keep this paper of manageable size, these results are not shown here. At this juncture, it wouid have been appropriate to obtain results using a finer grid, but with fimited sup~~computer resources, it was decided to use the (444 x 101) grid with 44844 mesh points in the caiculations. It should be noted that, at one point away from the body in the normal direction, a third-order accurate forward scheme is used for the vorticity, rather than the CDS used by Blodgett et al. (1990). With the grid size of (444 X 1011, the numerically permissible value of At was 0.05. The minimum grid spacings in the streamwise and normal directions near the leading edge (LEl were 0.003 42 and 0.000 838, respectively. To resolve the temporal scales of this unsteady separated flow problem and to achieve consistency of the solution, a value of At = 0.001 was chosen based on the numerical experiments.
Detailed results are obtained for two cases, fi) with Re = 1~~00 and (ii) for Re = 45flOO. For both of these cases, pitch-up rate Q, = 0.2 was used and the pitch axis was located at the quarter-chord point. The overall fully implicit methodology used for generating the results is highly vectorizable and currently achieves a computational index of 7 ps per time step per mesh point, using a single processor on the CRAY Y-MP 8/864 at the Ohio Supercomputer Center. Vectorization has been achieved satisfactorily, with the code currently performing in the range of 130-180 Mflops. The resufts for the two cases studied are discussed next. 5.3. Case ii): Re = 10 000
Comprehensive results for this case are presented in figs. 3-6. The inertia1 vorticity shown in figs. ?a--3h depicts the formation and subsequent shedding of the dynamic stall vortex, whose evolution in time and space wilt now be discussed. Additional rest&s avaiIabIe in the form of animation for this pitch-up motion are also employed in describing this dynamic stall event. The NACA 0015 airfoil used is started impulsively from rest at 19= 0” and calculations
K. N. Ghta et al. / Forced unsteady flow on NACA 0015 airfoil
0 BU
II 8” C’
Fig.
m
. OM
3. Instantaneou<
OX
100
150
2.m
z.91
.
. 50
tm
. LldO
1.00
150
2.m
vorticity for NACA 0015 airfoil undergoing constant d pitch-up motion: H(deg) = 0.0 (a). 11.4 (h). 14.8 Cc). 22.X Cd). 27.4 CC), 334.3 (f). 36.6 (g), 54.9 (h).
L.SO
Re = 10000.
Fig. 4. Instantaneous
stream
function
for NACA
0015 airfoil
undergoing
constant
i pitch-up
B(deg) = 0.0 (a), 11.4 (b), 14.8 cc), 22.8 cd), 27.4 (e), 34.3 (f), 36.6 (g), 54.9 (h).
motion;
Re = 10000.
(a)
(0 \
l-7
(cl
(h)
Fig. 5. Instantaneous
C, distribution
for NACA
0015 airfoil undergoing
constant
H pitch-up
motion:
Re = IOOO~I.
361
K N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
THETA-4.502
THETA- 18.241
-1
(9
_‘~:___$L_-_ -2mo.0
-.zaoD
-
THETA-7.928
wmo-
-4um.o -I -0.10
(4
O.&l
THETA-21.678
THEM- 11.36 mao-
do
Pi0
0.b
o.ia da CHORD
Fig. 6. Wall vorticity
0.43
for NACA
o.jo
0.b
o.hl
1.60
--a
0015 airfoil undergoing
THETA-25.1 16
(h)
a10
O.&a do
constant
oh
ah
b pitch-up
O.&l oh CHORD
motion:
ah0
da
0.b
Re = 10000.
ah
1.60
362
K. N. Ghiu rt ul. / Forced unsteady ,flow on NACA 0015 airfoil
are carried out until an asymptotic unsteady state is reached. The unsteady shedding process is set up quickly and it is observed that, between 20 and 25 characteristic times, the flow field results do not show any change, there by establishing that an asymptotic state has been achieved. The pitch-up motion is then initiated; separation occurs near the trailing edge (TE) at 6’ = 1.25” and propagates upstream towards the QCP by /3 = 7.93”. An extremely thin separation bubble is simultaneously observed at 0 = 7.93” near QCP; this is barely visible in the inertial surface vorticity plot of fig. 6. As the pitch-up motion continues to 6’ = 12.52”, the LE shear layer develops a kink and the energetic clockwise spinning fluid sets up a pressure gradient, which results in formation of a counterclockwise spinning vortex within the boundary layer and its strength continues to intensify near 0 = 15.95”. At H = 17.1”, a tertiary vortex with clockwise spinning fluid is formed underneath the secondary vortex. As 0 increases further, all of the vortices gain strength and, around 0 = 22.80”, a quaternary vortex emerges near the surface directly underneath the tertiary and secondary vortices. Also, these vortices slowly displace towards QCP; this is also evident in fig. 5 for C,,. It was observed by Smith and Walker (1992) and Walker et al. (1987) that, as these vortices intensify to their maximum strength, strong updrafts begin to develop and, at 0 = 28.56”, the unsteady boundary layer erupts, leading to unsteady inviscid-viscid interactions and the boundary-layer fluid is ejected into the inviscid flow, with some of this fluid being entrained in the formation of the dynamic stall vortex. As 19increases further. The stall vortex slowly convects over the airfoil surface, as shown in fig. 4g; the complete evolution of the stall vortex and its passage on the suction surface is seen in figs. 4c-4h. The various inviscid-viscid interactions near TE are not discussed here. It should be noted that, at higher values of 0, the maximum positive vorticity occurs near TE. 5.4. Case (ii): Re = 4.5 000 This case is simulated to compare the results of the present study with the experimental data of Walker et al. (1985) and, as such, the configuration parameters Re = 45 000, n,, = 0.2, with pitch axis at QCP are used in this study. For this case also, the flow was initiated from impulsive start at 19= 0” and an asymptotic state was attained with the results not varying between characteristic time 20 and 25. A key observation from the results not shown here is that the vorticity field did show the shedding process at early time, but as the time increased the shedding process ceased. Also, at 0 = O”, the vorticity field near TE showed a symmetric needle-type structure with opposite-sign vortices embedded in the flow field. An effort is underway to critically examine these features. The pitch-up motion is then initiated starting with 6’ = 0” and the vorticity contours for this case are depicted in figs. 7a-7h. The discussion that follows is also partly inferred from the animation sequence available to the authors. At 0 = 6.78”, separation is seen on the suction surface, extending between 0.6 and 0.75 c (c = chord length) and at 0 = 7.93” this range extends further as well as a separation bubble first appears just downstream of the LE on the suction surface. As the pitch-up motion continues, the leading-edge separation bubble embedded in the boundary layer intensifies and the upper-surface boundary layer beyond 0.3 c is lifted off up to 0.85 c, as seen from surface vorticity predictions in fig. log. At 0 = 14.80”, a tertiary vortex with clockwise rotating fluid is formed underneath the secondary vortex at about 0.1 c from the LE. With further increase in the value of 0. the secondary bubble grows abruptly in height around 0 = 17.1” and, by t) = 19.39”, the height decreases as the secondary vortex is being squeezed by a second clockwise spinning vortex formed just upstream of it. As the pitch-up motion continues, a quarternary vortex emerges at 0 = 21.68” and the eruption of the secondary vortex is initiated. The process is not as simple as in the case of Re = 10000, due to the presence of a second clockwise spinning vortex. The eruption process appears to
Fig. 7. Instantaneous
vorticity for NACA 0015 airfoil undergoing constant 4 pitch-up motion: (a), 10.2 (b), 12.5 Cc), 20.5 td), 22.8 (e), 34.3 (f), 40.0 (g), 44.6 (h).
Re = 45000. e(deg) = 0.0
Fig. 8. Instantaneous
stream function for NACA 0015 airfoil undergoing constant d pitch-up motion; fKdeg) = 0 (a), 10.2 (b), 12.5 cc). 20.5 cd), 22.8 (e), 04.3 (f), 40.0 (g), 44.6 (h).
Re = 45000.
KN. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
365
(h)
Fig, 9. Instantaneous
C, distribution for NACA 0015 airfoil undergoing constant i pitch-up motion; Re = 45000. 8(deg)= 0.0 (a), 10.2 (b), 12.5 CC), 20.5 cd), 22.8 (e), 34.3 (f), 40.0 (g), 44.6 (h).
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
366
THETA-14.603
CHORO
THETA-18.241
META=
(4
Fig. 10. Wall vorticity
,676
THETA-11.36
for NACA
0015 airfoil undergoing
constant
0 pitch-up
motion;
Re = 45 000.
K.N. Gizia ei al. / Forced unsteady flow on NACA 001.7 airfoil
367
36X
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
have been completed by 8 = 26.27” when the LE shear layer is pinched off by the secondary vortex and subsequently roils into a dynamic stall vortex. The subsequent passage of the dynamic stall vortex over the suction surface is depicted well in the stream function contours in figs. 8f-8h. The Cp distribution for this case is shown in fig. 9, which also shows the movement of the complex vortex structure towards QCP before the secondary vortex erupts. Finally, for this case, the streaklines were compared with the experimental data of Walker et al. (198.5) and the agreement as seen in fig. 11 is quite satisfactory. Rumsey and Anderson (1988a, bl have used a thin-layer compressible Navier-Stokes formulation to give results for this same NACA 0015 airfoil. Their results also showed qualitative agreement with the same experimenta data of Walker et ai. (198.5). In addition, they had also obtained results for a compressible turbulent flow configuration using Baldwin-urns as well as Johnson-~ng models; the latter model compared better with the experimental data used.
6. Conclusions
The unsteady NS analysis of Osswald et al. (1990) is modified suitably to treat higher-Re, constant R,-pitch-up motion for a NACA 0015 airfoil. The third-order accurate biased UDS used appears to be satisfactory for obtaining results up to Re = 45 000. A consistent secondorder accurate wall vorticity boundary condition is used to provide an accurate vorticity field. The numerical method used is a fully vectorized ADI-BGE technique of the authors. From the numerical simulation results for Re = 10000 and 45 000, the structure of the boundary layer is depicted and the separation process is partly quantified. It appears that both cases exhibit a leading-edge separation bubble embedded in the bounda~ layer that eventually leads to a complex vortex structure with four vortices. When the secondary vortex structure erupts, it sets up strong unsteady inviscid-viscid interactions and the stall vortex is created from the LE shear layer fluid. It appears that, if the separation bubble is controlled through suction appropriately between LE and QCP, the dynamic stall may be delayed and/or prevented, depending on the parameters of the problem. In the immediate future, it is planned to carry out numerical experiments to control the dynamic staI1 vortex. Also, software has been developed to visualize streaklines as well as the structure of the boundary layer; this should aid in comparing the unsteady computed results with the time-dependent flow visualization data more readily. Currently, effort is also under way to extend this work to study dynamic stall on wings.
Acknowledgement
The authors are sincerely indebted to Thomas Rohling for creating the animation for this study, to Keith Blodgett for providing help in various phases of the work and to Paul Melde, Christopher No11 and Hugh Thornburg in the preparation of the manuscript using LaTeX.
References Beam, R.M. and R.F. Warming (197X) An implicit factored scheme for the compressible Navier-Stokes equations, AI“54 L 16, 393401. Blodgett, K.J. (1990) Unsteady separated flow past an eliiptic cylinder using the ho-dimensional incompressible Navier-Stokes equations, nit..% ‘2%&s, University of Cincinnati. BIodgett, K., K.N. Ghia, G.A. Osswald and U. Ghia (19901 Analysis of unsteady separated flow past an elliptic cylinder, BUN. Am. Phys. Sot. 35. 2244.
K.N. Ghia et al. / Forced unsteady flow on NACA 0015 airfoil
369
Carr, L.W. (1988) Progress in analysis and prediction of dynamic stall, J. Aircraft 25, 6-17. Douglas, J. and J.E. Gunn (1964) A general formulation of alternating direction methods, Part I - parabolic and hyperbolic problems, Num. Math. 6, 428-453. Ghia, K.N., G.A. Osswald and U. Ghia (1985) Analysis of two-dimensional incompressible flow past airfoils using unsteady Navier-Stokes equations, in: Numerical and Physical Aspects of Aerodynamic Flows, Vol. III, ed. T. Cebeci (Springer, Berlin) pp. 318-338. Ghia, K.N., U. Ghia and G.A. Osswald (1992) Characterization of dynamic stall phenomenon using two-dimensional unsteady Navier-Stokes equations, Proc. NASA /AFOSR / AR0 Workshop Physics of Forced Unsteady Separation, Moffett Field, California, April; to appear as NASA CP3144, ed. L. Carr, pp. 129-147. Helin, H.E. (1989) The relevance of unsteady aerodynamics for highly maneuverable and agile aircraft, Proc. 4th Symp. Numerical and Physical Aspects of Aerodynamic Flows, Long Beach, CA (January 1989). Osswald, G.A., K.N. Ghia and U. Ghia (1985) An implicit time-marching method for studying unsteady flow with massive separation, AIAA CP 854, pp. 25-37. Osswald, G.A., K.N. Ghia and U. Ghia (1987) A direct algorithm for solution of incompressible three-dimensional unsteady Navier-Stokes equations, AIAA CP 874, pp. 408-421. Osswald, G.A., K.N. Ghia and U. Ghia (1989) Analysis of potential and viscous flows past arbitrary two-dimensional bodies with sharp trailing edges, AZAA CP 895, pp. 668-677. Osswald, G.A., K.N. Ghia and U. Ghia (1990) Simulation of dynamic stall phenomenon using unsteady Navier-Stokes equations, Proc. IMACS Isf Int. Conf: Comp. Physics (June 1990) Boulder CO, pp. 252-261; Comput. Phys. Commurr. 65 (1991) 209-218. Rumsey, CL. and W.L. Anderson, (1988a) Some numerical and physical aspects of unsteady Navier-Stokes computations over airfoils using dynamic meshes, AIAA Paper No. 88-0329. Rumsey, C.L. and W.L. Anderson (1988b) Parametric study of grid size, time step and turbulence modelling on Navier-Stokes computations over airfoils, AGARD CP 437, Validation of Computational Fluid Dynamics Vol. I, Paper No. 5. Smith, C.R., J.D.A. Walker (1992) Some aspects of unsteady separation, Proc. NASA/AFOSR/ARO Workshop on Physics of Forced Unsteady Separation, Moffett Field, CA (April 1990), NASA CP3144, ed. L. Carr, pp. 203-208. Smith F.T. (1990) Finite time singularities in unsteady separation, Proc. Workshop Analytical Methods in Unsteady Separation, Columbus, OH (January 1990) pp. l-22. Visbal, M.R. and J.S. Shang (1989) Investigation of the flow structure around a rapidly pitching airfoil, AU.4 J. 27, 104441055; Proc. AFSOR /FJSRL /DFAM/ UC Workshop II on Unsteady Separated Flows, FJSRL-TR-88-084, pp. 91-108. Walker, J.M., H.E. Helin and J.H. Strickland (1985) An experimental investigation of an airfoil undergoing large amplitude pitching motions, AHA J. 23, 1141-t 142; AIAA Paper No. 85-0039. Walker, J.D.A., C.R. Smith, T.L. Doligalski and A.W. Cerra (1987) Impact of a vortex ring on a wall, J. Fluid Mech. 1x1. 90- 140.