clear all
clc

%RK IV

f=@(x,y) exp(1).^(-1*x) -2*x;
y0 = -1;

x1(1) = 0;
y1(1) = y0;
n1 = 1;
h1 = 0.01;
ilosc_dzialan_dla_h1 = 0;
tic
while x1<1
    %tic;
    x1(n1+1) = x1(n1) + h1;
    k11 = h1*f(y1(n1),x1(n1));
    k21 = h1*f(y1(n1)+(k11/2),x1(n1)+(h1/2));
    k31 = h1*f(y1(n1)+(k21/2),x1(n1)+(h1/2));
    k41 = h1*f(y1(n1)+k31,x1(n1)+h1);
    phi1 = (1/6)*(k11+2*k21+2*k31+k41);
    y1(n1+1) = y1(n1)+phi1;
    n1 = n1+1;
    ilosc_dzialan_dla_h1 = ilosc_dzialan_dla_h1 + 1;
    %toc;
end
ilosc_dzialan_dla_h1
y1
Time1 = toc
plot(x1,y1,'r-');

rown1=@(x,y) exp(1).^(-x)-2*x;
tspan=0:0.01:1;
[x,y]=ode45(rown1, tspan, -1);
plot(x,y,'r+')

hold on



x2(1) = 0;
y2(1) = y0;
n2 = 1;
h2 = 0.001;
ilosc_dzialan_dla_h2 = 0;
tic
while x2<1
    %tic;
    x2(n2+1) = x2(n2) + h2;
    k12 = h2*f(y2(n2),x2(n2));
    k22= h2*f(y2(n2)+(k12/2),x2(n2)+(h2/2));
    k32 = h2*f(y2(n2)+(k22/2),x2(n2)+(h2/2));
    k42 = h2*f(y2(n2)+k32,x2(n2)+h2);
    phi2 = (1/6)*(k12+2*k22+2*k32+k42);
    y2(n2+1) = y2(n2)+phi2;
    n2 = n2+1;
    ilosc_dzialan_dla_h2 = ilosc_dzialan_dla_h2 + 1;
    %toc;
end

ilosc_dzialan_dla_h2
y2
Time2 = toc
plot(x2,y2,'b-');

rown2=@(x,y) exp(1).^(-x)-2*x;
tspan2=0:0.001:1;
[x,y]=ode45(rown2, tspan2, -1);
plot(x,y,'y-o')

hold on


x3(1) = 0;
y3(1) = y0;
n3 = 1;
h3 = 0.0001;
ilosc_dzialan_dla_h3 = 0;
tic
while x3<1
    %tic;
    x3(n3+1) = x3(n3) + h3;
    k13 = h3*f(y3(n3),x3(n3));
    k23 = h3*f(y3(n3)+(k13/2),x3(n3)+(h3/2));
    k33 = h3*f(y3(n3)+(k23/2),x3(n3)+(h3/2));
    k43 = h3*f(y3(n3)+k33,x3(n3)+h3);
    phi3 = (1/6)*(k13+2*k23+2*k33+k43);
    y3(n3+1) = y3(n3)+phi3;
    n3 = n3+1;
    ilosc_dzialan_dla_h3 = ilosc_dzialan_dla_h3 + 1;
    %toc;
end

ilosc_dzialan_dla_h3
y3
Time3 = toc
plot(x3,y3,'g-');
hold on

rown1=@(x,y) exp(1).^(-x)-2*x;
tspan3=0:0.0001:1;
[x,y]=ode45(rown1, tspan3, -1);
plot(x,y,'-o')


%legend('Wykres dla kroku : h1','Wykres dla kroku : h2','Wykres dla kroku : h3','Dokładny wykres ode23');