NASA-TR-R-287-1968 Classical fifth- sixth- seventh- and eighth-order Runge-Kutta formulas with stepsize control《带有步长控制的经典五阶 六阶 七阶和八阶Runge-Kutta公式》.pdf

1、NASA TR R-287CLASSICAL FIFTH-, SIXTH-, SEVENTH-, ANDEIGHTH- ORDER RUNGE- KUTTA FORMULASWITH STEPSIZE CONTROLBy Erwin FehlbergGeorge C. Marshall Space Flight CenterHuntsville, Ala.NATIONAL AERONAUTICS AND SPACE ADMINISTRATIONFor sale by the Clearinghouse for Federal Scientific and Technical Informati

2、onSpringfield, Virginia 22151 - CFSTI price $3.00Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-TABLE OF CONTENTSINTRODUCTION .PART I. FIFTH-ORDER FORM

3、ULASSection I. The Equations of Condition for the Runge-KuttaCoefficients 4Section II. A Solution of the Equations of Condition for theRunge-Kutta C oefficients . 14Section III. The Leading Term of the Local Truncation Error . . 21Section IV. Example for a Fifth-Order Runge-Kutta Formula . . 24Secti

4、on V. Numerical Comparison with Other Fifth-OrderRunge-Kutta Formulas 27Page1PART II. SIXTH-ORDER FORMULASSection VI, The Equations of Condition for the Runge-KuttaCoefficients 31Section VII. A Solution of the Equations of Condition for theRunge-Kutta Coefficients . 37Section VIII. The Leading Term

5、of the Local Truncation Error . . 46Section IX. Example for a Sixth-Order Runge-Kutta Formula . . 48Section X. Numerical Comparison with Other Sixth-OrderRunge-Kutta Formulas 50PART III. SEVENTH-ORDER FORMULASSection XI. The Equations of Condition for the Runge-KuttaCoefficients . 52Section XII. A S

6、olution of the Equations of Condition for theRunge-Kutta C oeffi cients 5 6Section XIII. The Leading Term of the Local Truncation Error . . 63iiiProvided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-TABLE OF CONTENTS (Continued)SectionXIV. Example for a Se

7、venth-Order Runge-KuttaFormula Section XV. Numerical Comparison with Other Seventh-Order Runge-Kutta Formulas Pa ge6466PART IV. EIGHTH-ORDER FORMULASSection XVI. The Equations of Condition for the Runge-Kutta Coefficients Section XVII. A Solution of the Equations of Condition for theRunge-Kutta Coef

8、ficients Section XVHI. Example for an Eighth-Order Runge-KuttaFormula .Section XIX. Numerical Comparison with Other Eighth-OrderRunge-Kutta Formulas 6872768OREFERENCES . 81i%:Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-LIST OF TABLESTABLE I.TABLE

9、 II.TABLE IN.TABLE IV.TABLE V.TABLE VI.TABLE VII.Equations of Condition for Sixth-Order FormulaRK 5(6) .Kutta 1 .Error Coefficientsfor Kutta 1 .Kutta 2 .Error Coefficients for Kutta 2 .Comparison of Fifth-Order Methods forExample (53)TABLE VIII. RK 6(7) .TABLE IX. Comparison of Sixth-Order Methods f

10、orExample (53) TABLE X. RK 7(8) .TABLE XI. Comparison of Seventh-Order Methods forExample (53) TABLE XII. RK 8(9) .TABLE XIII. Comparison of Eighth-Order Methods forExample (53) .Page526272828293O495165667780vProvided by IHSNot for ResaleNo reproduction or networking permitted without license from I

11、HS-,-,-Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-CLASSICAL FIFTH-,SIXTH-, SEVENTH-,ANDElGHTH-ORDER RUNGE-KUTTA FORMULASW ITHSTEPS IZECONTROLINTRODUCTIONIn two earlier papers 1 , 2, the author has described a Runge-Kuttaprocedure which provides

12、a stepsize control by one or two additionalevaluations of the differential equations. This earlier procedure, re-quiring an m-fold differentiation and a suitable transformation of thedifferential equations, yielded (m+4)-th order Runge-Kutta formulas -as well as (m+5)-th order formulas for the purpo

