A method for the numerical evaluation of an integral: Occuring in the expression for the wave drag of slender bodies

Co A RsDort No. ICO

Kanaalitraal 10 - DELFT

- 3 AUG. 1956

THE COLLEGE OF A E R O N A U T I C S

C R A N F I E L D

A METHOD FOR THE NUMERICAL EVALUATION

OF AN INTEGRAL

(occurring in the expression for the wave drag of slender bodies)

by

REPORT MD. 100 mBCE, 1956. T H E C O L L E G E O F A E R O N A U T I C S C R A N F I E L D A Method f o r Numerical E v a l u a t i o n of t h e I n t e g r a l 01 01 F ( x , y ) l o g Ü _{o <-' o} 1 x - y dx dy - b y - . T, Nonv/eiler, B . S c , SUIIIARY

This report presents formulae ard data for the num-erical evaluation of the double integi-al named in the title in the foiTa ^ ^ Zli ^ y ^ ( ^ »s ) , where numerical values for the Vveighting coefficients c and the lattice points (s ,3 ) are given. The method depends on a double Fourier series representation of F(x,y) in terms 6 = cos (l-2x) and 0 = cos (l-2y); the lattice points s , s are in fact . a* equally spaced intervals in 6 and 0,

The method is particularly applicable to integration vvhere F has half-order singularities (or zeros) along the borders of the region of integrationj but it mcy also be

adequate for the treatment of continuous functions P which are finite at these borders, (an accuracy of \ïithin less than 1 per cent being then expected from the use of a 1 1 x 1 1 lattice). Comments on applications of the formula to the evEiluation of wave drag are given.

-2-CONTENTS

List of Symbols 1• Introduction

2 , Formulae for Numerical Integration

2.1 The function F(x,y) bounded over the complete range of integration

2.2 The function F(x,y) v/ith half order singularities at the botindaries of the region of integration 3» I.Iodification v/here integration limits are altered 4. Simplification resulting from particular symmetries in

P(x,y)

5 , Applications in Aerodynamics Ackn wledgement

References

Appendix? Derivation of Integration Formula

LIST OF SYlfflOLS P = P ( x , y ) f u n c t i o n of x , y i n t h e d o u b l e i n t e g r a l t o be e v a l u a t e d F c o e f f i c i e n t d e f i n e d by e q u a t i o n ( l ) of Appendix I ()Ö)- i n t e g r a l d e f i n e d by e q u a t i o n ( 8 ) of Appendix I t h e double i n t e g r a l t o b e e v a l u a t e d S(x) c r o s s - s e c t i o n a l a r e a of body a t f r a c t i o n x of i t s l e n g t h from n o s e b = 1 - - i f 6 + 6 , N ") ^ 2 \^ nfi ( - n ) n ;

c weighting coefficients of integration formulae ^ tabulated in section 2.1 for various |i,v and n, c . defined in section 3(m)

m ^iv -'\ /

d weighting coefficients of integration formulae ^ tabulated in section 2.2

n number specifying lattice spacing in integration

formula

