1- Bisection method:
function root = bisectnew(func,xl,xu,Ead)
% bisectnew(xl,xu,es,maxit):
% uses bisection method to find the root of a function
% with a fixed number of iterations to attain
% a prespecified tolerance
% input:
% func = name of function
% xl, xu = lower and upper guesses
% Ead = (optional) desired tolerance (default = 0.000001)
% output:
% root = real root
if func(xl)*func(xu)>0 %if guesses do not bracket a sign change
error('no bracket') %display an error message and terminate
end
% if necessary, assign default values
if nargin<4, Ead = 0.000001; end %if Ead blank set to 0.000001
% bisection
xr = xl;
% compute n and round up to next highest integer
n = round(1 + log2((xu - xl)/Ead) + 0.5);
for i = 1:n
xrold = xr;
xr = (xl + xu)/2;
if xr ~= 0, ea = abs((xr - xrold)/xr) * 100; end
test = func(xl)*func(xr);
if test < 0
xu = xr;
elseif test > 0
xl = xr;
else
ea = 0;
end
end
root = xr;
2- False-position method:
function root = falsepos(func,xl,xu,es,maxit)
% falsepos(xl,xu,es,maxit):
% uses the false position method to find the root
% of the function func
% input:
% func = name of function
% xl, xu = lower and upper guesses
% es = (optional) stopping criterion (%) (default = 0.001)
% maxit = (optional) maximum allowable iterations (default = 50)
% output:
% root = real root
if func(xl)*func(xu)>0 %if guesses do not bracket a sign change
error('no bracket') %display an error message and terminate
end
% default values
if nargin<5, maxit=50; end
if nargin<4, es=0.001; end
% false position
iter = 0;
xr = xl;
while (1)
xrold = xr;
xr = xu - func(xu)*(xl - xu)/(func(xl) - func(xu));
iter = iter + 1;
if xr ~= 0, ea = abs((xr - xrold)/xr) * 100; end
test = func(xl)*func(xr);
if test < 0
xu = xr;
elseif test > 0
xl = xr;
else
ea = 0;
end
if ea <= es | iter >= maxit, break, end
end
root = xr;
3- Secant method
function root = secant(func,xrold,xr,es,maxit)
% secant(func,xrold,xr,es,maxit):
% uses secant method to find the root of a function
% input:
% func = name of function
% xrold, xr = initial guesses
% es = (optional) stopping criterion (%)
% maxit = (optional) maximum allowable iterations
% output:
% root = real root
% if necessary, assign default values
if nargin<5, maxit=50; end %if maxit blank set to 50
if nargin<4, es=0.001; end %if es blank set to 0.001
% Secant method
iter = 0;
while (1)
xrn = xr - func(xr)*(xrold - xr)/(func(xrold) - func(xr));
iter = iter + 1;
if xrn ~= 0, ea = abs((xrn - xr)/xrn) * 100; end
if ea <= es | iter >= maxit, break, end
xrold = xr;
xr = xrn;
end
root = xrn;
4- Modified secant method
function root = modsec(func,xr,delta,es,maxit)
% secant(func,xrold,xr,es,maxit):
% uses the modified secant method
% to find the root of a function
% input:
% func = name of function
% xr = initial guess
% delta = perturbation fraction
% es = (optional) stopping criterion (%)
% maxit = (optional) maximum allowable iterations
% output:
% root = real root
% if necessary, assign default values
if nargin<5, maxit=50; end %if maxit blank set to 50
if nargin<4, es=0.001; end %if es blank set to 0.001
if nargin<3, delta=1E-5; end %if delta blank set to 0.00001
% Secant method
iter = 0;
while (1)
xrold = xr;
xr = xr - delta*xr*func(xr)/(func(xr+delta*xr)-func(xr));
iter = iter + 1;
if xr ~= 0, ea = abs((xr - xrold)/xr) * 100; end
if ea <= es | iter >= maxit, break, end
end
root = xr;
5- Newton-Raphson method
function [x,fx,xx] = newton(f,df,x0,TolX,MaxIter)
%newton.m to solve f(x) = 0 by using Newton method.
%input: f = ftn to be given as a string ’f’ if defined in an M-file
% df = df(x)/dx (If not given, numerical derivative is used.)
% x0 = the initial guess of the solution
% TolX = the upper limit of |x(k) - x(k-1)|
% MaxIter = the maximum # of iteration
%output: x = the point which the algorithm has reached
% fx = f(x(last)), xx = the history of x
h = 1e-4; h2 = 2*h; TolFun=eps;
if nargin == 4 & isnumeric(df), MaxIter = TolX; TolX = x0; x0 = df; end
xx(1) = x0; fx = feval(f,x0);
for k = 1: MaxIter
if ~isnumeric(df), dfdx = feval(df,xx(k)); %derivative function
else dfdx = (feval(f,xx(k) + h)-feval(f,xx(k) - h))/h2; %numerical drv
end
dx = -fx/dfdx;
xx(k+1) = xx(k)+dx; %Eq.(4.4.2)
fx = feval(f,xx(k + 1));
if abs(fx)<TolFun | abs(dx) < TolX, break; end
end
x = xx(k + 1);
if k == MaxIter, fprintf(’The best in %d iterations\n’,MaxIter), end