ABSTRACT
A. zeroin Given the values, x0 and y0 say, of the end points of an interval assumed to contain a zero of the function f(x), zeroin attempts to find, by use of a combination of linear interpolation, extrapolation, and bisection, values, xn and yn say, of the end points of a smaller interval containing the zero. The method used is described in detail in [BusD75]. At successful exit f(xn)f(yn)≤0, f(xn)≤f(yn), and xn-yn≤2*tol(xn) where tol is a tolerance function prescribed by the user. The value of tol(x) may be changed during the course of the computation or left constant, as the user requires. For example tol(x) may be the function x*re+ae (where re, the permissible relative error, should be greater than the machine precision, and ae, the permissible absolute error, should not be set equal to zero if the prescribed interval contains the origin). The algorithm used results in the construction of five sequences of real numbers xi, ai, bi, ci, and di, and two integer sequences j(i), k(i) (i=1,...,n) as follows. Initially, with x1=y0, if f(x1)≤f(x0) then b1=x1, a1=c1=x0, otherwise b1=x0, a1=c1=x1. Then, for i=2,...,n: (a) j(i) is taken to be the largest integer in the range 1<j(i)<i for which
bj(i)-cj(i)≤½bj(i)-1-cj(i)-1 with j(i)=1 if this condition is never satisfied; (b) with h(b,c) = b + sign(c-b)tol(b) m(b,c) = (b+c)/2 l(b,a) = b - f(b)(b-a)/(f(b)-f(a)) (f(a)≠f(b)) l(b,a) = ∞ (f(a)=f(b)≠0) l(b,a) = b (f(a)=f(b)=0) and w(l,b,c) = l if l is between h(b,c) and m(b,c); w(l,b,c) = h(b,c) if l-b≤tol(b) and l does not lie outside the interval
bounded by b and m(b,c); w(l,b,c) = m(b,c) otherwise and r(b,a,d) determined by the condition that (x-r(b,a,d))/(px+q) = f(x) when x=a,b,d, xi is determined by use of the formula xi = w(l(bi-1,ai-1), bi-1,ci-1) if j(i) ≥ i-2 xi = w(r(bi-1,ai-1,di-1), bi-1,ci-1) if j(i) = i-3 xi = m(bi-1,ci-1) otherwise
(c) k(i) is taken to be the largest integer in the range 0≤k(i)<i for which f(xk(i))f(xi)≤0, and if f(xi)≤f(xk(i), then bi=xi, ci=xk(i), ai=bi-1, otherwise bi=xk(i), ai=ci=xi; also, if bi=xi or bi=bi-1 then di=ai-1, otherwise di=bi-1;
(d) n is the smallest integer for which bn-cn≤2*tol(bn), and then xn=bn, yn=cn. The number of evaluations of fx and tolx is at most 4*log2(x-y)/tau, where x and y are the argument values given upon entry, tau is the minimum of the tolerance function tolx on the initial interval. If upon entry x and y satisfy f(x)*f(y)≤0, then convergence is guaranteed, and the asymptotic order of convergence is 1.618 for simple zeros. Procedure parameters: boolean zeroin (x,y,method) zeroin: on exit zeroin is given the value true when a sufficiently small
subinterval of J containing a zero of the function f(x) has been found; otherwise zeroin is given the value false;
x: double x[0:0]; entry: one endpoint of interval J in which a zero is searched for; exit: a value approximating the zero within the tolerance 2*tol(x) when
zeroin has the value true, and a presumably worthless argument value otherwise;
y: double y[0:0]; entry: the other endpoint of interval J in which a zero is searched for;
upon entry x < y as well as y < x is allowed; exit: the other straddling approximation of the zero, i.e., upon exit the
values of y and x satisfy (a) f(x)*f(y)≤0, (b) x-y≤2*tol(x), and (c) f(x)≤f(y) when zeroin has a value true; method: a class that defines two procedures fx, and tolx, this class must
implement the AP_zeroin_methods interface; double fx(double x[ ]) defines function f as a function depending on the actual parameter
for x; double tolx(double x[ ]) defines the tolerance function tol which may depend on the actual
parameter for x; one should choose tolx positive and never smaller than the precision of the machine's arithmetic at x, i.e., in this arithmetic x+tolx and x-tolx should always yield values distinct from x; otherwise the procedure may get into a loop.
