Numerical Methods
Ordinary Differential Equations:
Initial Value Problems (IVP)
By
Norhayati Rosli
Faculty of Industrial Sciences & Technology
norhayati@[Link]
Numerical Methods
by Norhayati Rosli
[Link]
Description
AIMS
This chapter is aimed to solve initial value problems of single ODE by using three different types of methods involving
Euler’s method, 2nd order Runge-Kutta method and 4th order Runge-Kutta method. In addition, for system of ODEs,
two types of methods are considered; Euler’s method and 4th order Runge-Kutta method. Steps by steps of solving
initial value problems for single ODE and system of ODEs are presented
EXPECTED OUTCOMES
1. Students should be able to solve initial value problems using Euler’s method, 2nd order Runge-Kutta method
and 4th order Runge-Kutta method.
2. Students should be able to solve system of ODEs using Euler’s method and 4th order Runge-Kutta method.
REFERENCES
1. Norhayati Rosli, Nadirah Mohd Nasir, Mohd Zuki Salleh, Rozieana Khairuddin, Nurfatihah Mohamad
Hanafi, Noraziah Adzhar. Numerical Methods, Second Edition, UMP, 2017 (Internal use)
2. Chapra, C. S. & Canale, R. P. Numerical Methods for Engineers, Sixth Edition, McGraw–Hill, 2010.
Numerical Methods
by Norhayati Rosli
[Link]
Content
1 Introduction to Ordinary Differential Equatios
2 Numerical Methods of ODEs (IVP)
2.1 Euler’s Method
2.2 Runge-Kutta (RK) Methods
2.2.1 Second Order Runge-Kutta (RK2)Method
2.2.2 Fourth Order Runge-Kutta (RK4) Method
3 System of ODEs
3.1 Euler’s Method
3.2 Fourth Order Runge-Kutta (RK4) Method
Numerical Methods
by Norhayati Rosli
[Link]
INTRODUCTION
ODEs refers as a rate of equation .
It expresses the rate of change of a variable as a function of variables and
parameters .
ODEs provide a tool for better understanding the behaviour of many
biological and physical systems around us.
It forms the basis of simulation experiments in a realm where the
experiments
are often impractical or unethical.
By solving the underlying ODEs, one can identify trends and make forecasts
about the future path of the process.
dy
General Form of ODEs f ( x, y )
dx
Numerical Methods
by Norhayati Rosli
[Link]
INTRODUCTION (Cont.)
Most of ODEs cannot be solved analytically.
Its due to the complexity form of the equations.
Numerical methods offer a viable option to solve ODEs.
Numerical Methods for Solving ODEs (IVP)
1 Euler’s Method
Three Types
of Methods
Fourth Order
Runge-Kutta Method
3 2 Second Order
Runge-Kutta Method
Numerical Methods
by Norhayati Rosli
[Link]
EULER’S METHOD
Euler’s method is a one step method and can be formulated in
general as
yi 1 yi h
where ℎ denotes a step size and Φ is a slope estimate.
A new value of 𝑦𝑖+1 is extrapolated from an old value of 𝑦𝑖 over a
distance, ℎ.
𝑑𝑦
For Euler’s method the first derivative |𝑥𝑖 = 𝑓 𝑥𝑖 , 𝑦𝑖 provide the
𝑑𝑥
slope estimate at 𝑥𝑖 such that
f xi , yi
Numerical Methods
by Norhayati Rosli
[Link]
EULER’s METHOD (Cont.)
Euler’s Method Formula
yi 1 yi f xi , yi h
xi 1 xi h
Figure 1: Graphical Illustration of Euler’s Method
Numerical Methods
by Norhayati Rosli
[Link]
EULER’s METHOD (Cont.)
Example 1
Use Euler’s method to solve the following ODE (IVP)
dy
exp( x) 2 y, y (0) 2
dx
for 0 ≤ 𝑥 ≤ 2 with a step size, ℎ = 0.5.
Solution
Identify the estimate slope, 𝑓 𝑥, 𝑦 = exp −𝑥 − 2𝑦,
Step 1
initial values, 𝑥0 = 0, 𝑦0 = 2 and ℎ = 0.5.
Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) over the interval 0 ≤ 𝑥 ≤ 2
Step 2
by using Euler’s method.
Numerical Methods
by Norhayati Rosli
[Link]
EULER’s METHOD (Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2 Second iteration: 𝑖 = 1, 𝑥1 = 0.5, 𝑦1 = 0.5
y1 y0 f x0 , y0 0.5 y2 y1 f x1 , y1 0.5
2 f (0, 2)(0.5) 0.5 f (0.5, 0.5)(0.5)
2 (exp(0) 2(2))(0.5) 0.3033
0.5 x2 x1 h
x1 x0 h 1.0
0 0.5 0.5 y2 y (1.0) 0.3033
y1 y (0.5) 0.5
Numerical Methods
by Norhayati Rosli
[Link]
EULER’s METHOD (Cont.)
Solution (Cont.)
Third iteration: 𝑖 = 2, 𝑥2 = 1.0, 𝑦2 = 0.3033 Fourth iteration: 𝑖 = 3, 𝑥3 = 1.5, 𝑦3 = 0.1839
y3 y2 f x2 , y2 0.5
y4 y3 f x3 , y3 0.5
0.3033 f (1.0, 0.3033)(0.5)
0.1839 f (1.5, 0.1839)(0.5)
0.1839
0.1116
x3 1.5
x4 2.0
y3 y (1.5) 0.1839
y4 y (2.0) 0.1116
Numerical Methods
by Norhayati Rosli
[Link]
EULER’s METHOD (Cont.)
Solution (Cont.)
The solution is summarised in the following table and figure. Figure below shows the
comparison of the approximate solutions for Example 1 and the exact solutions.
𝒊 𝒙𝒊 𝒚𝒊
0 0 2
1 0.5 0.5
2 1.0 0.3033
3 1.5 0.1839
4 2.0 0.1116
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA
METHODS
General Form of Second Order Runge-Kutta Methods
yi 1 yi (a1k1 a2 k2 )h
k1 f ( x i , yi )
k2 f ( xi p1h, yi q11k1h)
General Formula of RK2
where Methods
a1 a2 1
1
a 2 p1
2
1
a2 q11
2 Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Based on the General Form of RK2, Three Types of Methods
are Developed.
1 Heun Method
Three Types
of RK2
Methods
Ralston Method
3 2 Midpoint Rule
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Heun Method with a Single Corrector (𝑎2 = 1/2)
1
yi 1 yi ( k1 k 2 )h
2
xi 1 xi h
Heun Method Formula
where
k1 f ( xi , yi )
k 2 f ( xi h, yi k1h )
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Midpoint Method (𝑎22 = 1)
yi 1 yi k2 h
xi 1 xi h
where
Midpoint Method Formula
k1 f ( xi , yi )
1 1
k2 f ( xi h, yi k1h)
2 2
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Ralston Method (𝑎2 = 2/3)
1
yi 1 yi (k1 2k2 )h
3
xi 1 xi h
Ralston Method Formula where
k1 f ( xi , yi )
3 3
k2 f ( xi h, yi k1h)
4 4
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Example 2
Use RK2 of Heun method to solve the following ODE (IVP)
dy
exp( x) 2 y, y (0) 2
dx
for 0 ≤ 𝑥 ≤ 2 with a step size, ℎ = 0.5.
Solution
Identify the estimate slope, 𝑓 𝑥, 𝑦 = exp −𝑥 − 2𝑦,
Step 1
initial values, 𝑥0 = 0, 𝑦0 = 2 and ℎ = 0.5.
Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) over the interval 0 ≤ 𝑥 ≤ 2
Step 2
by using Heun method.
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2 Second iteration: 𝑖 = 1, 𝑥1 = 0.5, 𝑦1 = 1.1516
k1 f x0 , y0 f (0, 2) 3 k1 f (0.5,1.516) 1.6967
k2 f x0 h, y0 k1h k2 f (1.0, 0.3033) 0.2387
f (0.5, 0.5) 0.3935 y2 0.6678
1
y1 y0 (k1 k2 )(h) x2 1.0
2 y1 y (1.0) 0.6678
1
2 (3 (0.3935)) 0.5
2
1.1516
x1 x0 h
0 0.5 0.5
Numerical Methods
y1 y (0.5) 1.1516 by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
Third iteration: 𝑖 = 2, 𝑥2 = 1.0, 𝑦2 = 0.6678 Fourth iteration: 𝑖 = 3, 𝑥3 = 1.5, 𝑦3 = 0.3897
k1 f (1.0, 0.6678) 0.9677 k1 f (1.0, 0.3897) 0.5562
k2 f (1.5, 0.1839) 0.1447 k2 f (2.0, 0.1116) 0.0878
y3 0.3897 y4 0.2287
x3 1.5 x4 2.0
y3 y (1.5) 0.3897 y4 y (2.0) 0.2287
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
The solution is summarised in the following table and figure. Figure shows the
comparison of the approximate solutions using Heun method, Euler method and the
exact solutions.
𝒊 𝒙𝒊 𝒚𝒊
0 0 2
1 0.5 1.1516
2 1.0 0.6678
3 1.5 0.3897
4 2.0 0.2287
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Example 3
Use RK2 of Midpoint method to solve the following ODE (IVP)
dy
exp( x) 2 y, y (0) 2
dx
for 0 ≤ 𝑥 ≤ 2 with a step size, ℎ = 0.5.
Solution
Identify the estimate slope, 𝑓 𝑥, 𝑦 = exp −𝑥 − 2𝑦,
Step 1
initial values, 𝑥0 = 0, 𝑦0 = 2 and ℎ = 0.5.
Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) over the interval 0 ≤ 𝑥 ≤ 2
Step 2
by using midpoint method.
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2 Second iteration: 𝑖 = 1, 𝑥1 = 0.5, 𝑦1 = 1.1394
k1 f x0 , y0 f (0, 2) 3 k1 f (0.5,1.1394) 1.6723
k2 f (0.75, 0.7213) 0.9703
k2 f x0 12 h, y0 12 k1h
y2 0.6543
f (0.25,1.25) 1.7212
x2 1.0
y1 y0 k2 (h)
y1 y (1.0) 0.6543
2 (1.7212) 0.5
1.1394
x1 x0 h
0 0.5 0.5
y1 y (0.5) 1.1394
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
Third iteration: 𝑖 = 2, 𝑥2 = 1.0, 𝑦2 = 0.6543 Fourth iteration: 𝑖 = 3, 𝑥3 = 1.5, 𝑦3 = 0.3784
k1 f (1.0, 0.6543) 0.9407 k1 f (1.5, 0.3784) 0.5337
k2 f (1.25, 0.4191) 0.5517 k2 f (1.75, 0.2450) 0.3162
y3 0.3784 y4 0.2203
x3 1.5 x4 2.0
y3 y (1.5) 0.3784 y4 y (2.0) 0.2203
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
The solution is summarised in the following table and figure. Figure below shows the
comparison of the approximate solutions using midpoint method, Euler method and the
exact solutions.
𝒊 𝒙𝒊 𝒚𝒊
0 0 2
1 0.5 1.1394
2 1.0 0.6543
3 1.5 0.3784
4 2.0 0.2203
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Example 4
Use RK2 of Ralston method to solve the following ODE (IVP)
dy
exp( x) 2 y, y (0) 2
dx
for 0 ≤ 𝑥 ≤ 2 with a step size, ℎ = 0.5.
Solution
Identify the estimate slope, 𝑓 𝑥, 𝑦 = exp −𝑥 − 2𝑦,
Step 1
initial values, 𝑥0 = 0, 𝑦0 = 2 and ℎ = 0.5.
Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) over the interval 0 ≤ 𝑥 ≤ 2
Step 2
by using Ralston method.
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2 Second iteration: 𝑖 = 1, 𝑥1 = 0.5, 𝑦1 = 1.1458
k1 f x0 , y0 f (0, 2) 3
k1 f (0.5,1.1458) 1.6850
k2 f x0 34 h, y0 34 k1h k2 f (0.875, 0.5139) 0.6109
f (0.375, 0.8755) 1.0627 y2 0.6613
1 x2 1.0
y1 y0 (k1 2k2 )(h)
3 y1 y (1.0) 0.6613
1
2 (3 2(1.0627)) 0.5
3
1.1458
x1 x0 h
0 0.5 0.5
y1 y (0.5) 1.1458
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
Third iteration: 𝑖 = 2, 𝑥2 = 1.0, 𝑦2 = 0.6613 Fourth iteration: 𝑖 = 3, 𝑥3 = 1.5, 𝑦3 = 0.3843
k1 f (1.0, 0.6613) 0.9547 k1 f (1.5, 0.3843) 0.5454
k2 f (1.375, 0.3033) 0.3537 k2 f (1.875, 0.1797) 0.2061
y3 0.3843 y4 0.2247
x3 1.5 x4 2.0
y3 y (1.5) 0.3843 y4 y (2.0) 0.2247
Numerical Methods
by Norhayati Rosli
[Link]
SECOND ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
The solution is summarised in the following table and figure. Figure below shows the
comparison of the approximate solutions using Ralston method, Euler method and the
exact solutions.
𝒊 𝒙𝒊 𝒚𝒊
0 0 2
1 0.5 1.1458
2 1.0 0.6613
3 1.5 0.3843
4 2.0 0.2247
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA
METHOD
The most popular Runge-Kutta method is often referred to as fourth order
Runge-Kutta (RK4).
It was developed around 1900 by the German Mathematicians C. Runge
and M. W. Kutta.
RK4 is normally known as classical fourth–order RK method.
Fourth Order Runge-Kutta Method Formula
The next value of 𝑦(𝑥𝑖+1 ) is determined by the sum of the current value of
𝑦(𝑥𝑖 ) and the weighted average of four increments. The terms 𝑘’𝑠 represent:
𝑘1 is the increment of the slope at the beginning of the interval, using 𝑦
𝑘2 is the increment of the slope at the midpoint of the interval, using 𝑘1
𝑘3 is the increment of the slope at the midpoint of the interval, using 𝑘2
𝑘4 is the increment of the slope at the end of the interval, using 𝑘3
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA
METHOD (Cont.)
Fourth Order Runge-Kutta Method Formula
h
yi 1 yi k1 2k2 2k3 k4
6
xi 1 xi h
where
k1 f xi , yi
Fourth Order Runge-Kutta
1 1
Method Formula k 2 f xi h, yi k1h
2 2
1 1
k3 f xi h, yi k 2 h
2 2
k 4 f xi h, yi k3 h
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA METHODS
(Cont.)
Example 5
Use RK4 method to solve the following ODE (IVP)
dy
exp( x) 2 y, y (0) 2
dx
for 0 ≤ 𝑥 ≤ 2 with a step size, ℎ = 0.5.
Solution
Identify the estimate slope, 𝑓 𝑥, 𝑦 = exp −𝑥 − 2𝑦,
Step 1
initial values, 𝑥0 = 0, 𝑦0 = 2 and ℎ = 0.5.
Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) over the interval 0 ≤ 𝑥 ≤ 2
Step 2
by using RK4 method.
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2
k1 f (0, 2) 3
k2 f 0 12 (0.5), 2 12 (3)(0.5)
f (0.25,1.25) 1.7212
k3 f 0 12 (0.5), 2 12 (1.7212)(0.5)
f (0.25,1.5697) 2.3606
k4 f 0.5, 2 (2.3606)(0.5)
f (0.5, 0.8197) 1.0329
0.5
y1 2 (3 2(1.7212) 2( 2.3606) ( 1.0329))
6
0.9836
x1 0.5
y1 y (0.5) 0.9836
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
Second iteration: 𝑖 = 1, 𝑥1 = 0.5, 𝑦1 = 0.9836 Third iteration: 𝑖 = 2, 𝑥2 = 1.0, 𝑦2 = 0.5106
k1 f (0.5, 0.9836) 1.3607 k1 f (1.0, 0.5106) 0.6532
k2 f (0.75, 0.6434) 0.8145 k2 f (1.25, 0.3473) 0.4080
k3 f (0.75, 0.7800) 1.0876 k3 f (1.25, 0.4086) 0.5306
k4 f (1.0, 0.4398) 0.5118 k4 f (1.5, 0.2453) 0.2674
y2 0.5106 y3 0.2774
x2 1.0 x3 1.5
y2 y (1.0) 0.5106 y3 y (1.5) 0.2774
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
Fourth iteration: 𝑖 = 3, 𝑥3 = 1.5, 𝑦3 = 0.1839
k1 f (1.5, 0.2774) 0.3317
k2 f (1.75, 0.6434) 0.2152
k3 f (1.75, 0.7800) 0.2734
k4 f (2.0, 0.4398) 0.1460
y4 0.1562
x4 2.0
y4 y (2.0) 0.1562
Numerical Methods
by Norhayati Rosli
[Link]
FOURTH ORDER RUNGE-KUTTA METHODS
(Cont.)
Solution (Cont.)
The solution is summarised in the following table and figure. Figure below shows the
comparison of the approximate solutions using RK4 method, Heun method, Euler
method and the exact solutions.
𝒊 𝒙𝒊 𝒚𝒊
0 0 2
1 0.5 0.9836
2 1.0 0.5106
3 1.5 0.2774
4 2.0 0.1562
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs
Many practical problems in science and engineering need to be modelled
in the form of a system of ODEs rather than single ODE.
In general, such system can be represented as
dy1
f1 ( x, y1 , , yn )
dx
dy2
f 2 ( x, y1 , , yn )
dx
dyn
f n ( x, y1 , , yn )
dx
which requires 𝑛 initial conditions at the starting values of 𝑥.
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
System of ODEs with 𝑛 = 2
dy
f1 ( x, y, z )
In general, a system of two
dx
first order ODEs with 𝑦 and 𝑧 dz
f 2 ( x, y , z )
referred to as dependent dx
variables and 𝑥 referred to as for the domain x0 x xn
independent variable has the
form with initial condition
y ( x0 ) y0 and z ( x0 ) z0
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Euler’s Method for System of ODEs
xi 1 xi h
Euler’s Method Formula yi 1 yi f1 ( xi , yi , zi )h
zi 1 zi f 2 ( xi , yi , zi )h
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
RK4 Method for System of ODEs
yi 1 yi
h
6
k y ,1 2k y ,2 2k y ,3 k y ,4
zi 1 zi k z ,1 2k z ,2 2k z ,3 kz,4
h
6
xi 1 xi h
RK4 Method
k y ,1 f1 xi , yi , zi , kz,1 f 2 xi , yi , zi
Formula
1 1 1
k y ,2 f1 xi h, yi k y ,1h, zi kz,1h
2 2 2
1 1 1
kz,2 f 2 xi h, yi k y ,1h, zi kz,1h
2 2 2
1 1 1
k y ,3 f1 xi h, yi k y ,2 h, zi kz,2 h
2 2 2
1 1 1
kz,3 f 2 xi h, yi k y ,2 h, zi kz,2 h
2 2 2
k y ,4 f1 xi h, yi k y ,3h, z i kz,3h
kz,4 f 2 xi h, yi k y ,3h, z i kz,3h
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Example 6
Use Euler’s method to solve the following system of ODEs
dy
2 y 4e x , y (0) 2
dx
dz yz 2
, z (0) 4
dx 3
for 0 ≤ 𝑥 ≤ 1 with a step size, ℎ = 0.5.
Solution
Identify 𝑓1 𝑥, 𝑦, 𝑧 and 𝑓2 (𝑥, 𝑦, 𝑧)
Step 1
dy
f1 ( x, y, z ) 2 y 4e x , y (0) 2
dx
dz yz 2
f 2 ( x, y , z ) , z (0) 4
dx 3
Step 2 Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) and 𝑧𝑖+1 = 𝑧(𝑥𝑖+1 ) over the
interval 0 ≤ 𝑥 ≤1 by using Euler’s method.
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2, 𝑧0 = 4
y1 2 f1 (0, 2, 4)(0.5)
2
z1 4 f 2 (0, 2, 4)(0.5)
1.3333
x1 0.5
y1 y (0.5) 2
z1 z (0.5) 1.3333
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
Second iteration: 𝑖 = 1, 𝑥1 = 0.5, 𝑦1 = 0.5, 𝑧1 = −1.3333
y2 2 f1 (0.5, 2, 1.3333)(0.5)
1.2131
z2 1.3333 f 2 (0.5, 2, 1.3333)(0.5)
1.9259
x2 1.0
y2 y (1.0) 1.2131
z2 z (1.0) 1.9259
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
The solution is summarised in the following
𝒊 𝒙𝒊 𝒚𝒊 𝒛𝒊
0 0 2 4
1 0.5 2 -1.3333
2 1.0 1.2131 -1.9259
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Example 7
Use RK4 method to solve the following system of ODEs
dy
2 y 4e x , y (0) 2
dx
dz yz 2
, z (0) 4
dx 3
for 0 ≤ 𝑥 ≤ 1 with a step size, ℎ = 0.5.
Solution
Identify 𝑓1 𝑥, 𝑦, 𝑧 and 𝑓2 (𝑥, 𝑦, 𝑧)
Step 1
dy
f1 ( x, y, z ) 2 y 4e x , y (0) 2
dx
dz yz 2
f 2 ( x, y , z ) , z (0) 4
dx 3
Step 2 Approximate iteratively 𝑦𝑖+1 = 𝑦(𝑥𝑖+1 ) and 𝑧𝑖+1 = 𝑧(𝑥𝑖+1 ) over the
interval 0 ≤ 𝑥 ≤1 by using RK4 method.
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2, 𝑧0 = 4
h k y ,1 k
k y ,1 f1 ( x0 , y0, z0 ) k y ,2 f1 ( x0 , y0 h, z0 z ,1 h)
2 2 2
f1 (0, 2, 4) f1 (0.25, 2,1.3334)
2(2) 4e (0) 2(2) 4e (0.25)
0.8848
0
h k y ,1 k z ,1
k z ,1 f 2 ( x0 , y0, z0 ) k z ,2 f 2 ( x0 , y0 h, z 0 h)
2 2 2
f 2 (0, 2, 4) f 2 (0.25, 2,1.3334)
2(4) 2 2(1.3334) 2
3 3
1.1853
10.6667
Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2, 𝑧0 = 4
h k y ,2 k z ,2
k y ,3 f1 ( x0 , y0 h, z 0 h)
2 2 2
f1 (0.25,1.7788, 3.7037)
2(1.7788) 4e (0.25)
0.4424
h k y ,2 k z ,2
k z ,3 f 2 ( x0 , y0 h, z 0 h)
2 2 2
f 2 (0.25,1.7788, 3.7037)
(1.7788)(3.7037) 2
3 Numerical Methods
8.1335
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2, 𝑧0 = 4
k y ,4 f1 ( x0 h, y0 k y ,3 h, z0 k z ,3 h)
f1 (0.5,1.7788, 0.0667)
2(1.7788) 4e (0.5)
1.1315
k z ,4 f 2 ( x0 h, y0 k y ,3 h, z0 k z ,3 h)
f 2 (0.5,1.7788, 0.0667)
(1.7788)( 0.0667) 2
3
2.6379 10 3 Numerical Methods
by Norhayati Rosli
[Link]
SYSTEM OF ODEs (Cont.)
Solution (Cont.)
First iteration: 𝑖 = 0, 𝑥0 = 0, 𝑦0 = 2, 𝑧0 = 4
h
y1 y0 (k y ,1 2k y ,2 2k y ,3 k y ,4 )
6
0.5
2 0 2(0.8848 0.4424) 1.1315
6
1.6845
y (0.5) y1 1.6845 x1 x0 0.5 0.5
h
z1 z0 (k z ,1 2k z ,2 2k z ,3 k z ,4 )
6
4
0.5
6
10.6667 2( 1.1853 8.1335) 2.6379 10 3
1.5578
z (0.5) z1 1.5578 x1 x0 0.5 0.5
Numerical Methods
Repeat the process for 𝑖 = 1 by Norhayati Rosli
[Link]
Conclusion
RK4 method has better order of convergence than Euler and RK2
methods.
However, the main computational effort in applying the RK4 method is one
needs to evaluate four functional evaluations per step.
For instance, in comparing with RK2 method, RK4 method requires twice
as many evaluations per step.
The approximation solution that is obtained by using RK4 method will
provide better approximate solution than Euler and RK2 methods.
Numerical Methods
by Norhayati Rosli
[Link]
Author Information
Norhayati Binti Rosli,
Senior Lecturer,
Faculty of Industrial Sciences & Technology (FIST),
Universiti Malaysia Pahang,
26360 Gambang, Pahang.
SCOPUS ID: 36603244300
UMPIR ID: 3449
Google
Scholars: [Link]
e-mail: norhayati@[Link]
Numerical Methods
by Norhayati Rosli
[Link]