Ma2002 1.9 rm

Numerical solution of ordinary and
partial dierential Equations
Module 9: Runge-Kutta methods - IV
Dr.rer.nat. Narni Nageswara Rao
£
August 2011
1 System of Equations
Consider the system of n equations
dy
dt
= f(t; y1; y2; ¡¡¡ ; yn)
y(t0) = (1)
where
y = [y1; y2; ¡¡¡ ; yn]T
; f = [f1; f2; ¡¡¡ ; fn]T
; = [1; 2; ¡¡¡ ; n]T
All the methods derived earlier can be used to solve the system of equations
(1), by writing the methods in vector form. Let us now apply some of these
methods to solve (1).
1.1 Taylor Series Method
We can write the Taylor series expansion in vector form as
yj+1 = yj + hyH
j +
h2
2!
yHH
j + ¡¡¡+
hp
p!
y(p)
j ; j = 0; 1; 2; ¡¡¡ ; N 1 (2)
where
y(k)
j =
2
6664
(y
(k)
1 )j
(y
(k)
2 )j
...
(y
(k)
n )j
3
7775 =
2
64
dk 1
dtk 1 [f1(t; y1; y2; ¡¡¡ ; yn)]j
...
dk 1
dtk 1 [fn(t; y1; y2; ¡¡¡ ; yn)]j
3
75
£nnrao maths@yahoo.co.in
1

1.2 Runge-Kutta method of Second order
Consider the Euler-Cauchy or the Heun method. We write it in vector form
as
yj+1 = yj +
1
2
(K1 + K2) (3)
where
Ki1 = hfi(tj; (y1)j; (y2)j; ¡¡¡ ; (yn)j)
Ki2 = hfi(tj + h; (y1)j + K11; (y2)j + K21; ¡¡¡ ; (yn)j + Kn1)
i = 1; 2; ¡¡¡ ; n
In an explicit form, (3) becomes
2
6664
y1
y2
...
yn
3
7775
j+1
=
2
6664
y1
y2
...
yn
3
7775
j
+
1
2
0
BBB@
2
6664
K11
K21
...
Kn1
3
7775 +
2
6664
K12
K22
...
Kn2
3
7775
1
CCCA; j = 0; 1; ¡¡¡ ; N 1
(4)
1.3 Runge-Kutta (Classical) Method of Fourth Order
The vector format of fourth order Ruge-Kutta method is
yj+1 = yj +
1
6
(K1 + 2K2 + 2K3 + K4) (5)
where
K1 =
2
6664
K11
K21
...
Kn1
3
7775; K2 =
2
6664
K12
K22
...
Kn2
3
7775; K3 =
2
6664
K13
K23
...
Kn3
3
7775; K4 =
2
6664
K14
K24
...
Kn4
3
7775
Ki1 = hfi(tj; (y1)j; (y2)j; ¡¡¡ ; (yn)j)
Ki2 = hfi

tj +
h
2
; (y1)j +
K11
2
; (y2)j +
K21
2
; ¡¡¡ ; (yn)j +
Kn1
2

Ki3 = hfi

tj +
h
2
; (y1)j +
K12
2
; (y2)j +
K22
2
; ¡¡¡ ; (yn)j +
Kn2
2

Ki4 = hfi (tj + h; (y1)j + K13; (y2)j + K23; ¡¡¡ ; (yn)j + Kn3)
i = 1(1)n
2

In explicit form, (5) becomes
2
6664
y1
y2
...
yn
3
7775
j+1
=
2
6664
y1
y2
...
yn
3
7775
j
+
1
6
0
BBB@
2
6664
K11
K21
...
Kn1
3
7775 + 2
2
6664
K12
K22
...
Kn2
3
7775 + 2
2
6664
K13
K23
...
Kn3
3
7775 +
2
6664
K14
K24
...
Kn4
3
7775
1
CCCA(6)
j = 1; 2; ¡¡¡ ; N 1
ExampleCompute an approximation to y(1), yH(1) and yHH(1) with the Taylor
series method of Second order and step length h = 1, for the initial value
problem
yHHH + 2yHH + yH y = cos t; 0 t 1
y(0) = 0; yH(0) = 1; yHH(0) = 2
after reducing it to a system of

