Category AERODYNAMICS

The representation of a continuous singularity distribution by a series of constant-strength elements results in a discontinuity of the singularity strength at the panel edges. To overcome this problem, a linearly varying strength singularity element can be used. The requirement that the strength of the singularity remains the same at the edge of two neighbor elements results in an additional equation. Therefore, with this type of element, for N collocation points 2N equations will be formed (see examples in Chapter 11).

Consider a doublet distribution along the x axis consisting of elements pointing in the z direction [ц = (0, /і)], as shown in Fig. 10.7. The influence at a point P(x, z) is an integral of the influences of the point elements between x,-*x2:

FIGURE 10.7 Constant-strength along the x axis.
doublet	distribution
(10.26)
(10.27)
(10.28)
Comparing this expression to the potential of a point vortex (Eq. (10.8)) indicates that this constant doublet distribution is equivalent to two point vortices with opposite sign at the panel edges such that Г = fi (see Fig. 10.8). Consequently, the velocity components are readily available by using Eqs. (10.9) and (10.10): u=±———- ?] (Ю.29)
(10.30)
FIGURE 10.8 Equivalence between a constant-strength doublet panel and two point vortices at the edge of the panel.

When the point P is on the element (z = 0±, хг<х <x2) then following Section 3.14

Ф(х, 0±) = T^

and the velocity components become

and the w velocity component is singular at the panel edges.

Details of the integral for the velocity potential appear in Appendix В (Eq. (B.10)) and in terms of the distances and angles of Eqs. (10.15) and

FIGURE 10.9

Constant strength vortex distribution along the x axis.

(10.16) (as shown in Fig. 10.6) it becomes

Ф=~2ж [^_ДСі)01_^_д:2^2 + ^1п^] (10.37)

and in terms of the x-z coordinates:

^ У Г, . Z. . _. z z, (x-x^f + z2!

Ф = ~2л l * _ tan 7=^-(Ar""2)tan ^ + 2,n0c-x2)2 + z2J

(10.38)

Following the formulation used for the constant source element, and observing that the и and w velocity components for the vortex distribution are the same as the corresponding w and u components of the source distribution, these equations become

The influence of the element on itself at z = 0± and (д:г <x <x2) can be found by approaching from above the x axis. In this case вг = 0, в2 = я and

Ф(*> 0+) = К* – *i)0 — (де — х2)я] – % (x – x2) (10.41a)

Similarly, when the element is approached from below, then

Ф(х,0-) = -^(х-х2) (10.416)

The x component of the velocity can be found by observing Eq. (10.24) for the source or by recalling Section 3.14,

n(x,0±) = ±^ (10.42)

and the w velocity component is similar to the и component of the source (Eq.

(10.23) ):

K4*;0±)=£lnji^ (10.43)

In most situations the influence is sought at the center of the element where rx = r2 and consequently w(panel – center, 0±) = 0.

Consider a source distribution along the x axis as shown in Fig. 10.5. It is assumed that the source strength per length is constant such that o(x) = a — const. The influence of this distribution at a point P is an integral of the influences of the point elements (described in the previous section) along the segment Xi~*x2.

The integral for the velocity potential (Eq. (10.11)) appears in Appendix B, Eq. (B.5), (note In r2 = 2 In r is used in the derivation) and in terms of the

• p

U. z)

Returning to x—z variables these equations become

Ф = 4jt ~ln “Xl)>2 + ~(x~*2),n К* – xz)2 + z2]

Of particular interest is the case when the point P is on the element (usually at the center). In this case z = 0± and the potential becomes

Ф(х, 0±) = [(x – Xi) ln (x – x,)2 – (x – x2) ln (x – x2)2] (10.22)

and at the center of the element it becomes

