See discussions, stats, and author profiles for this publication at: https://www.researchgate.

net/publication/320124650

A Case Study on Runge Kutta 4 th Order Differential Equations and Its

Application

Article  in  "Imperial Journal of Interdisciplinary Research (IJIR) · January 2017

CITATIONS READS

0 4,077

1 author:

Gowri Palanisamy
Sri Krishna Arts and Science College, Coimatore
7 PUBLICATIONS   6 CITATIONS   

SEE PROFILE

All content following this page was uploaded by Gowri Palanisamy on 30 September 2017.

The user has requested enhancement of the downloaded file.

Imperial Journal of Interdisciplinary Research (IJIR)
Vol-3, Issue-2, 2017
ISSN: 2454-1362, http://www.onlinejournal.in

A Case Study on Runge Kutta 4th Order

Differential Equations and Its Application
Gowri.P1 , Priyadharsini.S2, Maheswari.T3
1
Asst.Prof Sri Krishna Arts and Science College
2,3
Sri Krishna Arts and Science College
Abstract: In this paper it is discussed about Runge Runge-Kutta 4th order method is a numerical
Kutta 4th order with differential equations and its technique used to solve ordinary differential
applications. The diferential equation problems equation of the form
has sloved by Runge Kutta 4th order method and dy/dx = f(x,y),y(0)=y0
application problems are discussed with Runge
Kutta 4th order and codeded in C programming. So only first order ordinary differential equations
can be solved by using the Runge-Kutta 4th order
1.Introduction method. Euler and Runge-Kutta methods are used
to solve higher order ordinary differential
First Order Runge-Kutta Method equations or coupled (simultaneous) differential
Consider the following case: we wish to equations
use a computer to approximate the solution of the
differential equation 2. Derivation
dy(t)dt+2y(t)=0dy(t)dt+2y(t)=0 The runge –kutta 4th order method is based on the
following
or dy(t)dt=−2y(t)
yi 1  yi  a1 k1  a2 k 2  a3 k 3  a4 k 4 )h
Second Order Runge-Kutta Method (1)
The first order Runge-Kutta method used Where knowing the value of y=yi at xi, we can find
the derivative at time t₀ (t₀=0 in the graph below) the value of y=yi+1 at xi+1 , and
to estimate the value of the function at one time
step in the future. If you are not familiar with it, h  xi 1  xi
you should read the section entitled: A First Order
Linear Differential Equation with No Input. We
Equation is equated to the first five terms of Taylor
repeat the central concept of generating a step
series
forward in time in the following text.
dy(t)dt+2y(t)=0ordy(t)dt=−2y(t)

Third Order Runge-Kutta Method

This method is a third order Runge-Kutta method (2)

for approximating the solution of the initial value
problem y'(x) = f(x,y); y(x0) = y0 which evaluates Knowing that dy/dx=f(x,y) and xi+1 –xi =h
the integrand, f(x,y), three times per step. For
step i+1, yi 1  yi  f ( xi , yi )h  1/ 2! f ' ( xi , yi )h2  (3)
yi+1 = yi + 1/6 ( k1 + 4 k2 + k3 ), 1/ 3! f ' ( xi , yi )h3  1/ 4! f "' ( xi , yi )h4
where
k1 = h f(xi, yi), Based on equating equation (2)and equation(3),one
k2 = h f(xi + h / 2, yi + k1 / 2 ), of the popular solutions used is
k3 = h f(xi + h, yi - k1 + 2 k2 ),
Yi+1=yi+1/6(k1+2k2+2k3+k4)h
and xi = x0 + i h.
k1=f (xi,yi)
Fourth Order Runge-Kutta Method
k2=f(xi+h/2,yi+k1/2 h)

Imperial Journal of Interdisciplinary Research (IJIR) Page 134

Imperial Journal of Interdisciplinary Research (IJIR)
Vol-3, Issue-2, 2017
ISSN: 2454-1362, http://www.onlinejournal.in

