UC-NRLF $B 531 D3D £Mnbuvob nDatbematical XTvacts INTERPOLATION AND NUMERICAL INTEGRATION A COURSE IN INTERPOLATION AND NUMERICAL INTEGRATION for the Mathematical Laboratory by DAVID GIBB, M.A., B.Sc, Lecturer in Mathematics in the University of Edinburgh Xonbon: G. BELL & SONS, LTD., YORK HOUSE, PORTUGAL STREET 1915 Digitized by the Internet Archive in 2008 with funding from IVIicrosoft Corporation http://www.archive.org/details/courseininterpolOOgibbriGh PREFACE nPHE present work is intended for the use of students who are learning the practice of Interpolation and Numerical Integration. The advantages of a practical kno wedge of this part of Mathematics ai-e so obvious that it is needless to insist on them here : and these subjects form an important part of the course in the modern Mathematical Laboratory. There are, however, so many claims on the time of students that the extent of this course, as of all others, must be kept within narrow limits : and it has therefore been necessary to restrict the treat- ment to the most central and indispensable theorems. A large number of numerical illustrations and examples has been given of a kind likely to occur in the applications of Mathematics. My thanks are due to Professor Whittaker and to my colleague, Mr E. M. Horsburgh, M.A., B.Sc, Assoc. M.Inst. C.E., for their valuable criticisms and suggestions. D. G. The Mathematical Laboratory, University of Edinburgh, July 1915. 325S65 CONTENTS PAGE Preliminary Note on Computation - - - i CHAPTER I. Some Theorems in the Calculus of P^inite Differences. Section 1. Differences ------- 3 2. Properties of the Operator A - - - - 5 3. Expansion of /(.v + z/tc) for Integral Values of;? - 7 4. Relations between the Differences and the Differential Coefficients of a Function . . . . g Miscellaneous Examples - - - - - 12 CHAPTER II. Formulae of Interpolation. 5. Introduction ------- 14 6. Lagrange's Formula of Interpolation - - -15 7. Periodic Interpolation - - - - - 18 8. Notation ------- 18 9. Newton's Formula of Interpolation - - - 19 10. Stirling's Formula of Interpolation - - - 24 11. Gauss' and Bessel's Formulae of Interpolation - - 26 12. Evaluation of the Derivatives of the Function - - 29 Miscellaneous Examples - - - - -33 CHAPTER III. Construction and Use of Mathematical Tables. 13. Interpolation by First Differences - - - - 3^ 14. Case of Irregular Differences - - - - 3^ 15. Construction of New Tables - - - - 39 16. Tabulation of a Polynomial - - - - 39 VUl CONTENTS Section page 17. Calculation of the Fundamental Values for a Table of Logarithms - - - - - - 40 18. Amplification of the Table by Subtabulation - - 43 19. Radix Method of Calculating Logarithms - - 46 20. Inverse Interpolation - - - - - 51 21. Various Applications of Inverse Interpolation - - 58 Miscellaneous Examples - - - - - 61 CHAPTER IV. Numerical Integration. 22. Introduction - - - - - - - 63 23. Evaluation of Integrals by the Integration of Infinite Series Term by Term - - - - - 64 24. A Formula of Integration based on Newton's Formula - 65 25. Case when the Upper Limit is not a Tabulated Value - 70 26. A Formula based on Bessel's Expansion - - 73 27. The Euler-Maclaurin Expansion - - - - 76 28. An Exceptional Case - - - - - 78 29. The Formulae of Woolhouse and Lubbock - - 79 30. Formulae of Approximate Quadrature - - - 81 31. Comparison of the Accuracy of these Special Formulae 85 32. Gauss' Method of Quadrature - - - - 86 Miscellaneous Examples - - - - - 89 PRELIMINARY NOTE ON COMPUTATION. DEFORE introducing the formulae required in interpolation and numerical integration, it may be advisable to mention a few points which will facilitate the work of computing. Where computation is performed to any considei'able extent, computer's desks will be found useful. Those used in the mathe- matical laboratory of the University of Edinburgh are 3' 0" wide, r 9" from front to back, and 2' 6|" high. They contain a locker, in which computing paper can be kept without being folded, and a cupboard for books, and are fitted with a strong adjustable l)ook-rest. Thus the computer can command a large space and utilize it for books, papers, drawing-board, arithmometer, or instruments. Each desk is supplied with a copy of Barlow's tables (which give the square, square root, cube, cube root, and the reciprocal of all numbers up to 10,000), a copy of Crelle's multi- plication table (which gives at sight the product of any two numbei's each less than 1000), and with tables giving the values of the trigonometric functions and logarithms. These may, of course, be supplemented by a slide rule, or any of the various calculating machines * now in use, and such books of tables as bear particularly on the subject in question. Success in computation depends partly on the proper choice of a formula and partly on a neat and methodical arrangement of the work. For the latter computing paper is essential. A convenient size for such a paper is 26" by 16"; this should be divided by faint ruling into ^" squares, each of which is capable of * For descriptions of these see Modern Instruments and Methods of Cal- cidation, edited by E. M. Horsburgh. London : G. Bell & Sons. 1914. e. 1 ,2 INTERPOLATION AND NUMERICAL INTEGRATION holding two digits. It will be found conducive towards accuracy and speed if, instead of taking down a number one digit at a time, the computer takes it down two digits at a time, e.g., instead of taking down the individual digits 2, 0, 4, 7, 6, 3, it will be found better to group them together, as 20, 47, 63. Every computation should be performed with ink in preference to pencil ; this not only ensures a much more lasting record of the work but also prevents eye-strain and fatigue. It need hardly be emphasized that seven-place accuracy cannot be obtained by the use of four-figure tables. A fairly safe rule is to make use of tables containing one digit more than the accuracy requii-ed ; to use more accurate tables than these is simply to increase the amount of labour without increasing the efficiency. For the same reason contracted methods of multiplication and division should be adopted when machines and tables are not at hand. CHAPTER I. SOME THEOREMS IN THE CALCULUS OF FINITE DIFFERENCES. 1. DiflFerences. Let \f{x) denote the increment of a function f{x) correspond- ing to a given constant increment w of the variable a-, so that ^f{x)=f{x + w)-f{x). The expression l^/{x), which is usually called the First Differ- ence of the function f{x) corresponding to the constant increment w of the variable x, will, in general, be another function of x, and as such will have a first difference which will be obtained by operating on A/ (a-) with \ Calling the result of this operation A^f{x), we have AV(x) = A{A/(x-)} = A{f{x + w)-/{x)} = {/{x + 2w) -f{x + w)}- {/{x + ic) -f{x)] =f{x + '2w)-2f{x + w)+f{x). This is termed the Second Difference of f{x) corresponding to the constant increment w of x. By repeating the process we obtain in turn ^^f{x), ^^f{x), ... i\'\f{x), which are called the 3rd, 4th, ... nth differences of /{x) corresponding to this constant increment 7V. Xote. — It is always possible by means of a suitable trans- formation to make the value of w unity. For since w is the increment of x in f{x), then unity is the increment of — . Replacing x in f(x) by wy, we obtain a new w function F {y), which is such that ^F{y) = F{y + \)-F{y). INTERPOLATION AND NUMERICAL INTEGRATION [ch. Example 1. — Tabulate and difference the successive values of f[x): from a; = 110 to a; = 118. A3 A* X 7? A A2 no 1331000 36631 111 1367631 37297 666 112 1404928 37969 672 113 1442897 38647 678 114 1481544 39331 684 115 1520875 40021 690 116 1560896 40717 696 117 1601613 41419 702 118 1643032 In the above table the column headed x^ was obtained from Barlow's Tables. The column A gives the first differences, being obtained by sub- tracting each entry of the column x^ from the entry following it, e.g. if x=114, a;+l = 115, A/(a:) = (a;-i- l)3-a;3= 1520875- 1481544 = 39331. This result is set on a line midway between the numbers corresponding to a;=114 and a;= 115. In the same way the 2nd differences, A-, are obtained from the 1st differences, the results being set on a line midway between the latter, and so on. It will be noticed that the 3rd differences of ar* are constant and that its 4th differences are zero. It will be shown later that the nth differences of a function of the nih. degree are constant and that the differences of higher order are all zero. Example 2. — To determine the successive differences of sin (ax + 6). Let /(a;) = sin(aa:-l-t). Then A/{a:) = sin {a x+w-vh) -sin [ax-vh] = 2 sin ^a IV cos (ax + b + ^aio) : 2 sin ha w sin (ax + b + haw + ir) A-'/{x) = A{ 2 sin I a w sin {ax + 6 + ^ a ?t' + tt)} = 2 sin ^ as = {2 sin haw}- cos {ax + b + ^2aw + v) ia w { sin {a x + io + b + ^aw + ir) -sin (a x + b + ^ a w + w)} ■{2smhaiv}"am (ax + b-\- ^2aiv + 2 ir). 1, 2] THEOREMS ON DIFFERENCES Similarly, it can be shown that A»/(a;) = { 2 sin ^awf sin{a x + h + ^3aw+^Tr) and, more generally, that A"f{x) = { 2 sin law}" sin {ax + b + ^n aw + 7nr). This can be proved by induction in the following way :- A"+i/(a;) = A [{ 2 sin | a w}" ain {ax +b + ^naw + n ir)'] : { 2 sin i a m; }" [sin (ax+w + b + ^naw + mr) - sin {ax + b + ^naw + mr)] = { 2 sin I a iv }"+i cos{ax + b + ^{n + l)aw + mr) = { 2 sin i a «}«+i sin {a x + b + ^ {n+ I) a lu + {n+ I) n). But the theorem has been proved for n—l, 2. Hence it is true always. Example 3. — Find the first differences of the functions cos {ax + b) ; tan— ^ ax ; «''•"=. Example 4. — Find the nth differences of the functions — ; a-' ; cos^(a.v + 6). X 2. Properties of the Operator ^. The operator A obeys certain of the fundamental laws of algebra, viz. — (i) The Law of Distribution. For d {fix) + ^ {x)) = {/(x + IV) + c/, (x- + w)) - {/(.«) + c/) ix)) = [fix + w) -f{x)] + {c^{x + w)- {x)] = A/(a;) + Ac/,(x). (ii) The Law of Com mutation if the coefficients of the various terms are constants. For, if a is a constant, A [a .f{x)] = af{:c + w) - af{x) = a{f{x^-w)-f{x)] = a.\f{x). (iii) The Law of Indices. For A'- . \^f{x) == {A . A . A ... (?• factors) .AAA ... (s factors)} /(a;) = {A . A . A ... (r + s factors) }/(x-) The indices r and a" are, of course, integers : the operator A" was defined in the previous section only for integral values of ii. 6 INTERPOLATION AND NUMERICAL INTEGRATION [ch. i Example 1. — If f{x) is a polynomial in x of the nih. degree, then its nth difference is constant and its differences of higher order than n are zero. Let f(x) — anX'^+an-\x^-'^+ ... +aix + ao in which the coefficients a„,an-i, . . ai, ao are constants. Then Af{x) = {a„ {x + lo)" + «"-i (x + w)"-^ + ... + ai (x + ir) + a^, } — {a„ X" + a,i-i x"-'^ + ... + ai a; + ao } = 7t an x"-^ ir + h„_.2 X"-- + ... +bi x + h„ where 6„_2, ..., ^i, 60 are constants. We thus see that A/(x) is a polynomial in x of degree (71 -1). Repeating the operation, we have A^f{x) = {n a„ [x + »')"~^ "' + t„_.j {x + *")"-- + ... + ?>i (x -^ »•) + ?>,j} -{nanX^'-'^ ii; + h„-2X"--+ ... +h.^x-\-h(^} = n(,n-\)a„x"'-->if- + Cn~'iX'^-^-\- ... +Cia; + Co where c„_3, ..., q, c^ are constants. Hence A-/(x') is a polynomial in x of degree (/i-2). It is evident that by continued repetition of this process the degree of the function will be reduced to zero, and that we shall obtain A«/(a;) = »i (vi - 1 ) (?j - 2) . . . 3 . 2 . 1 a„ »•" , i.e., the ?ith difference of the function is a constant. Obviously the next and all differences of higher order will be zero. This result is of great import- ance ; probably no other theorem will be applied as frequently as this. For instance, it enables us to complete a mathematical table in which there are gaps. This is illustrated in the next example. Example 2. — If /(a:) = 0-5563, /(x f l)-0-p682, /(x + 3) = 0-5911, /(a; + 5) =^0-6128, obtain /(a; + 2) aad/(a; + 4) on the assumption that /(x) is a polynomial of the third degree. It has been shown in § 1 that ^""/{oc) =/(.T + 2 w) - 2f[x + w) +f(x). By repeating the process we get A'f{x)=f{x + 4.ir)-4:f(x + 3w) + 6/ix + 2,r)-4/(x + >r)+f(x). Bat the given function is of the third degree. Hence for it we have A-'/(.x-) = 0. Noting that w is unity we obtain f{x + 4) -■i/(x + 3) + &/{x + 2)-i/{x+l)+/{x) = or, denoting /(x + Ji.) by/,,, /4-4/; + 6/,-4/.-f-/o = 0. Similarly /o - 4/, + 6/ - 4/, 4-/1 - 0. 2, 3] THEOREMS ON DIFFERENCES 7 Substituting the given values in these two equations, we get /4 + 6/3 = 4 -0809 /4+ y;- 1-1819 whence /, = 0-5798 , /j = 0-6021 , i.e., /(X- + 2) = 0-5798, /(x-r 4) = 0-6021. Example 3. — If when a; = 0, 1, 3, 4, 6, the corresponding values of f {x) are 1, 1, 79, '253, 1291, 6nd approximate values for f{x) when x = 2 and when .1=0. Under what conditions will the result be accurate ? 3. Expansion of f{x + mv) for Integral Values of n. We shall now establish certain formulae which express J{x-\-mv) (where n is an integer) in terms of differences of various orders. These formulae will be shown in Chapter II, to be valid under certain conditions even when n is not an integer. It will be found convenient to introduce another operator E defined by the equation E=i+x y This new operator will obviously obey the same laws as A. Now EJ (x) = (1 + ^)/{x) = f{x) + \f{x) ^/{x + w). Thus the effect of the operator E on the function f{.v) is to increase its argument x by tlie given constant quantity w. Repeating the process we get E"J{x)=f{x + nw) and since the operator sS. may under the conditions of § 2 be treated as an algebraical quantity, this last equation may be written f{x + nw) = {l + ^Yf{x) =f{x)+,fi,^f{x) + „c,^''f{x)+ ... +^'\f{x) (1) where ,/.v denotes the coefficient of a'' in the binomial expansion of (1+a)". The function /(a; + 7iM/-) is here expressed in terms of f{x) and its successive differences. It is evident that, since the effect of the operator E on f{x) is to increase x by w, the result of operating on f{x) with E~'^, where E~^ is defined by the equation EE'^^^1, will be to subtract w from X, i.e., E-^/{x)=^/{x-w). 8 INTERPOLATION AND NUMERICAL INTEGRATION [cu. i Hence the above formula (1 ) may be transformed as follows : J {x + nw) =/{x) + ,.ci ^f^x) + {l + ^) E-' {,c, A- fix) + „c, ^y(x-) + ... +S"/{x)} =/(•«'•) + ,fii ^J (a) + (1 + Aj {„c, A- f{x - iv) + ^^Cs^'/{x-7J))+ ... +1"/{X-W)} + (l+A)^-M„+iC,Ay(x-«;)+ ... } + (l+A){„^,c,Ay(x-2u;)+ ... } + ,,+1^4 AV'(« - 2^) + „+oCg A'f{x - 2m;) + ... Operating next with (1 + A)^-' on the 7th and subsequent terms, then on those beginning with the 9th term, and so on, we obtain f{x + nw) =/{x) + „ci ^/(x) + „c., /\'/{x -w) + „+iC3 Ay (a; - w) + „+,c, A V(«-' - 2*«) + „+.c, Ay (X -2w)+... (2) the formula ending at the term .X-"~^/(x - n - I w). In the same way it may be shown that f{x + n w) ^f{x) + „c, Ay'(A' - w) + ,,^^c. A-/(^; - lo) + „^,c, AV(-*; - '^^) + ,,+.^4 A-* ;\x - 210) + „^.,C5 A'/'(x - otv) + ... {3} the formula ending at the term A-'"y"(;f - n w). Changing x into x + w and rt into n-\ in (3) we have J\x + 11 w) =/{x i-w) + „_iCi \f{x) + ,fi., A-y"(a;) + „c. A"/ (a; - w) + „+jC, AVXa; -w) + „+,c, A' fix -'2w) + ... (4) the formula ending at the term A-"/(x -n-2 w). J, 4] THEOREMS ON DIFFERENCES Adding (2) and (3) we obtain f{x + n lo) =f{x) + ,fi, ,^ ^ ' + — ,fi, A-7 (,x - w) A\f{x-w) + ^'/(x-2w) + ^„^,c,i\Y{x-2w) + ... (5). Similarly, addition of (2) and (4) gives Formulae (1), (2), (3), (4), (5), (6) may be compared with the expansions introduced in Chapter II. under the names Newton's, Stirling's, etc., Formulae of Interpolation. Again, from the definition of the operator E, we see that Ay\x) = {E-\r/{x). Expanding the operator in this last relation we obtain A"/{x) = E''/{x) - ,fi, ^-VX*-) + .^. E^'-'fi^) ••• + ( - 1)V(^-) =f{x -\-nw) - „Ci/(x- + n-lw)+ „Co/(.r + w - 2 w) + ... +{-\rf{x). In this result we have the 7ith difference of the function f{x) expressed explicitly in terms of f{x) and its successive values. 4. Relations between the Differences and the Dif- ferential Coefficients of a Function. The second example of g 1 may ha\ e reminded the reader of a somewhat similar result in the diflerential calculus, viz., -— sin (a x + h) = a" sin (« x i-b + hn tt). dx This result may be derived from the example in question. For if we put iv = .\x, where Ax' denotes a small inci-ement in x, and, 10 INTERPOLATION AND NUMERICAL INTEGRATION [ch. i after dividing both sides of the equation by (A x)", make A x tend to zero, we get d'\f {^) T (-2 sin ia . Ax^" . , . , . , .. ^d^ = lu['""tx~ ) sm{a*' + 6 + l7.(a.Aa; + 7r)}. Aa;->0 = a" sin {ax + b + J n tt). This naturally suggests that there may be some relations connecting the successive differences of fix) and the successive differential coefficients of the same function. That such a relation does exist we now proceed to show. Applying Taylors Theorem* to the function /(a; + w) we obtain f{:x + w) =/{x) + wf (a;) + ^y"' (.r-) + ... + ^./■'"" {^^ + • ■ • or, since f{x + id) -J\x) = \f{x), A/{x) = w/'{x)+ '^^/"{x)+ ... + ^f-'(x) + ... which is a relation of the kind sought for. Again, by Taylor's Theorem, ./(^- + "')=/ (-»•)+ wj'{x)+ — _/'(y/)+ ... + --/'""(.'/;)+... 'J> ! in\ f(x+2to)=/{x) + 2w/'{x)+^^^/"{x)+ .. +l^'/<™)(a;) + ... and generally /{x + n 2v) =/ (*) + n to/' (x) + ^' /" (x) + ...+ ^' f"" (o:) + ... Hence '^"/(•^-') =/(''' + '* '^') - »Ci /{-^ -t- w - 1 tv) + „c,/{x + 11-2 iv) - ... +{-\)y(x) = 2 -^ {"'" - „Ci {n - 1)'" + ..Co [n - 2)'" - . . . }/'"" (x). * It would be inopportune at thi.s stage to introduce a discussion of the questions of convergence whicli arise in connection with the application of Taylor's Theorem. We shall assume a convergence sufficient for our purpose. 4] THEOREMS ON DIFFERENCES 11 But ^"' ~ «Ci (w - 1 )'" + ,fi2 ('«■ - 2)'" + ... = if m<7i* = nl if ni = 71 = i ri (w + 1 ) ! if m = ?i + 1 = — ?i(3« + 1) (?i + -2) ! ifm = 7i + 2 etc. Making use of these results in the expi'ession for A"/(x) we find that 15 7r (?i + 1) „ . ,„ + ^^ 'w"^'J^'^+-{x) + ... Corollary 1. — Since the rzth differential coefficient of a function _/'(x') of the nth degree is a constant, and the derivatives of higher ordei- ai-e zero, we see again from this expression for A"/{x) that the ?ith difference of a polynomial of degree ?*. is a constant, and the differences of higher order zero. Corollary 2. — By reversing the series for Z\"/"(.^'), we obtain an expansion for the nth differential coefficient of /{x) in terms of its differences, viz. lon{n + 2) {n + 3) .... . ^ 6! - "-^^' CoruUary 3. — If the ?«th differences of a function are constant, then the function is a polynomial of the nth degree. For since the nth differences are constant, the higher differences are zero, and therefore the nth. differential coefficient is constant, Chrystal's Algebra, Fart ii., Chap, xxvii., § 9. 12 INTERPOLATION AND NUMERICAL INTEGRATION [ch. i and the differential coefficients of higher order zero. Thus the function, when expanded by Taylor's Theorem, takes the form /(,«) =f(X) ^(x- x)f (X) + '•' ^ f^\ r (X) + ... .'-^^■,/ (.r, " n ! and this is of the ?ith degree in x. Corollary i. — Since e^ = 1 + j-: + — + ... -\ . + •■• we see that the series which represents /{x + w), viz. f{x) + ivf'(^x)+ y,/"(a^) + ... may be symbolically represented as d e"~^f{x). Hence Ef{x) = e'^f{x). MISCELLANEOUS EXAMPLES. (Asstime in the followiny exercises that w is unity, unless otherwise indicated). 1. Find the first differences of the functions tail a X ; log x ; ( a-^ - a* ) / (« - 1 )• 2. Prove that sin A tan a; = eosx d cos(a,-+ \)d' 3. Prove that 1 1 e ^ X"'"^^ taii^;^ . 2^ tan— 4. Find the u'-'' differences of the functions x+\ ^ . , —- _ — . x™+" ; sin ax cos ox. X- \Zx ' 5. Find the »ith difference of each of the factorials a;(a; + a) (x + 2a)... and 1 a;(a; + a)(x + 2a). there being m factors in each factorial and the increment 4] THEOREMS ON DIFFERENCES 13 of X being a: and show how to develop any polynomial in x in a series of factorials. In particular show that x^ =x^lx{x-\)^^x{x-\){x-'l)-Vx{x-Y){x~1){x-'>,). 6. If f(x) = a?,^ + 6 2^- prove that AV(a;)-3A/(x') + 2/(x) = 0. 7. Prove that A X" - 2A--^ x" + 3 A3 a;" - = (a- ^ 1 )" - (.r - 2)" . 8. Establish the formula for tJ e Jith difference of the product of two functions A« {f.x) . 4> {X) } = A«/(.T) • 0(-i-) + n A«-i/(.r + 1 ) . A ( .r) 7i(n-l) + ^yy- A«- V(a' + 2) . A'-^,^(a-) + . . . 9. Show that (E) a\f{x) = a^

