Roots of Equations - Numerical Methods in Engineering with MATLAB - page 157

Newton-Raphsoniterationdoes not stay within the brackets, it is disregarded and
replacedwith bisection. Since newtonRaphson uses the function f(x) as well as its
derivative, functionroutines forboth (denotedby func and dfunc in the listing) must
be providedbythe user.
function root = newtonRaphson(func,dfunc,a,b,tol)
% Newton-Raphson method combined with bisection for
%findingarootoff(x)=0.
% USAGE: root = newtonRaphson(func,dfunc,a,b,tol)
% INPUT:
% func = handle of function that returns f(x).
% dfunc = handle of function that returns f'(x).
% a,b = brackets (limits) of the root.
% tol = error tolerance (default is 1.0e6*eps).
% OUTPUT:
% root = zero of f(x) (root = NaN if no convergence).
if nargin < 5; tol = 1.0e6*eps; end
fa = feval(func,a); fb = feval(func,b);
if fa == 0; root = a; return; end
if fb == 0; root = b; return; end
if fa*fb > 0.0
error('Root is not bracketed in (a,b)')
end
x=(a+b)/2.0;
fori=1:30
fx = feval(func,x);
if abs(fx) < tol; root = x; return; end
% Tighten brackets on the root
iffa*fx<0.0;b=x;
else;a=x;
end
% Try Newton--Raphson step
dfx = feval(dfunc,x);
ifabs(dfx)==0;dx=b-a;
else; dx = -fx/dfx;
end
x=x+dx;
% If x not in bracket, use bisection
if(b-x)*(x-a)<0.0