Transcription of Runge-Kutta-Fehlberg Method (RKF45)
1 S EC . RUNGE -K UTTA M ETHODS 497. Runge-Kutta-Fehlberg Method ( rkf45 ). One way to guarantee accuracy in the solution of an is to solve the problem twice using step sizes h and h/2 and compare answers at the mesh points corresponding to the larger step size. But this requires a significant amount of computation for the smaller step size and must be repeated if it is determined that the agreement is not good enough. The Runge-Kutta-Fehlberg Method (denoted rkf45 ) is one way to try to resolve this problem. It has a procedure to determine if the proper step size h is being used. At each step, two different approximations for the solution are made and compared. If the two answers are in close agreement, the approximation is accepted. If the two answers do not agree to a specified accuracy, the step size is reduced. If the answers agree to more significant digits than required, the step size is increased.
2 Each step requires the use of the following six values: k1 = h f (tk , yk ), 1 1. k2 = h f tk + h, yk + k1 , 4 4. 3 3 9. k3 = h f tk + h, yk + k1 + k2 , 8 32 32. (28) 12 1932 7200 7296. k4 = h f tk + h, yk + k1 k2 + k3 , 13 2197 2197 2197. 439 3680 845. k5 = h f tk + h, yk + k1 8k2 + k3 k4 , 216 513 4104. 1 8 3544 1859 11. k6 = h f tk + h, yk k1 + 2k2 k3 + k4 k5 . 2 27 2565 4104 40. Then an approximation to the solution of the is made using a Runge-Kutta Method of order 4: 25 1408 2197 1. (29) yk+1 = yk + k1 + k3 + k4 k5 , 216 2565 4101 5. where the four function values f 1 , f 3 , f 4 , and f 5 are used. Notice that f 2 is not used in formula (29). A better value for the solution is determined using a Runge-Kutta Method of order 5: 16 6656 28,561 9 2. (30) z k+1 = yk + k1 + k3 + k4 k5 + k6 . 135 12,825 56,430 50 55. The optimal step size sh can be determined by multiplying the scalar s times the current step size h.
3 The scalar s is 1/4 1/4. tol h tol h (31) s= 2|z k+1 yk+1 | |z k+1 yk+1 |. 498 C HAP. 9 S OLUTION OF D IFFERENTIAL E QUATIONS. Table rkf45 Solution to y = 1 + y 2 , y(0) = 0. RK45 approximation True solution, Error k tk yk y(tk ) = tan(tk ) y(tk ) yk 0 1 2 3 4 5 6 7 8 9 10 where tol is the specified error control tolerance. The derivation of formula (31) can be found in advanced books on numerical anal- ysis. It is important to learn that a fixed step size is not the best strategy even though it would give a nicer-appearing table of values. If values are needed that are not in the table, polynomial interpolation should be used. Example Compare rkf45 and RK4 solutions to the y = 1 + y2 with y(0) = 0 on [0, ]. An rkf45 program was used with the value tol = 2 10 5 for the error control tolerance. It changed the step size automatically and generated the 10 approximations to the solution in Table An RK4 program was used with the a priori step size of h = , which required the computer to generate 14 approximations at the equally spaced points in Table The approximations at the right endpoint are y( ) y10 = and y( ) y14 = and the errors are E 10 = and E 14 = for the rkf45 and RK4 methods, respectively.
4 The rkf45 Method has the smaller error.. S EC . RUNGE -K UTTA M ETHODS 499. Table RK4 Solution to y = 1 + y 2 , y(0) = 0. RK4 approximation True solution, Error k tk yk y(tk ) = tan(tk ) y(tk ) yk 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 Numerical Methods Using Matlab, 4th Edition, 2004. John H. Mathews and Kurtis K. Fink ISBN: 0-13-065248-2. Prentice-Hall Inc. Upper Saddle River, New Jersey, USA.