(x) does not become infinite when x = a^, o.,, . . . a„. If, however, /(x) is a function of which nothing is known except that it takes certain definite values corresponding to certain definite values a^, a.., ... a„ of the variable, then the simplest assumption which we can make regarding its form is, that it is the function of degree (n - 1) which sa,tisfies the given conditions ; and this function is given by equation (2), The expression on the right-hand side of (2) may therefore be taken as the expression for /{x) for all values of x in the interval considered. It may be written in the form fix) = Z -, , , , J (a.) ,.^1 (x- a,.) (b (a^) where (/> (x) = {x - a^) {x - a^) ... (x - a„). This is known as Layrange^s Formula of Interpolation. The advantage which it has over other formulae for the same purpose is that it is in a form suitable for logarithmic computation. This latter form of Lagrange's formula may be obtained directly from the formula previously given. It can also be obtained very readily by the Method of Partial Fractions, /(£l {X) {x) = [x - a,) [x - ff.,) ... (.T - a„ ). For consider the proper fractional function , , . where 6] FORMULAE OF INTERPOLATION 17 Decomposing it into partial fractions we have

{^) is the polynomial of lowest degree which passes through the points. Suppose, now, that the values of a function are given not merely for a finite set of values of its argument, but for the infinite set x = a + riv {r = Q,\,'2,'c>, ... ad inf.). Then we can, by Newton's formula, construct a polynomial of degree {k - \) which has the given values for the first k of these values of the argument, where k is any positive integer : and it may happen that, as k increases indefinitely, the series on the right in Newton's formula becomes (for some range of values of a;) a con- vergent series, having for its sum a function f{oi). The function /(*■) will now have the given values for each of the values a + rw (r = 0, 1, 2, 3, ... ad inf.) of the argument. There are, however, an infinite number of functions which satisfy this condition and /(a-) is that one of them which is the limit, when ^ -> x , of the polynomial of lowest degree which has k of the given values. An expression for the diff'erence between the value of the function and the sum of the first (»" + l) terms of the limiting form of this polynomial will now be obtained. Let n{n-\) ... (ti- r ■\-\) . , „ t\a + mv)=S, + n^J. + ... + '—^^ 67^ + R. Then " r/ ^ r s/ m{m-\) ... {m-r+\) f{a + m9v) -J^-niSji^ - ... — — 6/^ - H Cf. § 3, Formula (1) INTERPOLATION AND NUMERICAL INTEGRATION [en. will be zero for ?n = 0, 1, 2, ... r and n: its derivate of order (r+ 1) will therefore be sero for some value n of 7n lying between the least and the greatest of the quantities 0, r, n, and hence the residue is given by E = „C,+,to'-+' /"■+'' (a + n' n-). It may be remarked that Newton discovered the binomial coefficients when studying the expansion in series of a function by interpolations. He first of all considered the expansions of the expressions (1 -a'-)", (l-.i-)^, (1 -.r-)'-^, ... and deduced from them the expansions of (1— ,<-)'^, (1 -x^fi , ... . Then, by analogy, he obtained the expression for the general term in the expansion of a binomial, and thus established the binomial theorem. determine the value oi/{3-5). X fix) 1 0-208460 Example 1. — From the following table of values of the function f{x) A3 •29242 29029 28789 28523 0-237702 0-266731 0-295520 0-324043 -213 -240 -266 27 In this example a=l, n = ^, w=\ /y =/(!) = 208460. Ph Hence /(3-5)=^ + i 5/, - I S'-^j = 0-208460 14621 27-2 = 0-223106. E.ramplc J.— Corresponding values of x and/(,t:) are given in the fuUowiiig table : 4 4-1 4-2 4-3 4-4 4-5 /(-■) 54-5982 60-3403 66-6863 73-6998 81-4509 90-0171 Obtain the value of /(4-05). 9] FORMULAE OF INTERPOLATION 23 A certain amount of discretion is necessary in the application of Newton's formula. Let us assume, for instance, that the f!;iven observations are finite — say seven — in number, and correspond to the values a - 3 w, a-2w, ... a + Ziv of the argument as in the scheme, § 8. If a + n w lies between a - 3 w and a-'lw, then Newton's formula fory(a + nw) will involve the quantities /_:„ S/_ s, S-/_._;, S\f_ ^, 8^/_i, 8^/_ ^ and 8^/1, all of which occur in the scheme. Thus a good approximation to the value of /(a + nw) may be obtained by means of the formula. But if a + 7iw lies between a + 2w and a + 3w, then in the Newton's formula for /{a + n w) the only known terms are the first two depending on /a and 5/g. In general, these are quite insufficient to determine the value of f{a + n w) to the necessary degree of accuracy. It is therefore expedient in this case to find a formula involving a larger number of known quantities. A suitable one may be obtained by replacing iv in the original formula by - w. Then, since we must replace t>j\ by Similarly, since /,-/-! or -S/_i. we must replace S'-yj by ./; - i/Li +/L, or 671i . In the same way we can determine the functions which must take the places of S*./#, S^/j, .... Newton's formula then becomes n(n-\) ^, , n(7i-l) hi-2) ^. , J{a-mv)=/,-ndJ^.^+ \^ ^ <^-./-i - -^ 37 ■^V-s+ ••• This expansion involves diiFerences along a line parallel to the lower side in the scheme, and it may therefore be regarded as a formula for Backward Interpolation. When a + n w lies between a + 2 w and a + 3 1<;, the formula will depend on the quantities fi-, 5/s, S-/^, S^/|, S^/i, 8f/^, and S^/o) all of which are known: it is therefore very convenient when the value of fix) is required for values of X near the end of the series of observations. 24 INTERPOLATION AND NUMERICAL INTEGRATION [ch. ii Example 3. — From the following table of values oi f{x) determine the values of /( 104 -25). X /(.f) A A--^ 101 1030198 30906 102 1061104 31518 612 103 1092622 32136 618 104 1124758 32760 624 105 1157518 In this example a — 105, n — 0'75, "■ = 1 /y=/(105) = 1157518 /(104-25) 0-75 (0-75 - 1) 0-75 (0-75- 1) (0-75 - 2) -/o-0-75o/.^ + ^^— ^5y_,- --^ 5y-«+... =/o-0-75 5/'_, - 0-09375 5-/_, - 0-0390625 sy. + ... = 1157518-0-75 (3-2760) -0-09375 (6-24) - 0-0390'625 (6) = 1157518-24570 58-5 0-234375 = 1132889-265625. Example 4- — Using the data of Example 2, find the value of/(4-45). 10. Stirling's. Formula of Interpolation. As Newton's formula involves differences -which lie along a slanting line it will frequently not be suitable for calculating the value of the function for a value of the argument near the middle of a series of observations. By a suitable transformation of this formula, however, we may derive another involving differences which lie on a horizontal line : such a formula is termed a Central- Difference Formula. Since 6J\ - 8/_^ = S^Jl we have Hence 8jl =ijlSj[, + }.8'-/o. ), lUj FORMULAE OF INTERPOLATION In the same way the higher differences which occur in Newton's formula may be determined in terms of the differences which occur on the same horizontal line Sis/,,. Substituting these in Newton's formula we obtain /{a + n w) - /; + 7i{ji 3 /; + }, 3;/;} + ^^^^^ — - { S-/o + /^ 5"_/o + -h ^Vo} n {n - 1 ) {to - 2) 3"! {fx oVo + § 8\/l + II S% + i S%} + . . . =Jo + n II oy,3 + — 6-/o + -^^-^ — jx b-J, + ^ SVo n {ii- - 1) (vr — 4j i^ 8% This is known as Stirling's Formula of Interpolation. It should not be used for determining the value of a function for values of the argument near the beginning or the end of a series of observations, as it would not then be possible to obtain a sufficiently high order of differences. Example i.— Using the data of §9, Example 1, determine the value of /(3-5). .r fU-) A A- A3 1 0-208460 29242 2 0-237702 29029 -213 -27 3 0-266731 (28909) 28789 -240 (-27) -26 4 "295520 (28656) 28523 -266 (-26) - 26 5 0-3-24043 28-231 -292 6 0-35-2274 In this case a = 3, IV = 1, 11 - --\ fo = /(3) = 0-26673L /(3-5) = /o + i^5/o + i- 5Vo - tV 5Vo = 0-266731 14455 2- 30 = '281158. Examj}/e ^^ — Using the data of § 9, Example 2, Hnd the value of /(4-21). 26 INTERPOLATION AND NUMERICAL INTEGRATION [cir. ii If n is greater than 0'5, it is better to use backward differences. The corresponding formula may be obtained by merely changing the sign of w in Stirling's Forward Difference Formula, viz. :- — /{a - n iv) = /; -nix of, + ^^ 8\tl - ^^ ^'^ ~ /x S\fl 71- (n- - 1 ) ^ , , 71 (71- - 1 ) («- - 4) .. . 4 ! Of course, either formula may be applied whether it is greater or less than 0-5 : thus, when a check is necessar)^, it will be found better to compute/(.r) by the formula not already used. Example 5.— Use this Bickward Difference Formula to check Example \. Here a = 4, c- = 1, n — ^ /o=/(4) = 0-295520 /(3-5) -/, - i M 5/o + I 5- A + -h 5Vo = 0-295520 - 14328 33 2 = 0-281157. The agreement with the former result is perfect when we take into consideration the fact that the last figure is forced in both cases. Example 4.— Using the data of § 9, Example 2, find the value of/(4-29). 11. Gauss' and Bessel's Formulae of Interpolation. Another central difference formula may be obtained by trans- forming Stirling's formula. Since /^S/, = 8y,- |8V; ,xSi^j^ = S^f^-i8% etc., we see that Stirling's formula may be written in the form J {a +71 IV) =j, + 71 8/^ + ^^, sv^ + ^ — Yr — - ^A (n+l)7i{7l-\) ( 71-2) ^ {71 +2) {71+ I) 71 (71- I) {71 -2) ^^^. ^ 5 ! • I ■■■ In this result, which is frequently known as Gauss Formula, the /x-terms do not occur. The differences which occur in this 10, 11] FORMULAE OF INTERPOLATION 27 formula are those which lie on two adjacent horizontal lines of the scheme in § 8. Again, since t\ +/o ^ - /Vl hav( ^f, = \^^f,-\^f, etc. w€ nave and hence f\a + n w) (^t + 2)(7^+ 1) 7i (« - 1) (?t - 2) ,. + g-^ 0/^+ ... , ■. . 71 (lb -\) ^, 71 (n — \) (n - h) ^.. . {n + l)n{ n-l){n--2) + - ^ ^-^ — /^^A + ;-; 6'yj^ + . . . This is frequently known as BesseVs Formula. In it the differ- ences are those wdiich occur on a horizontal line mid-way between the lines on which /J, and/j stand. The formula most suitable for the construction of mathematical tables, when differences of higher order than the fourth are neglected, is obtained from the first of these on putting When this substitution is made we get ./ (a + 71 iv) = /o + n 6/, + \^ ^ ' 6% \ ' ^J\ W (1 - 71-) (2 - 7l) ,, . 28 INTERPOLATION AND NUMERICAL INTEGRATION [cii. ii Bessel's formula will be found particularly suitable when n is equal to |-, for in this case the formula becomes /(»+j„j.,4+ i^.ov,. (i±Mli^ll(izlVsv>... in which only even differences occur. As with Stirling's formula, none of these results should be employed when the value of the function is required for values of the argument near the beginning or the end of a series of observations. Stirling's formula will be found more convenient when 71 is very small, or when the number of observations given is odd : Bessel's when n is nearly equal to a half, or the number of observations even. The backward difference formula corresponding to Bessel's forward difference formula may be obtained in the same way as in Newton's and Stii-ling's cases. But, except as regards sign, the coefficients in the two formulae will be the same, and therefore the one cannot be used as a check on the other. When a check is required to a result computed from Bessel's formula, it is better to compute it again using Stirling's formula. Example /. —The values of a function for the values -2, -1, u, 1, 2, of the argument are respectively 0-12569, 0-17882, 0-23004, 0-27974, 0-32823. Using Gauss' formula, show that when the argument has the value 0-4, the value of the function is 0-25008. X fix) 1 X' A3 -2 0-1-2569 5313 - 1 0-17S82 - 191 5122 39 0-23004 -152 4^ ,^~ 31 1 0-27974 - 121 4849 2 0-328-23 In this example a = 0, n = O'i, ic = 1 /;, =/(0) = 0-23004. Hence /(0-4) = /„ + 4 5/^ - 0-12 o'Vo 0-056 5^;.. = 0-23004 1988 18-2 = 0-25008. 11, 12] FORMULAE OF I?^TERPOLATION 29 Example 2. — Using the data of §9, Example 2, find the values of/(4"21) and /(4 -29). 12. Evaluation of the Derivatives of the Function. When, as often happens, the values of the differential co- efficients of a function which is given by a series of observations are required, use may be made of the symbolic equality established in § 4, Cor. 4. For, since e"'~^f{:c) = Ef{x) we have w — .J\x) = log E .f{x). The method will be evident from the following example : — Examph 1. — The values of a function corresponding to the values I'Ol, 1-02, 1-03, 104, 1 05 of the argument are 1030301, 1-061208, 1092727, ri24864, 1*157625 respectively. Find the value of the first differential coefficient of the function when the argument has the value 1 "03. Here ,r = O'Ol. Let /o = /(I -01). d Then *'• j^/(l-03) = log £". /(1'03) = log E.f, = log E.E^f, = log(l+A).(l + A)2./„ / A2 A^ X = 1-^ - T+ T -...). (1 +2 A + A--'j./o. Neglecting differences of higher order tlian the 4th, this gives 0-01-^^,/(l-03) = (A + |A2+iA^'-AA^)/„ = {j\-iE+^E'-j\E')/, = i(/3-/i)-TV(/4-/«) = 0-042437 -0-010610 = 0-0318-27 /'(l-03) = 3-1827. Examplt 2. — The values of a function corresponding to the values 1-5, 2, 2-5, 3, 3-5, 4, of the argument are 4-625, 2-000, -0125, -1, 0-125, 4000, respectively. Find the value of the first differential coefficient of the function at the points where the argument has the values 2-5 and 3-0. 30 INTERPOLATION AND NUMERICAL INTEGRATION [en. ii The method is perfectly general. The same result may be obtained by diiierentiating any of the formulae of interpolation already established. For example, if we differentiate with respect to n the formula of Stirling, f{a + n w) = /; + nix 8jl + J^ S\f„ + ^ '^~ /^ ^Vo 7r{n--l) + 4!"" -^"^ - we have wf {a + nw) = tJ. 8/0 + n 8' f, + ^{Sn--l)ix S\fl + ^\,{2rv'-n)S\f,+ ... w-f" {a f n w) = S-_/q -\-nii S-'/q + J-.- (6 rr - 1 ) S-*/', + . . . Putting ?i = these become wf (a) = IX 8/1 -}ix 8% + gV H- ^Vo - • • • w"-/" (a) = S\/l - ^\ 8\fl + ^ 8% -... These formulae are of use in determining the turning values of a function, the values of the argument corresponding to given values of the differential coefficients, the points of inflexion on a curve of which isolated points are given, and also in determining whether a series of observations is periodic. Example 3. — The consecutive daily observations of a function being 0-099833, 0-208460, 0-314566, 0-416871, 0-514136, 0-605186, 0-688921, '764329, show that the function is periodic and determine its period. From the given observations we obtain the table : — X 1 y=f{x) 0-099833 A 108627 A" A" A* 2 0-208460 106106 -2521 -1280 3 0-314566 102305 -3801 -1239 41 4 0-416871 97-265 -5040 -1175 64 5 0-514136 91050 - 6215 -1100 75 6 0-605186 83735 -7315 -1012 88 7 0-688921 75408 -8327 8 0-7643-29 12] FORMULAE OF INTERPOLATION 31 Taking a = 3, /« -/(a) =./"(3), we have, since tv =\, /"(3) = svo- w5*/;-. = -3801-iV41 = -3804 1 ^, 3804 •■• /(3) •/"(^) = - 314566 = - 0-0121. In the same way we find that (i) When a = 4, /, -/(4) /"(4)= -5045, j;^/"(4)= -0-0121 (ii) When a = 5, /„ = /{5) f"(o)= -6-221, jr^f"(b)= -0-0121 (iii) When a = 6, ./'o =./'(6) /"(6)= -732-2, ^-:^-/"(6)= -0-0121. 1 Hence in all cases 77T\/"('t') - -0-0121. Putting y = J {x), this gives the differential equation 1 (Py - -r.,:= -0-0121 y dx- d-y whence -^., + 0-0121 3/ = y - A cos 0*11 a- + 5sin011.T. This result shows that y or f(x) is a periodic function of x, and that its period is 2 tt/ Oil, or 57-12 days. Example 4. — The following table represents one of the Bessel func- tions, / (x) : — 1-0 0-765198 1-1 0-719622 1-2 0-671133 1-3 0-620086 1-4 0-566855 1-5 0-5118-28 1-6 0-455402 dy d-y Find the values of -j-, and -j-,, for x — 1 -3, and hence, by substituting (with four-place accuracy) in the equation ah^ \ dy / 71^ \ d^^ + lc d^ + {^-lc^)y = ^' determine the value of n. INTERPOLATION AND NUMERICAL INTEGRATION [ch. ii Forming the table of differences we have Here Hence 0-765198 0-7 19622 0-671133 0-620086 0-566855 0-511828 0-455402 A - 45576 - 48489 -51047 (-52139) - 53231 - 55027 - 56426 2913 ■2558 2184 1796 1399 A^* 355 374 (381) 388 397 w = 0-1, rt = 1-3, /o =/(a) =/(l-3 0-l/'(l-3) = -0-052139 -i(0-0003Sl) /'(1-3) = -0-52203 •01/" (1-3) = -0-002184 -iV (0-000014) /"(I -3) = -0-2185. Substituting these in we have t4 + - T^ + (1 - — )y = 0, dx-^ X d X \ x^ } " ' •0-2185 + j-:^( -0-5220) + (l - Tf^.) (0-6201) = -0-2185 - 0-4016 + (l - TfWi) (0-6-201) = -0-6201 + (l - TyW,) (0-6201) = 1 + 1 (1-3)-^ ti = 0. Hence the Bessel function in question is /(,(»'). Example 'j. — Show that the following observations are these of a periodic function, and find its period : — X (years) ,/■(.»■) 1 198669 2 0-295520 3 0-389418 4 0-479425 5 0-564642 6 0-644217 7 0-717356 8 0-783327 12] FORMULAE OF INTERPOLATION 33 Example tf. — Observations were taken of the cooling of a boiler, and the results are given in the following table : - t (time in hours) 1-2 19 26 3-3 40 e (temperature °F) 182-3 1752 ITl'l 1643 16U0 (le Determine as carefully as you can the value of -^ wlien < = 2-6 hours. MISCELLANEOUS EXAMPLES. 1. Find a polynomial of the third degree in x which assumes the values 37, II, L3, 29 when x is equal to 2, 1, -1. -2 respectively. Also find the minimum value of the function. 2. A steam electric generator on three long trials, each with a different point of cut off on steady load, is found to use the following amounts of steam per hour for the following amounts of power : — Lb. of steam per hour 40-20 66o0 10800 Indicated horse-power 210 480 706 Kilowatts produced 114 290 435 Find the indicated horsepower and the weight of steam used per hour when 330 kilowatts are being produced. 3. Using first differences only, find the value of /(4817'36) if .<■ 4816 4817 4818 f{x) 69.3974063 694046108 694118146 To what extent may the result be relied on ? 4. Interpolate for .>; = 35-2967 from the following table and state the accuracy of the result :^ X /(-f) 35-28 162021-40139 35-29 162043-22198 35-30 16-2065-04777 5. If the sun's apparent longitude at O. M.N. on 1914 Aug. lat was 128° 2-2' 25"-8 ; on Aug. 3rd 130" 17' 15" -2 ; on Aug. 5th 13-2" 12' 7"-8 ; find the longitude at G.M.N, on Aug. 2nd. 6. At the following draughts in sea water a particular vessel has the following displacements : — Draught in feet - - 15 12 9 6 Displacement in tons - '2098 1512 1018 536 What are the probable displacements when the draughts are 11 and 13 feet respectively ? 7. The following table gives the expectation of life for males of the ages mentioned : — Present age .... 5 10 15 25 35 Expectation of life in years - 50-87 47-60 43-41 35*68 28-64 Find the expectations of 3 men whose ages are 18, 20, 30 years respectively. 34 INTERPOLATION AND NUMERICAL INTEGRATION [ch. ii 8. The HM table of the Institute of Actuaries gives the following values of h :- a: l^ 10 100000 12 99113 14 98496 16 97942 18 97245 20 96223 Obtain the values of Ix when x has the values 11, 13, 15, 17, 19. 9. The sun's apparent right ascension at mean noon is given in the follow- ing table : — June 1 4 34 0-69 6 4 54 31-75 11 5 15 10-64 16 5 35 55-34 21 5 56 43-18 26 6 17 30-99 Determine the apparent right ascension for 1914 June 15d. 16h. 10. Let the number of recruits in the United States army whose height did not exceed x inches be denoted by y. Determine the number whose height exceeded 66 inches, but did not exceed 66| inches, if 62 618 63 1855 64 3802 65 6821 66 10296 67 14350 68 17981 69 21114 70 23189 71 24674 m's apparent declination 3 at mean noon ^ 314, May 1 + 14° 54' 19" -6 3 + 15° 30' 25" -0 5 + 16° 5' 28" -4 7 + 16° 39' 27" -6 9 + 17° 12' 20" -3 11 + 17° 44' 4" -3 Determine the hourly difference in ? for 1914 May 5 and 7. MISC. Exs.] FORMULAE OF INTERPOLATION 35 12. In the following table s denotes the space described at the time t by a point in a mechanism. Tabulate the values of the velocity and the acceler- ation, and obtain their values when t = OAo. t (time in sees.) O'l 0-:2 03 0-4 0-5 O'O O'T O'S 0-9 1-0 1-1 « (space in feet) O-OS-' 0-12 0-2S1 0-J2j 0-804 1-200 1-725 2-229 2-753 3-2S2 3-Sll 13. There is a piece of mechanism whose weight is 200 lbs. The following values of s in feet show the distance of its centre of gravity from some point in its straight path at the time t seconds from some era of reckoning. Find its acceleration at the time < = 2-05, and the force in pounds which is giving this acceleration to it. t 200 2-02 2-04 2-06 2-OS 2-10 s 0-3090 0-4931 0-6799 0-S701 10643 1-2631 14. The following values of p (lbs. per sq. foot) and d"C are for saturated steam. Calculate with as much accuracy as the numbers will allow the value dp of ^for ^=115T. e 105 110 115 120 125 p 2524 2994 3534 " 4152 4854 15. From the following table of the secular variations of the magnetic declination determine the chief period : — Year. Declination. 1540-5 - 7° -463 1580-5 - 8° -500 1620-5 - 6°-250 1660-5 - 0°-500 1700-5 + 8° -228 1740-5 +15° -755 1780-5 +20° -832 1820-5 +22° -380 1860-5 +19° -364 1900-5 +14° -833 [CH. Ill CHAPTER III THE CONSTRUCTION AND USE OF MATHEMATICAL TABLES 13. Interpolation by First Differences. All mathematical tables are alike in that they give numerical values of the function for certain values of the argument, which are so chosen that intermediate values of the function may be derived by interpolation. The interval of the argument varies in the different mathematical tables, its choice in each case being determined by the rapidity of variation of the function. In most of the published elementary tables, which are intended for the use of students, the interval of the argument is so chosen that the second difference is smaller than one unit in the last place of decimals retained : on this account it is not necessary to consider second differences in the interpolation, and the function can be calculated for any value of the argument intermediate between two of the tabulated values by a simple proportion in the following manner : — Consider a function f{x) which is tabulated for the values .To, Xq + w, Xq + 2w,.... of the argument. Then, since the function increases uniformly, when digits beyond a certain place of decimals are neglected, we have /(xq + w') -/(Xq) _ «;' f{Xo + w)-f{Xa) w where 0„2 1 o \ ' jo \ y^ -X" + —J 27- The Euler-Maclaurin Expansion. In the foregoing the value of the integral is expressed as a series involving differences. In many cases in which the function is known analytically, it will be found more convenient to use what is known as the Euler-Maclaurin expansion which employs 'f differential coefl&cients instead of differences. The rapidity of ' convergence is of the same order as that given by Bessel's formula. It may be obtained as follows : — Since = ^[J{a+w)-j\a-w)] W' = h {/(«) + "^/' («) + y/" («) + ••• -fia) + wf (a) - ^/" (a) + . . . } = w/'{a)+ '^/"'(a)+... and similarly /x Sy, = tvj ' {a + r to) + — y '" {a + r to) + w'' Ij. B^f,. = vff" (a + r iv) + —./"*'■' (a + r w) 4 we have by substituting in the formula If*'" fix) dx=\j\ +j, +y; + . . . +./;._, + y. 26, 27] NUMP:RICAL INTEGRATION 77 the result ^j"^' " ./■ (^^ dx = 1 /; + /; + .. . + /;_, + 1 /;. -TroiTo"'-M/.'^'-/o'") + ... This is the EuIerMacIauriii* expansion. The coefficients of the various items in it are proportional to the well-known Bernoullian numbers. For if f we put f{x) = e^-" in the formula and at the same time chansfe the limits to a ~ w and a + tv we iret e-" dx = ifo +/] + y: + A w {fj -/;) + b «- (/;'" -./;'") where A, B, ... are constants. But — I e" " dx ^ — I W Ja-w W J-„. dx — (e'" - e-'"). w Hence — (e'" - e-"') = .i/o + /; + 1/, + ^ 7r (/,' -./;/) + B >r (./:'" - /;,'") + . . . = I e-'" 4-1+1- e'" + A w (e"' - e'"') + B ir^ (e'" - e~'"} + ... from which it is easy to deduce that ^; ::. :, = I - A iv - B >r* - ... w to , , „ r. ■ T— . cot -T—. = \ - Aw- - Bijcr - ... '1% 2 1 But the Bernoullian numbers are the quantities ^j, B.., ... in the expansion and therefore .4, ^, ... are proportional to B^., B.,, ... * For a rigorous discussion, c/". Whittaker and Watson's Modern Analysis, p. 128. t Poisson, Mem. de I'Acad. des Sc, lb23. 1 1 78 INTERPOLATION AND NUMERICAL INTEGRATION [ch. iv Example i. — Calculate — to 8 places of decimals, and check the J 95 '^ result by computing it as - loge(l - ygy) by the logarithmic series. Here /o = tfV, ^"=^, r = ^- Then Jp5 ^ 2 ■ 95 ' 96 97 98 99 2 100 T^ ^ 100- 95-7 ^ - " V " '*»' 95V w -O-OOOlOO0OO\i/ -0 000000010 ^ - 0-005-63158 1^ + 0-000110803/ ^'-"\ +0-000000012/ 010416667 0-010309278 0-010204082 0-010101010 0-005000000 = 0-051294195 -0-000000900 = 0-051293295 Again , /, 5 \= ^^. (0-05;^ (0-05)3 (005)^ (0-05)3 (o-05)« -i°g^(i-ioo; o'^«+-T-+-^-+— r-^~5-^-6~~ = 0-050000000 0-001250000 0-000041667 0-000001563 0-000000063 0-000000003 = 0-051293296 The two results will be seen to agree to 8 places of decimals. fl05 d.r Example 2. — Evaluate I -7- to 8 places of decimals. IT Example 3. — Approximate to tlie value of 1 co^-xdx. 28- An Exceptional Case- Legendre* remarked that if the odd differential coefficients above a certain order take the same value at both limits, the formula fails to give an accurate value of the definite integral. For instance, if the integral under * Fonctions Elliptiques, Vol. 11., p. 57. 28, 29] NUMERICAL INTEGRATION 79 discussion were the elliptic integral of the second kind taken between limits and 7r/2, viz. rTT's I J l—JS^ sin^ X dx, Jo all the odd differential coefficients become zero at both limits, and thus the value of the integral would appear to be given by i/o +/i +/2 +•• ■+./;-! +i/;-. But the value given by this expression is not in close agreement with the true value of the integral. The reason of this is that the numerical coefficients introduced in the successive differentiations increase without limit, so that each term takes the form oc x 0, which is, of course, inde- terminate. Better methods of evaluating the three kinds of elliptic integrals will be given in another tract in the present series. 29. The Formulae of Woolhouse and Lubbock. We shall next derive from the Euler-]\Iaclaui-in expansion a formula which is of use when it is required to evaluate the sum of a large number of terms. If the interval between consecutive values of the argument in the Euler-Maclaurin formula is w, we have and if this interval is w j m the formula becomes - f{x) dx= ./; + /■ , +/ , + . . . +./;. - 1 (/o +/,.) n- J f, mm -A^(./;.'-y:')+Tk|(./;--y;"')+... Subtracting the latter from m times the former we deduce /o+/i +J I + .-• +J>m{f,+J\+... +/.) - '"-^ (/;+/.) 12 m^'^' '^'■> m* - 1 tv' ^ .,„ ^.,„^ 80 INTERPOLATION AND NUMERICAL INTEGRATION [ch. iv Tliis is Woolhouse's Formula. OV)viou,sly by using the right- hand expression instead of the left-hand one, the labour involved in determining the sum ^ 4-/2 +./" 1 + ■■■ + A- is enormously reduced. A somewhat similar formula, in which differences are employed instead of differential coefficients, is due to Lubbock. Expressing the derived functions in Woolhouse's formula in terms of the successive differences of /q, J\, ■■■/,■ we obtain the formula in question, m — \ 12m '"5 2 24 ?/t (m-^-l)(19m^-l) I 20 7H' ' ■' ^ (m--l){9m-- 1) 480 m" {8V>-2 + ^V2) 3j0 1 E.rample 1. — To determine 2 — • '/! = 300 " Let w = 10, ?o=10. Then from Woolhouse's formvila we have 350 J_ M +J_ +J_ +J_ +J_ + J_\ rMw ^ "^^^300 310 320 330 340 350/ 2 VSOO 350/ 12 V 350- 300V = 018512761-0-02785714 2432 = 0-15724615 Losing Lubbock's formula we get 2 IT = 10 (/o +/i + • • ■ +^5) - s (./'o +./;) - -/-vo (5 n - 5/a ) - i?^ roV\ + s% > n = 300 720000 ^ -'-V •'4' 480000 ^"-'-i ^ "-'2^ ^0-18512761 -0027S5714 193S 487 3 3 = 15724616 29, 30] NUMERICAL IXTEGRATION 81 To test the accuracy of these answers, the reciprocals of the numbers 300, 301, 302, ..., 350 were taken from Barlow's tables, and summed by means of a comptometer. The result was found to be 0'15724616. The closeness of the approximation is remarkable. 175 1 ^.rrt??ip/e ^.—Obtain the value of 2 — • n = ioo " Examplf 3. — By giving n increasing positive values, show that T 1 30. Formulae of Approximate Quadrature The expansions which have been obtained in the preceding articles furnish the value of an integral to any required degree of accuracy, provided a sufficient number of terms is taken. We shall now consider certain formulae which are in theii' nature only aj>proxiinate, and which, though frequently useful, cannot be used where great accuracy is required. Suppose that it is required to integrate a f unction /(.x) between limits, which may, without loss of generality, be taken to be - 1 and + 1. We shall endeavour to obtain expressions of the type ^o/(^'u) + H,f(K) + ...+ iiJiK) which represent the integral as closely as possible, where Aq, h^, ... h„ are (n+l) values of x within the range of integration and h^, h^, ... h,„ Bq, //j, ... Zr„ do not depend on the function /(.r). Let us first choose B„, B^, ... B,„ so as to make the formula strictly accurate, so long as f{x) is a polynomial of degree less than n + \. Let F{x) = (.X - /g {x-h,) ... (x - h„). F(x) Then -— is a polynomial of degree less than n+l, so that the X - h,. •' ° F (x) formula must be strictly accurate when for /(xj we put — ^ . " X - h,. This gives F(x) B.L BrF'{K). 82 INTERPOLATION AND NUMERICAL INTEGRATION [cii. iv Hence, whatever hg, //j, ... h„ may be, in order tliat the above con- dition may be satisfied, the values of Hq^ IIj, ... H^ must be given by In practice it is convenient to make the intervals between successive values of the argument equal to each other. Suppose, then, that A,i, ^i, ... A„ are chosen so as to divide the whole range of integi'ation into n equal parts. Then A^ = - 1 1 1 K = -i + — n /i., = -1 + 1 n H, '■ {h, -K) ... (h, - /i,,„, ) {h, - A,+, ) . . . (A, - A J X (x- Ao) • • ■ (•>' - K-i)(^' - h,.+i) ■ . . (x - h„) i r 2(r-l) 2- -2 -4 _ 2 (n - r) ?i 7J n n ti n ( - 1)"-^ ("— Y\-! (ri-r)! Now^ let t = -^{x-\-l). Then //.= ^~^^" '^^ , f%(^-l)... (<-r+l)(<-r-l) ... (t-n)dt. n . r\ [n - r) I Jo Hence, finally, if we put 1*0 = a, ... h^ = a + r w 30] NUMERICAL INTEGRATION wo obtain \y{x) dx^j: II, f {a + r w) J a )-=o whei'i //,= ^-^— ^ -^ rV<-l) ... (t-r+l){t-~r-l) ... {t-n)dL r I {71 - r) ! Jo This is known as Cotes Formula of Quadrature. By means of it we may calculate the value of the integral without first of all forming a table of diifei'ences. Special Cases. (i.) Let?i=l. Tlien f 4 w f{x) d X = nj{a) + II, /{a + n-) — iv r^ where H, = — — - {t-l)dt = h rr ! I ! J Hence [^'\f{x)dx = ^ {f(^a)+f{a + iv)]. This is known as the Trapezoidal Rule. (ii.) Letw = 2. Then where r^' "'./•(.r) d X = HJ{a) + HJ (a + w) + H.J(a + 2 w) Hence 1+2 19 ?0 ,, 4 7t' ./ {x) dx =^ — y («) + ^y (a + «;) + y y (a + 2 w). This is called Simpson'' s * i??i/e or the Parabolic Rule. * It was, however, first given by James Gregory. 84 INTERPOLATION AND NUMERICAL INTEGRATION [ch. iv Corollary. — If we divide the whole range of integration, say {a, b), into 2 n equal parts, and apply Simpson's Rule to each of them, the value of the integral may be expressed as / = -g^ {yo + 2/2n+2(?/o + y4+ ...+?/.2„-:)+4(t/,+2/:;+..-+?/2„_l)} where i/q, y^, ... y..,, are the values oi /(x) "when x has the values h - a 2 n form. The rule may be stated thus : — Divide the area into an even number of strips by equidistant ordinates y^„ t/j, ... v/.j„ : to the sum of the extreme orJinates add ttvice the sum of the other ordinates with even suffixes, and four times the sum of the ordinates with odd suffixes, and multiply this total by one-third of the common distance betiveen the ordinates.* (iii) Let n = 3. Then fa+3w) f{x) dx = HJ{a) + II J {a + ic) + HJ{a + 2iv) + HJ{a + 3 to) where H, = -^^^^^{t-\){t -•2)(t -?,) dt = %w Hence T"^ VH dx ^%w {f(a) + 3/(« + tw) + 3 /'(a + 2 w) + /{a + 3 w)], which is often termed the Three- Eighths Ride. (iv) Let n = 6. Tlien it may be shown, as in the above special cases, that f(x)dx = lU {Uf(a) + 2\Gf{a + w) + 27 f{a -{- 2 tv) + 272/(a + 3 m;) + 27 f {a + 4 w) + 216 f {a + 5 w) + 41/(a + 6M7)}. * The reader is warned that the notation varies in the different text- books, so that in all cases aa examination of the notation used should first of all be made. i: 30,31] NUMERICAL INTEGRATION - 85 Adding ih S!/*(« + 3 w) = -^Ijj {/{a) - 6/(a + w) + 15/(ffl + 2 w) - 20,/"(« + 3 w) + 1 bj\a + 4 w) - 6/(a + 5 w;) +/(a + 6M;)} we obtain Ttt + Ow /(,«) ci X- + ^i^ S7 (a + 3 w;) = ^io {^V\<^) + ^lOy^a + t^) ^^ + 42/(a + 2 w) + 252/(« + 3 m;) + 42/(a + 4 w;) + 210/(a + 5 w) + 42/(a + 6to)} = TO {/■(«) + 5/(« + ^<^) +/(« + - «*^") + 6/(a + 3 Mj) +/(« + 4 ztf) + 5/(a + 5 w) +/(« + 6 i«) } . Neglecting the term -^\-^h^ f {a-^?> lo), which will be fairly small, we get Weddl&s Rule of Quadrature for 6 intervals. J{x)dx = j% { /(a) + 5/ {a + lo) +/(a + 2 w) + 6/ (a + ow) "^ +/{a + 4 w) + of {a + 5 iv) +f{a + 'oiv)]. 31. Comparison of the Accuracy of these Special Formulae. We may compare the accuracy of these formulae by evaluating, by each n fix of these methods, the definite integral I ~^^, , the correct value of which is log^7 = 1-94591. (i) Using the Trapezoidal Rule six times, i.e., dividing the region into six equal portions and applying the rule to each of them, we get f7 dx — = \f{(«) denotes any polynomial of degree <(n+ 1), and where F{x) = {x - h,} {x - h,) (x -h.^ ... {x - h„) we see that the condition to be satisfied may be written f F{x) {■'>') dx = and, conversely, this condition is sufficient to determine the Legendre polynomials. Hence we have (save for a multiplicative constant) F{x) = P„_,,(x) and therefore the quantities Z/^, /tj, ... h„ must be the roots of the equation ^,,-1-: (^■) - 0. The foUuwing are the values of the constants h^, h^, ... h , //„, H^, ... H in the simplest cases. (i) Whenn^l. K = 1 "v^3 K- 1 (ii) When7i = 2. /in = - ^'k li, = /l, = •Ji H, = 1 H, - 1 H, _ 5 — y H, = 1! H., _ 5 — y * Cf. Whittaker and Watson's Modern Analysis or Byerly'a Fourier's Series and Spherical Harmonics. 88 INTERPOLATION AND NUiVlERICAL INTEGRATION [ch. (iii) When?t = 3. K= - s^( 'i + is v^30) H, = h- rk ^/30 '>i= - V ( ? - i's v'30) Hi = i + sV v/30 ''i = >/{■?- Tft V30) ^2 = i + gV x'30 '^.i = s'( f + y^ V^30) i^:, - i - 3^6 V30 (iv) When w = 4. K= - - s'({; + ^n'70) K= - • N'(^-i^V70) h.,= h,= xAa-Av^70) ^4 = N^(H6lTV^70) ^^ ^ 322 + 13 v/70 900 900 128 225 322+13\/70 900 322-13V70 ^^= 900 ^' 900 Example I. — Determine the val (i) Using 3 urdinates r' ^ = 5 _J__ J-] 3+^ 9 -3- ^1-4 in J-l 3+« 9 ■ 3 9 • 3+ x'f = 0-69312165. (ii) Using 5 (jrdinates 3-22- 13v/70 1 3-22+ 13 V 900 1 3+ V(|- 70 6-3 1 900 1-28 + 225 322- 0-693147157. 3- N^(^ . 1 3 3 + -13V70 900 + 6^V70) ' 22+13v70 900 L 3 - s!{^ - h sHO] v/70) 3+ v'(| + ff-ov70) The correct value is 0-69314718, so that using 3 ordinates the error is in the 5th decimal place : using 5 ordinates, in the 8th decimal place. ^xa»ip/e i^.— Determine the values of (i) I -r4> ' MISC. Exs.] NUMERICAL INTECJRATION MISCELLANEOUS EXAMPLES. L The values of a function /(.r) for the values 1-050, 1060, 1-070, 1-080, 1-090, 1-100 of the argument are 1-25386, 1-26996, 1-28619, 1-30254, 131903, 1-33565. Use Newton's formula of quadrature to show that the value of fi-inu .fix)dx J 1-050 is 0-06 472. 2. The ordinates of a plane curve are of the following lengths, 3-5, 4-7, 58, 68, 7-6, 8-1, 80, 7-2 feet, and they are 4 feet apart. What is the area included between the two end ordinates? and what is the distance of the centre of gravity of that area (i) from the extreme left-hand ordinate, and (ii) from the base line ? 3. Approximate to the values of the following definite integrals : — (i) r dx.^{S-.v + x^) (ii) c/a i: IT f? dx 4. Use (i) Lubbock's formula, (ii) Woolhouse's formula, to determine the value of an annuity-certain for 36 years, interest being at the rate of 3%. 5. Verify that the area of the curve y~A + Bx + Cx- + Dx^ between the limits x = h and x= -h is equal to the product of h into the sum of the ordinates at ,r = A/^'3 and x= -h/^/S. In the case of the curve i/=f(x) = A + Bx+ Cx'^ + Dx^ + Ex^ ■\- Fx-' verify in like manner that the area between x = h and x= -h is equal to { 5f(h Jl ) + 8/ (0) + 5/( - h J J ) } /t/9. 6. Corresponding values of .r and y are given in the following table, the unit being the inch. Find the volume of the solid of revolution formed when this curve revolves about the .r-axis. Find also the centre of gravity of this uniform solid. 4 8 12 16 20 0-50 2-58 4-41 5-40 5-10 0-50 90 INTERPOLATION AND NUMERICAL INTEGRATION [ch. iv 7. A body of 200 lbs. is acted on by a variable force F. No other force acts on the body. When the body has passed through the distance x feet the force in pounds is as follows : — 0-1 0-2 3 0-4 Oo 0-6 07 0-8 0-9 1 20 21 21 20 19 18 5 ISO 1.3-5 9 4-5 Determine the work done upon the body from .r = to x = 0'4. If the velocity is when x = 0, what is its value when a; = 0'4 ? 8. The semi-ordinates in feet of the deck plan of a ship are respectively 3, 16-6, 25-5, 28-6, 298, 30, 29-8, 29-5, 28-5, 24-2, and 6-8, the common interval between them being 28 feet. Find the area of the deck. 9. The "half" ordinates in feet of the mid ship section of a vessel are 12-5, 12-8, 12-9, 13, 13, 128, 124, US, 104, 68, Oo respectively, and the common interval is 2 feet. Find the area of the whole section and the position of its centre of gravity. 10. The area of a ship's load-water plane is 8000 sq. ft. ; the body below the load-water plane is divided into six portions by equidistant horizontal sections, 3 feet apart, whose areas in sq. ft. are respectively 7600, 7000, 6000, 4500, 2800, and 250. Find (a) the displacement in tons, and (b) the number of tons which must be taken out of the ship to lighten her 4 inches. Printed by Lindsay & Co., 17 Blackfriars Street, Edinburgh 14 DAY USE RETURN TO DESK FROM WHICH BORROWED LOAN DEPT. This book is due on the last date stamped below, or on the date to which renewed. Renewed books are subject to immediate recall. "aTFe^saWJ REC'D LD T-^m- HECEIVED MAKixms- 'VMdyWpW DEC 2 7 '66 -4 PM LOH DEPT. RECrP \ TX TfWHM96e- ^0^^ar'6St0 ^REOU-LO wty-'S^^i LD 21A-50m-9,'£ (6889sl0)476B General Library University of Calif. Berkeley