ф(£іті-г’0±)-й(Іг-^1п(і1?;і)2 (Щ22л)

The x component of the velocity at z = 0 becomes

“(j:’0±)’5lni(f^)i (ШЗ>

which is zero at the panel center and infinite at the panel edges.

For evaluating the w component of the velocity, it is important to distinguish between the conditions when the panel is approached from its upper or from its lower side. For the case of P being above the panel 0, = 0, while в2 = л. These conditions are reversed on the lower side and therefore

w(x, 0±) = ±^

This is the same result obtained in Section 3.14 for the source distribution.

The discretization of the source, doublet, or vortex distributions in the previous section led to discrete singularity elements that are clearly not a continuous surface representation. A more refined representation of these singularity element distributions can be obtained by dividing the solid surface boundary into elements (panels). One such element is shown schematically in Fig. 10.4, and both the surface shape and the shape of the singularity strength distribution are approximated by a polynomial. In this section, for the surface representation, a straight line will be used. For the singularity strength, only the constant, linearly varying, and quadratically varying strength cases are presented, but the methodology of this section can be applied to higher-order elements.

In this section, too, three examples will be presented (source, doublet, and vortex) for evaluating the influence of the generic panel of Fig. 10.4 at an arbitrary point P. For simplicity, the formulation is derived in a panel-attached coordinate system, and the results need to be transformed back into the global coordinate system of the problem.

Consider a point vortex with the strength Г located at (x0, z0). Again using the definitions of the points, as in Fig. 10.2, and the results of Section 3-Я the increment to the velocity potential at a point P is

(10.8)

and the increments in the velocity components are

__T_____ Z ~Zq

u 2л(х – x0)2 + (z – 2o)2 – Г X – x0

Ш —- .1. —P —и

2я (x – x0)2 + (z – Zo)2

Note that all these point elements fulfill the requirements presented in Fig. 10.1. That is, the increments of the velocity components and potential at P depend on the geometry (x, z, x0, Zq) and the strength of the element.

These elements are probably the simplest and easiest to use and also the most efficient in terms of computational effort. Consequently, even when higher – order elements are used, if the point of interest is considered to be far from the element, then point elements can be used to describe the “far-field” effect. The three point elements that will be discussed are source, doublet, and vortex, and their formulation is given in the following sections.

10.1.1 Two-Dimensional Point Source

Consider a point source singularity at (x0, z0), with a strength a, as shown in Fig. 10.2. The increment to the velocity potential at a point P (following Section 3.7) is then

ф(*> z) = TIln V(* ~ *<>)* + (z – z0)2 (10.1)

and after differentiation of the potential, the velocity component increments are

ЭФ a x — x0

Эх 2л (x – x0)2 + (z – z0)2 ЭФ a z — Zq

dz 2л (x – x0)2 + (z – Zq)2

10.1.2 Two-Dimensional Point Doublet

Consider a doublet that is oriented in the z direction [ц = (0, /і)] as in Section 3.7. If the doublet is located at the point (x0, z0), then its incremental influence on the velocity potential at point P (Fig. 10.2) is

FIGURE 10.2 FIGURE 10.3

The influence of a point singularity element Transformation from panel to global coord-

at point P. inate system.

and the velocity component increments are

__дф _jl (x ~ Xp)(z – Zq)

Эх n[(x-x0)2 + (z-Zo)2]2 _ Эф (x ~ Xpf – (z – Zo)2

dz 2я [(x – x0)2 + (z – z0)2]2

In the case when the basic singularity element is given in a system (x*, z*), which is rotated by the argle /3 relative to the (x, z) system, as shown in Fig. 10.3, then the velocity components can be found by the transformation

It was demonstrated in the previous chapters that the solution of potential flow problems over bodies and wings can be obtained by the distribution of elementary solutions. The strengths of these elementary solutions of Laplace’s equation are obtained by enforcing the zero normal flow condition on the solid boundaries. The steps toward a numerical solution of this boundary-value problem are described schematically in Section 9.7. In general, when increas­ing the complexity of the method, mostly the “element’s influence” calculation becomes more elaborate. Therefore, in this chapter, emphasis is placed on presenting some of the typical numerical elements upon which some numerical solutions are based (the list is not complete and an infinite number of elements can be developed). A generic element is shown schematically in Fig. 10.1, and it requires the information on the element geometry, and strength of singularity, in order to calculate the induced potential and velocity increments at an arbitrary point P (xP, yP, zP).

xP, yP, zP Panel geometry Singularity strength

FIGURE 10.1

Schematic description of a generic panel influence coefficient calculation.

For simplicity, the symbol A is dropped in the following description of the singularity elements. However, it must be clear that the values of the velocity potential and velocity components are incremental values and can be added up according to the principle of superposition.

In the following sections some two-dimensional elements will be pre­sented, whose derivation is rather simple. Three-dimensional elements will be presented later and their complexity increases with the order of the polynomial approximation of the singularity strength. Also, the formulation is derived in the panel frame of reference and when these formulas are used in any other “global coordinate system,” the corresponding coordinate transformations must be used (for rotations and translations).

The thin wake behind lifting wings tends to follow the local velocity (Eq. (9.18)) which results in a wake rollup due to the lifting surface and the wake-induced velocity. This condition causes the shape of the wake to be unknown when the boundary conditions for the potential flow are established. Traditionally, the shape of the wake is assumed to be known (e. g., planar vortex sheet) and after the solution is obtained the validity of the initial wake shape can be rechecked. In Chapter 14 two methods used by two panel codes will be presented to calculate the wake shape (VSAERO-wake relaxation,92 and PMARC-time stepping9 6’9’7).

Since the wake is modeled by a doublet/vortex distribution, it is possible to extend this method for modeling jets and even shear layers of separated flows. Details about models to represent the effect of jets can be found in Refs 9.2 and 9.6, while models to treat some effects of flow separation will be presented in Chapter 14.

REFERENCES

9.1. Lamb, H., Hydrodynamics, Dover Publications, 1945.

9.2. Maskew, B., “Program VSAERO Theory Document,” NASA CR 4023, Sept. 1987.

9.3. Johnson, F. T., “A General Panel Method for the Analysis and Design of Arbitrary Configurations in Incompressible Flows,” NASA CR 3079, May, 1980.

9.4. Gerald, F. C., Applied Numerical Analysis, 2d edn., Addison-Wesley, Reading, Mass., 1980.

9.5. Bristow, D. R., and Hawk, J. D., “Subsonic Panel Method for Designing Wing Surfaces from Pressure Distributions,” NASA CR-3713, July 1983.

9.6. Ashby, D., Dudley, M., and Iguchi, S., “Development and Validation of an Advanced Low-Order Panel Method,” NASA TM 101024, 1988.

9.7. Ashby, D., Dudley, M., Iguchi, S., Browne, L., and Katz, J., “Potential Flow Theory and Operation Guide for the Panel Code PMARC,” NASA TM 102851, March 1990.

PROBLEMS

9.1. Solve the problem of a flat plate at an angle of attack of a using the lumped-vortex element. Divide the chord to five equal panels of length a, as shown in Fig. 9.19.

(a) Calculate the influence coefficient matrix a, t. Is this is a diagonally dominant matrix?

(b)

Calculate the lift and moment coefficients. How do these compare with the analytical results of Chapter 5.

FIGURE 9.19

Discrete vortex model for the flat plate at angle of attack.

9.2. Calculate the lift and moment coefficient (about the origin, x = 0) of the two flat plates shown in Fig. 9.20. Use a single-element lumped-vortex model to represent each plate and investigate the effect of the distance between the two plates on their lift (repeat your calculation with gap values of c/2, c, 2c, 4c).

Z(l

Г r2

The potential flow model presented in this chapter results in a very simple mathematical model that can be transformed into a very efficient and economical numerical solution. This resulted in the development of three – dimensional “panel codes” for arbitrary geometries and, naturally, modifica­tions were sought to improve these methods beyond the limits of incompres­sible inviscid flows. Some of these modifications are listed here.

Effects of Compressibility

The first and most straightforward modification to an incompressible potential – flow based method is to incorporate the effects of “low-speed compressibility” (e. g., for M <0.6). This can be obtained by using the Prandtl-Glauert rule, as developed in Section 4.8. Thus, small-disturbance flow is assumed, and a compressibility factor f} can be defined as

j8 = Vl – M2 (9.44)

If the free stream is parallel to the x coordinate then the x coordinate is being stretched with increased Mach number while the у and z coordinates remain unchanged. Consequently an equivalent incompressible potential ФА/=0 can be defined such that

(9.45)

Once the x coordinate is transformed, the equivalent incompressible potential problem is solved as described previously. This procedure results in an increase in the aerodynamic forces (as noted in Section 4.8) and

Effects of Thin Boundary Layers

When analyzing high Reynolds number flows in Section 1.8, it was assumed that the boundary layer is thin, and the boundary conditions are specified on the actual surface of the body. However, by neglecting the viscosity terms in the momentum equation, the information for calculating the viscous surface friction drag is lost too.

It is possible to account for the viscosity effects such as displacement thickness and friction drag by using a boundary-layer solution that can be matched with the potential flow solution. Two of the most common methods for combining these two solutions are as follows.

1. The first approach is to use a boundary-layer solution, usually a two – dimensional model along a streamline, which will work for simple wings and bodies. The solution begins by solving the inviscid potential flow, which results in the velocity field and the pressure distribution. This data is fed into two-dimensional boundary-layer solutions that provide the local wall friction coefficient and the boundary-layer thickness. The friction coefficient can be then integrated over the body surface for computing the friction drag. If the displacement thickness effect is sought, then a second iteration of the potential flow computation is needed, but now with modified surface geometry. This
modification can be obtained by displacing the body panels according to the local boundary-layer displacement, and the procedure can be reiterated until a satisfactory solution is obtained. Some of the principles of a computer program (e. g., the MCAIR panel code) that uses this method are presented in Ref. 9.5.

2. The second approach to incorporate boundary layer solutions into panel codes is to follow the procedure described above, but to account for the displacement effects by a modification of the boundary conditions instead of a change of the surface geometry. In this case, at each panel the normal flow is given a certain blowing value that accounts for the local boundary-layer displacement thickness d*. The formulation can be derived, using the results of Section 4.4, as follows:

Д Oi=^{qb*) (9.48)

ds

where q is the local streamwise velocity component of the potential flow (outside the boundary layer) and the differentiation takes place along a streamline s. The transpiration velocity then becomes qn. = Да,. For more details on this approach see Ref. 9.2.

As a first example for demonstrating the principle of numerical solutions, let’s consider the solution for the symmetric, thin airfoil. Because the airfoil is thin, no sources will be used, while the doublet distribution will be approximated by two constant-strength doublet elements (/r1( ju2, pointing in the – z direction). This element is equivalent to two concentrated vortices at the panel edges (see Fig. 9.16). However, the geometry of the “lumped-vortex” model was developed in Chapter 5, and by placing the vortex at the quarter chord and the collocation point at the three-quarter chord point of the panel the Kutta condition is automatically satisfied. Using this knowledge the equivalent discrete vortex model (with only two elements) for the airfoil is shown in Fig. 9.17. Also for the thin lifting surface only the Neumann (velocity) boundary condition can be used, because of the zero thickness of the airfoil. (Note that the doublet representation in Fig. 9.16 clearly indicates the existence of a

starting vortex at a large distance behind the airfoil which is shown in Fig. 9.17 too.)

Selection of singularity element. For this very simple example the lumped – vortex element is selected and its influence is derived in terms of the geometry involved. Such an element is described in Fig. 9.18л, and it consists of a concentrated vortex at the panel quarter chord, a collocation point, and a normal vector n at the three-quarter chord. It is important to remember that this element is a simplification of the two-dimensional continuous solution and therefore accounts for the Kutta condition at the trailing edge of the airfoil.

If the vortex of circulation Г of the element is at (jc0, Zo), then at an arbitrary point P(x, z), according to Section 3.8, the velocity induced by this element will be

Г (г-2b)

U 2л (x — x0)2 + (z — Zo)2 (лг-лср)

W 2л (x – x0f + (z – Zo)2

and in a matrix form that is more useful for computations

where

This can be programmed as an influence coefficient subroutine in the manner shown in Fig. 9.186. Let us call this routine VOR2D and an algorithm based on Eq. (9.31) will have the form

(и, w) = VOR2D (Г, x, z, *o, Zq) (9.32)

Discretization of geometry and grid generation. For this example, the thin airfoil case is being solved (Fig. 9.17). For simplicity, only two elements will be used so that no computations are necessary. At this phase the geometrical information on the grid has to be derived. This can be automated by computer routines for more complex situations, but for this case the vortex point locations are

(•*oi. *o,) = (c/8, 0) (*o2, Z02) = (5c/8, 0)

and the collocation points are

(*ci, Zci) = (3c/8, 0) (*e2, zc2) = (7c/8, 0)

The normal vectors n must be evaluated at the collocation points, and for an arbitrary element і are
where the angle /3, is defined in Fig. 9.18a. In this particular case, when the airfoil has no camber and is placed on the z = 0 plane, both normals are identical:

ІІ! = (0, 1) n2 = (0,1)

Influence coefficients. Here the condition requiring zero velocity normal to the airfoil will be enforced. This boundary condition, according to Eq. (9.4), is

(q + Q«) • n = 0 (9.34)

The velocity q is induced by the unknown vortices, while the free-stream normal component can be calculated directly, and hence is moved to the right-hand side of the equation:

q • n = Q, • n (9.34a)

Because, in this case, the airfoil was divided into two elements with two unknown vortices of circulation Г1( T2, two equations based on the zero flow across the airfoil boundary condition will be derived at the collocation points. Defining as positive Г a clockwise rotation, the velocity induced by a unit strength vortex at point 1 on collocation point 1 is calculated with Eq. (9.32):

(«11, vt’n) = VOR2D (1.0, *cl, zcl, xm, z01) = (o, – ^ ^—j

and the velocity induced at collocation point 1, by a unit vortex at point 2, is

(«12, и’и) = VOR2D (1.0, xcl, zcl, x^, zm) = (o,

The velocity induced at collocation point 2, by a unit vortex at point 1 is

(«21, w2i) = VOR2D(1.0,xc2, zc2,x0l, Zoi) = (о, ~2д. .зс/4)

and the velocity induced at collocation point 2, by a unit vortex at point 2 is

(«22, *v22) = VOR2D (1.0, xc2, zc2, jc02, z02) = (o, – — *

The influence coefficients a,; are really the normal component of the flow velocity induced by a unit strength vortex element Г; at collocation point і:

a,/ = Ч, ДГ; = 1) • n, (9.35)

For the current problem, Eq. (9.35) is applied to collocation point 1 and to vortex point 1, thus

-2

Similarly for the second vortex,

0l2 = (Ml2, W’izXO, 1) = —

and for the second collocation point,

2

az = (M2i. И’гіХО» 1) = —

2

a22 = (U22> И>22)(0) 1) =

ЛС

Note that the left-hand side of Eq. (9.34a) can be described now as

q • n = 2 «оГ, for i = 1, 2 (9.36)

;=i

Establish RHS. The solution is based on enforcing the boundary condition of Eq. (9.34a) at the collocation points. Since the product Q„ • n is known it is transferred to the right-hand side of the equation:

q • n = – Q* • n = RHS (9.37)

It is useful to express the component of the free stream in vector form to allow easy vector operations; for this particular case the right-hand side (RHS) is

RHS, = -(l/., W.) • n, (9.38)

where (t/„, W„) = Q„(cos a, sin a). Computing the RHS vector for the two collocation points results in

RHSj = – gooSin a

RHS2= – Q„sin a

Solve linear set of equations. The results of the previous computations can be summarized as

І a, jTi = RHS,

/=і

and explicitly, for this particular case,

/ —2/яс 2/яс /ГЛ -2/Зяс —2/яс/Г2/

which can be solved by standard matrix methods

Secondary computations: pressures, loads, velocities, etc. The solution of the

above set of algebraic equations is

The resulting pressures and loads can be computed by using the Kutta – Joukowski theorem (Section 3.11):

A U = pQSi

and by assuming a constant pressure distribution along the element, the pressure difference becomes

Ap, = pQSja

where a is the panel length. The lift and moment are then

2

L = 2 A A, = рОІ, лс sin a (9.40)

i=i

2 c2

M = — 2 A L/Х/ = – pj22jr—sin a (9.41)

i = l **

These results are similar to those for a zero-thickness symmetrical airfoil (Section 5.4) and equal to the exact flat plate solution (Section 6.5.1) and the method can easily be expanded to various camberline shapes and even multielement lifting airfoils.

Description of more complex methods for solving the potential flow problem will be presented in the following chapters.