function [z,varargout] = newton(func,z0,tol,maxiter,varargin)
hopt = 2*sqrt(eps);
z = z0;
value = feval(func,z,varargin{:});
iter = 1;
while abs(value) >=tol
deriv = (feval(func,z+hopt,varargin{:}) - value)/hopt;
z = z - value/deriv;
value = feval(func,z,varargin{:});
iter = iter+1;
if iter >= maxiter
fprintf('maxiter exceeded, no convergence')
break
end
end
if nargout >= 2
varargout{1}=value;
if nargout == 3
varargout{2} = iter;
end
end