Application of the method of integral relations (MIR) to transonic airfoil problems. Part I: Inviscid supercritical flow over symmetric airfoils at zero angle of atack

September 1970

DEPARTMENT OF THE NAVY

NAVAL SHIP RESEARCH AND DEVELOPMENT CENTER

WASHINGTON. o. o 20034

APPLICATION OF THE ÌTHOD OF INTEGRAL RELATIONS (nR)

TO TRANSONIC AIRFOIL PROBLEMS. PART I - INVISCID

SUPERCRITICAL FLOW OVER SYTRIC AIRFOILS AT ZERO ANGLE OF ATTACK

by

Tsze C. Tai

This document has been approved for public release and sale; its distribution is unlimited.

NSRDC Report 3424 Aero Report 1176

TABLE 0F CONTENTS

Page

SUMMARY y

INTRODUCTION I

ANALYSIS 2

BASIC FLOW EQUATIONS 2

METHOD OF INTEGRAL RELATIONS 3

NUMERICAL PROCEDURES FOR SY1'TRIC AIRFOILS AT

ZERO ANGLE OF ATTACK 4

APPLICATION OF FREE-STREAM BOUNDARY CONDITIONS 4

EVALUATION OF INITIAL VALUES 6

EFFECT OF BLUNTNESS 7

TREATMENT OF SONIC POINT 7

LOCATION OF IMBEDDED SHOCK 8

RESULTS OF EXPLORATORY WORK 9

CONCLUSIONS 11

ACKNOWLEDGEMENT 12

APPENDIX A - ORDINARY DIFFERENTIAL EQUATIONS REDUCED WITH A

TWO-STRIP APPROXIMATION IN CARTESIAN COORDINATES 13

APPENDIX B - EQUATIONS FOR EVALUATION OF INITIAL VALUES USING

TWO-STRIP APPROXIMATIONS 18

APPENDIX C - TWO-DIMENSIONAL INVISCID FLOW EQUATIONS IN A BODY

ORIENTED ORTHOGONAL CURVILINEAR COORDINATE SYSTEM 21

REFERENCES 32

LIST OF FIGURES

Page

Figure 1 - Coordinate Systems and Strip Boundaries 22

Figure 2 - Extension of Free-Stream Boundary to "Infinity" . 23 Figure 3 - Accuracy of Different Approximations for the

Method of Integral Relations 24

Figure 4 - Velocity Gradients Along Surface of an NACA 0012 Airfoil

at M = 0.75 25

oe

Figure Sa - Downstream Flow Field of an NACA 0012 Airfoil

at M = 0.75, for a "CorrectT' Shock Location 26

oe

Figure 5b - Downstream Flow Field of an NACA 0012 Airfoil

at M = 0.75, for a "Wrong" Shock Location 27

Figure 6 - Iterative Processes for a Supercritical Flow

Over Symmetric Airfoils 28

Figure 7 - Pressure Distribution on a Symmetric 8.4 Percent

Circular-Arc Airfoil at Moe = 0.85 29

Figure 8 - Pressure Distribution on an NACA 0012 Airfoil

atM

0.75 30

oe

Figure 9 - Mach Number Contours for NACA 0012 Airfoil

atM =0.75

31

SYOL S

A,B,Q Functions in Equation (6)

ak Constants in Equation (7)

C C

1 ± 2/[(y-1)]

c Chord length of an airfoil

D Function defined by Equation (A-11) E1 thru E4 Functions defined by Equations (A-12) F1 thru F5 Functions defined by Equations (A-13)

G Function defined by Equation (A-14)

K K l/(YN)

M

Mach number

N Number of strips

P Static pressure normalized by its free-stream value

q1,q2 Functions defined by Equations (A-15)

R Local body radius of curvature normalized by the chord length

s,n Orthogonal curvilinear coordinates measured along and normal

to the airfoil surface normalized by the chord length

u,v Velocity components in cartesian coordinates, normalized

by the free-stream velocity

Velocity components along and normal to the airfoil surface, normalized by the free-stream velocity

x,y Cartesian coordinates normalized by the chord length

Y Specific heat ratio

ô Normal distance between the airfoil surface and system boundary normalized by the chord length

9 Surface inclination angle with respect to the direction of free-stream velocity

p Static density normalized by its free-stream value

x Angle between shock tangent and normal to the airfoil axis of symmetry

S UB S C PTS Conditions at airfoil surface

c,x,y Evaluated by using continuity, x-momentum and y-momentum equations, respectively

ô Conditions at system boundary

0, 1, 2 Along lines 0, 1, and 2 in a two-strip approximation

Conditions at free-stream

S tTMMA RY

The feasibility of applying the method of integral relations _(MIR) to transonic flows over symmetric airfoils has been studied. In order to take account of the severe transversal flow variation and still retain simplicity in computation, the method is modified so that the number of strips used may be considerably higher than the order of the polynomial which approximates the integrand. _{An important feature of}

this modification, however, is its capability to extend the free-stream boundary to "infinity". _{Using one-strip and two-strip approximations,} flow equations are reduced to a set of ordinary differential equations in a cartesian coordinate system. _{Numerical procedures, including the} treatment of the sonic point and the determination of shock location, are also formulated. _{The fourth-order Runge-Kutta method was used in} numerical computation.

The advantage of small computer capacity and time required by the method is evidenced by exploratory calculations for a symmetric circular-arc airfoil and an NACA 0012 airfoil, traveling at supercritical speeds. The calculated results compare fairly well with those obtained by the unsteady finite difference scheme for the case of a circular-arc airfoil and very well with experimental data for that of an NACA 0012 airfoil.

