function xNew = stephensen(f,xOld,tol)
fprintf(' Stephesens root-solving method\n\n');
fprintf(' n |     x(n)     |  |x(n)-x(n-1)|  |     f(x(n)) \n');
fprintf('_____________________________________________________\n');


for n = 1 : 1000
    xNew = xOld - f(xOld)^2 / (f(xOld + f(xOld)) - f(xOld));
    fprintf('%2d | %12.10g | %15.10g | %15.10g \n', n, xNew, abs(xNew-xOld), f(xNew));
    if (abs(xNew - xOld) < tol)
        return
    end
    xOld = xNew;
end
fprintf('Method did not converge after 1000 steps\n');
