c05azf locates a simple zero of a continuous function in a given interval by using Brent's method, which is a combination of nonlinear interpolation, linear extrapolation and bisection. It uses reverse communication for evaluating the function.

2 Specification

Fortran Interface

Subroutine c05azf (

x, y, fx, tolx, ir, c, ind, ifail)

Integer, Intent (In)	::	ir
Integer, Intent (Inout)	::	ind, ifail
Real (Kind=nag_wp), Intent (In)	::	fx, tolx
Real (Kind=nag_wp), Intent (Inout)	::	x, y, c(17)

C Header Interface

#include <nag.h>

void	c05azf_ (double x, double y, const double fx, const double tolx, const Integer ir, double c[], Integer ind, Integer *ifail)

The routine may be called by the names c05azf or nagf_roots_contfn_brent_rcomm.

3 Description

You must supply x and y to define an initial interval

[a, b]

containing a simple zero of the function

f (x)

(the choice of x and y must be such that

f (x) \times f (y) \leq 0.0

). The routine combines the methods of bisection, nonlinear interpolation and linear extrapolation (see Dahlquist and Björck (1974)), to find a sequence of sub-intervals of the initial interval such that the final interval

[x, y]

contains the zero and

| x - y |

is less than some tolerance specified by tolx and ir (see Section 5). In fact, since the intermediate intervals

[x, y]

are determined only so that

f (x) \times f (y) \leq 0.0

, it is possible that the final interval may contain a discontinuity or a pole of

f

(violating the requirement that

f

be continuous). c05azf checks if the sign change is likely to correspond to a pole of

f

and gives an error return in this case.

A feature of the algorithm used by this routine is that unlike some other methods it guarantees convergence within about

{(\log_{2} [(b - a) / δ])}^{2}

function evaluations, where

δ

is related to the argument tolx. See Brent (1973) for more details.

c05azf returns to the calling program for each evaluation of

f (x)

. On each return you should set

fx = f (x)

and call c05azf again.

The routine is a modified version of procedure ‘zeroin’ given by Brent (1973).

4 References

Brent R P (1973) Algorithms for Minimization Without Derivatives Prentice–Hall

Bus J C P and Dekker T J (1975) Two efficient algorithms with guaranteed convergence for finding a zero of a function ACM Trans. Math. Software 1 330–345

Dahlquist G and Björck Å (1974) Numerical Methods Prentice–Hall

5 Arguments

Note: this routine uses reverse communication. Its use involves an initial entry, intermediate exits and re-entries, and a final exit, as indicated by the argument ind. Between intermediate exits and re-entries, all arguments other than fx must remain unchanged.

1: $x$ – Real (Kind=nag_wp) Input/Output

2: $y$ – Real (Kind=nag_wp) Input/Output

On initial entry: x and y must define an initial interval

[a, b]

containing the zero, such that

f (x) \times f (y) \leq 0.0

. It is not necessary that

x < y

On intermediate exit: x contains the point at which

f

must be evaluated before re-entry to the routine.

On final exit: x and y define a smaller interval containing the zero, such that

f (x) \times f (y) \leq 0.0

, and

| x - y |

satisfies the accuracy specified by tolx and ir, unless an error has occurred. If

ifail = 4

, x and y generally contain very good approximations to a pole; if

ifail = 5

, x and y generally contain very good approximations to the zero (see Section 6). If a point x is found such that

f (x) = 0.0

, on final exit

x = y

(in this case there is no guarantee that x is a simple zero). In all cases, the value returned in x is the better approximation to the zero.

3: $fx$ – Real (Kind=nag_wp) Input

On initial entry: if

ind = 1

, fx need not be set.

ind = −1

, fx must contain

f (x)

for the initial value of x.

On intermediate re-entry: must contain

f (x)

for the current value of x.

4: $tolx$ – Real (Kind=nag_wp) Input

On initial entry: the accuracy to which the zero is required. The type of error test is specified by ir.

Constraint:

tolx > 0.0

5: $ir$ – Integer Input

On initial entry: indicates the type of error test.

$ir = 0$: The test is: $| x - y | \leq 2.0 \times tolx \times \max (1.0, | x |)$ .
$ir = 1$: The test is: $| x - y | \leq 2.0 \times tolx$ .
$ir = 2$: The test is: $| x - y | \leq 2.0 \times tolx \times | x |$ .

Suggested value:

ir = 0

Constraint:

ir = 0

1

2

6: $c (17)$ – Real (Kind=nag_wp) array Input/Output

