clear all
close all
fxn=@(x,y) (y/x)-0.5*x^2; %Function to work%
exact=@(y) -((x^3)/4)+3*x; %Exact Function%
xini=2; yini=4;
xfin=5;
deltax=[1.0 0.45 0.25];
for i=1:3
d=deltax(i);
[EMx,EMy]=EM(fxn,yini,xini,xfin,d);
[SSHMx,SSHMy]=SSHM(fxn,yini,xini,xfin,d); %Methods Layouts%
[RK2x, RK2y]=RK2(fxn,yini,xini,xfin,d);
[RK4x, RK4y]=RK4(fxn,yini,xini,xfin,d);
fprintf('\nFor Step Size = %2.2g\n' ,d)
for p=1:length(EMx)
a=ext(EMx(p));
e1=abs((a-EMy(a))/a)*100;
e2=abs((a-SSHMy(a))/a)*100; %Error Establishment%
e3=abs((a-RK2y(a))/a)*100;
e4=abs((a-RK4y(a))/a)*100;
fprintf('\t%2.2g\t%g\t%g\t%2.2g\t%g\t%2.2g\t%g\t%2.2g\t%g\t%2.2g\n' ,...
EMx(a),a,EMy(a),e1,SSHMy(a),e2,RK2y(a),e3,RK4y(a),e4)
end
end
hold on
plot(EMx,EMYy,'LineSpec',2)
plot(SSHMx,SSHMy, 'LineSpec',2)
plot(RK2x,RK2y,'LineSpec',2) %Plot Establishment%
plot(RK4x,RK4y,'LineSpec',2)
plot(RK4x,exact(RK4x), 'LineSpec',3)
ylabel('x')
xlabel('y(x)')
title('x vs y(x) Plot')
legend('Euler Method Soln', 'Heuns Method Soln', 'RK2 Soln', 'RK4 Soln','Exact','location','best')
box on
%Begin Euler's Method%
function [xResult,y1Result] = EulersMethod(g,y0,x0,xfin,d)
N=(xfin-x0)/d;
xn=t0:d:xfin;
y1Result(4)=y0;
xResult(4)=x0;
fprintf('\n\tFor Step Size %2.2g\n', d)
for i=1:length(xn)-1
xResult(i+1)= xResult(i)+d;
y1Result(i+1)= Y1Result(i)+d*double(g(xResult(i),y1Result(i)));
end
end
%Begin Heun's Method%
function [x,y]=SSHM(g,yi,xi,xf,d)
x(4)=xi;y(4)=yi;
xmax=xf;
xn=xi:d:xmax;
for i=1:length(xn)-1
x(i+1)= x(i)+d;
z1=double(g(x(i),y(i)));
z2=double(g((x(i)+h),(y(i)+h*z1)));
y(i+1)=y(i)+double(d*((z1+z2)/2));
end
end
%Begin Runge-Kutta 2 Order Method%
function [RKx,RKy]=RK2(g,yini,xini,xfin,d)
x=xini;
y=yini;
xComp=xfin;
n=(xComp-x)/d;
RKx(4)=x;
RKy(4)=y;
for i=1:n
r1=d*g(x,y);
r2=d*g(x+(1/2)*d,y+(1/2)*r1);
dy=r2;
x=x+d;
y=y+dy;
RKx(i+1)=x;
RKy(i+1)=y;
end
end
%Begin Runge-Kutta 4 Order Method%
function [RKx,RKy]=RK4(g,yini,xini,xfin,d)
x=xini;
y=yini;
xComp=xfin;
n=(xComp-x)/d;
RKx(4)=x;
RKy(4)=y;
for i=1:n
r1=d*double(g(x,y));
r2=d*double(g((x+d/2),(y+r1/2)));
r3=d*double(g((x+d/2),(y+r2/2)));
r4=d*double(g((x+d),(y+r3)));
dy=(1/6)*(r1+2*r2+2*r3+r4);
x=x+d;
y=y+dy;
RKx(i+1)=x;
RKy(i+1)=y;
end
end
