Initial Value Problems for
Ordinary Differential Equations
INTRODUCTION
• Equations which are composed of an unknown
function and its derivatives are called differential
equations.
• It becomes an initial value problem when the original
solution is to be found, given the initial condition.
2
Differential equations play a fundamental role in
engineering because many physical phenomena are
best formulated mathematically in terms of their rate of
change.
• When a function involves one dependent variable, the
equation is called an ordinary differential equation (or
ODE). A partial differential equation (or PDE)
involves two or more independent variables.
• Differential equations are also classified as to their
order.
– A first order equation includes a first derivative as
its highest derivative.
– A second order equation includes a second
derivative, as its highest derivative.
• Higher order equations can be reduced to a system of
first order equations, by redefining a variable.
Single-step method: calculation at current
step only uses value at last step (yn yn+1).
Multi-step method: calculation at current step
uses values at several previous steps.
Types of Errors :
Absolute error : is the numerical difference between
the true value and approximate value of a quantity.
Truncation errors: They are caused by using
approximate results or on replacing an infinite process
by a finite one.
• Round off errors : They arise when a calculated
figure is rounded off to a fixed number of digits to
represent the exact figure, the difference between
such rounded figure and exact figure is called
round off error.
Taylor Series Method
The method gives a straightforward adaptation of
classic calculus to develop the solution of the infinite
series. It is a powerful single step method used to find
the successive derivative easily. It is useful for finding
the starting values for powerful method like Runge-
kutta method, Milne's method .
Example 1:
Solve the initial value problem y' = -2xy2, y(0) = 1 for
y at x = 1 with step length 0.2 using Taylor series
method of order .
Solution:
Given y' = f(x,y) = -2xy2 , y'' = -2y2 - 4xyy'
y''' = -8yy' - 4xy'2 - 4xyy''
yiv = -12y'2 - 12yy'' - 12xy'y'' - 4xyy'''
yv = -48y'y'' - 16yy''' -12xy''2 - 16xy'y''' - 4xyyiv
The fourth order Taylor's formula is
y(xi+h) = y(xi) + h y'(xi, yi) + h2 y''(xi, yi)/2!
+ h3 y'''(xi, yi)/3! + h4 yiv(xi, yi)/4! + . . .
given at x=0, y=1 and with h =0 .2 we have
y' = -2(0)(1)2 = 0.0
y'' = -2(1)2 - 4(0)(1)(0) = -2
y''' = -8(1)(0) - 4(0)(0)2 - 4(0)(1)(-2) = 0.0
yiv = -12(0)2 - 12(1)(-2) - 12(0)(0)(-2) - 4(0)(1)(0) = 24
y(0.2) = 1 + .2 (0) + .22 (-2)/2! + 0 + .24 (24)/4! = .9615
now at x = .2 we have y = .9615
y' = -0.3699, y'' = -1.5648, y''' = 3.9397
and yiv = 11.9953 then
y(0.4) = 1 + .2 (-.3699) + .22 (-1.5648)/2! +
.23(3.9397)/3! + .24 (11.9953)/4! = 0.8624
y(0.6) = 1 + .2 (-.5950) + .22 (-0.6665)/2! + .23
(4.4579)/3! + .24 (-5.4051)/4! = 0.7356
y(0.8) = 1 + .2 (-.6494) + .22 (-0.0642)/2! + .23
(2.6963)/3! + .24 (-10.0879)/4! = 0.6100
y(1.0) = 1 + .2 (-.5953) + .22 (-0.4178)/2! + .23
(0.9553)/3! + .24 (-6.7878)/4! = 0.5001
now at x = 1.0 we have y = .5001,which is the required
solution
Euler’s Method and Modified Euler’s Method
The Taylor –series method may be difficult to
apply if the derivatives become complicated and in
this case, it is difficult to determine the error.
The error in a Taylor series will be small if the step
size h is small. In fact, if we make h small enough,
we may only need a few terms of Taylors series
expansion for good accuracy.
In Euler’s method if h is small, the method is too
slow and if h is large, it gives in accurate value .
It is one of the oldest method.
Euler method Formula: yi+1 = yi+ h f(xi, yi)
Modified Euler Formula :
yi+1 = yi + ½ h [f (xi, yi) + f (xi + h, yi + hf (xi, yi)) ]
Example 2: Find y(1.0) accurate up to four
decimal places using Modified Euler's method by
solving the IVP y' = -2xy2, y(0) = 1 with step
length 0.2.
Solution: f(x, y) = -2xy2
y' = -2×x×y×y, y[0.0] = 1.0 with h = 0.2
Given
y[0.0] = 1.0
Euler Solution: y(1) = y(0) + h× (-2×x×y×y)(1)
y[0.20] = 1.0
Modified Euler iterations:
y(1) = y(0) + 0.5 × h ×((-2 × x × y × y)(0)
+ (-2 × x × y × y)(1)
y[0.20] = 1.0 y[0.20] = 0.95999
y[0.20] = 0.96313 y[0.20] = 0.96289
y[0.20] = 0.96291
Euler Solution:
y(2) = y(1) + h × (-2 × x × y × y)(2)
y[0.40] = 0.88873
Modified Euler iterations:
y(2) = y(1) +0.5 × h × ((-2 × x × y × y)(1)
+ (-2 × x × y × y)(2)
y[0.40] = 0.88873 y[0.40] = 0.86263
y[0.40] = 0.86629 y[0.40] = 0.86578
y[0.40] = 0.86585
Euler Solution:
y(3) = y(2) + h ×(-2 × x × y × y)(3)
y[0.60] = 0.74589
Modified Euler iterations:
y(3) = y(2) + .5 × h × ((-2 × x × y × y)(2)
+ (-2 × x × y × y)(3)
y[0.60] = 0.74589 y[0.60] = 0.73910
y[0.60] = 0.74031 y[0.60] = 0.74010
y[0.60] = 0.74014
Euler Solution:
y(4) = y(3) + h × (-2 × x × y × y)(4)
y[0.80] = 0.60866
Modified Euler iterations:
y(4) = y(3) + .5 × h ×((-2 × x × y × y)(3)
+ (-2 × x × y × y)(4)
y[0.80] = 0.60866 y[0.80] = 0.61512
y[0.80] = 0.61385 y[0.80] = 0.61410
y[0.80] = 0.61405
Euler Solution:
y(5) = y(4) + h ×(-2 × x × y × y)(5)
y[1.00] = 0.49340
Modified Euler iterations:
y(5) = y(4) + .5 × h ×((-2 × x × y × y)(4)
+ (-2 × x × y × y)(5)
y[1.00] = 0.49340
y[1.00] = 0.50504
y[1.00] = 0.50272
y[1.00] = 0.50318
y[1.00] = 0.50309
y[1.00] = 0.50311
Runge -Kutta method
• Runge -Kutta methods generate an accurate solution
without the need to calculate high order derivatives.
• Second order Runge-Kutta method have local
truncation error of order O(h3) and global truncation
error of order O(h2).
• Higher order Runge-Kutta method have better local
and global truncation errors.
4th Order Runge-Kutta Formula:
k1 = f ( xi , yi )
h 1
k2 = f ( xi + , yi + k1h )
2 2
h 1
k3 = f ( xi + , yi + k2h )
2 2
k4 = f ( xi + h, yi + k3h )
yi +1 = yi + (k1 + 2k2 + 2k3 + k4 )
h
6
5 4
Local error is O ( h ) and global error is O ( h )
Example 3 :
dy
= 1 + y + x 2 , y( 0 ) = 0.5
dx
Use Runge − Kutta method to find y( 0.2 ), y( 0.4 )
Solution:
Step 1:
h = 0.2, f(x,y) =1 + y + x 2 , x0 =0, y0 =0.5
k1 = f ( x0 , y0 ) = ( 1 + y0 + x0 ) = 1.5
2
1 1
k 2 = f ( x0 + h , y0 + k1h )
2 2
= 1 + ( y0 + 0.15) + ( x0 + 0.1) = 1.64
2
1 1
k 3 = f ( x0 + h , y 0 + k 2 h )
2 2
= 1 + ( y0 + 0.164 ) + ( x0 + 0.1) = 1.654
2
k 4 = f ( x0 + h , y 0 + k 3 h )
= 1+ (y + 0.16545) + ( x + 0.2 ) = 1.7908
2
0 0
h
y1 = y0 + (k1 + 2k 2 + 2k 3 + k 4 ) = 0.8293
6
Step 2:
h = 0.2, f ( x , y ) = 1 + y + x , x1 = 0.2 , y1 = 0.8293
2
k1 = f ( x1 , y1 ) = 1.7893
1 1
k 2 = f ( x1 + h , y1 + k1 h ) = 1.9182
2 2
1 1
k 3 = f ( x1 + h , y1 + k 2 h ) = 1.9311
2 2
k 4 = f ( x1 + h , y1 + k 3 h ) = 2.0555
0.2
y 2 = y1 + (k1 + 2k 2 + 2k3 + k4 ) = 1.2141
6
Summary of the solution :
xi yi
0.0 0.5
0.2 0.8293
0.4 1.2141
Predictor-Corrector Method
In numerical analysis, a predictor–corrector method is an
algorithm that proceeds in two steps. First, the prediction
step calculates a rough approximation of the desired
quantity. Second, the corrector step refines the initial
approximation using another means.
There are two such methods namely Milne’s method and
Adam’s method
Milne’s Method :
It is also a multi-step method where we assume that
the solution to the given IVP is known at the past four
equally spaced points t0, t1, t2 and t3.
4h
Pr edictor : yn +1 = yn −3 + (2 yn′ − 2 − yn′ −1 + 2 yn′ )
3
h
Corrector : yn +1 = yn −1 + ( yn′ −1 + 4 yn′ + yn′ +1 )
3
Example 4:
Tabulate the solution of
dy
=
t + y, y (0) =
1
dt
in the interval [0, 0.4] with h = 0.1, using Milne’s
Predictor-Corrector method.
Solution :
To find the solution at first four points t0, t1, t2 and
t3, as Milne’s P-C method is not a self – starting
method, we shall use R-K method of fourth order to
get the required solution and then switch over to
Milne’s P – C method.
Thus, taking t0 = 0, t1 = 0.1, t2 = 0.2, t3 = 0.3
we get the corresponding y values using R–K method
of 4th order; that is y0 = 1, y1 = 1.1103, y2 = 1.2428
and y3 = 1.3997
y1′ =t1 + y1 =0.1 + 1.1103 =1.2103
y2′ =t2 + y2 =0.2 + 1.2428 =1.4428
y3′ =t3 + y3 =0.3 + 1.3997 =1.6997
Using Milne’s predictor formula :
4h
P : y4 = y0 + (2 y1′ − y2′ + 2 y3′ )
3
4(0.5)
= 1+ [ 2(1.21103) − 1.4428 + 2(1.69971)]
3
= 1.58363
Before using corrector formula, we compute
y4′= t4 + y4 ( predicted value)
=
0.4 + 1.5836 =
1.9836
Finally, using Milne’s corrector formula, we compute
h
C : y4 =y2 + ( y4′ + 4 y3′ + y2′ )
3
0.1
= 1.2428 + (1.9836 + 6.7988 + 1.4428)
3
= 1.5836
The required solution is:
t 0 0.1 0.2 0.3 0.4
y 1 1.110 1.242 1.399 1.583
Adms-Bashforth method
Adams- predictor formula :
h
y *
n +1
= y n + ( 55 y ′n − 59 y ′n −1 + 37′n −2 − 9 y ′n −3 ) ,
24
Adams- corrector formula :
h
y n +1 = y n + ( 9 y n′+1 + 19 y n′ − 5 y n′−1 + y ′n −2 )
24
Example 5:
Use the Adams-Bashforth- Method with h = 0.2 to
obtain y(0.8) for the solution of y′ = x + y − 1 , y ( 0) = 1
Solution:
With h = 0.2, y(0.8) will be approximated by y4. To get
started, we use the RK4 method with x0 = 0, y0 = 1, h =
0.2 to obtain
y1 = 1.02140000, y2 = 1.09181796,
y3 = 1.22210646
Now with x0 = 0, x1 = 0.2, x3 = 0.4, x4 = 0.6, and
f (x, y) = x + y – 1, we find
y0′ = f ( x0 , y0 ) = ( 0 ) + ( 1 ) − 1 = 0
y1′ = f ( x1 , y1 ) = ( 0.2 ) + ( 1.02140000 ) − 1 = 0.22140000
y′2 = f ( x2 , y 2 ) = ( 0.4 ) + ( 1.09181796 ) − 1 = 0.49181796
y3′ = f ( x3 , y3 ) = ( 0.6 ) + ( 1.22210646 ) − 1 = 0.82210646
Then predictor formula gives
0.2
y 4 = y3 +
*
(55 y3′ − 59 y2′ + 37 y1′ − 9 y0′ ) = 1.42535975
24
To use the corrector formula, we need
′
y4 = f ( x4 , y4 ) = 0.8 + 1.42535975 − 1 = 1.22535975
*
0.2
y 4 = y3 + (9 y4′ + 19 y3′ − 5 y2′ + y1′ ) = 1.42552788
24
Applications
Runge - Kutta method is applied in many
chemical and biological settings.
In hydrology , equations of this type model the
transport and fate of adsorbing contaminants and
microbe-nutrient systems in groundwater.
In chemistry , this equation can stimulate the air
pollution in environmental cases.
Adams–Bashforth scheme used in modeling the
population dynamics in erbium-doped fiber
amplifiers and lasers based on highly doped fibers.
Euler method is used to find the salt
concentration in chemical Engineering and to find
the voltage across the capacitor in electronic
Engineering.
