**Modified Euler’s Method** is a popular method of numerical analysis for integration of initial value problem with the best accuracy and reliability. It solves ordinary differential equations (ODE) by approximating in an interval with slope as an arithmetic average. This method is a simple improvement on Euler’s method in function evaluation per step but leads to yield a second order method.

In earlier tutorials, we’ve already gone through a C program and algorithm/flowchart for this method. Here, we’re going to write a program for **Modified Euler’s method in Matlab**, and discuss its mathematical derivation and a numerical example.

## Derivation of Modified Euler’s Method:

Consider y which is the function of t and solution to an ordinary differential equation. y can be graphically represented as follows:

Let h be height of small intervals taken during the process of approximation.

h = (x_{n} – x_{0})/n

We are provided y(t) as an initial value and we are to evaluate y(t+h)

By using Taylor’s Series:

y(t ) = y( t+ h ) – h [dy ( t + h )]/dt + h^{2}/2 [d^{2}y ( t + h )]/dt^{2 } – O(h^{3})

Rearranging the terms:

y( t+ h ) = y (t ) + h [dy ( t + h )]/dt – h^{2}/2 [d^{2}y ( t + h )]/dt^{2 } + O(h^{3})

Solving this equation:

y( t + h ) = y (t) + h/2 [ dy(t)/dt + dy (t +h ) / dt ] + O(h^{3})

The values of O(h^{3}) and other higher terms are very small and generally neglected. Based on this derivation, we’ll work out the code for Modified Euler’s method in Matlab.

## Modified Euler’s Method in MATLAB:

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 |
function a= eular( df ) %df = @(x, y) -2xy^2 % for calculating value of dy/dx at some particular pt using modified euler's method % asking intial conditions x0 = input('Enter initial value of x : '); y0 = input ('Enter initial value of y : '); x1 = input( 'Enter thevalue of x at which y is to be calculated : '); tol = input( 'Enter desired level of accuracy in the final result : '); % calaulating the value of h n =ceil( (x1-x0)/sqrt(tol)); h = ( x1 - x0)/n %loop for calculating values for k = 1 : n X(1,1) = x0; Y (1,1) = y0; X( 1, k+1) = X(1,k) + h; y_t = Y(1,k) + h* feval( df , X(1,k) , Y(1,k));% Eular's formula Y(1 ,k+1) = Y(1,k) + h/2* (feval( df , X(1,k) , Y(1,k)) + feval( df , X(1,k) + h , y_t ) ) ; %improving results obtained by modified Eular's formula while abs( Y(1,k+1) - y_t ) > h y_t = Y(1,k) + h*feval( df , X(1,k) , Y(1,k+1)); Y( 1 ,k+1) = Y(1,k) + h/2* (feval( df , X(1,k) , Y(1,k)) + feval( df , X(1,k) + h , y_t ) ); end end %displaying results fprintf( 'for \t x = %g \n \ty = %g \n' , x1,Y(1,n+1)) %displaying graph x = 1:n+1; y = Y(1,n+1)*ones(1,n+1) - Y(1,:); plot(x,y,'b') xlabel = (' no of intarval '); ylabel = ( ' Error '); |

The above source code for **Modified Euler’s Method in Matlab** is written for solving ordinary differential equation: **y’ = -2xy ^{2 }**with the initial value condition that, x

_{0}= 0 and y

_{0}= 1. The program can be modified to solve any equation by changing the value of ‘df’ in the code.

This code is a four-parameter input program: it needs initial value of x, initial value of y, the final value of x at which y is to be evaluated and allowed error as the inputs. When all the inputs are given to the program, it displays the value of step height, h and the final value of y at the desired x as output.

The sample output of the program for the function and initial value condition mentioned above is given below:

## Numerical Example of Modified Euler’s Method:

Now, let’s analyze mathematically the above program code of Modified Euler’s method in Matlab. The questions here is:

*Find the value of y(1) by solving the following initial value problem. (Use Modified Euler’s Method.)*

* y’ = -2xy ^{2 }, y(o) = 1 and h = 0.2*

**Solution:**

Given function is, f(x, y) = -2xy^{2 }

Initial value condition, x_{0} = 0 and y_{0} = 1

Step height, h = 0.2

**1.** Approximation for y(t +h) i.e. y( 0 + 0.2) = y ( 0.2 )

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.9599999988079071 * y[0.20] = 0.9631359989929199

y[0.20] = 0.9628947607919341 * y[0.20] = 0.9629133460803093

**2.** Approximation for y(t +h) i.e. y( 0.2 + 0.2) = y ( 0.4 )

y(2) = y(1) + 0 .5 * h * ((-2*x*y*y) (1) + (-2*x*y*y) (2)

y[0.40] = 0.8887359638083165 y[0.40] = 0.8626358081578545

y[0.40] = 0.8662926943348495 y[0.40] = 0.8657868947404332

y[0.40] = 0.865856981554814

**3.** Approximation for y(t +h) i.e. y( 0.4 + 0.2) = y ( 0.6 )

y(3) = y(2) + 0.5 * h* ((-2*x*y*y)(2) + (-2*x*y*y) (3)

y[0.60] = 0.7458966289094106 y[0.60] = 0.7391085349039348

y[0.60] = 0.7403181774980547 y[0.60] = 0.7401034281837107

y[0.60] = 0.7401415785278189

**4.** Approximation for y(t +h) i.e. y( 0.6 + 0.2) = y ( 0.8 )

y(4) = y(3) + 0 .5* h * ((-2*x*y*y) (3) + (-2* x* y* y) (4)

y[0.80] = 0.6086629119889084 y[0.80] = 0.6151235687114084

y[0.80] = 0.6138585343771569 y[0.80] = 0.6141072871136244

y[0.80] = 0.6140584135348263

**5.** Approximation for y(t +h) i.e. y( 0.8 + 0.2) = y ( 1.0 )

y(5) = y(4) + 0.5 * h * ((-2*x*y*y)(4) + (-2*x*y*y)(5)

y[1.00] = 0.49340256427369866 y[1.00] = 0.5050460713552334

y[1.00] = 0.5027209825340415 y[1.00] = 0.5031896121302805

y[1.00] = 0.5030953322323046 y[1.00] = 0.503114306721248

which is the required answer. It is same as that obtained from the program of Modified Euler’s method in Matlab.

If you’ve any questions regarding Modified Euler’s method, its MATLAB program, or its mathematical derivation, bring them up from the comments section. You can find more Numerical Methods tutorials using Matlab here.