13、se of stepsize control.The stepsize control was based on a complete coverage of the leadinglocal truncation error term. The procedure required altogether sixevaluations of the differential equations, regardless of m.Since for m=0 no differentiations have to be performed, our earlier for-mulas repres

14、ent, in this special case, classical Runge-Kutta formulasof the fourth order, requiring six evaluations of the differential equationsand including a complete coverage of the leading truncation error termfor the purpose of stepsize control.o In this paper we shall derive classical Runge-Kutta formula

15、s of the fifth,sixth, seventh, and eighth order including a stepsize control procedurewhich is again based on a complete coverage of the leading local trunca-tion error term. Naturally, these new formulas require more evaluationsper step of the differential equations than the known classical Runge-K

16、uttaformulas without stepsize control procedure.Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-However, they require fewer evaluations per step than the knownclassicalformulas if Richardsons extrapolation to the limit is used for such for-mulas as a

17、 stepsize control device. Since the application of Ric_hardsonsextrapolation to the limit meanspractically a doublingof the computationaleffort for the benefit of the stepsize control only, it is worthwhile to lookfor a less expensive stepsize control procedure.Less expensiveprocedures have been sug

18、gestedby different authors.However, these procedures generally do not make any effort to cover thetruncation error, but rather try to estimate the truncation error from thelast or the last few considered terms of the Taylor series. Since such aprocedure has no mathematical base, these estim;xtesare

19、rather unreli-able. Generally, since the terms in a convergent Taylor series are de-creasing with increasing order, the last considered term will be largerthan the leading truncation error term. Therefore, a stepsize controlprocedure based onthe last or last few considered terms of the Taylorseries

20、will, in general, largely underestimate the permissible stepsize,thereby wasting computer time and building up unnecessarily largeround-off errors,. The new formulas of our paper contain one or m_re free parameters. Bya proper choice of these parameters the leading term of the local trunca-tion erro

21、r reduces substantially. This, in general, results in an increaseof the permissible stepsize. This increase, together with the smallernumber of evaluations per step, accounts for the superiority of our newformulas compared with the known Runge-Kutta formulas operated withRichardsons principle as ste

22、psize control procedure. Naturally, the new classical Runge-Kutta formulas of this paper, being ofthe_ eighth or lower order, are in general less economical than our earlierRunge-Kutta transformation formulas 1, 2 which represent Runge-Kutta formulas of any desired order. However, our new formulas h

23、avecertain advantages, since they require no preparatory work (like differen-tiation of the differential equations) by the programmer. The new formulasincluding the stepsize control procedure can easily be written as a sub-routine.Provided by IHSNot for ResaleNo reproduction or networking permitted

24、without license from IHS-,-,-Onemight try to further raise the order of our new classical Runge-Kuttaformulas hoping to make them still more economical. However, we be-lieve that we have m_re or less reached the optimum with our eighth-orderformula. The examplesthat we ran show that the gain of our

25、eighth-orderformula RKS(9)if comparedwith our seventh-order formula RK7(8) is notvery substantial any more.We madethe same experience with a ninth- and a tenth-order Runge-Kuttaformula that E.B. SHANKShas developedrecently. These (not yet published)new formulas of SHANKSdo not bring a substantial ga

26、in any more comparedwith SHANKSseighth-order Runge-Kutta formula.Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-SECTIONPART I. FIFTH-ORDER FORMULASI. THE EQUATIONS OF CONDITION FOR THERUNGE-KUITA COEFFICIENTS5. Although all formulas of this paper ho

27、ld for systems of differentialequations, we state them-for the sake of brevity-for a single differentialequationy = f(x,y) (i)only.In the customary way, we set:fo = f(xo, Yo)K-1f =f(Xo+_Kh, Yo +h _ fl_)f)K)_=0and we require:5y=yo+h_ c fKKg=07 2y=yo+h _ fKKg=O(K =1,2,3,.7)+ O(he) i+ O(h v)(2)(3)with

28、h as integration step size.Equations (3) mean that we try to determine the coefficients _ , flg g)_Ac , c in such a way, that the first formula (3) represents a fifth-order andE Kthe second formula (3) a sixth-order Runge-Kutta formula. According to (2)Provided by IHSNot for ResaleNo reproduction or