rst order equations.
Solution
Set y = v1; vH
1 = v2; vH
2 = v3
Note that v2 = vH
1 = yH
v3 = vH
2 = yHH and
vH
3 = vHH
2 = yHHH
The system of equations is
vH
1 = v2; v1(0) = 0
vH
2 = v3; v2(0) = 1
vH
3 = cos t 2v3 v2 + v1; v3(0) = 2
or
vH =
2
4
v1
v2
v3
3
5
H
=
2
4
v2
v3
cos t 2v3 v2 + v1
3
5; v(0) =
2
4
0
1
2
3
5
The Taylor series method of second order is
v(t0 + h) = v0 + hvH
0 +
h2
2!
vHH
0
= v0 + vH
0 +
1
2
vHH
0
3

Since h = 1, we have
vH
0 =
2
4
v2(0)
v3(0)
1 2v3(0) v2(0) + v1(0)
3
5 =
2
4
1
2
4
3
5
Now,
vHH =
2
4
vH
2
vH
3
sin t 2vH
3 vH
2 + vH
1
3
5and vHH
0 =
2
4
2
4
7
3
5
Hence,
v(1) = v0 + vH
0 +
1
2
vHH
0 =
2
4
0
1
2
3
5 +
2
4
1
2
4
3
5 +
1
2
2
4
2
4
7
3
5 =
2
4
2
1
3=2
3
5
) y(1) = 2; yH(1) = 1; yHH(1) =
3
2
Example Solve the system of equations
xH = 3x + 2y; x(0) = 0
yH = 3x 4y; y(0) = 0:5
with h = 0:2 on the interval [0; 0:4]. Use the
(i) Euler-Cauchy method and
(ii) Classical Runge-Kutta fourth order method.
Solution (i) The Euler-Cauchy method is given by
xj+1 = xj +
1
2
(K1 + K2); j = 0; 1
K1 = hf(tj; xj)
K2 = hf(tj + h; xj + K1)
For j = 0, we have
t0 = 0; x0 = 0; y0 = 0:5
4

