Transcript Runge-Kutta 2nd Order Method for Solving Ordinary

Runge 2nd Order Method
Chemical Engineering Majors
Authors: Autar Kaw, Charlie Barker
http://numericalmethods.eng.usf.edu
Transforming Numerical Methods Education for STEM
Undergraduates
7/21/2015
http://numericalmethods.eng.usf.edu
1
Runge-Kutta 2nd Order Method
http://numericalmethods.eng.usf.edu
Runge-Kutta 2nd Order Method
For
dy
 f ( x, y ), y (0)  y0
dx
Runge Kutta 2nd order method is given by
yi 1  yi  a1k1  a2 k2 h
where
k1  f xi , yi 
k2  f xi  p1h, yi  q11k1h
3
http://numericalmethods.eng.usf.edu
Heun’s Method
Heun’s method
Slope f xi  h, yi  k1h 
y
Here a2=1/2 is chosen
1
a1 
2
yi+1, predicted
Slope f xi , yi 
p1  1
q11  1
resulting in
1 
1
yi 1  yi   k1  k2 h
2 
2
where
k1  f xi , yi 
Average Slope 
yi
xi
1
 f xi  h, yi  k1h   f xi , yi 
2
xi+1
x
Figure 1 Runge-Kutta 2nd order method (Heun’s method)
k 2  f xi  h, yi  k1h
4
http://numericalmethods.eng.usf.edu
Midpoint Method
Here a2  1 is chosen, giving
a1  0
p1 
1
2
q11 
1
2
resulting in
yi 1  yi  k2h
where
k1  f xi , yi 
1
1


k 2  f  xi  h, yi  k1h 
2
2


5
http://numericalmethods.eng.usf.edu
Ralston’s Method
Here a 2  2 is chosen, giving
3
1
3
3
p1 
4
3
q11 
4
resulting in
a1 
2 
1
yi 1  yi   k1  k2 h
3 
3
where
k1  f xi , yi 
3
3


k 2  f  xi  h, yi  k1h 
4
4


6
http://numericalmethods.eng.usf.edu
How to write Ordinary Differential
Equation
How does one write a first order differential equation in the form of
dy
 f  x, y 
dx
Example
dy
 2 y  1.3e  x , y 0   5
dx
is rewritten as
dy
 1.3e  x  2 y, y 0  5
dx
In this case
f x, y   1.3e  x  2 y
7
http://numericalmethods.eng.usf.edu
Example
The concentration of salt, x in a home made soap maker is
given as a function of time by
dx
 37.5  3.5 x
dt
At the initial time, t = 0, the salt concentration in the tank is
50g/L. Using Euler’s method and a step size of h=1.5 min,
what is the salt concentration after 3 minutes.
dx
 37.5  3.5 x
dt
f t , x  37.5  3.5x
8
1 
1
xi 1  xi   k1  k 2 h
2 
2
http://numericalmethods.eng.usf.edu
Solution
Step 1: i  0, t0  0, x0  50
k1  f t0 , xo   f 0,50  37.5  3.550  137.50
k2  f t0  h, x0  k1h   f 0  1.5,50   137.501.5  f 1.5,156.25
 37.5  3.5 156.25  584.38
1 
1
x1  x0   k1  k 2 h
2 
2
1
1

 50    137.50  584.381.5
2
2

 50  223.441.5
x1 is the approximate concentration of
 385.16g / L
salt at
t  t1  t0  h  0  1.5  1.5 min
x1.5  x1  385.16g/L
9
http://numericalmethods.eng.usf.edu
Solution Cont
Step 2: i  1, t1  t 0  h  0  1.5  1.5, x1  385.16g / L
k1  f t1 , x1   f 1.5,385.16  37.5  3.5385.16  1310.5
k 2 f t1  h, x1  k1h   f 1.5  1.5,385.16   1310.51.5  f 3,1580.6
 37.5  3.5 1580.6  5569.8
1 
1
x2  x1   k1  k 2 h
2 
2
1
1

 385.16    1310.5  5569.81.5
2
2

 385.16  2129.6 1.5
x1 is the approximate concentration of
 3579.7 g / L
salt at
t  t2  t1  h  1.5  1.5  3 min
x3  x1  3579.7 g/L
10
http://numericalmethods.eng.usf.edu
Solution Cont
The exact solution of the ordinary differential equation is
given by
x(t )  10.714 39.286e 3.5 x
The solution to this nonlinear equation at t=3 minutes is
x3  10.715 g/L
11
http://numericalmethods.eng.usf.edu
Comparison with exact results
Figure 2. Heun’s method results for different step sizes
12
http://numericalmethods.eng.usf.edu
Effect of step size
Table 1. Effect of step size for Heun’s method
Step size, h
3
1.5
0.75
0.375
0.1875
x3
Et
|t | %
1803.1
−1792.4
16727
3579.6
−3568.9
33306
442.05
−431.34
4025.4
11.038 −0.32231
3.0079
10.718 −0.0024979 0.023311
x(3)  10.715 (exact)
13
http://numericalmethods.eng.usf.edu
Effects of step size on Heun’s
Method
Figure 3. Effect of step size in Heun’s method
14
http://numericalmethods.eng.usf.edu
Comparison of Euler and RungeKutta 2nd Order Methods
Table 2. Comparison of Euler and the Runge-Kutta methods
Step size,
h
3
1.5
0.75
0.375
0.1875
x(3)
Euler
Heun
Midpoint
Ralston
−362.50
720.31
284.65
10.718
10.714
1803.1
3579.6
442.05
11.038
10.718
1803.1
3579.6
442.05
11.038
10.718
1803.1
3579.6
442.05
11.038
10.718
x(3)  10.715 (exact)
15
http://numericalmethods.eng.usf.edu
Comparison of Euler and RungeKutta 2nd Order Methods
Table 2. Comparison of Euler and the Runge-Kutta methods
Step size,
h
3
1.5
0.75
0.375
0.1875
t %
Euler
Heun
Midpoint
Ralston
3483.0
6622.2
2556.5
0.023249
0.010082
16727
33306
4025.4
3.0079
0.023311
16727
33306
4025.4
3.0079
0.023311
16727
33306
4025.4
3.0079
0.023311
x(3)  10.715
16
(exact)
http://numericalmethods.eng.usf.edu
Comparison of Euler and RungeKutta 2nd Order Methods
Figure 4. Comparison of Euler and Runge Kutta 2nd order
methods with exact results.
17
http://numericalmethods.eng.usf.edu
Additional Resources
For all resources on this topic such as digital audiovisual
lectures, primers, textbook chapters, multiple-choice
tests, worksheets in MATLAB, MATHEMATICA, MathCad
and MAPLE, blogs, related physical problems, please
visit
http://numericalmethods.eng.usf.edu/topics/runge_kutt
a_2nd_method.html
THE END
http://numericalmethods.eng.usf.edu