function [x,k,xk,fxk] = newton2(f,df,x,tol,maxit)
%newton2 Finds root
%https://en.wikipedia.org/wiki/Newton's_method
    k=1;
    fx = f(x);
    xk(k) = x;
    fxk(k) = fx;
    while (k == 1 || norm(x-xk(end-1)) >= tol) && k < maxit
        x = x - fx/df(x);
        fx = f(x);
        k = k+1;
        xk(k) = x;
        fxk(k) = fx;
    end
end
