% This MATLAB code implements the Newton method to solve for roots of the
% input equation myfun with derivative myfunprime.
% Call as newton(initapprox,'myfun','myfunprime',tol)
function[approx]=newton(initapprox,myfun,myfunprime,tol)
% Set up iteration control
maxiter=10000; i = 1;
% Set up initial approx
approx = initapprox;
nextadjust=feval(myfun,approx)/feval(myfunprime,approx);
while abs(nextadjust)>=tol && i<=maxiter;
    % Implement Newton's formula to update the approximation
    approx=approx-nextadjust;
    % Update iteration count
    i=i+1;
    % Print a message if maximum iteraton is reached:
    if i==maxiter;
        fprintf('Max number of %d iterations reached; root may not\n',maxiter);
        fprintf('have been found to within tolerance\n');
    end;
    nextadjust=feval(myfun,approx)/feval(myfunprime,approx);
end;
fprintf('%d iterations',i);