29、 networking permitted without license from IHS-,-,-and (3) this means that the coefficients a and flKk have to be the same in both Kformulas for the first six evaluations of f. This restriction explains why weshall need eight evaluations (instead of seven) for our sixth-order formula.6. Since the eq

30、uations of condition for the Runge-Kutta coefficients, result-ing from Taylor expansions of (2) and (3), are well known in the literature, werestrict ourselves to stating these equations. For Runge-Kutta formulas up tothe eighth order these equations - in condensed form - are listed in a paper byJ.C

31、. BUTCHER (3, Table 1). For the convenience of the reader, we listthese equations for a sixth-order formula such as the second formula (3) in thecustomary summation form. However, since these equations then become some-what lengthy, we introduce the following abbreviation:k kfiKl c_l +ilK2 c_2 + +fl

32、KK-laK-1 = PKk (K= 2,3,.,7;k= 1,2,3 (4)The 37 equations of condition for our sixth-order formula, listed in the sameorder as in BUTCHERs paper, then read:TABLE I. EQUATIONS OF CONDITION FOR SIXTH-ORDER FORMULA7(i,1) Z - 1 : 0Kg=07A 1(II, 1) 2J c a - 0g K 2g=l7 A iK=2(HI,2) 2 K K - 6K=I1K flKxPkl - 2

33、4 - 0g=31 x-,7 1(IV, 2)_.z - -0cKPK2 24_2 _-5Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-TABLE I. (Continued)(IV, 3) _Kc_KPK1 - _ - 0g=2(IV,4) _ K=I K K-24 -0K=4 K Lk=3 21 7 1CK flK_.B_2 - 120 - 0(V, 2) 2 g=3 =2(V, 3) flKXaXPXl - 40 - 0K=3 K k= 2

34、1 _7 ,_ P 1(V, 4) _K=2 cK K3- 120 - 0(v, 5) o _ /_ PKxPx -3-o-K=3 K g k=21 7 I(v,o)Z_ - _o= K K K2 301 _? 1 -0(V, 7) 2K=2_KP2K1 -40(V, 8) _K=2 g K KI -20-0(V, 9)(VI, 1)(VI, 2)1 -012072O6Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-TABLE I. (Contin

35、ued)(VI, 3) _ fiKk fl)_#olp%K=4 g LX=3 2(VI, 6)(VI,4) _ _3 _K flKkPk - 720 - 0K k=2(vI, 5) _ _Kxax _xPK=4 KLt=3 2eK _KxaxPx - 180 - 02 K=3 _X=2K-1(VI, 7)(vi,8)- 240 - 01(VI, 9)( VI, 10)(vl, li)K=3 Cg_K-t )_ 2p _1 2 _ _Kx% xl 12oK=3 K =21 _ cAKPK4 12-4 K=2 - 720 - 01CKC_K flK_ flk#P_K=4 Lk=3 2-0e_ _X

36、PX7 -144 =02K= 3 K K 2E A /K-I )VI, 12) c o_ l_ 1(VI,13) 6K-2 K K K3-144 =0(VI, 14) _ _1: K flKxlUM - 7-_ =0K=3 1 2K=2 PK1P/2 - _ = 0(VI, 15)1- - 02401m= 01801_ =0144Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-TABLE I. (Concluded)K-1 1/ 11 _ _,_

