@prefix vivo: . @prefix edm: . @prefix ns0: . @prefix dcterms: . @prefix skos: . vivo:departmentOrSchool "Science, Faculty of"@en, "Mathematics, Department of"@en ; edm:dataProvider "DSpace"@en ; ns0:degreeCampus "UBCV"@en ; dcterms:creator "Creemer, Albert Lee"@en ; dcterms:issued "2011-11-30T16:36:09Z"@en, "1962"@en ; vivo:relatedDegree "Master of Arts - MA"@en ; ns0:degreeGrantor "University of British Columbia"@en ; dcterms:description """The purpose of this thesis is to study the factors involved in determining a most efficient method for the numerical integration of the differential equation x' = f(t,x) . By "a most efficient method" we mean a method requiring a minimum of computation to obtain a solution within prescribed error bounds. We outline two computational procedures and derive estimates for the propagated error of a general multi-step method when based on either procedure. These estimates, lead us to conclude that a stable single-iterate procedure, involving one evaluation of f at each step, will determine a solution most efficiently. In particular, this procedure based on Adams formulas is recommended. Experimental results support our conclusion in all stable cases. However, these results also indicate that the role of stability in the choice of a most efficient procedure is in need of further investigation."""@en ; edm:aggregatedCHO "https://circle.library.ubc.ca/rest/handle/2429/39377?expand=metadata"@en ; skos:note "EFFICIENT .METHODS. FOR THE NUMERICAL INTEGRATION OF ORDINARY DIFFERENTIAL EQUATIONS \"by ALBERT LEE CREEMER B.A., U n i v e r s i t y of B r i t i s h Columbia, 1956 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF ARTS i n the Department of MATHEMATICS We accept, t h i s ' t h e s i s as conforming ,.to:.-the' r e q u i r e d standard THE UNIVERSITY OF BRITISH COLUMBIA April,.1962 In presenting t h i s t h e s i s i n p a r t i a l f u l f i l m e n t of the requirements f o r an advanced degree at the U n i v e r s i t y of B r i t i s h ' C o l u m b i a , I agree t h a t the L i b r a r y s h a l l make i t f r e e l y a v a i l a b l e f o r reference and study. I f u r t h e r agree that permission f o r extensive copying of t h i s t h e s i s f o r s c h o l a r l y purposes may be granted by the Head of my Department o r by h i s r e p r e s e n t a t i v e s . I t i s understood t h a t copying or p u b l i c a t i o n of t h i s t h e s i s f o r f i n a n c i a l g a i n s h a l l not be alloived without my w r i t t e n permission. Department of Mathematics The U n i v e r s i t y of B r i t i s h Columbia, Vancouver 8, Canada. Date_ A p r i l 13, 1962 ABSTRACT The purpose of t h i s t h e s i s i s t o study the fac t o r s , i n v o l v e d : i n determining a most e f f i c i e n t method f o r the numerical i n t e g r a t i o n of the d i f f e r e n t i a l equation x' = f ( t , x ) . By \"a most e f f i c i e n t method\" we mean a method.requiring a minimum of computation to o b t a i n a s o l u t i o n w i t h i n p r e s c r i b e d e r r o r bounds. We o u t l i n e two computational procedures and de r i v e estimates f o r the propagated e r r o r of a general m u l t i - s t e p method when based on e i t h e r procedure. These estimates, l e a d us to conclude t h a t .a s t a b l e s i n g l e - i t e r a t e procedure, i n v o l v i n g one e v a l u a t i o n of f at each step, w i l l determine a s o l u t i o n most e f f i c i e n t l y . In p a r t i c u l a r , t h i s procedure based on Adams formulas, i s recommended. Experimental r e s u l t s support our c o n c l u s i o n i n a l l s t a b l e cases. However, these r e s u l t s a l s o i n d i c a t e t h a t the r o l e of s t a b i l i t y i n the choice of a most e f f i c i e n t procedure i s i n need of f u r t h e r i n v e s t i g a t i o n . I hereby c e r t i f y t h a t t h i s a b s t r a c t i s s a t i s f a c t o r y . i i i . . ACKNOWLEDGEMENTS The author wishes, to acknowledge his•indebtedness to Dr. T-E. H u l l f o r advice and guidance r e c e i v e d i n f o r m u l a t i n g and pre p a r i n g t h i s t h e s i s . Thanks are a l s o due t o the s t a f f of the Computing Centre of the U n i v e r s i t y of B r i t i s h Columbia. T h e i r e f f o r t s made p o s s i b l e the p r a c t i c a l v e r i f i c a t i o n of t h e o r e t i c a l r e s u l t s . i v . TABLE OF CONTENTS CHAPTER I I n t r o d u c t i o n and Summary 1 CHAPTER I I Estimates of the Propagated e r r o r 3 P ( E C ) m Procedure, Flow Chart '. h P E ( C E ) m Procedure, Flow Chart 6 E r r o r Estimate, PE(CE) m Procedure - 7 E r r o r Estimate, P ( E C ) m Procedure . . 12 CHAPTER I I I Comparison and Choice of Methods . . . . . . . . . . 15 CHAPTER IV Numerical R e s u l t s 20 Comparison of S i n g l e - I t e r a t e Procedures . 21 Conclusion 22 \"BIBLIOGRAPHY 23 CHAPTER I INTRODUCTION AND SUMMARY An important c l a s s of methods f o r the numerical i n t e g r a t i o n of the i n i t i a l value problem: x' = f ( t , x ) , x ( t o ) = X Q i s the c l a s s of i n t e g r a t i o n formulas of the p r e d i c t o r - c o r r e c t o r type. The purpose of t h i s study i s to i n v e s t i g a t e c e r t a i n p r o p e r t i e s of t h i s wide c l a s s , w i t h the u l t i m a t e g o a l of choosing the method r e q u i r i n g the l e a s t amount of computing f o r the s o l u t i o n of ( l ) to a p r e s c r i b e d accuracy. Our main concern i s f o r problems i n v o l v i n g a complicated f u n c t i o n f. The computing time f o r the numerical s o l u t i o n of such problems i s almost e n t i r e l y determined by the number of e v a l u a t i o n s of f r e q u i r e d by the i n t e g r a t i o n technique over a s p e c i f i e d i n t e r v a l of t . Under these circumstances, the many c a l c u l a t i o n s of f r e q u i r e d by the Runge-Kutta method would be p r o h i b i t i v e and f o r t h i s reason, i n searching f o r a most e f f i c i e n t means of i n t e g r a t i o n , we c o nsider only the c l a s s of m u l t i - s t e p p r e d i c t o r - c o r r e c t o r formulas. To achieve our aim, we begin by d e r i v i n g expressions d e s c r i b i n g the behaviour of the e r r o r introduced by the p r e d i c t o r - c o r r e c t o r procedure. These expressions w i l l be employed i n a study of the c h a r a c t e r i s t i c s of such methods and i n p a r t i c u l a r , t o determine c r i t e r i a f o r the s e l e c t i o n of a most e f f i c i e n t i n t e g r a t i o n procedure. These c r i t e r i a l e a d us to consider m u l t i - s t e p i n t e g r a t i o n techniques i n v o l v i n g a s i n g l e e v a l u a t i o n of f at each i n t e g r a t i o n step. At l e a s t as long as there i s no d i f f i c u l t y w i t h s t a b i l i t y , these techniques are found t o be best i n the sense th a t over a s p e c i f i e d i n t e r v a l of t , r e l a t i v e l y few such c a l c u l a t i o n s are r e q u i r e d , i n order t o meet p r e s c r i b e d requirements of p r e c i s i o n . In p a r t i c u l a r , techniques of t h i s type based on Adams formulas of h i g h order, prove to be h i g h l y . e f f i c i e n t and convenient. S i m i l a r p r e d i c t o r - c o r r e c t o r procedures f o r the numerical s o l u t i o n of ( l ) have r e c e n t l y been proposed by Hamming [l] , Milne and^Reynolds [V], Nordsieck [5j , and Ra l s t o n [b] . These methods r e q u i r e only a s i n g l e a p p l i c a t i o n of the c o r r e c t o r formula but u s u a l l y two e v a l u a t i o n s of f at each step of i n t e g r a t i o n . To some extent, these methods a l s o r e q u i r e s p e c i a l techniques which are not g e n e r a l l y a p p l i c a b l e . In the next chapter, c e r t a i n p r e d i c t o r - c o r r e c t o r , techniques are descr i b e d and.estimates are d e r i v e d . f o r the propagated e r r o r . In Chapter I I I , observations are made on these estimates, to draw.conclusions as to the r e l a t i v e c h a r a c t e r i s t i c s of s e v e r a l p r e d i c t o r - c o r r e c t o r methods. F i n a l l y , i n Chapter XV, experimental evidence- i s o f f e r e d i n support of the accuracy of the e r r o r estimates and,the v a l i d i t y of a s s o c i a t e d c o n c l u s i o n s . These r e s u l t s a l s o i n d i c a t e the need f o r a more complete i n v e s t i g a t i o n of the r o l e of s t a b i l i t y i n the choice of a most e f f i c i e n t method. CHAPTER I I ESTIMATES OF THE PROPAGATED ERROR 3-As mentioned p r e v i o u s l y , should the f u n c t i o n f on the right-hand .side of ( l ) be d i f f i c u l t t o evaluate, the number of such e v a l u a t i o n s becomes the f o r e m o s t . f a c t o r i n determining the c o s t of numerical i n t e g r a t i o n . An e f f i c i e n t m u l t i - s t e p method would be one which r e q u i r e d r e l a t i v e l y few ev a l u a t i o n s of f i n i n t e g r a t i n g ( l ) to w i t h i n a p r e s c r i b e d accuracy,.across a s p e c i f i e d i n t e r v a l of t . To determine i f a p a r t i c u l a r m u l t i - s t e p formula w i l l meet p r e s c r i b e d requirements on p r e c i s i o n , we need an estimate of e r r o r , an expression g i v i n g a reasonable approximation .to the e r r o r propagated through numerical c a l c u l a t i o n s . However, our purpose, i s t o s e l e c t an e f f i c i e n t method f o r the i n t e g r a t i o n of ( l ) . For t h i s reason, we s h a l l r e s t r i c t our use of e r r o r estimates to q u a l i t a t i v e a p p r a i s a l s of one method r e l a t i v e to other methods r a t h e r than as a means f o r p r e c i s e q u a n t i t a t i v e measures of accuracy. In d e r i v i n g an e r r o r estimate, we c a r e f u l l y f o l l o w the computational procedure o f . a p p l y i n g a p r e d i c t o r - c o r r e c t o r p a i r . We s h a l l consider two such procedures, each l e a d i n g t o a d i s t i n c t form of the e r r o r . I t w i l l be seen th a t a f t e r many.iterations of the c o r r e c t o r formula,.these procedures become equi v a l e n t . P ( E C ) m Procedure ( P r e d i c t , Evaluate,.Correct, m i t e r a t i o n s of the c o r r e c t o r ) . A l o g i c a l f l o w diagram best describes the c a l c u l a t i o n and f o r t h i s procedure i s shown i n Figure 1. PREVIOUS STEP PREDICT CORRECT J L TEST P(EG)™,Procedure Figure 1. 5-The 'PREDICT1 and 'CORRECT' stages represent the c a l c u l a t i o n s i n v o l v e d i n a p p l y i n g the r e s p e c t i v e formulas. At the 'EVALUATE f stage, recourse i s made t o some subprogram f o r the c a l c u l a t i o n o f f ( t , y ) , where y i s the approximation t o the t r u e s o l u t i o n x ( t ) , obtained at the immediately preceding 'PREDICT'.or 'CORRECT' stage. The 'TEST', a p p l i e d at the f i n a l stage of the procedure, may determine i f two successive i t e r a t e s , y^-''\"^ and y ^ ^ , s a t i s f y some c r i t e r i o n of cl o s e n e s s , or the t e s t may merely represent a count of i t e r a t i o n s t o some f i x e d number, m. Throughout t h i s d i s c u s s i o n , we s h a l l assume f o r s i m p l i c i t y , t h a t t h e • l a t t e r i s the case. I f m does vary from.step to s t e p , . t h i s assumption o f f e r s a reasonable approximation. P E ( C E ) m Procedure. A fl o w diagram f o r t h i s procedure i s shown i n Figure 2 . The stages of c a l c u l a t i o n shown i n Figure 2 may be rearranged i n var i o u s ways to give e q u i v a l e n t c o n f i g u r a t i o n s . In any case, the P E ( C E ) m procedure, as compared to a P ( E C ) m procedure w i t h the same c r i t e r i o n of 'TEST', i n v o l v e s an e x t r a e v a l u a t i o n of f . In d e r i v i n g an estimate f o r the e r r o r propagated by a m u l t i - s t e p p r e d i c t o r - c o r r e c t o r method, we consider any f i x e d number m of i t e r a t i o n s at e v e r y . i n t e g r a t i o n step. However, our main concern w i l l be the case m=l. As m in c r e a s e s , the accuracy and e f f i c i e n c y of the P ( E C ) m and P E ( C E ) m techniques w i l l d i f f e r only s l i g h t l y . However, f o r small m, the d i s t i n c t i o n becomes important. C o n s i d e r i n g both accuracy and e f f i c i e n c y , i t w i l l be shown t h a t provided both are s t a b l e , the PEC technique i s s u p e r i o r t o PECE. For t h i s purpose, we begin by d e r i v i n g an estimate of the e r r o r propagated by a p r e d i c t o r -c o r r e c t o r method when c a l c u l a t i o n s are extended by the P E ( C E ) m procedure. Such an estimate i s a l s o i n s t r u c t i v e , i n t h a t i t o f f e r s an i n s i g h t i n t o the behaviour of the approximation as m grows l a r g e . PREVIOUS STEP 5_ PREDICT ( C E ) m Procedure Figure 2. 7-E r r o r estimate,-PE(CE) m procedure. Let be approximations to the tr u e s o l u t i o n , obtained from ;the mth . i t e r a t i o n of the c o r r e c t o r formula.at the p o i n t s Then at t = t Q + (ft+k)h, the p r e d i c t e d o r d i n a t e i s k - i - »« (2) I f f i s evaluated .using the ( q . - l ) s ^ i t e r a t e of the corr e c t o r , formula, the . q^1 i t e r a t e i s given by (3) C = ~2_ *« + ^ X 1 3 ' ^ < x\"' y \" < : 7 > ) ^ ^ w i t h W « W By f o r m a l l y . s u b s t i t u t i n g the tr u e s o l u t i o n i n t o (2) and (3)> we o b t a i n the - t r u n c a t i o n e r r o r of each formula, and , r e s p e c t i v e l y : However, any numerical e v a l u a t i o n of (2) and (3) r e s u l t s i n numbers &ntv and ^*+h > % = 1^2, . . ., m, such t h a t 8. and (my) w i t h zrn+fe. = n^*-!, . The q u a n t i t i e s n^w, Av • i n equations (6) and (7) have been determined at the k+1 previous steps, whereas X>+k. and Xm denote the e r r o r due to•round-off i n e v a l u a t i n g the p r e d i c t o r and qth . i t e r a t e of the c o r r e c t o r at t = t + (n+k)h. We define and :and subtract (1+) from (6) and (.5) from (7) to o b t a i n the e r r o r equations xte-i K CO In equations (8) and (9)> ( 1 0 ) = 4(t«.x„j - J7*..,, _„i), 1,-so t h a t iL.^ i s some value of > . at l e a s t when the d e r i v a t i v e e x i s t s , at ~t»v+j , and on an open., i n t e r v a l between and v^»ti . I f i n equation (9)> we set q = m and r e c u r s i v e l y s u b s t i t u t e to (1'0 corresponding, r e l a t i o n s f o r , q = m, m-l> . . ., 1, then f i n a l l y , we o b t a i n •*— o+i^ i n terms of C;n<.| , i = 1,2, n+k-1. I f f s a t i s f i e s a L i p s c h i t z c o n d i t i o n and h i s s u f f i c i e n t l y s m a l l , as m in c r e a s e s , the r e s u l t a n t equation w i l l g i ve an expression f o r the e r r o r had the c o r r e c t o r been i t e r a t e d i n d e f i n i t e l y . We-now.make the assumption t h a t a constant f o r q=0,l,2,..., m and n .= 1,2, Then equation.(9) becomes The above-.is a .l i n e a r , d i f f e r e n c e equation w i t h constant c o e f f i c i e n t s . I f we take q = m i n the corresponding homogeneous equation, a f t e r r e c u r s i v e s u b s t i t u t i o n s , we o b t a i n (11) „(»*»>) + e *»0 (2 + 4 L ? In equation ( l l ) , © ~ h ^ t , , and i n p r a c t i c e , j@| < I I f i n equation ( l l ) i s re p l a c e d by ^ , the general s o l u t i s determined by.the r o o t s o f the polynomial, i o n (12) .ft) We now assume t h a t the q u a n t i t i e s 'n«-h , 'n+b , > and VJ^ k are constants i n n and q , denoted by T c , Tp , r c , and rp r e s p e c t i v e l y . Then a p a r t i c u l a r s o l u t i o n of the.;.nonhomogeneous equation-corresponding t o ( l l ) when q.= m,.is (13) .1-© ^ fegj - ^ Pi) ^ ( £ < - h £ ^ 10. At t h i s p o i n t we a s s e r t that our i n t e r e s t i s con f i n e d t o c o r r e c t o r formulas which are c o n s i s t e n t i n the sense t h a t the tr u e s o l u t i o n w i l l s a t i s f y (3) t o w i t h i n terms t h a t are o(h). For our purposes, the p r e d i c t o r formula should at l e a s t supply an exact s o l u t i o n to ( l ) when f = 0. I f CK^~ ^t,^ ' ' t h e n l t i s necessary (as w e l l as s u f f i c i e n t ) f o r the c o e f f i c i e n t s of any p r e d i c t o r - c o r r e c t o r p a i r to s a t i s f y : i n order to be c o n s i s t e n t i n the sense described ( H u l l , Luxemburg [2] ). I f we apply the f i r s t of r e l a t i o n s (lh) to (13), we o b t a i n , a f t e r regrouping terms: (15; U>- — j /._ ; Z _ r — ' I f the k+1 r o o t s of (12) are denoted by which f o r s i m p l i c i t y , are assumed d i s t i n c t f o r a l l m, the complete s o l u t i o n to (10) when q = m, w i l l have the form (16) e T - A ~ . , X 1 , + A-xA-t + ••• + A ^ , / w , \" ° ^ • The use of a k^*1 order m u l t i - s t e p formula i n v o l v e s r e p l a c i n g the f i r s t order d i f f e r e n t i a l equation ( l ) , by a k**1 order d i f f e r e n c e equation. This process introduces k-1 extraneous s o l u t i o n s . However, the coupled p r e d i c t o r -c o r r e c t o r p a i r introduces an a d d i t i o n a l \" i t e r a t i v e r o o t \" whose c o n t r i b u t i o n to the e r r o r w i l l vanish as the e f f e c t of the p r e d i c t o r i s suppressed. That i s , as m i s increased. To show t h i s to be the case, we note t h a t e x a c t l y one root of (12), \\ I- A - o. say f o r def i n i t e n e s s , has the property t h a t 11. I t remains t o show t h a t i n (l6), the m u l t i p l i e r of t h i s root remains bounded wi t h i n c r e a s i n g m. Let £ Q , £, , . .., G-b, be the e r r o r s i n . t h e k+1 o r d i n a t e s used t o s t a r t the i n t e g r a t i o n procedure. Also,.define For any m, the m u l t i p l i e r s , i n (l6) are determined by a system of equations of the form: r — A _. _\\ -4- • ' • • - + • _V E =A~,X~, +A~_X,*a + • • • + A - . f c t l X ^ s sA^X +A«Xx + - • + A~^ A - - * . w i t h determinant I I • - . 1 x«,( X ^ i • * • « • * I f — i s expanded by minors of the ( k + l ) s t column, then •I. A =• I» -ji? l inn. \"^!*>,n +| I I ivyi ^ We consider A-By use of equation (15), we f i n d t h a t i f h i s s u f f i c i e n t l y small A Thus remains bounded as m . increases. 12. We conclude t h a t the component of the e r r o r due t o one root of (12) w i l l . d i m i n i s h w i t h r e p e a t e d . i t e r a t i o n of the c o r r e c t o r formula. F u r t h e r , the remaining roots, tend.towards those of the c h a r a c t e r i s t i c polynomial of the c o r r e c t o r a c t i n g alone. I f we l i m i t our choice to s t a b l e methods which are c o n s i s t e n t i n the sense t h a t r e l a t i o n s (lh) h o l d , then f o r any m, equation (12) w i l l have one 0<3h root which i s approximately . The s t a b i l i t y p r o p e r t y i n s u r e s t h a t the components of the e r r o r due to the remaining r o o t s w i l l be n e g l i g i b l e . For t h i s choice of i n t e g r a t i o n formulas, equation (l6) gives ( i 7 ) = £ e s r t\" ' +-' \\m approximately, as the propagated e r r o r of the PE(CE) procedure. In d e r i v i n g (l7)> the e r r o r s i n the s t a r t i n g values are assumed equal and denoted by € .-In summary, equations (12) and (17) i n d i c a t e t h a t repeated i t e r a t i o n has the e f f e c t of suppressing ;the i n f l u e n c e of the. p r e d i c t o r on the s t a b i l i t y and accuracy c h a r a c t e r i s t i c s of the method. As m.increases, equations (12) and (17) become r e s p e c t i v e l y , the c h a r a c t e r i s t i c polynomial and propagated e r r o r estimate d e r i v e d by H u l l and Newbery Q[) f o r s t a b l e and c o n s i s t e n t c o r r e c t o r formulas. E r r o r estimate, P ( E C ) m procedure. We begin by w r i t i n g expressions f o r the approximations to the tr u e s o l u t i o n obtained from the p r e d i c t o r and each subsequent i t e r a t i o n of the c o r r e c t o r up t o the m^. Since the PEC procedure does not i n c l u d e an e v a l u a t i o n (a) of f u s i n g the f i n a l approximation at each step, the equations f o r V^ *^ } »^o,i,...,/m w i l l be 13-H-l Procgeding as i n the PE ( C E ) m case, we introduce the tr u e s o l u t i o n , i n t o each equation. From the r e s u l t i n g r e l a t i o n s , we s u b t r a c t the corresponding equations w i t h J n + j r e p l a c e d by the computed q u a n t i t i e s , j q = 0,1,2,..., We assume t h a t , given by equation (10), as w e l l as the e r r o r s due to t r u n c a t i o n and round-off are constant i n n and q. Then u s i n g the same n o t a t i o n as b e f o r e , we are l e d t o the f o l l o w i n g system of equations f o r the e r r o r s at m. each i t e r a t i o n : K-i K (18) 1=-o K-i I f we r e c u r s i v e l y s u b s t i t u t e i n t o the expression f o r W , each preceding equation of ( l 8 ) , then the system reduces t o the p a i r of equations: k-i C—-0 —r— e (19a) + h3> P'^ ' + _ k-i -e Ik. (19b) e „ t R where 0 = ^ ( 3 ^ and \\o\\ < I I f we e l i m i n a t e from (19),.the r e s u l t a n t equation f o r ^-rw^ admits the p a r t i c u l a r s o l u t i o n ^ i - e ~ 0 ^ Z _ f t | + CTP+ er- h 9 7 ft assuming the method i s c o n s i s t e n t i n the sense of (lk). Under the l a t t e r assumption, one of the ro o t s of the c h a r a c t e r i s t i c polynomial of system (19) w i l l , d i f f e r from by terms t h a t are o(h) as n->0 I f we assume th a t the errors, i n the s t a r t i n g values are constant and denoted by £ , then the propagated e r r o r of a s t a b l e , c o n s i s t e n t method based on a P ( E C ) m procedure w i l l be approximately, w i t h W^m^ given by equation (20).. We may v e r i f y from the preceding d i s c u s s i o n , t h a t i f jfc^l^. I , the s t a b i l i t y and,accuracy c h a r a c t e r i s t i c s of the P ( E C ) m and P E ( C E ) m procedures become i n d i s t i n g u i s h a b l e w i t h repeated i t e r a t i o n . For, as m i n c r e a s e s , the P ( E C ) m e r r o r , given by.the homogeneous equations corresponding to (l9)> w i l l s a t i s f y the same r e l a t i o n as the P E ( C E ) m e r r o r given by the l i m i t case of ( l l ) This i n d i c a t e s t h a t the c h a r a c t e r i s t i c polynomials of the procedures w i l l c o i n c i d e a f t e r many i t e r a t i o n s . F u r t h e r , both -U^m) andU>( m) become - ^ w i t h i n c r e a s i n g m. CHAPTER I I I COMPARISON AND CHOICE OF METHODS Our purpose i s t o s e l e c t from among the c l a s s of m u l t i - s t e p methods, the most e f f i c i e n t f o r the i n t e g r a t i o n of ( l ) i n circumstances i n v o l v i n g a complicated f u n c t i o n f . In such circumstances, the time taken at any i n t e g r a t i o n step to evaluate the p r e d i c t o r - c o r r e c t o r formulas may be regarded as n e g l i g i b l e when compared t o even a s i n g l e c a l c u l a t i o n of f. Thus the s e l e c t e d method and technique of a p p l i c a t i o n w i l l be c h a r a c t e r i z e d by a need f o r r e l a t i v e l y few eva l u a t i o n s of f , when determining a s o l u t i o n over a s p e c i f i e d i n t e r v a l of t , to w i t h i n p r e s c r i b e d e r r o r bounds. As a guide i n making our choice, the e r r o r estimates d e r i v e d i n the previous chapter w i l l be used to provide c r i t e r i a f o r comparing v a r i o u s i n t e g r a t i o n methods. These estimates w i l l be employed only i n t h i s q u a l i t a t i v e sense, r a t h e r than as exact measures of p r e c i s i o n . However, i n the next chapter, some numerical evidence i s o f f e r e d i n support of the accuracy of these estimates. When usin g the estimates, we assume the word l e n g t h adequate t o make the e f f e c t s of round-off s m a l l , compared t o the t r u n c a t i o n e r r o r of any formula considered. In the d i s c u s s i o n t h a t f o l l o w s , a comparison i s made of v a r i o u s p r e d i c t o r - c o r r e c t o r methods from the standpoint of e f f i c i e n c y . To advance t h i s d i s c u s s i o n , we assume t h a t a l l methods considered are both c o n s i s t e n t and s t a b l e . I t w i l l be seen i n the next chapter that s t a b i l i t y w i l l l i m i t the a p p l i c a b i l i t y of methods w i t h high e f f i c i e n c y c h a r a c t e r i s t i c s . However, these methods prove u s e f u l f o r problems solved by 1+th or 5th order formulas and, i n a r e s t r i c t e d sense, when solved by 6th order formulas. We f i r s t compare the e f f i c i e n c y of the P ( E C ) m technique, t a k i n g m = N>1 i t e r a t i o n s per step w i t h the same process having m=l. Provided the methods are s t a b l e , the behaviour of the propagated e r r o r i s de s c r i b e d by equation (21) with m=N or m=l. I f we assume the e r r o r s i n the s t a r t i n g values are small and neglect 1 \\N round-off, the propagated e r r o r f o r P(ECj w i l l be determined by 16. l - e ff 7mo i-o S i m i l a r l y , i f the s i n g l e - i t e r a t e technique remains s t a b l e w i t h step s i z e h' , the propagated e r r o r is. determined by (22) U ) -In equation.(22), the terms \"Te and \"Tp , are expressions of the form X (\"t) C^) ' '^le c o e^i c-'- e n\" t : > °^ e i t h e r case, i s determined by the choice of i n t e g r a t i o n formula. For most methods of i n t e r e s t , ^ /?}, and are approximately C*0 one*, and .usually |@1 « I • Comparing these expressions f o r U> and oJ\" , we see t h a t the accuracy of a s t a b l e PEC technique w i t h h 7 = h would.be comparable t o t h a t of P ( E C ) N . F u r t h e r , such a method would be much l e s s c o s t l y in terms of e v a l u a t i o n s of f. If we set h' = h/N ,,the f a c t o r ( i ) ^ ^ appearing i n (22), would, i n st cases, make a s t a b l e PEC process s u b s t a n t i a l l y more accurate than P(EC)^ w i t h mo step s i z e h . For t h i s choice of h' ; the cost of both techniques would be the same. Between these extremes of accuracy and c o s t , the PEC procedure can be made more e f f i c i e n t than P(EC)' N i f we can choose a step s i z e h'> h/N , f o r which the former i s s t a b l e . b I B -* P r e d i c t o r s and c o r r e c t o r s of the Adams type have ^ ^ — ^* ~ I » 2-o •17. We conclude from.the above•discussion,. t h a t , provided.the method remains, s t a b l e , there i s a range of h f o r which the propagated e r r o r of the s i n g l e - i t e r a t e procedure remains w i t h i n prescribed'bounds, and f o r this choice of h,. the method i s more e f f i c i e n t than a s i m i l a r procedure w i t h m >.l . A s i m i l a r comparison of PE ( C E ) m procedures would l e a d to the same c o n c l u s i o n : provided the method remains s t a b l e , the s i n g l e - i t e r a t e procedure w i l l determine a s o l u t i o n w i t h i n p r e s c r i b e d e r r o r bounds, i n the most e f f i c i e n t manner. These c o n c l u s i o n s l e a d us to compare the PECE and PEC i n t e g r a t i o n methods w i t h step s i z e h and h/2 r e s p e c t i v e l y . N e g l e c t i n g round-off, the •propagated e r r o r i s determined by Tc + eTP when.the c a l c u l a t i o n s are extended by the PECE procedure, and.by the case N^2 in.equa t i o n (22), f o r PEC. I f h i s s u f f i c i e n t l y s m a ll to keep the e r r o r of the former method w i t h i n p r e s c r i b e d bounds, the same would be tru e of the PEC procedure unless simultaneously, » __. i _ _ ( |eTP| » y X \\ I f the l a t t e r were t r u e , the f i r s t i n e q u a l i t y would be u n l i k e l y t o h o l d . f o r most methods of i n t e r e s t , since JG| = |^^^»»| • The second i n e q u a l i t y i s by i t s e l f i n v a l i d f o r Adams' formulas of a l l orders. Methods of t h i s type,.of k t h order, have A — —— \\ I/LL^^CZU*-^) ••• (.1** +*fc-0 d • Consequently, I ~ k_ ^ P, - 1 and ( 6 h f l > k+l Pk ' The f o r e g o i n g leads us t o conclude t h a t p rovided the PEC method remains s t a b l e , the f a c t o r (•g-)^ '*' i n u/'' > would permit a s u b s t a n t i a l . increase i n h w i t h consequent improvement i n e f f i c i e n c y . Therefore, we are l e d to the PEC procedure 18. as the most' e f f i c i e n t means of a p p l y i n g a p a i r of p r e d i c t o r - c o r r e c t o r - formulas for. the s o l u t i o n , of ( l ) . That i s , i f the method.remains s t a b l e , t h i s technique w i l l p rovide a s o l u t i o n of d e s i r e d accuracy and do so w i t h a minimum number of e v a l u a t i o n s of f. Having e s t a b l i s h e d .the most e f f i c i e n t means of extending c a l c u l a t i o n s , we must now adopt a s u i t a b l e p a i r of m u l t i - s t e p formulas. The c r i t e r i a governing t h i s choice are accuracy and s t a b i l i t y . To meet the accuracy requirement, we may r e s o r t to h i g h order formulas, since the time taken to evaluate the p r e d i c t o r and c o r r e c t o r i s n e g l i g i b l e when f i s complicated. Thus the s t a b i l i t y requirement remains as the formost f a c t o r governing our choice. Formulas w i t h good s t a b i l i t y p r o p e r t i e s are those of the Adams type; having extraneous s o l u t i o n s only at the o r i g i n .of the complex plane when h=o. Besides t h i s d e s i r a b l e property, these methods are simple t o program f o r computation, r e q u i r i n g at each step only the o r d i n a t e determined.at the immediately previous step i n a d d i t i o n t o the d e r i v a t i v e s In summary, the estimates of propagated error, have l e d us to the PEC procedure as the most e f f i c i e n t f o r the- i n t e g r a t i o n of ( l ) . In the next chapter, a comparison i s made between PEC and PECE procedures'using Adams formulas. CHAPTER IV EXPERIMENTAL RESULTS The t h e o r e t i c a l d i s c u s s i o n of the f o r e g o i n g s e c t i o n s has l e d us to s e l e c t the PEC procedure \"based on Adams formulas, as a most e f f i c i e n t and convenient means of s o l v i n g ( l ) , to w i t h i n a p r e s c r i b e d t o l e r a n c e of e r r o r . In t h i s chapter, we o f f e r experimental evidence i n support of our c h o i c e . S p e c i f i c a l l y , our aims are to e s t a b l i s h the accuracy of the e r r o r estimate d e r i v e d i n Chapter I I , and to demonstrate the h i g h e f f i c i e n c y of the s e l e c t e d methods. The d i f f e r e n t i a l equation used i n a l l experiments was (23) w i t h o * t 6 S r a d i a n s . The exponential component of the s o l u t i o n may be e l i m i n a t e d by s u i t a b l e choice of i n i t i a l c o n d i t i o n s , so t h a t OtSiryb't. + b co% b't X(-t) -z. Various values were assigned to a and b^ and the problem was solved u s i n g Adams formulas w i t h k=3,^ ,5> o r 6. In each case, we used a p r e d i c t o r and c o r r e c t o r of the same order. The t r u n c a t i o n e r r o r s of these formulas may be expressed as w i t h fcl\" and given i n Table 1 f o r formulas of both types. Table 1 Truncation E r r o r C o e f f i c i e n t s , Adams Methods. k 3 1+ 5 6 H k 251 95 19O87 5257 P r e d i c t o r 720 258 60I+80 17280 H k -19 -3 -863 -275 r\\. Corrector 720 16*0 60480 24192 .20. In Table 2, the row labeled'\"Max. Abs. E r r o r \" shows the magnitudes of the maximum e r r o r observed u s i n g the PEC procedure. The a p r i o r i bounds have been c a l c u l a t e d from equation (2l).with m=l, g=a, and \\C\\t [T\"p| t jv^j < Icf ' t In a l l cases we have taken a=b and h=0.04. . Table 2 PEC, m=l; a=b, h=0.0U, 0 i S t - 8 . k . 3 h 5 6 a 1 -1 l -1 1 -1 1 -0.625 Max. Abs. E r r o r x,10' O.53 3-8 0.018 •0.26 0.0006 0.25 0.000U A p r i o r i . Bound x 10? . 658 0-73 81 3-025 •5 0.0017 .5 0.003 R a t i o 0.08 0.7 O.05 0.7 0.05 O.k 0.05 0.1 In experiments w i t h ;5th and 6th order formulas, the e r r o r was predominantly due to round-off. Consequently, i n c a l c u l a t i n g the e r r o r bounds f o r these methods, we have allowed f o r the randomness of the e r r o r by co n s i d e r i n g the mean e r r o r to be zero, and re c o r d i n g the standard d e v i a t i o n , which i s .(computed bound) V 200 very•approximately, for.-200 i n t e g r a t i o n steps. With f i x e d a and b i n Table 2,,the r a t i o of e r r o r t o bound remains approximately constant i n :k, at l e a s t r e l a t i v e to the wide v a r i a t i o n i n the bounds themselves. This, i n d i c a t e s t h a t when comparing methods, the estimate w i l l at l e a s t e x h i b i t a q u a l i t a t i v e behaviour s i m i l a r t o t h a t of the e r r o r . A study, of the e r r o r estimate f o r the PECE procedure gives r e s u l t s very s i m i l a r to the PEC case, and.are shown.in Table 3-21. Table 3 k 3 k 5 6 a 1 . -1 1 - l 1 -1 .1 -0.625 Max. Abs. E r r o r x 1CK 16 O.O78 3-2 0:001*9 0.021 O.COO30 0.0090 0.000003 A p r i o r i Bound x 105 150 .11 60 0.007 33 D .000^ 0.36 0.000:015 R a t i o 0.1 0 . 7 0.05 0 .7 0.06 0 . 7 0 .03 0 . 2 Comparison of s i n g l e - i t e r a t e procedures. We f i r s t d efine the cost of i n t e g r a t i o n per u n i t l e n g t h of i n t e g r a t i o n i n t e r v a l as ^ s , where h i s the step s i z e and ^ \" the number of eva l u a t i o n s of f at each step. For the PEC and PECE procedures, V i s one and two r e s p e c t i v e l y . To compare the r e l a t i v e e f f i c i e n c y of these procedures, we have•integrated.equation ( 2 l ) u s i n g both, on an equal cost b a s i s . That i s , we i n t e g r a t e d the - d i f f e r e n t i a l equation u s i n g Adams formulas of v a r i o u s orders w i t h step h . f o r the PECE procedure and h/2 f o r PEC Tables k and 5 show the r a t i o of the magnitude of the maximum observed e r r o r i n a PECE i n t e g r a t i o n to t h a t of PEC, f o r the same case. In a l l cases a=b, and the Adams p r e d i c t o r and c o r r e c t o r formulas are of the same order. Table k Max. PECE e r r o r a k +1 -1 -2 -3 3 35 15 Ik k 67 27 27 •* 5 8 52 •* .6 k * * * PEC procedure unstable. 22. Max. PECE e r r o r J Max. PEC e r r o r c=20, h=0.10> \\ \\ a k +i - l -2 -3 3 35 12 Ik ih •h 170 2.5 * •* 5 50 1* * -* 1 6 7 •* •* * PEC procedure unstable. In each.stable case,.the r a t i o of the maximum PECE e r r o r to the maximum-PEC error, i s c o n s i d e r a b l e g r e a t e r than. one. These . r a t i o s give an i n d i c a t i o n of the improvement i n e f f i c i e n c y t o be obtained by i n c r e a s i n g the step s i z e used f o r the PEC i n t e g r a t i o n . Conclusion Both t h e o r e t i c a l and experimental r e s u l t s have i n d i c a t e d t h a t when s t a b l e , the PEC procedure based on Adams formulas i s a most e f f i c i e n t means of determining a s o l u t i o n t o w i t h i n e r r o r t o l e r a n c e s . The experimental r e s u l t s have a l s o exphasized the•importance of s t a b i l i t y when choosing a most e f f i c i e n t method. Our o r i g i n a l assumption had been t h a t i f a m u l t i - s t e p formula based.on the PECE procedure was s t a b l e w i t h step s i z e h, then.the corresponding PEC procedure would a l s o be s t a b l e w i t h h/2 I f t h i s was t r u e , the c r i t e r i o n of choice would be-the r e l a t i v e s i z e of the propagated e r r o r . However, our numerical experiments have shown otherwise and consequently, a new p o s s i b i l i t y has a r i s e n . To achieve a p r e s c r i b e d accuracy, the s u p e r i o r s t a b i l i t y p r o p e r t i e s of PECE co u l d permit the use of'higher order formulas at lower i n t e g r a t i o n c ost than a s t a b l e PEC procedure. The i n v e s t i g a t of t h i s p o s s i b i l i t y remains as a problem f o r the f u t u r e . BIBLIOGRAPHY 23-Hamming, R.W. Stable P r e d i c t o r - C o r r e c t o r Methods f o r Ordinary D i f f e r e n t i a l Equations, J.Assoc.Comput.Mach..vol.6, 1959, PP-37-^7-Hull,-T.E. and Luxemburg, W.A.J. Numerical Methods and Existence Theorems, f o r Ordinary D i f f e r e n t i a l Equations, Numerische Math, v o l . 2 , I960, pp.30-41. H u l l , T.E. and Newbery, A.C.R. I n t e g r a t i o n Procedures which Minimize Propagated E r r o r s , J.Soc.Indust. Appl.Math. vol-9> 19^1, pp.31-47. M i l n e , W.E. and Reynolds, R.R. F i f t h - O r d e r Methods f o r the Numerical S o l u t i o n of Ordinary D i f f e r e n t i a l Equations, J.Assoc.Comput.Mach. v o l . 9 , 1962, pp.64-70. Nordsieck, Arnold. On Numerical I n t e g r a t i o n of Ordinary D i f f e r e n t i a l Equations, Math, of Comp..vol.16, 1962,.pp.22-49-R a l s t o n , A. Some T h e o r e t i c a l and Computational Matters R e l a t i n g to P r e d i c t o r - C o r r e c t o r Methods of Numerical I n t e g r a t i o n , Computer J . vol.4, 1961, pp.64-67. "@en ; edm:hasType "Thesis/Dissertation"@en ; edm:isShownAt "10.14288/1.0080563"@en ; dcterms:language "eng"@en ; ns0:degreeDiscipline "Mathematics"@en ; edm:provider "Vancouver : University of British Columbia Library"@en ; dcterms:publisher "University of British Columbia"@en ; dcterms:rights "For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use."@en ; ns0:scholarLevel "Graduate"@en ; dcterms:title "Efficient methods for the numerical integration of ordinary differential equations"@en ; dcterms:type "Text"@en ; ns0:identifierURI "http://hdl.handle.net/2429/39377"@en .