In principle, the results may be improved by increasing the number of strips and extending the free-stream boundary condition sufficiently far away from the airfoil. It is also found that in transonic airfoil problems using a cartesian coordinate system, the order of polynomial approximating the integrand in the method of integral relations should be at least second-order or higher.

INTRODUCTION

The determination of the transonic flow over an arbitrary two-dimensional airfoil has long been one of the most difficult problems in aerodynamics. The essential nonlinearity in the flow governing equations prevents an exact analytical solution even for an inviscid

case. _{A simplified solution for thin, sharp nosed airfoils was}

develop-ed by Spreiter and Alksne (Reference 1) in which an iterative process is

required. Other works such as the method of parametric differentiation suggested by Rubbert and Landahi (Reference 2) , and the integral method

employed by Crown (Reference 3) are generally concerned with sharp nosed thin symmetric airfoils. _{Problems associated with the arbitrary lifting} airfoils involve more difficulty for analytical treatment, and there has been no alternative approach other than numerical procedures.

The numerical method of finite difference scheme was applied by Magnus, Gallaher and Yoshihara (Reference 4) to calculate the inviscid

supercritical airfoil flows. The approach is an unsteady one where the desired steady flow is obtained as an asymptotic limit for large times from a given initial flow. Although there is no restriction of the method on the shape of airfoils, its practical usefulness is limited by the large computer time and capacity required. It will be even more difficult then to consider viscous and three-dimensional effects.

In order to meet the urgent need for the design of an efficient transonic airfoil, one must overcome the above difficulty by selecting a method which requires relatively small computer time and capacity, and thus may be able to handle the general lifting airfoil problem including the boundary layer effects. _{Among various approaches, the} method of integral relations (KIR), originally developed by Dorodnitsyn (Reference 5), appears to satisfy this requirement. The purpose of this study is to apply the method of integral relations to the transonic airfoil problems.

Part I covers exploratory calculations for symmetric airfoils traveling at supercritical speeds at zero angle of attack, using first- and second-order polynomial approximations. Cases for arbitrary lifting airfoils along with a third-order polynomial approximation will be presented in

Part II. Further extension to viscous flows will be investigated and presented in Part III.

-1-ANALYSIS BASIC FLOW EQUATIONS

In Cartesian coordinates, the governing equations of a steady, adiabatic, inviscid, non-heat conducting, and perfect gas flow are as

follows: Continuity O(pu) O(pv) -ox x-Momentum o 2 0

-

(KP + pu ) + - (puy) = O

ox

y-Momentum O (puy) +. (KP+ pv2) = O Energy p,/ (4)

where lengths are normalized with respect to the chord length, and pressure, density and velocity components with respect to their free-stream values. The symbol K represents

where for air, y = = 7/5. Note that for flows of interest, an isentropic relation has been used in place of the energy equation.

The boundary conditions are: at the airfoil surface, the normal velocity component is equal to zero, i.e.,

V = u sin e + e = o n

where e is the surface inclination angle with respect to the direction of free-stream velocity. At infinity, the flow is undisturbed, i.e.,

(5)

-3-p= -3-p=u= 1

v=O

(5h)

METHOD OF INTEGRAL RELATIONS

The method of integral relations, originally proposed by Dorodnitsyn (Reference 5), has been employed as a general method of numerical solution for nonlinear aero-hydrodynamic problems. It is the purpose of this study to apply the method to the transonic airfoil problems.

In applying the method of integral relations, the system of flow equations must be written in divergence form:

A(x,y,u,...) + B(x,y,u,...) = Q(x,y,u,...) (6)

Note that Equations (1), (2) and (3) are already in divergence form. Unlike the procedure in the hypersonic blunt body problem, some typical

streamlines are used here as strip boundaries to divide the flow field into strips, as shown in Figure 1. The reason for doing this is twofold, i.e., to facilitate the numerical procedure and to locate the streamlines. The divergence form of Equations (1) thru (3) may then be integrated outward from the airfoil surface (but not necessarily normal to the surface) to each strip boundary successively at some constant value of

x. This procedure reduces the partial differential equations (with independent variables of x and y) to ordinary ones (with independent variable x). In performing the integration, the variation of integrand along y must be known. A general approach is to approximate the inte-grands by interpolation polynomials, for example A by

A

=

ak (x)(y y0)k (7)

where N is the number of strips and ak, constants evaluated at strip boundaries. In principle, the actual variation may be represented more closely with an increasing number of strips.

To examine whether the method of integral relations may be applied to the transonic airfoil problem, one- and two-strip approximations are

used for exploratory computation. The complete set of ordinary differential equations derived by using two-strip approximation are presented in Appendix A.

Since the method is a general one, it can be applied to any shape of airfoils at zero or certain angle of attack. However, the essential difficulty in applying the method to the present transonic airfoil problem lies in the numerical procedures. The procedure for handling cases of symmetric airfoils are formulated and outlined in the subse-quent sections.

NUMERICAL PROCEDURES FOR SYMMETRIC AIRFOILS

AT ZERO ANGLE OF ATTACK

With the partial differential equations reduced to a set of ordinary differential equations by the method of integral relations, the numerical integration may be carried out along the longitudinal axis x, using a

standard scheme such as the Runge-Kutta method. In the course of numerical calculation, however, difficulty arises in the application of free-stream boundary condition, the evaluation of initial values for numerical inte-gration, treatment of the sonic point, location of the imbedded shock, and the effect of bluntness of the airfoil. These problem areas must be considered and overcome for a complete solution to exist.

APPLICATION OF FREE-STREAM BOUNDARY CONDITIONS