37、( _ _F_ - _ - o(VI, 16) 2K=3 K KX=2_ Ii _7_ _2F_2 _ - o(VI, 17) _K=2 K1 _7 _ c_ p2 1(VI, 18) _=2 K K K1 -48 -0(VI, 19) _K_2 K KPKI -72 -01 T i(vi,2o) _ _ _ -_ -o120 _ K 720K-1The Roman numerals in front of the equations in Table I indicate the order ofthe terms in the Taylor expansion.A similar tabl

38、e holds for a fifth-order formula such as the first formula(3). We obtain this table from Table I by omitting the sixth-order equations(VI, 1) through (VI, 20), and replacing, in the remaining equations,_ with cK Kand the upper limit 7 of the _-sums with 5.Naturally, all equations of this new table

39、for a fifth-order formula aswell as those of Table I have to be satisfied simultaneously.7. The equations of Table I which represent necessary and sufficient con-ditions for a sixth-order Runge-Kutta formula, can, however, be replaced bya much sfmpler system of sufficient condition.First, we make th

40、e assumptions: _6 cA7_5 =z7 = 1, ol6=0;_1=cl =0, cA2-c2, c a =c 3, c4-c 4, c5=0, = =c 6 (5)Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-which greatly facilitate the simultaneous solution of both systems of equations(for fifth-order and sixth-order

41、 Runge-Kutta formulas). Then, we introducethe following assumptions (A, (B), (C) which are well-known to reduce ourequations of condition to a great extent:1 1 2 1 (“:) P21 =_ ol_ (A2) (*) 1)51 =_-ors =_- (A5)i _ (A3) 1 2 (A6)(A) (*) P31=E P6_=Eo_ = o1 _ (A4)(*) P41=_ 1 2 1 (A7)PTi - 2 0t7 = -2-l_sa

42、 i = c4( 1 - o_4)t (B4)c5 (_5 +_75) = c5( t - ozs) = o (B5)c_? 6 = c5 (B6)(c)(,-)c2_zB21+e3_3f131+c4_4f141+csl_57_ I =0 (CI)c2_b21+ c3“b31 + c4“b41 + ecru = 0C30_3_32_21% C40_4(_42_21 + _43_31) + C5(_2fl21 + _73_31 + _74_41+_i +_i) =o ,(C2)(C3)Provided by IHSNot for ResaleNo reproduction or networki

43、ng permitted without license from IHS-,-,-The asterisk (*) in front of some of the equationsindicates that theseequationshold for the fifth-order formulas as well as for the sixth-orderformulas. If the equations marked by an asterisk are split in two lines, the topline holds for the fifth-order form

44、ulas and the bottom line for the sixth-orderformulas. The equations without an asterisk are required for the sixth-orderformulas only.8. Inserting the assumptions (A) into Table I, one immediately finds thefollowing identities:(III, 1) = (III, 2); (IV, 3) - 3(IV, 4); (V, 7) - 3(V, 9);(V, 8) -6(V,9);

45、 (VI, 8)= 2(VI, 7); (VI, 14)-(VI, 16);(VI, 15) - (VI, 17); (VI, lS) = 15(VI, 20); (VI, 19)- 10(VI, 20). I (6)Therefore, the equations on the left-hand sides of the identities (6) can beomitted from Table I.Using also assumption (B1), three more identies are obtained:(IV, i) -(IV,2); (V, 3) -3(V, 4);

46、 (VI, 7) -3(VI, 9); (7)thus eliminating equations (VI, 1), (V, 3), (VI, 7) from Table I,The assumptions (B) lead to the following identities:(IV, 2) = (III, 2) - 3(IV, 4);(V,2) - (IV,2) - (V, 6);(VI, 1) = (V,1) - (VI, IO);(VI, 4)- (V,4) - (VI, 13);(w, 6)_ C:,G)(v, 1) - (Iv, 1) - (v, 5);(v,4) =- (IV,

47、4) - 4(v, 9);(vI, 2) = (v, 2) - (vI, 11);(vi, 5)- (v, 5) - 2(vI, 16);- 2(vI, 17); (vI, 9)- (v, 9) - 5(vI, 20);(8)thus eliminating ten more equations from Table I.10Provided by IHSNot for ResaleNo reproduction or networking permitted without license from IHS-,-,-Finally, we make use also of the assum

48、ptions (C and obtain the followingidentities:(V,5) = (V,6); (VI, 3)-3(V,4) - 3(VI, 13);(VI, 10)- (VI, 11); (VI, 12) - 3(VI, 13);(VI, 16)- (VI, 17) .(9)9. Cancelling all equations listed on the left-hand sides of the identities (6),(7), (8), (9), Table I reduces to the following ten equations:(I, 1), (II, l) (III, 2), (IV,4), (V, 6), (V, 9), (VI, 11), (VI, 13),%)