k3=f(xi+h/2,yi+k2/2) =0.1360
Δy=1/6[0.1+0.2310+0.2344+0.1366]
k4=f(xi+h,yi+k3 h =0.1169
Y1=y0+ Δy
3. Problems =1+0.1169
=1.1169
1. Using fourth order runge-kutta method. y(0.1)=1.1169
Evaluate the value of y when x=1.1 given that For the second approximation we have x1=0.1
dy/dx+y/x=1/x2 y(1) =1 K1=0.1[0.1*(1.1169)+(1.1169)2 ]
=0.1359
Solution : K2=0.1[0.15(1.1849)+(1.1849)2]
The formula for the fourth order Runge-kutta =0.1582
method of the differential equation dy/dx=f(x,y) is K3=0.1[0.15(1.196)+(1.196)2]
given by =0.1610
k1=h *f (x0,y0) K4=0.1[0.2(1.2779)+(1.2779)2]
k2=h*f(x0+h/2,y0+k1/2) =0.1889
k3=h*f(x0+h/2,y0+k2/2) Δy=1/6[0.1359+0.3164+0.3220+0.1889]
k4=h*f(x0+h,y0+k3) =0.1605
Δy=1/6(k1+2k2+2k3+k4) Y2=y1+ Δy
Where h is the interval of differencing and (x0,y0) is =1.1169+0.605
the initial value =1.2774
Hence f (x,y)=1/x2-y/x; x0 =1 y0=1 h=0.1 Y(0.2)=1.2774
K1=(0.1)*(1/12-1/1) 3. Use the classical RK method to estimate
=0
y(0.4) when y' ( x)  x  y 2 with y(0)=0
2
K2=(0.1) [(1/(x0+h/2)2-(y0+k1/2/x0+h/2))]
=(0.1)[(1/(1+0.1/2)2)-(1/1+0.1/2)] assume h=0.2
=-0.00454 SOLUTION :
K3=(0.1)[(0.9070)-(1+(-0.00454/2)/1.05)] f ( x, y)  x 2  y 2
=0.1(0.9070-0.9502)
=-0.00432 m1  f ( x0 , y 0 )  0
K4=(0.1)[(1/(1.1)2)-(1-0.00432/1.1) h mh
=(0.1)(0.8264-0.9052) m2  f ( x 0  , y0  1 )
=-0.00788 2 2
Δy=1/6(0-0.00908-0.00864-0.00788)  f (0.1,0)  0.01
=-0.0042667
h mh
Y1=y(1.1) m3  f ( x0  , y 0  2 )
=y0+ Δy 2 2
=1+(-0.0042667) 0.2 0.01  0.2
=0.9957 = f ( , )
2. Compute y(0.1) and y(0.2) by runge kutta 2 2
method of fourth order differential equation =0.01
dy/dx=xy+y2 ,y(0) =1 m4  f ( x0  h, y 0  m3 h)
Solution:
 f (0.2,0.01  0.2)
The formula for the fourth order runge kutta
method are =0.04
k1=h *f (x0,y0) 0  2  0.01  2  0.01  0.04
k2=h*f(x0+h/2,y0+k1/2)
y(0.2)  0  0.2
6
k3=h*f(x0+h/2,y0+k2/2)
k4=h*f(x0+h,y0+k3) =0.002667
Δy=1/6(k1+2k2+2k3+k4) ITERATION 2
Where h is the interval of differencing and (x0,y0) is X1 =0.2
the initial value Y1 =0.002667
Hence f (x,y)=xy/y2; x0 =1 y0=1 h=0.1
K1=(0.1)(0+1)=0.1
4. Application Problem
K2=(0.1)[(0.05(1.05)+(1.05)2)]
1. In this program for Runge Kutta method in
=0.1155
C, a function f(x,y) is defined to calculate slope
K3=(0.1)[0.05(1.05775)+(1.05775)2
whenever it is called. f(x,y) = (x-y)/(x+y)
=0.1172
SOLUTION:
K4=(0.1)[(0.1)(1.1172)+(1.1172)2]
#include<studio.h>

Imperial Journal of Interdisciplinary Research (IJIR) Page 135

Imperial Journal of Interdisciplinary Research (IJIR)
Vol-3, Issue-2, 2017
ISSN: 2454-1362, http://www.onlinejournal.in

#include<math.h> the bacteria will reduce as fresh water enters the

float f(float x,float y); lake .Find the concentration of the pollutant
void main() after 7 weeks.
{
float x0,y0,m1,m2,m3,m4,m,y,x,h,xn; dC
printf("Enter x0,y0,xn,h:");  0.06c  0
dt
scanf("%f %f %f %f",&x0,&y0,&xn,&h);
x=x0;
y=y0; The differential equation that governs the
printf("\n\nX\t\tY\n"); concentration c of the pollution as a function of
while(x<xn) time (in week) is given by
{
m1=f(x0,y0); dc
m2=f((x0+h/2.0),(y0+m1*h/2.0));  0.06c  0, c(0)  10 6
m3=f((x0+h/2.0),(y0+m2*h/2.0)); dt
m4=f((x0+h),(y0+m3*h));
dc
m=((m1+2*m2+2*m3+m4)/6);  0.06c
y=y+m*h; dt
x=x+h;
printf("%f\t%f\n",x,y); f (t , c)  0.06c
}
getch();
} ci 1  c1  1 / 6(k1  2k 2  2k3  k 4 )h
float f(float x,float y)
{ For i=0,t0 =0,c0 =107
float m;
m=(x-y)/(x+y); k1  f (t 0 , c0 )
return m;
}
= f (0,10 7 )
OUTPUT:
Enter x0,y0,xn,h :
=  0.06(10
7
X= )
Y=
=  600000

2. Using fourth order runge-kutta method. k 2  f (t 0  1 / 2  h, c0  1 / 2k1h)

Evaluate the value of y when x=1.1 given that
dy/dx+y/x=1/x2 y(1) =1
Solution:  f (0  1 / 2  h, c0  1 / 2k1h)
Function[]=runge (f,x0,y0,xn,n)
x=x0, y1=y0;  f (0  1 / 2  3.5,10 7  1 / 2(600000)3.5)
h=xn-x0/n;
disp(‘x rk4’)
= f (1.75,8950000)
for i=1:n
k1=h*f(x,y1);
k2=h*f(x+0.5*h,y1+0.5*k1);  0.06(8950000)
k3=h*f(x+0.5*h,y1+0.5*k2);
k4=h*f(x+h,y1+k3); =  537000
Y1=y1+1/6(k1+2*k2+2*k3+k4);
X=x+h; k3  f (t 0  1 / 2  h, c0  1 / 2k 2 h)
disp ([x y1])
end
disp (‘RK 4th order method:’)  f (0  1 / 2  3.5,10 7  1 / 2(537000)3.5)
disp ([y1])
OUTPUT:  f (1.75,9060300)
f=@(x,y)1+y*y;
Runge(f,1,1, , 5)  0.06(9060300)
3. A polluted lake has an initial concentration of
a bacteria of 107parts/m3 .The concentration of

Imperial Journal of Interdisciplinary Research (IJIR) Page 136

Imperial Journal of Interdisciplinary Research (IJIR)
Vol-3, Issue-2, 2017
ISSN: 2454-1362, http://www.onlinejournal.in

 543620  f (3.5  1 / 2  3.5,8105900  1 / 2(440648)3.5)

k 4  f (t 0  h, c0  k3 h)
 f (5.25,7344100)
 f (0  3.5,10  (543620)3.5)
7

 0.06(7444100)
 f (3.5,8097300)
 440648
 0.06(8097300)
k 4  f (t1  h, c1  k 3 h)
 485840
 f (3.5  3.5,8105900  (440648)3.5)
c1  c0  1/ 6(k1  2k 2  2k3  k 4 )h
= f (7,6563600)
 10 7  1 / 6(600000  2(537000)  2(543620) 
 0.06(6563600)
(485840))3.5
 10 7  1 / 6(3247100)3.5  393820

 8.1059  10 6 parts / m 3 c2  c1  1 / 6(k1  2k 2  2k3  k 4 )h

C1 is the approximate concentration of bacteria at  8105900  1/ 6(486350  2  (435290) 

t=t1=t0+h=0+3.5=3.5 parts/m3
2  (440648)  (393820))  3.5

C(3.5)=C1  8.1059  10 6 parts / m 3

 8105900  1 / 6(2632000)  3.5
For i=1, t1=3.5,C1  8.1059  10
6
 6.5705  10 6 parts / m 3
k1  f (t1 , c1 )

 f (3.5,8.1059  10 6 ) C2 is the approximate concentration of bacteria at

 0.06( 8.1059  10 6 ) t 2  t1  h  3.5  3.5  7weeks C

 486350 C(7)= C2  6.5705  10 6 parts / m3

k 2  f (t1  1 / 2  h, c1  1 / 2k1h) PROGRAM:

#include<conio.h>
#include<iostream.h>
 f (3.5  1 / 2  3.5,8105900  1 / 2(486350)3.5) Void main()
{
float c[10], f ,t[10],h , n;
 f (5.25,7254800) cout<<”Enter the initial values of concentration of
bacteria:”;
int i;
 0.06(7254800) cin>>c[0];
cout<<”Enter the initial value of time:”;cin>>t[0];
 435290 cout<<”Enter the value of time in weeks at which
we want to see the concentration:”;cin>>f;
k3  f (t1  1 / 2  h, c1  1 / 2k 2 h) cout<<”Enter the difference: “;cin>>h;
n=(f-c[0])/h;
for(i=1;i<=n;i++)

Imperial Journal of Interdisciplinary Research (IJIR) Page 137

Imperial Journal of Interdisciplinary Research (IJIR)
Vol-3, Issue-2, 2017
ISSN: 2454-1362, http://www.onlinejournal.in

{  2.2067  10 12 (653.05 4  81 108 )

a=0
k1=-(.06*c[a]); =-0.38347
k2=-h*(.06*(c[a]+k1/2));
k3=-h*(.06*(c[a]+k2/2)); k 3  f (t 0  1 / 2h, 0  1 / 2k 2 h)
k4=-h*(.06*(c[a]+k3));
k=(k1+2*k2+2*k3+k4)/6;
=
y[i]=y[a]+1;
a++; f (0  1 / 2(240),1200  1 / 2(0.38347)  240)
}
Cout<<”/n table “;
For (i-0;i<n;i++) = f (120,1154.0)
{
Cout<<”\t x=%f\ty=%f”,val[i][0],val[i][1]; =  2.2067  10 12 (1154.0 4  81 108 )
Cout<<”\n”;
Getch();
} =  3.8954
Output:
Enter the initial values of concentration of bacteria: k 4  f (t 0  h, 0  k 3 h)
Enter the initial value of time:
Enter the value of time in weeks at which we want = f (0  240,1200  (3.894)  240)
to see the concentration:
Enter the difference:
4. A ball at 1200 K is allowed to cool down
 f (240,265.10)
in air at an ambient temperature of 300 K.
Assuming heat is lost only due to radiation,  2.2067  10 12 (265.10 4  81 108 )
the differential equation for the temperature
of the ball is given by  0.0069750
1   0  1 / 6(k1  2k 2  2k3  k 4 )h
d
 2.2067  10 12 ( 4  81  108 ), (0)  1200k
dt =1200+1/6(-4.5579+2(-0.38954)+0.069750))240

=1200+(-2.1848)  240
f(t,  )=  2.2067  10 12 ( 4  81 108 )
=675.65K
 i 1   i  1 / 6(k1  2k 2  2k3  k 4 )h
 1 is the approximate temperature at
For i=0, t0=0,  0 =1200k
t = t1

k1  f (t 0 , 0 ) =t0+h

=f(0,1200) =0+240

=240
 2.2067  10 12 (1200 4  81 108 )

=-4.5579  1 =  (240)

k 2  f (t 0  1/ 2h, 0  1/ 2k1h) ~675.65K

For i=1,t1=240,  1 =675.65K

=
f (0  1 / 2(240),1200  1 / 2 4.5579)  240
k1  f (t1 ,1 )

=f(120,653.05) = f (240,675.65)

Imperial Journal of Interdisciplinary Research (IJIR) Page 138

 Imperial Journal of Interdisciplinary Research (IJIR)
Vol-3, Issue-2, 2017
ISSN: 2454-1362, http://www.onlinejournal.in

=  2.2067  10 12 (675.65 4  81 108 ) =240+240

=-0.44199 =480

k 2  f (t1  1 / 2h,1  1 / 2k1h)  2 =  (480)

= =594.91K
f (240  1 / 2(240),675.65  1 / 2(0.44199)240)
5. References
= f (360,622.61) [1]Atkinson, Kendall A. (1989), An Introduction to
Numerical Analysis (2nd ed.), New York: John
=  2.2067  10 12 (622.614  81 108 ) Wiley & Sons, ISBN 978-0-471-50023-0.

=-0.31372 [2]Butcher, John C. (May 1963), study of Runge-

Kutta integration processes 3 (2), pp. 185–
k 3  f (t1  1 / 2h,1  1 / 2k 2 h) 201, doi:10.1017/S1446788700027932.

= [3]Butcher, John C. (1975), "A stability property of

f (240  1/ 2(240),675.65  1/ 2(0.31372)  240) implicit Runge-Kutta methods", BIT 15: 358–
361, doi:10.1007/bf01931672.

= f (360,638.00) [4]Butcher, John C. (2008), Numerical Methods for

Ordinary Differential Equations, New York: John
=  2.2067  10 12 (638.00 4  81 108 ) Wiley & Sons, ISBN 978-0-470-72335-7.

=  0.34775
[5] John.H.Mathews , KurtisD.Fink , Numerical
Method Using MATLAB, ISBN 978-81-203-2765-
8.
k 4  f (t1  h,1  k 3 h)
= f (240  240,675.65  (0.34775)  240) [6] Laurene V Fausett , Applied Numerical
Analysis Using MATLAB ,2009 ,ISBN 978-81-
= f (480,592.19) 317-2853-6

= [7]Pallab Ghosh , Numerical Methods with

 2.2067  10 12
(592.19  81 10 )
4 8 computer programs in c++ , ISBN 81-203-2987-2

=-0.25351 [8] Vedamurthy. V N Dr N Ch S Iyengar,

Numerical Methods ,Vellore , ISBN 978-81259-
 2  1  1 / 6(k1  2k 2  2k3  k 4 )h 0630-8.

=
675.65  1/ 6(0.31372)  2(0.34775)  (0.25351))  240
=675.65+1/6(-2.0184)  240

=594.91K

 2 is the approximate temperature at

t  t2

=t1+h

Imperial Journal of Interdisciplinary Research (IJIR) Page 139

View publication stats