1 (, . a%'\

x.,y

variables of integration

6 = 1 i f | i = v ; = 0 i f | i ? ^ v . The Kronecker delta

^(x) non-dimensional vdjig section ordinate at fraction

X of its chord from nose

Ö = cos""'' (l-2x)

|i,v nuiïibers characterising point of evaluation in

integration or interpolation scheme.

0 = cos""^ (l-2y)

'fV - 2 + 2n •

1• Introduction

A double integral of the type named in the title

appears for instance in expressions for the wave drag of

slender bodies and of svrept wings and presents some difficulty

in its evaluation, as the function F(x,y) is often only

known either by its numerical values at discrete points, or

else by a closed algebraic expression rendering formal

evalua-tion tedious, if indeed possible. For this reason numeidcal

quadrature is often a convenience, but the well-knovm formulae.

for the evaluation of double integrals (such as Vfeddle's Rule)

are inapplicable owing to the logarithmic singularity in the

integrand» ¥e present here a method for its approximate

numerical evaluation which depends essentially on the

rep-resentation of the function F(x,y) in terms of a finite

double Fourier series, which can be shovm to be equivalent

to its expression as a polynomial in x and y divided by

-4-This form of representation is evidently particularly suitable where P(x,y) has half-order singularities (or zero^ at the boundaries of the region of integration, which is some-times so, but admittedly V70uld introduce unnecessary errors where such conditions are not met. For this reason, the

quadrature formulae quoted are tested for a particular example of the latter category, where the function F is a constant over the range of integration. This probably represents the most severe test possible, implying as it does the represent-ation of (sin e) over a half period 0 < 0 < 7t by a finite cosine series,

Til ri\

0 o ly o v/here

The integration formulae are quoted in the form n-1 n-1

F(x,y) In dx dy -— j \ '

x-y I ^ |i=-n+1 v>=-n+1 %v ^^%>%^ 1 + sin ^ )

Using an 11 x 11 -point lattice (i,e, n = 6 ) , and putting P equal to a constant, the eri'Or is less than one per cent of the exact value. The numerical value is less than the tmae one? the precise error ;vith this and other lattice spacings is tabulated below.

n

1

2

3

4

5

6

l a t t i c e p o i n t s

1 X 1

3 x 3

5x5

7 x 7

9 X 9

11 X 11

e r r o r as

percentage of exact

value

32.71

6.ii5

2,91

1,66

1,07

0.75

For the reasons stated above, it is to be expected that accuracy in relation to functions F vïhich' have half-order

singularities or zeros on the borders of the region of

integration would be considerably better,'

Full details of the method of calculation and the

values of the weighting coefficients c are given in the

next section. The mathematical derivation of the formulae

is treated in the Appendix, The third section concerns the

evaluation of such integrals where the limits are not betv/een

zero and unity, and the next section concerns seme

simplific-ations resulting from symmetries in the function P(x,y).

Finally seme comments are made concerning the application of

the formulae to the evaluation of wave drag, by way of exaniple,

2. Formulae for Numerical Integration

The analysis of Appendix I leads us to distinguish

between two conditions n the behaviour of F(x,y), which is

assumed bounded everyvrhere inside the region of integration.

The distinctions arise from the presence or absence of

singu-larities in P(x,y) on the boundaries of the region of

integration.

2,1, The function F(x,y) bounded over the complete region

of integration

Equations (10) and (12) of Appendix I

show

that

JA1

r-1 , n-1 n-1

P(x,y)ln - ^ dx

dyzh. ..-^'^

. --'' v c F(s , s )

I I ^ i ' ^ ' x-y» ^ — -—--•. * :; liv ^ LL* v'

.) ; ^ |i=-n+1 v=-n+1

^ ^

I.: o -> o ^

where s =•=•11 + s i n ^ \ • The c o e f f i c i e n t s c

[Ji 2 \^ 2 n / p

are givei

for n = 2 ( l ) 6 . I t w i l l be observed t h a t

(together with numericEil values of s , s ) are given below

c . = c

-6-so that out of each array of (2n-l) coefficients, arranged 2

in a square, only n are different, and these appear in ea.ch of the four triangular parts of the array formed by and

including the diagonal elements,

It will also be seen that the lattice points s , s . are at equally spaced in terms of 6 and 0 where

0 = cos'''"(l-2x), 0 = cos"\l-2y), (O :$ d,0 ^ %) ;

so that it is convenient to suppose that F(x,y) = f(Q,0), say, in which event

F(s s ) - f/'^ + ^ ^ + ^ ^

^^%'^v^ - ï( 2 + 2n » 2 + 2n^

The lattice points corresponding to the values (i = +^ n, or V = +1^ n lie on the boundaries of the region of integration, ajid the values of P (and so of f) at these points are not used in the integration schemes. Given F in terms of numerical data at discrete points Tstó-ch do not coincide ïri.th the lattice points, its values required in the integration scheme would have to be found by interpolation? Y/hether or not this can more easily be accomplished using Ö and 0 as

by

independent variables, or/onterpolating Vidth respect to x and y, depends of covirse upon the nature of the derivation of the

C o e f f i c i e n t s c f o r n=2 "^v % 0.8536 0 . 5 0,1464 ^ 1 0 - 1 0.1464 -1 0.0137 0.1375 0.2193 0 . 5 0.8536 0 1 0.1375 0,2193 0.3873 0,1375 0,1375 0.0137 C o e f f i c i e n t s c f o r n=3 ^ 0,9330 0,75 0 , 5 0 . 2 5 0.0670 s p. 2 1 0 - 1 - 2 0.0670 - 2 0 , 0 0 0 8 O.OWi-0.0247 0.0555 0,0625 0 . 2 5 - 1 0.01 Vf 0.0297 0.0922 0.1599 0.0555 0 . 5 0 0.0247 0.0922 0.2035 0.0922 0.0247 0 . 7 5 1 0.0555 0.1599 0.0922 0.0297 0»0144 0.9330 2 0.0625 0.0555 0.0247 0.012^4 0,0008 C o e f f i c i e n t s c f o r n=4 % 0.9619 0.8536 0 . 6 9 1 3 0 . 5 O.3O87 0 , 1 4 6 4 O.038I s

3

2 1 0 -1 - 2

-3

O.038I

-3

0.00016 0.00266 0.00507 0.01227 0.01682 0.02452 0.02396 0.1464 - 2 0,00266 0,00566 0,01674 0,02643 0,04888 0,07028 0.02452 0.3087 0 . 5 - 1 0,00507 0.01674 0.02933 0.06520 0.11128 0.04888 0.01682 0 0,01227 0.02643 0.06320 0.12734 0.06520 0,02643 0.01227 0.6913 [ 0 , 8 5 3 6 1 0.01682 0.04888 0.11128 0.06320 0.02933 0.01674 0.00507 1

1

2 0.02452 0.07028 0.04888 0.02643 0.01674 0.00566 0.00266 0.9619

3

0.02396 0.02452 0.01682 0,01227 0,00507 0,00266 0.00016

8 -I UN Ü. ^

eg

> j l to ^ o • H Ü ' M

8

Ü UN UN ON • o •^ -d-O CTN . O CTN t o , ON r-• o ^ -d-UN MO • o UN . o i n • ^ -d-_N~\ . O 'T-* MO O CM . O UN UN ON O . o 1 ^ -d-CsJ O o o i co j—^ _ ^ l O , CM i ; — o T ~ 1 CM 1 ^n -d-1 i / / '^ / / •>• / 0^ / y ' ^ T * O . O 1 ^ CM T -o • o NN CyJ O T -O . o

3

co o o • o ^

'S

o . o 1 ^ UN l O , O • o

3

T —

g

, o (^ MO O Q , O

§

o o , o -d-UN U N r^ ON • O NN CM O , o CO !> J -m o • o MO ON MO CM o • o CM UN CO V O • o CO UN NN o • o MO MO r>-o o , o ON

8

• o ^-MO T -o o _e o ^ o O o • o NN ^ -d-_O ON • , o ^ lOs CM O O • o MO ON MO CM O , O CO r-o MO o , o o CM ON NN o . o ON UN NN CM O . O ^— 0 0 UN o . O MO

5

8

. o ON O o . o % ^— o o • o CM ON NN ON r-. o % CO

8

• o CM UN 0 0 V " o . o o CM ON NN O • o ^_ -d-_O co o • o CM o ^ O . o MO CM UN CM o • O T ~ cxo _UN

5

• o MO MO f ^ O O • O NN UN NN O O . O "> -d-UN MO . o ! >

8

. o 0 0 UN N ^ V -o . o ON U N NN CM o . O CM O ^ O • O r^ MO 0 0 O . O CM O Q , O ON UN NN

8

# o CO UN NN r -o • o

s

_Q O , o o UN • O NN UN NN

8

• o MO MO r-o o . o .^ 0 0 U-N T -o . o MO CM UN CM O . o CM ^ Q • O V --d" O co o . o o CM ON NN O . O CM U N CO •\— o • o

eg

cO

8

• o T -UN ^ -d-_NN . „ Q

-3

T —

8

. o CTN Q O • o M>

5

_o o . o T -co LTV T -O • O ON UN NN CM O . o O CM CTN NN O . o CO r^ o o . o MO ON MO CM O • O NN CM O O • O CM 1

Cö

_o CM . O MO o

8

• O T -MD T " O O • o CTN

:5

8

• o MO MO r^

8

. o co UN NN o « o CM "^ CO %*~ O • o MO ON MO

s

. o 0 0 r>--d-NN o , o NN CM «^ o , o NN 1 UN U N CTN O • O Po 1 O

8

, o f ^ M5

8

O 1 _{, I} O

3

• ' S " 1

8

, o NN UN NN

8

_, o •^

©

o , o

3

co o o . o •*•***"'' NN CM O

5

« o T — NN CM ^f~ O • o CTN • ^ ^ T -O • o -d-1 ^ -d-CM O . O

r \ j

t

s

i ^ 4 ^

S V \ . 1 V i ^^

1 0.9850

0.9530

0.8536

0.75

0,6294

iO.5

0.3706

1 0.25

0.12,^4

, 0.0670

0,0170

i

5

4

3

2

1

0

-1

-2

- 3

- 4

"5

!

0.0170 1 0.0670

-5

0.00017

0.00005

0.00070

0.00113

0.0025e

0,0052f4

0.00555

0,00606

0,00590

io.00717

,0.00561

- 4

0.00005

0,00042

0,00174

0,0022+S

0.00458

0,00674

0,01091

0.01209

0,01607

0.01885

0,00717

0,1464

- 3

0,00070

0.00174

0.00257

0.0058O

0.00391

0,01327

0.02075

0,02504

0.05475

0.01607

0,00590

i

0.25

-2

0,00115

0.0022^8

0.00580

0.00826

0,01526

0,01958

0,05143

0.04954

O0O2504

0,01209

0,00606

0.3706

-1

0.00258

0.00458

0.00591

0.01526

0.03003

0.03579

0.05059

0.05145

0.02075

0.01091

0.00555

0.5 10.6294 0.75 iO.3556 10.9550 0.9850

1 ! •. i l

! i i

0 1 1 1 2 1 5

0.0052*4

0.00674

0.01527

0.01958

0.05579

0.06559

0.05579

0.01958

0.01527

0.00674

0.0052*4

0.00555

0.01091

0.02073

0.03145

0.05039

0,03579

0.03003

0.01526

0.00391

0,00453

0.00258

0.00606

0.01209

0.02504

0.04954

0.03145

0.01958

0.01526

0.00826

O.CO58O

0.0022f8

0.00113

0.00590

0.01607

0.05475

0.02504

0,02075

0.01327

0,00391

O.OO58O

0.00257

0.00174

0.00070

4

0.00717

0.01885

0,01607

0.01209

0.01091

0.00674

0.00458

0.0022,5

0.00174;

0.002*2

0.00005

5

0.00561

0.00717

0.00590

0.006061

0.00355 i

0.0052*4!

0.00258I

0.00113;

0.00070!

0.000051

0.00017

-10-2,2 The function F(x,y) v/ith half-order singularities

at the boundaries of the region of integration, but

bounded over the interior of the repjon

E question (11) of Appendix I shows that, under the

stated conditions of the heading,

1^ i'1

i P(x,y)

In

Jo

1 i «sr^^-v

^ I d x d y ^ ^ p - ^

n

'^o

c P(s ,s .)

li=-n v=-h

where s and the values of c for -n <* u, v < +n are

as given in the previous section. But if (i or v = + n, we

must interpret

°uv ^(^u'^) = ^uv ^ ^

' ^

/:=<y(1-x)(l-y) P(x,y)

^ ^ ^ x-*s^y-*s^— -J

and the c o e f f i c i e n t s d are tabulated below for n = 6,

UV

Again it is found that a symmetry exists for the coefficients

d . such that

UV

and

UV

[XV

= d,

= d,

V|I

(-u)(-v)

so that it is necessary only to quote d for - n ^ v ^ nj

evidently

un

= d

nu

^-n)v = \(-v)

and d / \ = d / \

n(-n; n(-u)

Table of values of d for n;=6

nv

V

d

nv

V

\ v

-6

0,12796

1

0,02026

-5

0,1502,5

2

0.00828

-4

0.00962

3

0.00700

-3

0.06804

4

0,00082

-2

0.04562

5

-0.00226

-1

0.03300

6

-0.00077

0

0.02210

3, Modification v/here integration liraits are altered

A simple substitution vd.ll show that, if b _> a,

(\h

pb _ i\1

-A

! -1

dx dy

a u a

G(x,y)ln l ^ ^ ! dx dy = • ^ ^ 2 j P(x,y)Zn | ^

Ü o •'...' o f.1 '.1 _ i n ( ^ 1 I p ( ^ , y ) a x d y ( b - a ) ^ I ^ ' J o J o

where F ( x , y ) = G j (b-a)x+a, (b-a)y4-a |

The first integral o n the right-hand-side m a y b e evaluated b y the method of the previous section. A n analysis similar to that u s e d i n the Appendix I shovra that the second double integral

r^ 01 2 4:2, „iïL n n

F(x,y)dx dy = - ^ ' ^

~"'>'\v ^^

lim |/xy(l-x)(l-y) P(x,y)_|

I 16n -'^..•w -dl—K ^ x-*s y-4 B "•""

o '-'o |i=-n v=-n ^ u V

where b = 4 if U and v / 4 - n ! b = 2 if either u or v = 4- nj

( J I V " ' ^ — ' | i V ^ — '

and b = 1 if b o t h fi and v = +^ n , (The limits i n this double sum only require interpretation if n or v = 4^ n ) , Use of this formula can b e convenient particularly if P has half-order singularities or zeros at the edges of the region of integration, as it involves data w h i c h is required to evaluate the double integral containing the logarithmic singularity. O t h e r , and more accjurate, methods m a y of course b e u s e d to evaluate the

1 2

-TECHNISCHE HOGESCHOOL

VLIEGTUIGBOUWXUN'DE Kanaalstraat 10 - DELFT

4 , Simplification r e s u l t i n g from p a r t i c u l a r symiMctries i n F(x,y)

( i ) I f P(x,y) = F(y,x) -ae see from the previous sections

t h a t as c = c ,

UV vn '

4-n

_{+n +n u}

' \ ' ' c P(s , s ) = "'™' •<:'',c P ( s , s )

v=-n u=~ïi v=-n

Y/here

and v/here

1°UV = ^2-Vv^°UV

8 = 1 i f |i = v , = 0 i f u / v i s the ICronecker d e l t a ,

Thus i n s t e a d of a square l a t t i c e of p o i n t s at v/hich the fui^ction

F(s , s ) has t o be evaluated, only the teruis i n a t r i a n g u l a r

l a t t i c e formed by and including the diagonal elements on u = v

need be used. Thus for n = 5» ^-i""! &• function F bounded

over the e n t i r e region of i n t e g r a t i o n , xfe find t h a t i n s t e a d of

the array of 5 x 5 terms given i n § 2 , 1 , the t r i a n g u l a r array

for the c o e f f i c i e n t s c t

UV

\ ^ u

V x..,^

2

1

0

-1

-2

-2

0,0625

-1

0,1599

0

0,2055 0.182^4 0.1110 1 0.0494

1

0.1599 0.182*4 0.0594

2

0,0625 0,1110 0,0494 0,0288 0.0288 i 0.0016

of 15 teixis suffices. Similar tabulations may easily be constructed for other values of n,

(ii) If P(x,y) = P(l-x,y), then as s = 1 - s, \

n

n

+n "'•^c .. P(s ,s ) = ^ "^7^ „c P(s ,s ) (i=-n v=-n (i=o v=-n

,x'

where „c = c +(l-6 ) c, \

Here the scjuare lattice is replaced by a rectangular lattice formed by ahd including the elements along |i = 0, Thus agaiji, for n = 5» £>J^cL a function P bounded over the entire range of integration, the rectangular array for „^ ^^

V ^ \ . 2 /[

1 °

' -1 -2 0 0.0247 0.0922 0.2055 0.0922 0.0247 4-1 0,0699 0.1896 0,1844 0,1896 0,0699 +2 0.0653 j 0.0699 0.0494 I 0.0699 i 0.0635

which is again of I5 instead of 25 terms.

(iii) If P(x,y) = P(x,1-y), then similarly 4-n +n n n " ^ J ^ ' " c P(s ,s ) = ^ ^

|i=-n v=-n li=-n v=o Vi?here ,c . = c 4- (1-5 )c / \

5 UV UV ov' u(~v)

and again the square lattice is replaced by a rectangular one divided now along the elements v=o,

--/ ,c P(s ,s )

(iv) If P(x,y) = P(l-x,y) = F(x,1-y), v/hich is a combination of the conditions (ii) and (iii) above, then

+n n ^.-. ^ , UV n=-n v=-n P(! U' V ) = n n

V ;.G„. F(s„,sJ

|a=o v=o n' V' where , c „ = I 1 4- (I-6 )(l-5 „) | c 4- (2-5 - 5 ,.) C/ N„ 4 UV !_. uo ov' j |j.v ^ no ov' (-u)v The lattice used is now the square formed by and including the elements on u = 0 and v = 0, Taking the same example as before, the nine coefficients , c for n = 5» and with a

' 4 UV ' function P bounded everywhere, are»

- 1 2 * ^

2

1

0

0

0,0494

0,1844

0,2055

1

0,1598

0,3792

0,1844

2

0,1266

0.1598

0,0494

Other tabulations for different n can quickly be assembled,

(v) If the conditions (i) and (iv) are both satisfied, further reduction is possible. For then

n n

«srirv ^^^— ^'c P(S , S ) = ' ^ ' ' ^ rC P ( s , S ) |i=-n v=-n |j,=o 'v=o

where c-c ,, = ( 2 - 5 „ ) , c 5 UV UV 4 I UV

Thus i n t h e p r e v i o u s example o n l y 6 c o e f f i c i e n t s , c a r e new n e e d e d

2

1

0

0

0,2055

1

0,5792

0.5688

2

0,1266

0,2796

0,0988

Taking n = 6, and the fxinction P boiinded over

the ran^e of i n t e g r a t i o n the square array of 1 1 x 1 1

i e n t s c can be reduced t o a t r i a n g u l a r array of 21 coeffic-

coeffic-i e n t s f-c subject t o the s t a t e d c o n d coeffic-i t coeffic-i o n s ! these Ecoeffic-ire

t a b u l a t e d belcw,

i ^ - ^ ' .

_[

5

4

3

2

1 ''

! 0

0

0.06559

1

0,16084

0,14516

2

0,11560

0,18676

3

0.07464

0.12336

0.09856

4

0.03854

0.07124

0,05828

0.06196

5 i

0.01156

0,02888

0.0262,0

0.02876

0.02452

0.07852 ; 0.05508 1 0.02696 ; 0.01576 1

(vi) If F(x,y) = - F(l-x,y), then

_n_ +n n. n

'^>' " ^ ? c F(s ,s ) ="^-5

"^'JrC

F(s ,s )

^.A -i:l.N l^v ^ |i' v^ ^ ^ ,

_{li=-n v=-n}

_{(i=1 v=-n}

^:^6

[iv ^ n' v'

where

_6°nv

= c

UV - c (-U)v

Here the scjuare lattice is replaced by a rectangular lattice formed by, but not including the elements along u = 0, Thus for n = 3 aiKl a function P bounded over the entire region of integration the non-vanishing coefficients /-c are 8 in number, and can be calculated ast

V =

2

y- =

i 1

2

0,0617

0.0411

,

1

0,0411

0,1302

0

0

0

-1

-0,0411

-0,1302

2 !

-0,0617

-0,0411

(vii) Similarly if P(x,y) =

n n

'^' ^ 7 c F(s ,s )

|ji=-n |i=-n

P(x,1-y), then

n __n

|i=-n v=1

•vï^ere 7 UV UV u(-v)

1 6

-( v i ü ) If conditions -(vi) and -(vii) both apply then as

V=°(-^)(-v)

^ ; ^

n n

yc„„p(s„,sj = "^^ ^ ' s v ^(^n'^)* ^«^

. ^ s . ^ ^ UV ^ tx' V

U=-n n=-n u=1 v=1

QC = 2 [^gC '\ , For the example quoted before the 5 x 5 square lattice of coefficients is now reduced simply to a 2 x 2 lattice.

(ix) If conditions (i) and (viii) both apply, thenj

è ± > "'V'%) = .S ^: 9"Vv ^f^'^)

li=-n v=-ri |i=1 v=1 where c . = 2(2-6 ) /-c . ,

nv ^ nv 6 nv

In the example for n=5 and a function P bounded every.vhere, the number of coefficients „c is nov/ only 3 given by

V - ^ ^ 2 1 1 0,2604 2 0,1234 0,162,4

Taking n = 6, and the function F bounded over the range of integration the square array of 11 x 11 coefficients c can be reduced to a triangular array of 15 coefficients _c subject to the stated conditions! these are tabulated

5. Applications in Aerodynamics

As mentioned in the introduction to the paper integrals of the type referred to in this note appear in formulae for the wave drag at supersonic speeds. Two particular applications can be quoted! vd.th

P(x,y) = S'(x)

Viy)

the integral appears in the expression for the wave drag of a swept v/ing of infinite span vriLth subsonic leading edges, and with a section whose thickness is t and v/hose ordinates are

A

z = •g-t^(x) at a fraction x of the chord from the nose. The present integration formulae are particularly suited to evalua-tion of such integrals where the secevalua-tion leading-edge and

trailing edge are either rounded or cusped, because the function z'(x) has then half-order singularities at x = 0 or 1• Where such conditions are satisfied it is found that the formulae given above yield exact values in ccsnparison "v/ith all simple exaniples vrorked exactly (by formal integration) merely because the form of the interpolating formula (used in deriving the present method) describes the section shape exactly, T/here such conditions are not met only a fev/ examples have been calculated by exact methods i one is for the biconvex wing section, which corresponds to

-18-In this example F(x,y) is finite (and non-sero) at the borders of the region of integration: use of the integration formula with n = 6 then provides a result too small by a little less than 1 per cent. Application to sections with round noses and angular trailing-edges (so that ^' has the half-order singularity at one end of the rajigo of integration only) has been found to provide answers correct v/ithin less than 0,1 per cent; the same order of accioracy is also attain-able for the treatment of sections v\dth one cjusped edge. The method is certainly not applicable to the evaluation of

integrals v/here F (although bounded) has discontinuities! thus, in the application under discussion, if

^•(x) = 2 sgn(l-2x),

corresponding to the double vredge section, the integration formula provides too small an ansv/er by some 10 per cent,

The second application which comes to mind is in the evaluation of the expression for the v/ave drag of certain

2

classes of slender bodies vdiere in our notation an integral with

P(x,y) = S"(x) S"(y)

appears, S(x) denoting the body cross-section area at a fraction x of the length fran the nose. The integration

formula quoted here are particularly applicable if S'(x) can be represented by a finite sine series in 6, where

6 = cos (l-2x), because then the function P(x,y) has ha.lf-order singularities or zeros along the bha.lf-orders of the region of integration (as is assumed in the metliod for numerical integration). Many shapes of body, v/hich have been derived to prcovide minimum drag under stated conditions, fall within this category, and likevdse those obtained closely approximating to them. On the other hand many body shapes - such as

ellipsoids, or those vra.th parabolic meridian sections - have a cross-sectional area distribution S(x) v/hich is represented

by a polynomial in x, and which is equivalent to the expression of S'(x) by a cosine series. Applied to such shapes the integration formula are accurate only vdthin about 1 per cent, even vd.th n = 6j this is in direct analogy to the results cjuoted above for vdng sections with angular edges, Likev/ise, too, the formula are inapplicable where S " (x) is discontinuous (i,e, v/here the body has discontinuities in curvature),

It is understood that work on similar lines to tliat \xnder present discussion has been undertaken at the Royal Aircraft Establishment, Pamborough, in relation to this application to the wave drag of slender bodies. The accuracy cjuoted from their preliminary examples does not seen as close as that obtained from the use of the above formulae! this is probably due to the fact that the double differentiation of S(x) to obtain S " (x) is implicitly accomplished in

their technique, as a Fourier series representation of S(x) is used to fit its values at stipulated points, and not of S" (x). Numerical differentiation by Fourier series inter-polation is knov/n to be a relatively inaccxnrate device. This source of error does not appear in the examioles cjaoted here, as no attention has been paid to the problem of derivation of S" (x) given only numerical data for S(x), which of course often arises in practice,

To sum up, in aeronautical applications where F(x,y) is separable as f'(x) f'(y)> s^> the integration formulae have been found quite satisfactory v/here f'(x) is k n a m oxid

has half-or der singularities (or zeros) at x = 0 and/or 1, Ihere f'(x) is finite but non-zero at these end points, an accuracy of about 1 per cent (if the formulae relevant to n = 6 are used) is all that can be anticipated,

-20-In the failing cases, which certainly include those where f' (x) is discontinuous at interior points of 0 -,: x <; 1, the v/ork of lagendre may be of particular use. In this he

reduces the double integral of the type under consideration to one with a bounded integrand, together v/ith a line integral (and other terms if f(x) is discontinuous)! this new double integTEil may be treated by ¥eddle's Rule, or another similar formula. The flexibility of this Rule, in allov/ing arbitrarily close spacing of lattice points at which the integrand is

evaluated, together v/ith the fact that it is f(x), and not f'(x), which appears in the finite integrand derived by legendre, ccmmend this method, ïfeddle's Rule is strictly speaking not in this connection applicable if f'(x) has

half-order singularities or zeros at the end points of the range of integration, but as has been pointed cut, the method herein described is then quite adequate and probably rather simpler to apply,

Acknov/ledgement

The author v/ishes to acknowledge the assistance of Mr, S, Ingham who performed the canputations of the results quoted here,

1, T, Nonv/eiler

2, G,N, Ward

5. R, Legendre

References

Theoretical Supersonic Drag of Non-lifting

I n f i n i t e - s p a n 17ings Swept b e h i n d t h e iiach l i n e s ,

R, and 11. No» 2 7 9 5 . (1950) S u p e r s o n i c flow p a s t S l e n d e r P o i n t e d

B o d i e s .

Quart .J.I.Iech. and App.Maths . V o l , I I , P t . 1 (1949) L i m i t e sonique de l a r e s i s t a n c e d ' o n d e s

d ' u n £ieronef.

Ccciptes R-endus, Tome 256, N o , 2 6 , p , 2 4 7 9 ( 1 9 5 5 ) . \^ ('1

1 f'

^'o Jo

^ ^ l i - .

* The result is that I \ f'(x)f'(y)ïn - ~ | dxdy = -^ j

f(x)-f(.Y)

r2

x-y i • f x-y

and other terms, so that in fact the value of f'(x) is recjuired for the value of the integrand along x = y,

AEFENDIX

Derivation of Integration Formula

Y/e suppose that the function F(x,y) is expressible in the forms

F(x,y) = " 2 ^ ^ P sin 2nQ sin 2n^

n^-n v=-n ^ (cos-f' ~cos 6)(cos'f -cos 0)

' (1)

v/here n, n» ^^^ v are integral, P constant, and

e = cos''' (l-2x), 0 = cos"\l-2y), 'f^ = f + U ^^^ . ( 2 ) •

Using the i d e n t i t y !

s i ^ 2n6 2 ( l ) " " ^ ^ ^ , ,

-= —^ ^ •-. b, cos k tl cos kö

II, Q • Q .^"'^ k - n "^n

cos V' - COS 6 s m 6 ~-• '^ / - \

^1^ k=o (3)

where b = 1 for u / ± n , and b = ^ for u = + ï^>

n u

i t v/ill be seen t h a t {_sin 6 s i n 0 P ( x , y ) j i s represented as a

f i n i t e Fourier double cosine s e r i e s over the half-periods

{Of%) i n 6 and 0, By t a k i n g n as \inbounded the

v/ell-known methods of Fourier analysis may be used t o show t h a t a

large c l a s s of functions F(x,y) can be expressed i n the

chosen form, Hov/ever i n v/hat follows v/e take n as bounded,

and i n general any a r b i t r a r y function P(x,y) can only be

expressed approximately by the expression ( l ) ,

Now i t can be shov/n t h a t , i f p i s any i n t e g e r

lim / _ 2 i ï L 2 n 9 ^ ^ 2(-l)"-'^' n cosec il' , ±f p = U / ± n|")

Ö - ^ f \^cosy. - COS 6/ '^ [

= 0, i f p ;i^ n . (

(4)

2 2

-s l n 2nQ \ 4 n

\J

\ COS \u -

w

cos

V

S i n

, as 9"-.y;^

_.(5)

I t follows, from (4) and (5) i n ( l ) , t h a t

P

Izii^*" ,„

1 y . . t —

^'" 4n2 e 5 ! ^ 6 ^ ! / / . I V v ^ ^ Ö s i n 1^ P(x,y)

_%\y

(6)

where, as before, b = 1 if n ;^ + n, and b = 4- if n = + n»

The limits are of cxoiirse easily v/ritten dov/n unless n o^ v

equals ^ n. It v/ill be seen that the advantage of the

representation (l) lies in the ease v/ith which the coefficients

F of the series are expressible in terms of the values of the

function it represents. If the value of P(x,y) is knov/n only

at the lattice points (s ,s ) , where n>v = -n, -n+1, ,,, n-1,n,

and

^U "" ^ V "^ ^"^ 2^) ' ^'^^

% ^

°°^"'' ^^"^^U^* "^^^

then the right-hand side of (1) may be regarded as an

interpola-tion formula for F, giving its value at all intermediate points,

Let us nov/ put

1 {0) ^\

sin 6 sin 2nQ ^ 2

^ \

cos

^ -

COS Ö cos Q - cos

0

% •

do (8)

Then substituting from (3) and integrating by parts

• i T C

2n

I (i2^) = 2 (-1

T^^ \ ^<'

b, cos k

^v

cos kG In ;

U ;

li:!.:

k-n Tn

•} o k=o

cos 6-cos

0

de

= 2(-1)^•^^^

s^

In

1-=.'o

jcos 0-cos

0\

I dO

-2n • "'^

\ ' 1 ^ , 1 ' sin kOsin e ,„

- y-rb, cos ki/, s v dS

^ k k-n TfJ. cos 6 - oos

p

k=1 Jo

i . e . 1^(0) = 27c(-l)"-^'^ l o g 2 + "^N f b. c o s k </' cos kJ2f

/ 1

_••0]

^-^^ k k - n i n

k=1 --J

(9)

But the double i n t e g r a l of v/hich we r e q u i r e t h e v a l u e i s

n1 T'l

1 = 1 I P(x,y) In — I dx ay

.(10)

'..' o '.' o

T/hich on substitution from (l) and (2) becanes

T 1 I 1 '-^^/ *^~"' -ni sin 2nQ sin 2n2f sin Q sin 0 ^

I = T j

•;>

^ P '

'^— log

J o J o ^ n f = - n ^^ ( ° ° s ' f ^ - ° ° s G)(cos ^-^-cos 0) cos6-cosj2f dGdj2^ and from ( 8 ) . n n _ I = r '^~->' ' ' \ ' P s i n 2n;2f s i n jZf ^ ^^^ ^ ^ iïc:'. i -s.-:iix ^^v j ,, g^ n n:=-n v=-n '-J o j V ^ Using t h e r e s u l t (9) i t f o l l o w s t h a t ^ " ^

I = f -^V ^ ' ( - 1 ) ' ' ^ ^ F

n=-n v=-n

f

2n l o g 2 + ' ^ ' r - b T- b. c o s k'/^ c o s kjZf k=1 V s i n 2n(^ s i n 0 ,-/ cosil' - c o s 0 l y T ^ .

-= f ^ > " S ' (-1)^'^ F I ! log 2,sin 2ni^ ^

n=-n v - - n 2n

v.'o

V 2k*\-n°°^ ^T/^ ^'^^ 2n+k0 4- sin 2n-kJ2f) j

]c=1 c o s y:' -COS ^ 2 n ^ _ — J . n 7i:=-n v=-n , J n 1 2 l o g 2 4- ""> r t», c o s k^ COS. k y *c:L.N k k - n / n ^v

_ k=1 _i

I

-24-Using (6), this equation can be v/ritten as -2 __n.. _ n I - — 2 8n^ V7here .. UV lim lim U=-n v=-n rn V 2n

b b sin 6 sin

0

F(x,y) ƒ

V (11) 1 .2

^^v = ^°S 2 4. ^ " ^ 1 bj_^ cos ky.- cos kj/. ^

k=1

3

In particular if P(x,y) is bounded at the edges of the region of integration, from (7),

^n-1 n-1

n=-n+1 v=-n4-1

r

.(12)