## Description

This method is one of the most widely used method and again it is also an iterative method which includes the routine of well known Euler's Method.
Suppose dy/dx = f(x,y) , y = y0 when x = x0 to find y for x = x1 = x0 + h
we have to compute
k1 = hf(x0, y0)
k2 = hf(x0 + h/2, y0 + k1/2)
k3 = hf(x0 + h/2, y0 + k2/2)
k= hf(x0 + h/2, y0 + k3)
k = 1/6( k1 + 2k2 + 2k3 + k4), Now the value of y1 of y for x1 = x0 +k is given by y1 = y0 + k
similarly we have to repeat for next value of x , x2 ,x... and so on

In General for xn
k= hf(xn-1, yn-1)
k2 = hf(xn-1 + h/2, yn-1 + k1/2)
k3 = hf(xn-1 + h/2, yn-1 + k2/2)
k4 = hf(xn-1 + h/2, yn-1 + k3)
k = 1/6( k1 + 2k2 + 2k+ k4),   yn = yn-1 + k

• Programming implementation of RK4 in C++
`#include<iostream>`
`using namespace std;`

`double f(double x,double y){`
`    return x*y + 1;`
`}`
`double f1(double x, double y){`
`    return 3*x + y/2;`
`}`
`int main(){`
`    //Suppose dy/dx = f(x,y) and y = y0 when x=x0 `
`    double x0,y0;`
`    double xn,yn;`
`    double k1,k2,k3,k4, k;`
`    cout<<"Provide initial (x0,y0) whose differential equ is defined:";`
`     cin>>x0>>y0;`
`     cout<<"Enter the value of x for which you need to find out the y=?:";`
`     cin>>xn;`

`     double h = (xn-x0)/5;`
`     yn = y0;`
`     for(double xi = x0 ; xi < xn ; xi += h){`
`         k1 = h*f(xi , yn);`
`         k2 = h*f(xi + h/2, yn + k1/2);`
`         k3 = h*f(xi + h/2, yn + k2/2);`
`         k4 = h*f(xi + h/2, yn + k3);`

`         k = (k1 + 2*k2 + 2*k3 + k4)/6;`

`         yn = yn+k;`
`     }`

`     printf("Value of y at x=%f is %0.4f when h: %f\n",xn , yn, h);`
`    return 0;`
`}`
• Code Execution of RK4 for dy/dx = xy +1 , y(0)=1 at x = 0.1
`Provide initial (x0,y0) whose differential equ is defined:0 1`
`Enter the value of x for which you need to find out the y=?:0.1`
`Value of y at x=0.100000 is 1.1052 when h: 0.020000`
Note: check the value of h, calculate the accuracy of the method and compare with other methods of finding the solution of the ODE.

## Description

Modified Euler's Method is just the extended version of Euler's Method to increase the accuracy of the approximate solution.
It's recommended for you to go through Euler's Method.
The approximate value of y1 of y for x = x1 = x0 +h by Euler's method is given by
y1 = y+ hf(x0, y0)
so in edition to attend more accuracy we repeat the below formula,
y1(1) = y0 + h/2[ f(x0, y0) + f(x1, y1) ]
y1(2) = y0 + h/2[ f(x0, y0) + f(x1, y1(1))]
...
in pan and paper process continues untill two value coincides now we have value of y1
y2 = y1 + hf(x1, y1)
repeat the above process for y2
y2(1) = y1 + h/2[ f(x1, y1) + f(x2, y2) ]
y2(2) = y1 + h/2[ f(x1, y1) + f(x2, y2(1))]
...
• Correctness

• Programming implementation of modified Euler's Method in C++
`#include<iostream>`
`#include<math.h>`
`using namespace std;`

`double f(double x,double y){`
`    return x+y;  // dy/dx = x + y`
`}`
`double f1(double x,double y){`
`    return log2(x+y);`
`}`
`int main(){`
`    // COde for solving Ordinary differential equations`
`    //     dy/dx = f(x,y) `
`    //Define your equation`

`    //Given y(x0) = y0`
`    //Need to find y(xn) = yn`
`    double x0,y0;`
`    double xn,yn;`