K11 = hf1(t0; x0; y0) = 0:2( 3x0 + 2y0)
= 0:2( 3 ¢0 + 2 ¢0:5) = 0:2
K21 = hf2(t0; x0; y0) = 0:2(3x0 4y0)
= 0:2[3 ¢0:4 ¢0:5] = 0:4
K12 = hf1(t0 + h; x0 + K11; y0 + K21)
= 0:2[ 3(x0 + K11) + 2(y0 + K21]
= 0:2[ 3(0 + 0:2) + 2(0:5 0:4)]
= 0:08
K22 = hf2(t0 + h; x0 + K11; y0 + K21)
= 0:2[3(x0 + K11) 4(y0 + K21)]
= 0:2[3(0 + 0:2) 4(0:5 0:4)]
= 0:04
x(0:2) % x1 = x0 +
1
2
(K11 + K12)
= 0 +
1
2
(0:2 0:08)
= 0:06
y(0:2) % y1 = y0 +
1
2
(K21 + K22)
= 0:5 +
1
2
( 0:4 + 0:04)
= 0:32
For j = 1, we have
t1 = 0:2; x1 = 0:06; y1 = 0:32
K11 = hf1(t1; x1; y1) = 0:2( 3x1 + 2y1)
= 0:2( 3 ¢0:06 + 2 ¢0:32) = 0:092
K21 = hf2(t1; x1; y1) = 0:2(3x1 4y1)
= 0:2(3 ¢0:06 4 ¢0:32) = 0:22
K12 = hf1(t1 + h; x1 + K11; y1 + K21)
= 0:2[ 3(x1 + K11) + 2(y1 + K21)]
= 0:2[ 3(0:06 + 0:092) + 2(0:32 0:22)]
= 0:0512
5

K22 = hf(t1 + h; x1 + K11; y1 + K21)
= 0:2[3(x1 + K11) 4(y1 + K21)]
= 0:2[3(0:06 + 0:092) 4(0:32 0:22)]
= 0:0112
x(0:4) % x2 = x1 +
1
2
(K11 + K12)
= 0:06 +
1
2
(0:092 0:0512)
= 0:0804
y(0:4) % y2 = y1 +
1
2
(K21 + K22)
= 0:32 +
1
2
( 0:22 + 0:0112)
= 0:2156
(ii) For the classical Runge-Kutta fourth order method, we obtain the fol-
lowing results
For j = 0, we have
t0 = 0; x0 = 0; y0 = 0:5
K11 = hf1(t0; x0; y0) = h( 3x0 + 2y0)
= 0:2( 3 ¢0 + 2 ¢0:5)
= 0:2
K21 = hf2(t0; x0; y0) = h(3x0 4y0)
= 0:2(3 ¢0 4 ¢0:5)
= 0:4
K12 = hf1(t0 +
h
2
; x0 +
K11
2
; y0 +
1
2
K21)
= h[ 3(x0 +
1
2
K11) + 2(y0 +
1
2
K21)]
= 0:2[ 3(0:1) + 2 ¢0:3]
= 0:06
K22 = hf2(t0 +
h
2
; x0 +
1
2
K11; y0 +
1
2
K21)
= h[3(x0 +
1
2
K11) 4(y0 +
1
2
K21)]
= 0:2[3 ¢0:1 4 ¢0:3]
= 0:18
6

K13 = hf1(t0 +
h
2
; x0 +
1
2
K12; y0 +
1
2
K21)
= h[ 3(x0 +
1
2
K12) + 2(y0 +
1
2
K22)]
= 0:2[ 3 ¢0:03 + 2 ¢0:41]
= 0:146
K23 = hf2(t0 +
h
2
; x0 +
1
2
K12; y0 +
1
2
K22)
= h[3(x0 +
1
2
K12) 4(y0 +
1
2
K22)]
= 0:2[3 ¢0:03 40:41]
= 0:31
K14 = hf1(t0 + h; x0 + K13; y0 + K23)
= h[ 3(x0 + K13) + 2(y0 + K23)]
= 0:2[ 3 ¢0:146 + 2 ¢0:19]
= 0:0116
K24 = hf2(t0 + h; x0 + K13; y0 + K23)
= h[3(x0 + K13) 4(y0 + K23)]
= 0:2[3 ¢0:146 + 2 ¢0:19]
= 0:0644
x(0:2) % x1 = x0 +
1
6
(K11 + 2K12 + 2K13 + K14)
= 0 +
1
6
(0:2 + 0:12 + 0:292 0:0116)
= 0:1001
y(0:2) % y1 = y0 +
1
6
(K21 + 2K22 + 2K23 + K24)
= 0:5 +
1
6
( 0:4 0:36 0:62 0:0644)
= 0:2593
For j = 1, we have
t1 = 0:2; x1 = 0:1001; y1 = 0:2593
K11 = hf1(t1; x1; y1) = h( 3x1 + 2y1)
= 0:2( 3 ¢0:1001 + 2 ¢0:2593)
= 0:0437
7

K21 = hf2(t1; x1; y1) = h(3x1 4y1)
= 0:2[3 ¢0:1001 4 ¢0:2593]
= 0:1474
K12 = hf1(t1 +
h
2
; x0 +
1
2
K11; y0 +
1
2
K21)
= h[ 3(x1 +
1
2
K11) + 2(y1 +
1
2
K21)]
= 0:2[ 3 ¢0:1220 + 2 ¢0:1856]
= 0:0010
K22 = hf2(t1 +
h
2
; x1 +
1
2
K11; y1 +
1
2
K21)
= h[3(x1 +
1
2
K11) 4(y1 +
1
2
K21)]
= 0:2[3 ¢0:1220 4 ¢0:1856] = 0:0753
K13 = hf1(t1 +
h
2
; x1 +
1
2
K12; y1 +
1
2
K22)
= h[ 3(x1 +
1
2
K12) + 2(y1 +
1
2
K22)]
= 0:2[ 3 ¢0:1006 + 2 ¢0:2217] = 0:0283
K23 = hf2(t1 +
h
2
; x1 +
1
2
K12; y1 +
1
2
K22)
= h[3(x1 +
1
2
K12) 4(y1 +
1
2
K22)]
= 0:2[3 ¢0:1006 4 ¢0:2214]
= 0:1170
K14 = hf1(t1 + h; x1 + K13; y1 + K23)
= h[ 3(x1 + K13) + 2(y1 + K23)]
= 0:2[ 3 ¢0:1284 + 2 ¢0:1423] = 0:0201
K24 = hf2(t1 + h; x1 + K13; y1 + K23)
= h[3(x1 + K13) 4(y1 + K23)]
= 0:2[3 ¢0:1284 4 ¢0:1423] = 0:0368
x(0:4) % x2 = x1 +
1
6
(K11 + 2K12 + 2K13 + K14)
= 0:1001 +
1
6
[0:0437 + 0:0020 + 0:0566 0:0201]
= 0:1138
8

Ma2002 1.9 rm

More Related Content

What's hot

Viewers also liked

Similar to Ma2002 1.9 rm

More from Ramakrishna Paduchuri

Recently uploaded

Ma2002 1.9 rm