function[sol]= heun(f,x0,y0,xn,h)
% This program estimates the solution of the first order differential
% equation of the form y'=f(x,y) at given point xend using Modified- Euler
% Method.
% Input	-> f  : function in the equation which is passed as a string.
%				x0  : initial x value.
%				y0  : initial y value.
%				h	: incremental step-size.
%				xn	: end point where corresponding value of y to be
%                     estiamted.
%               h   : Step size length
% Output	-> The computed table for x and y values are displayed.
% clc;
% clear all;
% close all;
tic;
% *************************************************************************
% f=@(x,y) y;
% x0=0;y0=1;xn=1;h=0.1;
% *************************************************************************
x=x0:h:xn;
n=length(x);
y(1)=y0;
iteration=0;
for i=1:n-1
    y_e=y(i)+h*f(x(i),y(i));
    y(i+1)=y(i)+(h/2)*(f(x(i),y(i))+f(x(i+1),y_e));  
    iteration=iteration+1;
end
sol=y;
e_time=toc;
error=norm(y-exp(x),inf);
M=[x;y];
fprintf('\n\n----------------Modified-Euler Mehod--------------------\n\n');
fprintf('This program was written by Mr. Harish Bhatt on 02,08/2013');
fprintf('\n\n\t\t x-values\t\t y-values\n');
fprintf('\t\t ---------\t\t ---------\n\n');
fprintf('\t\t %f \t\t%f \n',M);
fprintf('\n\n\t\tThe estimated value of y at %f is %f\n',x(end),y(end));
fprintf('\n\n\t\tThis method took %3d step to get the solution\n',iteration+1);
fprintf('\n\n\t\t The elasped time is %f\n',e_time);
fprintf('\n\n\t\t The maximum error is %f\n',error);
% % ------------------------Graphical Representation-------------------------
% % Plot of the numerical solution VS analytical solution y=e^x
% plot(x,y,'^-r',x,exp(x),'-');
% title('Numerical solution with Modified-Euler VS Analytical Solution')
% xlabel('x-values')
% ylabel('y-values')
% legend('Numerical sol','Analytical sol')