As stated previously, the boundary conditions of the present problem are that at the airfoil surface, the normal velocity component must vanish and at infinity, the flow should be undisturbed. The application of free-stream condition at infinity requires a sufficiently large domain of

integration, that in turn, requires a large number of strips to be used. Since the number of strips is normally equal to the order of a polynomial which approximates the integrands, the actual number of strips must then be compromised by the numerical difficulties involved in solving the

coefficient of the integrand polynomial. The contradiction between the physical requirement and the mathematical limitation constitutes the difficulty in applying the method of integral relations to the transonic airfoil problem.

In order to apply the free-stream boundary condition at "infinity" without involving the above numerical problem, it is advisable to treat the whole integration domain as a series of different effective regions; within each region a small number of strips may be used. The idea for doing this is illustrated in Figure 2.

For simplicity, only the upper part of a symmetric airfoil is shown in Figure 2 along with two-strip approximations in the method of integral relations. Let lines (0, 1, 2) be the strip boundaries associated with a two-strip approximation. The solution may be integrated if the flow properties along line "2" are given. If line "2" coincides with 'b" of set (a, b, c), then the boundary conditions along line "2" are determined by integration along line "b" with another two-strip approximation for lines (a, b, c). Similarly, the values along line "c" are provided by results along line "Tr' obtained from integration of lines (,

_fl,

Ç) . The

integration of all sets should be carried out simultaneously with single-valued solution along each line. In principle, the procedure may be

repeated on and on to extend the final boundary at "infinity". In fact,

line "Ç", the outer strip boundary of the third set is already located sufficiently far away (by three to four chord lengths) from the airfoil and thus the flow there could be eventually undisturbed, or at most, very slightly disturbed. It is then permissible to impose the free-stream boundary condition along line "Ç" for practical computation.

Associated with the technique developed here are terminologies used later on for different integration schemes. For example, in Figure 2 three sets of two-strip approximations working together will be called as three-two-strip integration. By the same token, two-three-strip is for

two sets of three-strip approximations used. Note that two-three-strip differs from three-two-strip integrations.

The advantage of the idea is best explained with Figure 3, in which three kinds of approximations, namely one-strip, three-strip and two-two-strip, for the method of integral relations are compared. For a given integrand distribution along the transversal direction of a physical flow at a fixed x-axis, the one-strip approximation, as shown in Figure 3, is very inaccurate for transonic airfoil flows. The two-two-strip scheme, which yields three effective strips, is seen being even superior to the

actual three-strip approximation. Of course the less labor work for the two-two-strip compared to the three-strip approximation offers another dividend of the present approach. Due to its ability of using a large number of strips, the present idea is believed to pave a way to the successful application of the method of integral relations to the tran-sonic airfoil problems.

EVALUATION OF INITIAL VALUES

With the free-stream boundary condition properly applied, the numerical integration of reduced ordinary differential equations may be carried out if the initial values are known at the "initial points". The initial points are those at the strip boundaries corresponding to a

constant x-axis where the numerical integration starts. The initial point for the body streamline is usually located a small distance away

from the stagnation point.

To evaluate these initial values, the method of integral relations is applied in the upstream region of the airfoil. The equations for this integration are reduced from Equations (1) thru (3) , letting

y = y = O along the stagnation streamline. They are given in Appendix B. Since the integration starts sufficiently far upstream of the airfoil where the flow is essentially undisturbed, the free-stream condition can be used as initial values for this upstream integration.

With free-stream properties used as the initial values and also as the boundary condition along the outermost strip, the numerical integra-tion may proceed, yielding a uniform flow. This is because there is nothing to vary then. The method of integral relations works until a nonzero term , which may be called as the transverse velocity

\ Y/o

gradient along the stagnation streamline, was identified, as shown in Appendix B. The existence of the transverse velocity gradient causes the velocity to decrease along the stagnation streamline, and accordingly, the properties along other lines begin to vary. Physically this term reflects the flow feedback from the airfoil to the upstream region. The

expression for the transverse velocity gradient is readily derived within the approximation of the method of the integral relations and is presented in Appendix B.

EFFECT OF BLUNTNESS

The nonlinearity of the flow equations in transonic _{flow becomes} more significant as the airfoil is blunted. _{According to Pearcey}

(Reference 6) , a peaky velocity may be obtained through proper

blunt-ness of the airfoil. _{Unlike sharp-nosed airfoils, the blunted airfoil} has large and drastically changing body slopes in its fore portion,

especially near the nose. Since the flow variation is directly influenced by the body shape, the number of strips in the present method was sub-stantially increased for the nose region of the airfoil. The use of large number of strips, however, needs not correspond to a high-order polynomial approximation, as discussed above. _{Hence, the basic computer program was} still a second-order polynomial approximation for the present tnethod.

