clear;
clc;

a=-2;
b=-0.8;
r=10;
n_max=100;
h=0.001;
metoda=input("podaj numer metody 1/2/3 ")

select metoda
case 1
    disp("wybrano metodę Newtona");
case 2
    disp("wybrano metodę Steffensena");
case 3
    disp("wybrano metodę Halleya");
end

function y=f(x)
    y=x*exp(sin(x))+cos(x);
endfunction

function y=fp(x,h)
    y=(f(x+h)-f(x))/h;
endfunction

function y=fpp(x,h)
    y=(fp(x+h,h)-fp(x,h))/h;
endfunction

function y=g(x)
    y=(f(x+f(x))-f(x))/f(x);
endfunction

n=0;
exit_code=0;

if(f(a)*f(b)<0 & fp(a,h)*fp(b,h)>0 & fpp(a,h)*fpp(b,h)>0) then
    if(f(a)*fpp(a,h)>0) then
        x(1)=a;
        n=1;
    else
        x(0)=b;
        n=1
    end

    for n=2:n_max
        select metoda
        case 1 then //metoda Newtona
            x(n)=x(n-1)-f(x(n-1))/g(x(n-1));
        case 2 then //metoda Steffensena
            x(n)=x(n-1)-(f(x(n-1))/g(x(n-1)));
        case 3 then //metoda Halleya
            l=f(x(n-1))*fp(x(n-1),h);
            m1=(fp(x(n-1),h))^2;
            m2=f(x(n-1))*fpp(x(n-1),h)/2;
            x(n)=x(n-1)-l/(m1-m2)
        end
        err=abs(x(n)-x(n-1))*abs(f(x(n)))/abs(f(x(n))-f(n-1))
        if err<10^(-r)
            exit_code=3
            break;
    //while(n<n_max)
        //n=n+1;
        //x(n)=x(n-1)-(f(x(n-1))/g(x(n-1)));
        //if(((abs(x(n)-x(n-1))/abs(f(x(n))-f(x(n-1))))))<10^-r) then
          //  exit_code=3;
            //break;
        //end
        //if(n>n_max)
          //  exit_code=2
            //break;
        end
    end
else
    exit_code=1;
end