In
numerical analysis, numerical differentiation
algorithms estimate the
derivative of a
mathematical function or function
subroutine using values of the function and perhaps other knowledge about the function.
Finite differences
The simplest method is to use finite difference approximations.
A simple two-point estimation is to compute the slope of a nearby
secant line through the points (''x'', ''f''(''x'')) and (''x'' + ''h'', ''f''(''x'' + ''h'')). Choosing a small number ''h'', ''h'' represents a small change in ''x'', and it can be either positive or negative. The slope of this line is
:
This expression is
Newton's
difference quotient (also known as a first-order
divided difference
In mathematics, divided differences is an algorithm, historically used for computing tables of logarithms and trigonometric functions. Charles Babbage's difference engine, an early mechanical calculator, was designed to use this algorithm in ...
).
The slope of this secant line differs from the slope of the tangent line by an amount that is approximately proportional to ''h''. As ''h'' approaches zero, the slope of the secant line approaches the slope of the tangent line. Therefore, the true derivative of ''f'' at ''x'' is the limit of the value of the difference quotient as the secant lines get closer and closer to being a tangent line:
:
Since immediately
substituting 0 for ''h'' results in
indeterminate form , calculating the derivative directly can be unintuitive.
Equivalently, the slope could be estimated by employing positions (''x'' − ''h'') and ''x''.
Another two-point formula is to compute the slope of a nearby secant line through the points (''x'' - ''h'', ''f''(''x'' − ''h'')) and (''x'' + ''h'', ''f''(''x'' + ''h'')). The slope of this line is
:
This formula is known as the
symmetric difference quotient. In this case the first-order errors cancel, so the slope of these secant lines differ from the slope of the tangent line by an amount that is approximately proportional to
. Hence for small values of ''h'' this is a more accurate approximation to the tangent line than the one-sided estimation. However, although the slope is being computed at ''x'', the value of the function at ''x'' is not involved.
The estimation error is given by
:
,
where
is some point between
and
.
This error does not include the
rounding error due to numbers being represented and calculations being performed in limited precision.
The symmetric difference quotient is employed as the method of approximating the derivative in a number of calculators, including
TI-82,
TI-83,
TI-84,
TI-85, all of which use this method with ''h'' = 0.001.
Step size
An important consideration in practice when the function is calculated using
floating-point arithmetic of finite precision is the choice of step size, ''h''. If chosen too small, the subtraction will yield a large
rounding error. In fact, all the finite-difference formulae are
ill-conditioned
In numerical analysis, the condition number of a function measures how much the output value of the function can change for a small change in the input argument. This is used to measure how sensitive a function is to changes or errors in the input ...
[Numerical Differentiation of Analytic Functions, B Fornberg – ACM Transactions on Mathematical Software (TOMS), 1981.] and due to cancellation will produce a value of zero if ''h'' is small enough.
[Using Complex Variables to Estimate Derivatives of Real Functions, W. Squire, G. Trapp – SIAM REVIEW, 1998.] If too large, the calculation of the slope of the secant line will be more accurately calculated, but the estimate of the slope of the tangent by using the secant could be worse.
For basic central differences, the optimal step is the cube-root of
machine epsilon.
For the numerical derivative formula evaluated at ''x'' and ''x'' + ''h'', a choice for ''h'' that is small without producing a large rounding error is
(though not when ''x'' = 0), where the
machine epsilon ''ε'' is typically of the order of 2.2 for
double precision. A formula for ''h'' that balances the rounding error against the secant error for optimum accuracy is
:
(though not when
), and to employ it will require knowledge of the function.
For computer calculations the problems are exacerbated because, although ''x'' necessarily holds a
representable floating-point number in some precision (32 or 64-bit, ''etc''.), ''x'' + ''h'' almost certainly will not be exactly representable in that precision. This means that ''x'' + ''h'' will be changed (by rounding or truncation) to a nearby machine-representable number, with the consequence that (''x'' + ''h'') − ''x'' will ''not'' equal ''h''; the two function evaluations will not be exactly ''h'' apart. In this regard, since most decimal fractions are recurring sequences in binary (just as 1/3 is in decimal) a seemingly round step such as ''h'' = 0.1 will not be a round number in binary; it is 0.000110011001100...
2 A possible approach is as follows:
h := sqrt(eps) * x;
xph := x + h;
dx := xph - x;
slope := (F(xph) - F(x)) / dx;
However, with computers,
compiler optimization facilities may fail to attend to the details of actual computer arithmetic and instead apply the axioms of mathematics to deduce that ''dx'' and ''h'' are the same. With
C and similar languages, a directive that ''xph'' is a
volatile variable will prevent this.
Other methods
Higher-order methods
Higher-order methods for approximating the derivative, as well as methods for higher derivatives, exist.
Given below is the five-point method for the first derivative (
five-point stencil in one dimension):
:
where