KIL1005: NUMERICAL METHODS FOR ENGINEERING

SEMESTER 2, SESSION 2022/2023

LECTURE 9: ORDINARY DIFFERENTIAL

EQUATION (ODE)

DR MAHAR DIANA HAMID

DEPARTMENT OF CHEMICAL ENGINEERING
DIFFERENTIAL EQUATION
• Definition: a mathematical equation that E.g.: Newton’s Second Law
consists of an unknown function and its
derivatives. dv c
g v
• It expresses the relation between a dt m
function and its derivatives.
• Differential equations play a fundamental
(1)
role in engineering because many physical
phenomena are best formulated Where:
mathematically in terms of their rate of 𝑣 = velocity
change. g = gravitational constant
t = time
• Examples: reaction rate, mass flow, heat m = mass
transfer rate, acceleration, velocity, etc. c = drag coefficient
𝑣- dependent variable
t- independent variable
2
DIFFERENTIAL EQUATION
• When a function involves one independent variable, the
equation is called an ordinary differential equation (or ODE).
• A partial differential equation (or PDE) involves two or more
independent variables.

3
 EXAMPLE: ODE IN CHEMICAL REACTION
PROBLEM
Hence,
A chemical reaction of the type of dCA
k1  k1C AC B  k 2CC
A B C dt
k2

(2)
dCB
 k1C AC B  k 2CC ODE
dt
takes place in a reactor, the material balance can be dCC
applied as:  k1C AC B  k 2CC
dt
Input + Generation =Output +Accumulation
(3)
For a batch reactor,
Input + Generation =Output +Accumulation
4
CLASSIFICATION OF ODES
ODEs are classified according to their order and their linearity.

The order of an ODE is the order of the highest derivative present in that equation.
1st order
dy
 y  kx
dx
d2y dy
2nd order 2
 y  kx
dx dx
2
d3y d 2 y  dy 
3rd order 3
 a 2  b   kx
dx dx  dx 

5
CLASSIFICATION OF ODES
ODEs also can be classified according to its linearity, i.e. linear and nonlinear ODE.

An ODE is nonlinear if it contains product of the dependent variable or its derivatives or

both.
dy
Linear  y  kx
dx
Nonlinear d 2 y dy
 y  kx
dx 2 dx
2
d3y d2y  dy 
Nonlinear 3
 a 2
 b    kx
dx dx  dx 
6
DEGREE OF ODES
The degree of an ODE is the largest power to which the highest order derivative is
raised.

5 2
 d3y  d 2 y  dy 
 3   a 2  b   kx
 dx  dx  dx  Has degree of 5

2 5 10
 d3y   d2y   dy 
 3   a 2   b   kx
 dx   dx   dx  Has degree of 2
7
RUNGE-KUTTA METHODS
• However not all ODEs can be solved using analytical methods of calculus.
• Need to use numerical methods.
• Methods used to solve ODEs using iterative methods.
• Developed by German mathematicians C. Runge and M.W. Kutta.

dy
 f ( x, y)
dx
8
RUNGE-KUTTA METHODS
• Euler’s method
• Heun’s method
• Midpoint method
• Higher Order Runge-Kutta (RK) method

9
RUNGE-KUTTA METHODS
• Solve ordinary differential equations of the form
dy
 f ( x, y)
dx

• Numerical method in the general form of

𝑁𝑒𝑤 𝑣𝑎𝑙𝑢𝑒 = 𝑜𝑙𝑑 𝑣𝑎𝑙𝑢𝑒 + 𝑠𝑙𝑜𝑝𝑒 × 𝑠𝑡𝑒𝑝 𝑠𝑖𝑧𝑒
or in mathematical terms
𝑦 = 𝑦 + 𝜙ℎ

• Applied step by step to compute out the value in the future (trajectory of the solution).

10
RUNGE-KUTTA METHODS
• The type of RK methods can be differentiated by the manner in which the slope is
estimated.
• The easiest one-step method: Euler’s method.
• However, Euler’s method has the least accurate predictions compared to other RK
methods.

11
EULER’S METHOD

• Named after Leonhard Euler (1707-1783).

• Swiss mathematician and physicist.
• He introduced several mathematical
notations such as f(x), e (also known as
Euler’s number), ∑ for summation and i.
• He also had introduce the Euler’s method to
solve ordinary differential equation.
• It is the simplest method to solve an ODE
with a given initial value, yo.

12
EULER’S METHOD
• The first derivative provides a direct estimate of
the slope at xi
  f ( xi , yi )
where f(xi,yi) is the differential equation evaluated
at xi and yi. This estimate can be substituted into
the equation:

yi 1  yi  f ( xi , yi )h Figure 25.2 Euler’s method

• A new value of y is predicted using the slope to

