clear all
clc
format long

syms x

z = 1:500;
fz = 1:500;
prompt = 'Function: ';
f = input(prompt);
prompt = 'x0: ';
z(1) = input(prompt);
fz(1) = subs(f,x,vpa(z(1)));
prompt = 'x1: ';
z(2) = input(prompt);
fz(2) = subs(f,x,vpa(z(2)));
prompt = 'x2: ';
z(3) = input(prompt);
fz(3) = subs(f,x,vpa(z(3)));
prompt = 'Error: ';
r = input(prompt);

i = 3;
error = 10000;

while  error > r
    z(i+1)=(fz(i-2)*(z(i-1)^2-z(i)^2)+fz(i-1)*(z(i)^2-z(i-2)^2)+fz(i)*(z(i-2)^2-z(i-1)^2))/(2*fz(i-2)*(z(i-1)-z(i))+2*fz(i-1)*(z(i)-z(i-2))+2*fz(i)*(z(i-2)-z(i-1)));
    fz(i+1) = subs(f,x,vpa(z(i+1)));
    error = abs(z(i)-z(i-1));
    X = sprintf('x%d: %0.10f | f(x%d): %0.10f | Error: %0.10f',i,z(i+1),i,fz(i+1),error);
    disp(X)
    i=i+1;
    
    if ( i > 500 )
        break
    end
end