On initial entry: if

ind = 1

, no elements of c need be set.

ind = −1

c (1)

must contain

f (y)

, other elements of c need not be set.

On final exit: is undefined.

7: $ind$ – Integer Input/Output

On initial entry: must be set to

1

−1

$ind = 1$: fx and $c (1)$ need not be set.
$ind = −1$: fx and $c (1)$ must contain $f (x)$ and $f (y)$ respectively.

On intermediate exit: contains

2

3

4

. The calling program must evaluate

f

at x, storing the result in fx, and re-enter c05azf with all other arguments unchanged.

On final exit: contains

0

Constraint: on entry

ind = −1

1

2

3

4

Note: any values you return to c05azf as part of the reverse communication procedure should not include floating-point NaN (Not a Number) or infinity values, since these are not handled by c05azf. If your code does inadvertently return any NaNs or infinities, c05azf is likely to produce unexpected results.

8: $ifail$ – Integer Input/Output

On initial entry: ifail must be set to

0

−1

1

to set behaviour on detection of an error; these values have no effect when no error is detected.

A value of

0

causes the printing of an error message and program execution will be halted; otherwise program execution continues. A value of

−1

means that an error message is printed while a value of

1

means that it is not.

If halting is not appropriate, the value

−1

1

is recommended. If message printing is undesirable, then the value

1

is recommended. Otherwise, the value

−1

is recommended since useful values can be provided in some output arguments even when

ifail \neq 0

on exit. When the value $- 1$ or $1$ is used it is essential to test the value of ifail on exit.

On final exit:

ifail = 0

unless the routine detects an error or a warning has been flagged (see Section 6).

6 Error Indicators and Warnings

If on entry

ifail = 0

−1

, explanatory error messages are output on the current error message unit (as defined by x04aaf).

Errors or warnings detected by the routine:

$ifail = 1$: On entry, $f (x)$ and $f (y)$ have the same sign with neither equalling $0.0$ : $f (x) = ⟨ value ⟩$ and $f (y) = ⟨ value ⟩$ .

$ifail = 2$: On entry, $ind = ⟨ value ⟩$ .
Constraint: $ind = −1$ , $1$ , $2$ , $3$ or $4$ .

$ifail = 3$: On entry, $ir = ⟨ value ⟩$ .
Constraint: $ir = 0$ , $1$ or $2$ .

On entry, $tolx = ⟨ value ⟩$ .
Constraint: $tolx > 0.0$ .

$ifail = 4$: The final interval may contain a pole rather than a zero. Note that this error exit is not completely reliable: it may be taken in extreme cases when $[x, y]$ contains a zero, or it may not be taken when $[x, y]$ contains a pole. Both these cases occur most frequently when tolx is large.

$ifail = 5$: The tolerance tolx has been set too small for the problem being solved. However, the values x and y returned may well be good approximations to the zero. $tolx = ⟨ value ⟩$ .

$ifail = - 99$: An unexpected error has been triggered by this routine. Please contact NAG.
See Section 7 in the Introduction to the NAG Library FL Interface for further information.

$ifail = - 399$: Your licence key may have expired or may not have been installed correctly.
See Section 8 in the Introduction to the NAG Library FL Interface for further information.

$ifail = - 999$: Dynamic memory allocation failed.
See Section 9 in the Introduction to the NAG Library FL Interface for further information.

7 Accuracy

The accuracy of the final value x as an approximation of the zero is determined by tolx and ir (see Section 5). A relative accuracy criterion (

ir = 2

) should not be used when the initial values x and y are of different orders of magnitude. In this case a change of origin of the independent variable may be appropriate. For example, if the initial interval

[x, y]

is transformed linearly to the interval

[1, 2]

, then the zero can be determined to a precise number of figures using an absolute (

ir = 1

) or relative (

ir = 2

) error test and the effect of the transformation back to the original interval can also be determined. Except for the accuracy check, such a transformation has no effect on the calculation of the zero.

8 Parallelism and Performance

Background information to multithreading can be found in the Multithreading documentation.

c05azf is not threaded in any implementation.

9 Further Comments

For most problems, the time taken on each call to c05azf will be negligible compared with the time spent evaluating

f (x)

between calls to c05azf.

If the calculation terminates because

f (x) = 0.0

, then on return y is set to x. (In fact,

y = x

on return only in this case and, possibly, when

ifail = 5

.) There is no guarantee that the value returned in x corresponds to a simple root and you should check whether it does. One way to check this is to compute the derivative of