extrapolate linearly over the step size h.

13
dy
 2 x3  12 x 2  20 x  8.5
dx

14
ERROR ANALYSIS OF THE RK METHOD
Numerical solutions of ODEs involves two types of error:
i) Truncation/discretisation errors
Caused by the nature of the techniques employed to approximate the value of y.

ii) Round off errors

Caused by the limited numbers of significant digits

15
ERROR ANALYSIS OF THE RK METHOD
Truncation errors:
i) Local truncation error that results from an application of the method in question
over single step.
ii) Propagated truncation error that results from the approximation produced
during the previous steps.

16
ERROR ANALYSIS OF THE EULER’S METHOD (CONT.)

Total error = Local truncation error + Propagated truncation error

Global truncation error

17
ERROR ANALYSIS OF THE RK METHOD (CONT.)
Conclusions:
i) The global error can be reduced by decreasing h.
ii) The method will provide error-free predictions if the underlying function is
linear.

18
IMPROVEMENTS OF EULER’S METHOD
• A fundamental source of error in Euler’s method is that the derivative at the
beginning of the interval is assumed to apply across the entire interval.
• Two simple modifications are available to circumvent this shortcoming:

• Heun’s Method

• The Midpoint (or Improved Polygon) Method

19
 Heun’s Method

• One method to improve the estimate of the slope involves the

determination of two derivatives for the interval:

• At the initial point averaged to obtain an improved estimate of

• At the end point the slope for the entire interval.

20
 DERIVATION OF HEUN’S METHOD

• The slope at the beginning of • It provides an estimate of yi+1 -

an interval: an estimated slope at the end of
𝑦′ = 𝑓(𝑥 , 𝑦 ) the interval
𝑦′ = 𝑓(𝑥 , 𝑦 )

• Hence the average slope for the

• Using Euler’s method, the
interval
slope is used to extrapolate
linearly to yi+1: 𝑦′ + 𝑦′
𝑦′ =
2
𝑦 = 𝑦 + 𝑓 𝑥 ,𝑦 ℎ
𝑓 𝑥 , 𝑦 + 𝑓(𝑥 ,𝑦 )
=
Predictor equation 2
21
DERIVATION OF HEUN’S METHOD

• Using the averaged slope to extrapolate linearly from yi to yi+1:

𝑓 𝑥 ,𝑦 + 𝑓(𝑥 ,𝑦 )
𝑦 =𝑦 + ℎ
2
Corrector equation

• Heun’s method – predictor-corrector approach (only method for

ODE).

22
PREDICTOR-CORRECTOR
APPROACH
Predictor : yi01  yi  f ( xi , yi )h
f ( xi , yi )  f ( xi 1 , yi01 )
Corrector : yi 1  yi  h
2

23
Predictor

Corrector

24
HEUN’S METHOD

• Since the corrector equation has yi+1 on both LHS & RHS, it can be applied
in iterative fashion (repeatedly) in order to get an improved estimation of
yi+1.
• A termination criterion of a Heun’s method iterative process:

yy''j ji i11yy''j j11i i11

aa  jj
100
100%%
yy'' i i11
where:
yj-1i+1 : result from the prior iteration
yj-1i+1 : result from the present iteration

25
 .
𝑦 = 4𝑒 − 0.5𝑦

26
27
 THE MID-POINT METHOD

• Also known as improved polygon method or the modified Euler’s

method.
• Uses Euler’s method to predict a value of y at the mid-point of an
interval
h
yi 1/ 2  yi  f ( xi , yi )
2

• The predicted value is used to calculate a slope at the mid-point:

y 'i 1/ 2  f ( xi 1/ 2 , yi 1/ 2 )

which is assumed to be a valid approximation of the average slope

for the entire interval. 28
THE MID-POINT METHOD

• The slope is used to extrapolate linearly from xi to xi+1

yi 1  yi  f ( xi 1/ 2 , yi 1/ 2 )h

• No yi+1 on the RHS => cannot be applied iteratively.

29
THE MID-POINT METHOD

• Midpoint method is superior from

Euler’s method as it utilises a slope
estimation at the midpoint of the
prediction interval.

30
EXAMPLE 3

Solve the following ODE problem using midpoint method over the
interval from x = 0 to x = 2 with a step size of 0.5. The initial condition
at x = 0 is y = 1.

dy
 yx 2  1.2 y
dx

31
 .
𝑦 = 4𝑒 − 0.5𝑦

32
SUMMARY OF IMPROVED EULER’S METHODS
• Required more computational effort to determine the slope.
• Improved accuracy by reducing errors.

1. Heun’s Method (without iteration)

2. Euler’s Method One step approaches
3. Mid-Point Method
Runge-Kutta Methods

33