`     cout<<"Provide initial (x0,y0) whose differential equ is defined:";`
`     cin>>x0>>y0;`
`     cout<<"Enter the value of x for which you need to find out the y=?:";`
`     cin>>xn;`

`     double h = (xn-x0)/5;`
`     yn = y0;`
`     for(double xi = x0;xi < xn; xi += h){`
`         cout<<xi<<" "<<yn<<endl;`
`         double YnPrevious = yn;`

`         yn = yn + h*f(xi , yn);`
`         for(int i=0;i<5;i++){`
`             yn = YnPrevious + h*( f(xi,YnPrevious) + f(xi+h,yn))/2;`
`         }`
`     //  Yn = Yn-1 + h*f(Xn-1,Yn-1) `
`     }`

`     printf("Value of y at x=%f is %0.4f when h: %f\n",xn , yn, h);`
`    return 0;`
`}`
• Code Execution for dy/dx = x + y ,y(0) =1 ,for x=0.3
`Provide initial (x0,y0) whose differential equ is defined:0 1`
`Enter the value of x for which you need to find out the y=?:0.3`
`0 1`
`0.06 1.06371`
`0.12 1.13507`
`0.18 1.21456`
`0.24 1.30268`
`Value of y at x=0.300000 is 1.4000 when h: 0.060000`
Note: In code, h is defined constant depending upon interval size, but you can change as h cause accuracy difference, try out the code for other functions.

## Description

Consider a defferential equation dy/dx = f(x,y) with initial condition y = y0 when x = x0, suppose we want to find y for x = xn.
First we divide the interval [x0,xn] into n equal parts. x0, x1 = x0 +h , x2 = x0 + 2h, .. , xn where h is small. Using the equation of tangent at point (x0, y0)
y - y0 = dy/dx * (x - x0)  where dy/dx = f(x0,y0)

In the small interval [x0, x1] or [x0, x0+h]  we approximate the solution curve by the tangent, So the approximate value of y1 at point y1 can be given by
y1 - y0 = f(x0, y0)( x0 + h - x0)
y1 = y0 + hf(x0, y0)
Similarly
y2 = y+ hf(x1, y1)
...
yn = yn-1 + hf(xn-1, yn-1)

In the Euler's method we approximate the solution curve by tangent in each subinterval by a sequence of short lines and therefore unless h is very small the error is bound to be quite significant.
• Programming implementation of Euler's Method in C++
`#include<iostream>`
`#include<math.h>`
`using namespace std;`

`double f(double x,double y){`
`    return x+y;  // dy/dx = x + y`
`}`
`double f1(double x,double y){`
`    return log2(x+y);`
`}`
`int main(){`
`    // COde for solving Ordinary differential equations`
`    //     dy/dx = f(x,y) `
`    //Define your equation`

`    //Given y(x0) = y0`
`    //Need to find y(xn) = yn`
`    double x0,y0;`
`    double xn,yn;`

`     cout<<"Provide initial (x0,y0) whose differential equ is defined:";`
`     cin>>x0>>y0;`
`     cout<<"Enter the value of x for which you need to find out the y=?:";`
`     cin>>xn;`

`     double h = (xn-x0)/5;`
`     yn = y0;`
`     for(double xi = x0;xi < xn; xi += h){`
`         cout<<yn<<endl;`
`         yn = yn + h*f(xi , yn);`
`     //  Yn = Yn-1 + h*f(Xn-1,Yn-1) `
`     }`

`     printf("Value of y at x=%f is %0.4f when h: %f\n",xn , yn, h);`
`    return 0;`
`}`
• Code execution for dy/dx = x + y ,y(0) = 1 ,finding y for x = 0.1
`Provide initial (x0,y0) whose differential equ is defined:0 1`
`Enter the value of x for which you need to find out the y=?:0.1`
`1`
`1.02`
`1.0408`
`1.06242`
`1.08486`
`Value of y at x=0.100000 is 1.1082 when h: 0.020000`

Note: In code their h is defined constant, change it if needed as accuracy completely dependent on the value of h.