Another difficulty in the case of a blunted airfoil _{is that the} surface derivative dy/dx becomes infinite at the leading edge. In this connection, it is helpful to use a body oriented orthogonal curvilinear coordinate system near the nose surface, while the outer region still remains in cartesian coordinates. _{It may be called as a 'mixed coordinate} system", as shown in Figure 1. The solution of the outer region provides the boundary condition for the inner region. _{Since the width of the inner} region can be fairly thin (normally up to O.2c), _{one- or two-strip} approxi-mation in orthogonal coordinates should be good enough for this region.

Equations for the one-strip approximation in this coordinate system are given in Appendix C.

TREATMENT OF SONIC POINT

Near the sonic point when u2 _{KyP/p, the denominator of Equation (A-1)} goes to zero which might cause the velocity gradient to blow _{up. However,} physically there should be a continuous flow through the sonic point.

A

continuous solution exists only if the numerator of that equation _becomes zero simultaneously as the denominator does so that the ratio O/O _still yields a finite quantity. _{To achieve this goal involves the adjustment} of the height of the outermost strip boundary. _{Note that this does not}

contradict the argument for imposing the free-stream boundary as far as possible, as mentioned previously. _{In fact, the solution is dictated by}

-7-the particular integration scheme used. For instance, the three-two-strip integration scheme allows the free-stream boundary to be imposed farther from the airfoil surface than the two-two-strip one. For a given integra-tion scheme, with specified strip spacing at the initial staintegra-tion, the exact height of the outermost strip boundary is solely determined by the condition of a continuous solution at the sonic point.

This phenomenon is illustrated in Figure 4 where the velocity gradients, du0/dx for various outermost strip boundary heights are plotted versus the distance from the airfoil leading edge, x. The "wrong" outermost strip boundary heights cause the velocity to diverge. While the "exact" value

for the outermost strip boundary height is completely insignificant itself, yet the determination of this value is necessary in order to obtain a

continuous solution at the sonic point. This procedure is quite similar to the hypersonic blunt body problem in which the adjustment of the shock standoff distance is required for a continuous solution to exist at the sonic point (for example, see Reference 7). Even with the "exact" height of the outermost strip boundary, the numerical integration is still very sensitive and unstable near the sonic point. It has to be replaced by a polynomial curve fit through the sonic region, as shown in Figure 4,,

LOCATION OF INBEDDED SHOCK

For transonic flow with supercritical free-stream Mach number, in general, there will be an imbedded shock on the airfoil. The location of the foot of the shock is determined by the condition such that the flow returns to its undisturbed state in the far downstream. For an inviscid flow, where the velocity must be tangent to the surface, the normal shock

relations may be applied along the surface streamline to give flow

quantities after the shock. The determination of the shock location also requires an iteration process so as to satisfy the free-stream condition in the downstream. This procedure is analogous to the principle of nozzle flow where the location of the shock is determined by matching the flow pressure at the nozzle exit (Reference 8) Figures 5a and 5b show the presently calculated flow field in the downstream of an NACA 0012 airfoil

at M = 0.75. The flow field returns to its free-stream condition in the far downstream if the assumed shock location is "correct", as shown in Figure 5a; otherwise it diverges, as in Figure 5b.

The iterative process involved here should not be confused with that for the adjustment of the outermost strip boundary for continuous solution at the sonic point. As shown in Figure 6, each iteration has its own sole purpose so that a unique solution can be obtained. These two kinds of

iteration can be treated independently. The only relationship between them is that the height of the outermost strip boundary determined by the Iteration I (see Figure 6) will also be used in the Iteration II.

The above considerations cover the numerical procedure to be employed for an inviscid flow over a symmetric airfoil at zero angle of attack. Procedures concerning the arbitrary lifting airfoil case will be discussed in Part II and the extension to viscous flows, in Part III.

RESULTS OF EXPLORATORY WORK

Exploratory calculations have been carried out for inviscid flows over a symmetric 8.4-percent circular-arc airfoil at Moe = 0.85, and an NACA 0012 airfoil at Moe = 0.75. Both cases are at zero angle of attack. The choice of these airfoil shapes and flow conditions is based on the availability of corresponding theoretical and/or experimental data for comparíson purposes. The numerical procedures outlined above have been used in the present calculations. 4ith necessary iterations done and converged to within five percent, the results are considered as unique solutions within the framework of the method.

The surface pressure distributions for the circular-arc case, along with those obtained by the unsteady finite difference scheme (Reference 4) are presented in Figure 7. One-strip and two-two-strip integration schemes are used in the present method. As shown in Figure 7, the calculated

pressures compare fairly well with those of Reference 4. The use of the two-two-strip scheme not only corresponds to that of refining the mesh size in the finite difference scheme, but also extends the integration domain much farther. Here, for one-strip solution, the free-stream boundary

(outermost strip boundary) was at approximately 1.5 chord lengths away from the axis of symmetry; for two-two-strip scheme, it increased to about 2.5 chord lengths.

The results also reveal two inherent discrepancies in using the one-strip integration scheme. _{First, the linear variation of the integrand} is a very poor approximation for transonic airfoil problems where the flow

-9-properties vary severely in the transverse direction and usually having values exceeding those at the free stream and at the body. The linear variation of the integrand associated with the one-strip integration fails

to describe this physical phenomenon even in an approxímate manner. Second,

in the one-strip integration where the free-stream boundary is placed at constant distance away from the axis of symmetry, the flow resembles that for a convergent-divergent nozzle flow. As a consequence, the sonic point will always occur at the minimum area, or corresponding to the maximum thickness of the airfoil. Physically this is not always true. Based on the above argument, the one-strip integration scheme will not be used in future calculations.

Calculated pressure distributions for the case of an NACA 0012 airfoil traveling at = 0.75 are presented in Figure 8. Also shown in the figure are experimental data of NACA TN 2174 (Reference 9) and NACA TN 1746

(Reference 10) for comparison purposes. The N-2-strip integration scheme was used in the present calculation, where N varies from 3 to 6. The flow

in the inner stagnation region (see Figure 1) was calculated with the aid of a body oriented orthogonal curvilinear coordinate system. One-strip approximation was employed in this small and thin region.

The present results are labeled with two shock locations, namely, "correct" and "wrong". The criterion for a shock location to be "correct" or "wrong" is based on the far downstream flow conditions. As mentioned

in the previous section, if the assumed shock location allows the far downstream flow field to converge to its free-stream state, it is con-sidered as "correct"; otherwise it is "wrong". Although the pressure distributions for both shock locations are fairly close to each other, their resulting downstream situations are quite different, as shown in Figures Sa and 5b. It is seen in Figure 8 that present results agree fairly close to the experimental data, yet better agreement is obtained with that of the "correct" shock location.

Finally, Figure 9 presents the Mach number contours in the vicinity

of the NACA 0012 airfoil at M = 0.75. The overall domain of computation

was actually set about 4 c from the upstream to the leading edge and ended downstream also about 4 c from the trailing edge, with the outer-most strip boundary at approximately 4.5 c. It is generally known that

the flow becomes practically undísturbed beyond _{four chord lengths away} from the airfoil in most transonic wind-tunnel _tests.

It is seen from the figure that the flow variation becomes more and more severe as the surface is approached. _{Also observed in Figure 9 are} the flow recotnpression in front of the shock and weak _{re-expansion aft the}

shock. _{These pehnomena have been established as typical flow patterns in}

a supercritical flow in previous analyses (see References 11 and _{12, for} instance).

The shock shown is a curved one. _{The reason for this is that since} the normal shock relations were applied only at its foot but nowhere else, the shock is still free to bend away from the surface. _{Of course this} becomes possible only if there is no strip boundary going through _the

shock wave. _{Fortunately, the shock length is too short to touch any} field-strip boundaries, since the imbedded shocks are generally weak ones in a transonic flow.

From the experience of the present calculation, it is found that the non-inflective nature of a second-order polynomial curve (two-strip approxi-mation) makes it inadequate for approximating the integrand distribution _in

certain regions of a transonic flow. _{This is because in certain regions} the flow variation becomes so severe that its distribution curve will be

inflected. It is believed that a third-order polynomial, having possibility of inflection point, should give a better approximation than the second-order one.

The computation time on an IBM 7090 computer for results using the

two-two-strip scheme is about 0.07 seconds for each step size. For a typical run with a step size varying from 0.001 to 0.02, the number _of strips varying from four to seven, and 20 iterations, it requires,

conservatively speaking, about 15 minutes on the IBM 7090 machine. The

computer time is approximately proportional to the number of strips _used.

CONCLUS IONS

As a result of present work, the following conclusions _{y be drawn:} 1. _{The application of the method of integral relations to traasonic} airfoil problems beeoLles feaiible through a modification by _{which the}

number of strips may be considerably higher than the order of the poly-nomial approximation.

-11-An important feature of the present modification of the method, however, is its capability of extending the free-stream boundary condition to "infinity'.

Using a basic second-order polynomial approximation, the calcu-lated results for cases of a circular-arc and an NACA 0012 airfoil compare favorably well with those obtained by the unsteady finite difference schene and experimental data, respectively. In principle, the results may be improved by increasing the number of strips and extending the free-stream boundary condition sufficiently far away from the airfoil.

In transonic airfoil problems using a cartesian coordinate system, the order of polynomial approximating the integrand in the method of

integral relations should be at least second-order or higher.

ACKNOWLEDGMENT

The author wishes to thank Dr. S. de los Santos, head of the

High-Speed Aerodynamics Division at the Naval Ship Research and Development Center, for his suggestion of this study and his encouragement throughout all phases of the present work.

APPENDIX A

ORDINARY DIFFERENTIAL EQUATIONS REDUCED WITH A TWO-STRIP APPROXIMkTION IN CARTESIAN COORDINATES

With the use of two-strip approximation in the method of integral relations, the ordinary differential equation associated with strip boundaries 0, 1 and 2 are as follows:

du o / P 0 2

\,KY+u

E1 dx f P n (KY 2.

-

u2 D

rO\

_p

0J

o dx pl -dx dy

y =u

-o

odx

2 2 _

-v

0(

P o o du (KY + E - u E4 F F

-v E

dv F F - 1 2

i

4

_3y

p uD

i i

C-1

-13--u E

0 3

(A-l)

(A-2)

F 2

(y1

yo)

y

P =r

_i

₁

dy1

V1

dx

-

u1

where

C = i +

(y-i)M2

D = F1 F5 - F2 F4

E =F F

_{1 2 SC}

-F F

S _3C

(A-i

E2=F4F

_3c

-F F

₁

_ec

E =F F

_{3 2 Gx}

-F F

_E

(A-L

E F F

-F F

4 4

3x

1

ex

F1 = (y1-y0)

[(Y1-Yo)(Y2-Yo)(Y2_Y1

(y2-yo)2

(yi-yo)2

+ 3

(y2-y1)

2

-

)

(YY0

Y1YO)

(A-(A-9

(A-1

(A-1

(A- L

(A-1

(A-1

) ) O)

1)

2a)

e)

d)

(A-1: b)

J

(y1 - y0)4 dA

T(A

-

A)(1

- y0)(y

- y0)(y

-F3=

_{6 d}

+ (y

y +

+ (y

_{- y0)2(A1 -}

A0)(Y2

- Yo

₃

j

dx

rq

(y1

- y2) q

-

(y1

- Y0)2

L+

3

1 c

+ (B1 - B0)(Y1 - Y0)(Y

- Y0)(Y2 - Y1)

(A-13c)

- y0)2

F 5

(y2

- y0)4

6

-15-y1-y0\

(ya-yo) (Y2-Yi

±

3

(A-13 d)

(A-13e)

y0)

(y1

_-

-r(y1

- y0)2 -

(y

Yo)L

2

-

y)2

[A2

- A0)(

_{- y0)}

-

- A0)(Y2 - YO)

- y0

dy

_2(y1

_{- y0)3(A}

_{- A )}

dy

O

3

(A1 -A)j

3

F4 = (y2

F 2 0 1 0

dy2

+ 2

_(y

- y

)3(A

- A)

dx

3

(y2 - y0)q2

+ 3

IG+(B-B)(y-y)(y-y)(y-y)

2 ₀

10

20

2

dy

o

G =-(y2 - y1)(y2 ± y1 - 2y)

dy1

+ (y2 - y)(y2 - 2y1 + y)

a-dy

2

+ (y1 - y)(2y2 - y1 - y0)

;-(A1 - A)(y2 - y)2 -

(A2 A(y1

-(y1 - y0)(y2 - y0)(y2 - y1)

Yo2

_{L 2}

(A-13f)

(A-14)

(A-15a)

6 =

-

-

- yo)2

(y1 - y0)

- y0)

dA

2

2 3

dx

±

-(A

- A)(y

- y0y2 - y0)(y

- y)

+ (y

-

y)2

[(A

- A)(y1 - y) - (A

- A)(y

- Yo)

'2

- y)(A - A)

y2 - y0 \ dy1

-

(y

2

-y)2(A -A)

o 2 0

l

-(

-o ₃

yo)

2 q 1

q

-(A - A )(y._o - y ) - (A - A )(y - y )

i o i o 2 o

(y1

- y0) (y - y0) (y

- y1)

The subscripts c, x and y appearing in F3 and _FB imply that F3 and

F are evaluated by using the continuity, x-momentum and y-momentum equations, respectively. Accordingly, A and B represent the following:

-17-(A - 15 b)

Equa t ion A B

continuity pu pv

x-momen turn KP + pu2 puy

dx V

=0

o 1 Ìc - u: Y

P =p

_{o o} APPENDIX B

EQUATIONS FOR EVALUATION OF INITIAL VALUES USING TWO-STRIP APPROXIMATION

The equations for evaluation of initial values using two-strip approximation are simílar to those in Appendix A except that here

y, y

and dy/dx are identically zero and the equation along line

"Ofl is replaced by the stagnation streamline equation, Equation (B-l). With strip boundaries O, 1 and 2 where line "O' corresponds to the

stagnation streamline, the equations are as follows:

2 du KY

+ U

₎ E2 - 1 E4 ' i P p Ky - u2 D

1\

P1

1'

dv F F

-F F

-v E

4 3)7 1 S\J 1 2 dx p u D (B-2) (B-5) (B-6)

or

p-2 2 C - u1 - V1 dy1 V1 dx u1

All the symbols here have the same meaning as those in Appendix _A. The term

(òvI&y)

appearing in Equations (B-1), (B-5) and (B-6) is the transverse velocity gradient along the stagnation streamline. It may be determined within the approximation of the method of integral relations. For a two-strip approximation, the

integrand of y-momentum equation is represented by

Aa +a(y-y)+a(y-y)2

o i o o puy = p u y + a1

(y - y ) +

a _(y 2 3 y

-

y3)

o o o

=a

y+a

y2

ly

2J

Differentiate with respect to y to obtain

(puy) a

+2a

_y 'y 2y (B-11) since (pu) (puv) -

_{+ pu}

- -19-(B-7) (B-12) = ply (B-8)

Combining Equations (B-11) and (B-12) and evaluating at the stagnation streamline line (where y = O and y = O for axisymmetric case), there

resu Its a ly _(B-13)

pu

₀₀

re a = 'y o 2 2

puvy

₁₁₁₂

- puvy

₂₂₂₁

y y (y -

V )

12

2

-i

(B-14)

APPENDIX C

TWO-DIÌNSIONAL INVISCID FLOW EQUATIONS IN A BODY ORIENTED ORTHOGONAL CURVILINEAR COORDINATE SYSTEM

PARTIAL DIFFERENTIAL EQUATIONS Continuity (pV5)

ò [(1

+ )

v]

= o (C-1)

+-S-Momentum PV V 2 pV5 + KP) +

₁₍₁

+ PV y s n R)

sn

R n-Momentum ò -

(pV5V) +

Energy P = p''

ORDINARY DIFFERENTIAL EQUATIONS REDUCED WITH AN ( + KP) n -21-2

PV +KP

S R (C-2)

ONE - STRIP APPROXIMATION

dV (PbVSb + p V )

âs

idó

d(öSô)

+ ds (C-5) ds b -dô ds - tan(>-e) (C-6) 2

C-Vs

i (C-7)

C-1

= _(C-8)

Outer stagnation region where cartesian coordinates are used as in other regions Streamline used as Outermost strip boundary (free-stream condition) -'-I

"1

_{'-r, /}

X

Figure 1 - Coordinate Systems

arid Strip Boundaries

Body oriented orthogonal coordinate system for inner stagnation region

z

Line Ç Free-stream boundary

This is a three-two-strip integration

scheme. The solution along line "Ti"

of set (, _fl, Ç) provides the boundary

condition for set (a, b, c);

that along line "b" of Set (a, b, c), for set (0, 1, 2).

Line _'fl, c

Line b, 2

Line i

Airfoil Upper Surface

Line , a, O

-23-X

,

/

Integrand distribution of a physical flow

/

/'

/ Three-strip approximation I / / / I / / / One-strip approximation

1/

Two-two-strip Approximation Sec tio ri A

Transverse Distance From Axis of Symmetry, y

Figure 3

- Accuracy of Different Approximations for the Method of Integral Relations

AirfoiL X

-7.0 6.0 5.0 2.0 1.0 o 0.04

Q

Q

8

Q

Velocity Components -25-NACA 0012 Airfoil Sonic Point

Figure 4 - Velocity Gradients Along Surface of an NACA 0012 Airfoil at M = 0.75

0.06 0.08 0.10 0.12 0.14

0.9 0.8

y=1.1

2.2 0.7 0.6 0.5 0.4 1.0 Trailing Ed ge

fi-Free-stream Mach number

-y =

O _{Normalized Distance From Airfoil Leading Edge, x}

The flow field in the far downstream returns to its free-stream condition if the assumed shock location is "correct".

Figure 5a - Downstream Flow Field of an NACA 0012 Airfoil at M

= 0.75,

for a "Correct" Shock Location

/

4.0

Airfoil

X

2.0

0.9 0.8 y = 1.1

-

---y2.2

, 0.7 0.6 0.5 0.4

Free-stream Mach number

X

The flow field in the far downstream does not return to its free-stream condition if the assumed shock location is "wrong.

Figure 5b - Downstream Flow Field of an NACA 0012 Airfoil at M

= 0.75,

for a "wrong" Shock Location

1.0

2.0

3.0

4.0

Trailing

Normalized Distance From Airfoil Leading Edge, x

Adjustment of the height of the outermost s trip boundary

Outermost strip boundary

(free-stream condition)

Iteration I

Iteration II

make s

Continuous solution at sonic point

make s

Adjustment of shock location

Downs t ream flow converged to free-stream condition

Figure 6

- Iterative Processes

for a Supercritical Flow

1.8 1.6 1.4 1.2 1.0 0.8 0.6 0.4 O I I I

Method of Integral Relations (Present Method)

O

One-strip

A

Two-two-strip

Unsteady Finite Difference Scheme (Magnus, Gallaher and Yoshihara)

Coarse mesh Fine mesh

I

\

L

a a

A

A

O

\

O

w

Airfoil -Pressure at sonic point O O

-A

_O

Normalized Distance From Airfoil Leading Edge, x

Figure 7 _{- Pressure Distribution on a Symmetric 8.4 Percent} Circular-Arc Airfoil at M = 0.85

I I I t

Method of Integral Relations (Present Method)

'Correct" Shock Location "Wrong" Shock Location

Experiment, Wind Tunnel (NACA TN 2174)

Experiment, Free Air, N = 0.74

£

(NACA TN 1746)

NACA 0012

Pressure at Sonic Point

Figure 8 - Pressure Distribution on an NACA 0012 Airfoil at M = 0.75

0 0.2 0.4

O6

08

1.0

Normalized Distance From Airfoil Leading Edge, x

1.8 1.6 1.4 1.2 Q) -o Q) N 1.0 o

z

0.8 0.6 0.4

1.0 0.8 0.6 0.4 0.2 1.2 1.0 0.8 0.2 0.4 0.6

Normalized Distance From Airfoil Leading Edge,

x

Figure 9

- Mach Number Contours for NACA 0012 Airfoil at

= 0.75

REFERENCES

Spreiter, John R. and Alberta Y. Alksne. _{Theoretical Predictions of} Pressure Distribution on Nonlifting Airfoils at High Subsonic Speeds. Wash., 1955. 42 p. mcl. illus. _{(National Advisory Committee for}

Aeronautics. Rpt 1217)

Rubbert, Paul E. and Marten T. Landahi. Solution of the Transonic Airfoil Problem Through Parametric Differentiation. ALAA Journal

(N.Y.), y. 5, Mar 1967. p. 470-479.

Crown, J. C. Calculation of Transonic Flow Over Thick Airfoils by Integral Method. AIAA Journal (N.Y.), y. 6, Mar 1968. p. 413-423 Magnus, R., W. Gallaher and H. Yoshihara. Inviscid Supercritical

Airfoil Theory. In Transonic Aerodynamics. _{[Paris, 1969] Paper 3.} (Advisory Group for Aerospace Research and Development. Conference Proceedings 35)

Dorodnitsyn, A. A. A Contribution to the Solution of Mixed Problems of Transonic Aerodynamics. In International Congress for the

Aeronautical Sciences. 1st, Madrid, 1958. Advances in Aeronautical Sciences, y. 2; Proceedings. N.Y., Pergamon Press, 1959.

Pearcey, H. H. The Aerodynamic Design of Section Shapes for Swept

Wings. In International Congress for the Aeronautical Sciences. 2d, Zurich, 1960. Proceedings. _{N.Y., Macmíllan [1962] p. 277-322}

Tai, T. C. Nonequilibrium Flow Over Blunt Bodies Using Method of Integral Relations. In 8th U.S. Navy Symposium on Aeroballistics, Corona, Calif., May 1969. Proceedings, Vol. 2. China Lake, Calif., NWC, Jun 1969. p. 267-309.

Liepmann, Hans Wolfgang and A. Roshko Elements of Gasdynamics. N.Y.,

Wiley [l957i 439 p.

Amick, James L. Comparison of the Experimental Pressure Distribution on an NACA 0012 Profile at High Speeds With That Calculated by the Relaxation Method. Wash., Aug 1950. 9 p. mcl. illus. (National Advisory Committee for Aeronautics. TN 2174)

Emtnons, Howard W. Flow of a Compressible Fluid Past a Symmetrical Airfoil in a Wind Tunnel and Free Air. Wash., Nov 1948. 31 p.

Holder, D. W. The Transonic Flow Past Two-Dimensional Aerofoils. Royal Aeronautical Society., Journal (London) y. 68, Aug 1964, p. 501-516.

Magnus, R. and H. Yoshihara. Inviscid Transonic Flow Over Airfoils. N.Y., Jan 1970. 8 p. mcl. illus. _{(American Institute of}

Aero-nautics and AstroAero-nautics. Paper 70-47)

-33-1 (O3A) 1 (03C) 1 (320) 1 (32OC) 1 (5108) 1 (5323) 1 (530.14) 1 (531) 2 (604)

1 Chief of Naval Operations i Naval Ordnance Systems Command i Office of Naval Research

1 Naval Air Development Center 1 Naval Air Engineering Center 1 Naval Air Test Center

1 Naval Avionics Facility Indianapolis, Indiana 1 Naval Missile Center

Point Magu, California i Naval Ordnance Laboratory 1 Naval Postgraduate School i Naval Research Laboratory i Naval Weapons Center 1 Naval Weapons Laboratory

20 Defense Documentation Center 1 Institute for Defense Analysis

1 Office, Chief of R & D (ABÌA)

1 Army Missile Command Redstone Arsenal i Hdq., USAF (AFMINDE)

1 Aerospace Research Laboratories Dayton, Ohio

I Hdq., Arnold Eng. Development Ctr.

1 Office of Research Analysis (OAR) Holloman Air Force Base

3 _{Hdq., Washington, D. C.} 1 Flight Research Center 1 Goddard Space Flight Center i Langley Research Center i Lewis Research Center

1 Marshall Space Flight Center 3 Clemson University

Clemson, South Carolina 1 Dept. of Mech. Engr. 1 Dr. T. Yang

1 Library

1 The John Hopkins Univer5 Applied Physics Laborat Silver Spring, Md.

i University of Maryland College Park, Maryland 2 _{North Carolina State Uni}

Raleigh, North Carolina i Dept. of Mech-Aerosar 1 Dr. F. R. DeJarnette 2 _{Virginia Polytechnic Ins}

Blacksburg, Va. i _{Dept. of Aerospace Er} i Library INITIAL DISTRIBUTION Copies Copies 10 NAVAIR 8 NASA ity o ry versity ce Engr. titute gr.

UNCLASSIFIED

iii rus (!.IssufR atti

D D

FORM 1473

NOV 65 I (PAGE 1) _UNCLASSIFIED

DOCUMENT CONTROL DATA - R & D

Surtir, t i 1.1551I taliOn of tutte. lmdy ci f nhtrrnct attui uuIdcxuuuc) annotatisln nu, be entered when t!Ie overa!! report is claSsified)

ONC!NA TINO AC liv! lv (C.rpuur,uiC Ouithiurl Department of Aerodynamics

Naval Ship Research and Development Center

Washington, D, C. 20034

2a, REPORT SECURITY CLASSIFICATION Unclassified 2h. GROUP

REPORT TITLE

APPLICATION OF T} THOD OF INTEGRAL RELATIONS (MIR) TO TRANSONIC AIRFOIL PROBLEMS. PART I - LNVISCID SUPERCR.ITICAL FLOW OVER SYfl'TRIC AIRFOILS

AT ZERO ANCLE OF ATTACK

A DISC RIPrIVE NOTES (Type of report and inclus iv. dates)

Research and Development Report

S AU T i-IOR,SI çFrsst name. middle Initia!, ta.! n8me) Tsze C. Tai

E REPORT DATI

September 1970

la. TOTAL NO OF PAGES lb. NO OF REES

12 6A CONTRACT OR GRANT NO WR 009.0201 b. PROJECT NO e. TASK 10202 d. NSRDC 645-277

55, OPIGINAORS REPORT NUMBER(S)

NSRDC Report 3424 Aero Report 1176

9. O To ER RE PORT ROIS! (Any other numbers that may be assigned

this report)

10 DISTRIBUTION STATEMENT

This document has been approved for public release and sale;

its distribution is unlimited.

li SUPPLEMENTARY NOTES I2 SPONSORING MILITARy ACT! VITV

Commander

Naval Air Systems Command

Washington, D, C, 20360

The feasibility of applying the method of integral relations (MIR) _to tran-sonic flows over symmetric airfoils has been studied. _{In order to take account of the} severe transversal flow variation and still retain simplicity in computation, the method is modified so that the number of strips used may be considerably _{higher than}

the order of the polynomial which approximates the integrand. _{An important feature of} thís modification, however, is its capability to extend the free-stream _{boundary to}

'infinity't. Using one-strip and two-strip approximations, flow equations _{are reduced} to a set of ordinary differential equations in a cartesian coordinate system.

Numeri-cal procedures, including the treatment of the sonic point and the _{determination of}

shock location, are also formulated. _{The fourth-order Runge-Kutta method was used in} numerical computation.

The advantage of small computer capacity and time required by the method is evidenced by exploratory calculations for a symmetric circular-arc airfoil and an NACA 0012 airfoil, traveling at supercritical speeds. _{The calculated results compare fairly} well with those obtained by the unsteady finite difference scheme for _{the case of a cir} cular-arc airfoil and very well with experimental data for that of _{an NACA 0012 airfoil}

In principle, the results may be improved by increasing the number of strips and extending the free-stream boundary condition sufficiently far _{away from the}

air-foil. _{It is also found that in transonic airfoil problems using a cartesian coordinate}

system, the order of polynomial approxinting the integrand in the method of integral relations should be at least second-order or higher.

JNCLASSIFID ',urttv Classi 1iatiin

DDFORM 1473 (BACK)

I NOV C5 UNClASSIFIED

KEY WORDS LINi< A IINK P I.. INK C

POI_E WT POLE WI POLE WI

Integral Relations Transonic Airfoils