There are other methods for finding numeric derivatives. For example, you can fit a polynomial (either
Lagrange or least-squares) through some sample points, then find the derivative of the polynomial.
These methods might not have any speed advantage over nder2(), because an accurate, high-order
polynomial requires considerable time to find the coefficients, although it is fast to find the derivative
once you have them.
In my library of numerical methods books I find little reference to finding numerical derivatives. This is
perhaps because it is relatively easy to find derivatives of most simple functions.
Numerical Recipes in Fortran, 2nd edition, William H. Press et al. Section 5.7, p181 describes various
issues and problems in calculating numerical derivatives. The expressions for the errors of avgRC()
and nDeriv() are also found here.
More comments on nder1() and Ridders' method
nder1() is a 'last resort' for finding numerical derivatives, to be used when the other alternatives have
failed. The other alternatives are 1) finding a symbolic derivative, and 2) using the 89/92+ built-in
numerical derivative functions avgRC() and nDeriv().
For example, there is no point in using nder1() for a simple function such as tan(x). The 89/92+ can
easily find the symbolic derivative, which can then be evaluated numerically. However, you may have
complex, programmed functions for which 89/92+ cannot find a symbolic derivative.
The error will increase dramatically near function asymptotes, where the function is changing very
rapidly.
It is possible to adjust h in nDeriv() to get reduce the error at a given point, but that is not very valuable
if you don't know what the answer should be! However, it is possible with a program like nDer1(), which
returns an error estimate, to improve the answer by adjusting h. For example, suppose we want to find
the derivative of tan(x) as accurately as possible at x = 1.4. The basic idea is to call nder() with
different values of h, and try to locate the value of h with the smallest error. The table below shows
some results.
6.97E-93.05E-90.005
6.66E-109.06E-100.01
1.44E-93.83E-100.02
2.71E-104.06E-100.04
4.28E-104.58E-100.06
2.57E-101.03E-90.08
5.78E-107.54E-100.10
1.06E-91.77E-90.13
Actual errorReported errorInterval size h
This example shows that the error is not too sensitive to the interval size, and the actual error is usually
less than the reported error bound.
nder1() uses Ridders' algorithm to estimate the derivative. The general idea is to extrapolate the
central-difference equation to h=0:
f
∏
(x)=
f
(
x+h
)
−f(x−h)
2$h
6 - 35