Topic 14.3: 4th-Order Runge Kutta's Method (Examples)

Contents Previous Chapter Start of Chapter Previous Topic Introduction Notes Theory HOWTO Examples Engineering Error Questions Matlab Maple Next Topic Next Chapter

Example 1

Given the IVP y(1)(t) = 1 - t y(t) with y(0) = 1, approximate y(1) with one step.

First, let t0 = 0, y0 = 1, and h = 1. Thus, we calculate

K0 = f(0, 1) = 1
K1 = f(0.5, 1.5) = 0.25
K2 = f(0.5, 1.125) = 0.4375
K3 = f(1, 1.4375) = -0.4375
(K0 + 2K1 + 2K2 + K3)/6 = 0.3229166667

Therefore, approximation is y1 = y0 + h 0.3229166667 = 1.3229166667. The actual value is 1.331309118 and therefore the absolute error is approximately 0.0084, significantly smaller than the error using Heun's method.

Example 2

Given the same IVP shown in Example 1, approximate y(0.5).

K0 = f(0, 1) = 1
K1 = f(0.25, 1.25) = 0.6875
K2 = f(0.25, 1.171875) = 0.70703125
K3 = f(0.5, 1.353515625) = 0.3232421875
(K0 + 2K1 + 2K2 + K3)/6 = 0.6853841146

Therefore, approximation is y1 = y0 + h = 1 + 0.5⋅0.6853851146 = 1.342692057. The actual value is 1.342841185 and therefore the absolute error is approximately 0.00015, or approximately 1/56th the error when we used h = 1. This is significantly better than the estimated 1/32.

Example 3

Repeat Examples 1 and 2 but with with the initial value y(0.5) = 2.5 and approximating y(1.5) and y(1.0).

First, let t0 = 0.5, y0 = 2.5, and h = 1. Therefore

K0 = f(0.5, 2.5) = -0.25
K1 = f(1, 2.375) = -1.375
K2 = f(1, 1.8125) = -0.8125
K3 = f(0.5, 1.6875) = -0.153125
(K0 + 2K1 + 2K2 + K3)/6 = -1.026041667

Therefore, our approximation is y1 = y0 + h(-1.026041667) = 1.473958333. The actual value is 1.502483616 and therefore the absolute error is approximately 0.0285.

Next, let t0 = 0.5, y0 = 2.5, and h = 0.5. Therefore

K0 = f(0.5, 2.5) = -0.25
K1 = f(1, 2.4375) = -0.828125
K2 = f(1, 2.29296875) = -0.719726562
K3 = f(0.5, 2.140136719) = -1.140136719
(K0 + 2K1 + 2K2 + K3)/6 = -0.7476399739

Therefore, our approximation is y1 = y0 + h(-0.7476399739) = 2.126180013. The actual value is 2.126611964 and therefore the absolute error is approximately 0.00043. Again, the error is approximately 1/66th that of the previous calculation where h = 1.

If we tabulate the errors for various values of h, as is shown in Table 1, we note that as h gets smaller, the error begins to drops by a factor 1/32 each time we divide h by two.

Table 1. Errors when approximating y(t0 + h) for decreasing values of h.

hApproximation
of y(0.5 + h)
Error
11.4739583330.0285
0.52.1261800130.00043
0.252.3686006860.0000089
0.1252.4508056500.00000023
0.06262.4798088990.0000000064
0.031252.4910366620.00000000019
0.0156252.495804920.0000000000058

You will note that with the last step, the error goes down by approximately a factor of 1/32.7.

Copyright ©2005 by Douglas Wilhelm Harder. All rights reserved.