Numerical differentiation and derivative error control

Let $U\subseteq \mathbb{C}$ be open, $x\in U\cap \mathbb{R}$, and $h>0$ with $x+ih\in U$. For a holomorphic $f:U\to \mathbb{C}$ that is real-valued on some real interval $I\subseteq U\cap \mathbb{R}$ containing $x$, the complex-step estimate is

$C_{h}f(x):=\frac{\mathrm{Im}f(x+ih)}{h}\mathrm{.}$

For an inexact complex evaluation $\widehat{f}(x+ih)$, define its imaginary-component error by ${\delta }_{\mathrm{i}\mathrm{m}}:=\mathrm{\mid }\mathrm{Im}\widehat{f}(x+ih)-\mathrm{Im}f(x+ih)\mathrm{\mid }$. The inexact estimate is ${\widehat{C}}_{h}f(x):=\mathrm{Im}\widehat{f}(x+ih)\mathrm{/}h$; its evaluation-error contribution is controlled by ${\delta }_{\mathrm{i}\mathrm{m}}\mathrm{/}h$. No claim is made that arbitrary code representing a real function has a holomorphic extension.

Let $I\subseteq \mathbb{R}$ be open, $x\in I$, and $h>0$. When $x+h\in I$, the forward difference of $f:I\to \mathbb{R}$ is

$D_h^{+}f(x):=\frac{f(x+h)-f(x)}{h}\mathrm{.}$

When $x-h,x+h\in I$, the central difference is

$D_h^{0}f(x):=\frac{f(x+h)-f(x-h)}{2h}\mathrm{.}$

An evaluation with absolute error bound $\delta \ge 0$ is a value $\widehat{f}(y)=f(y)+e_y$ with $\mathrm{\mid }{e}_y\mathrm{\mid }\le \delta$ at each requested point $y$. The corresponding inexact quotients are obtained by replacing each $f(y)$ in the displayed stencil by $\widehat{f}(y)$. The bound is absolute and applies separately to each call; no correlation or relative-error assumption is included.

Let $m,n\ge 1$, $U\subseteq {\mathbb{R}}^n$ be open, $x\in U$, $F:U\to {\mathbb{R}}^m$, $f:U\to \mathbb{R}$, and let $v\in {\mathbb{R}}^n$. Use the Jacobian $JF(x)$ and gradient $\mathrm{\nabla }f(x)$ of The Jacobian matrix of partial derivatives and the gradient in the scalar-valued case, the Hessian $H_f(x)$ of The Hessian matrix and critical points of a scalar field, and the published Euclidean vector norm $\mathrm{\parallel }\cdot {\mathrm{\parallel }}_2$ of The $p$-norms $\parallel x{\parallel }_p$ for rational $p\ge 1$, and $\parallel x{\parallel }_{\mathrm{\infty }}$. The Jacobian-vector products below require $JF(x)$ to exist, and the Hessian-vector products require $H_f(x)$ to exist. Each gradient-difference Hessian-vector estimate also requires $\mathrm{\nabla }f$ to exist at every point in its stencil. The exact products are $JF(x)v$ and $H_f(x)v$, respectively.

For $h>0$ with the required points in $U$, define the forward and central Jacobian-vector estimates

$J_h^{+}(F,x,v):=\frac{F(x+hv)-F(x)}{h},J_h^0(F,x,v):=\frac{F(x+hv)-F(x-hv)}{2h},$

and the gradient-difference Hessian-vector estimates

$H_h^{+}(f,x,v):=\frac{\mathrm{\nabla }f(x+hv)-\mathrm{\nabla }f(x)}{h},H_h^0(f,x,v):=\frac{\mathrm{\nabla }f(x+hv)-\mathrm{\nabla }f(x-hv)}{2h}\mathrm{.}$

For each requested vector value $z(y)$, write an inexact evaluation as $\widehat{z}(y)=z(y)+e$. When an error bound is specified, it is per call: $\mathrm{\parallel }e{\mathrm{\parallel }}_2\le \delta$ for that call's $\delta \ge 0$. Errors may be dependent or random; a stated bound on a random error is understood to hold almost surely. Define the inexact estimates by substitution:

${\widehat{J}}_h^{+}(F,x,v):=\frac{\widehat{F}(x+hv)-\widehat{F}(x)}{h},{\widehat{J}}_h^0(F,x,v):=\frac{\widehat{F}(x+hv)-\widehat{F}(x-hv)}{2h},$

${\widehat{H}}_h^{+}(f,x,v):=\frac{\widehat{\mathrm{\nabla }f}(x+hv)-\widehat{\mathrm{\nabla }f}(x)}{h},{\widehat{H}}_h^0(f,x,v):=\frac{\widehat{\mathrm{\nabla }f}(x+hv)-\widehat{\mathrm{\nabla }f}(x-hv)}{2h}\mathrm{.}$

Let $m,n\ge 1$, $U\subseteq {\mathbb{R}}^n$ be open, $F:U\to {\mathbb{R}}^m$ be differentiable at $x\in U$, and $J:=JF(x)$ be its Jacobian as defined in The Jacobian matrix of partial derivatives and the gradient in the scalar-valued case. A structural row pattern is a list of sets $S_i\subseteq \{0,\dots ,n-1\}$, one for each row $0\le i<m$, such that $J_{ij}=0$ whenever $j\mathrm{\notin }{S}_i$. It is a safe pattern only when every actual nonzero lies in the declared support.

The column-interference graph has vertices $0,\dots ,n-1$ and an edge between distinct columns $j,k$ exactly when some row support $S_i$ contains both. A proper coloring is a map $\chi :\{0,\dots ,n-1\}\to C$ to a set of color labels such that adjacent columns have different labels. For a color $c\in C$, its seed vector $s^c\in {\mathbb{R}}^n$ has coordinate $s_j^c=1$ when $\chi (j)=c$ and $0$ otherwise. The compressed product is $J{s}^c$.

Let $I\subseteq \mathbb{R}$ be open, let $x\in I$, and let $h>0$ with $[x-h,x+h]\subset I$. Suppose $f\in C^5$ on a neighborhood of this interval and $\mathrm{\mid }{f}^{(5)}(y)\mathrm{\mid }\le M_5$ throughout it. Then for some real $r_h$,

$D_h^{0}f(x)=f^{\mathrm{\prime }}(x)+c_2{h}^2+r_h,c_2:=\frac{f^{\mathrm{\prime }\mathrm{\prime }\mathrm{\prime }}(x)}{6},\mathrm{\mid }{r}_h\mathrm{\mid }\le \frac{M_5{h}^4}{120}\mathrm{.}$

The same coefficient $c_2$ works for every admissible $h$ at the fixed point $x$.

Let $U\subseteq \mathbb{C}$ be open, let $x\in U\cap \mathbb{R}$, and suppose $f:U\to \mathbb{C}$ is holomorphic, real-valued on a real interval $I\subseteq U\cap \mathbb{R}$ around $x$, and $\{x+it:0\le t\le h\}\subset U$ for $h>0$. Suppose $\mathrm{\mid }\mathrm{Re}{f}^{\mathrm{\prime }\mathrm{\prime }\mathrm{\prime }}(x+it)\mathrm{\mid }\le M_3$ for $0\le t\le h$, with $M_3\ge 0$. Then

$\mathrm{\mid }{C}_{h}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{M_3{h}^2}{6}\mathrm{.}$

If the evaluated imaginary component has error at most ${\delta }_{\mathrm{i}\mathrm{m}}$, then

$\mathrm{\mid }{\widehat{C}}_{h}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{M_3{h}^2}{6}+\frac{{\delta }_{\mathrm{i}\mathrm{m}}}{h}\mathrm{.}$

Holomorphy and real-axis reality are essential to the formula: real differentiability of the restriction alone is insufficient.

Let $I\subseteq \mathbb{R}$ be open, $x\in I$, and $h>0$. For the forward bound assume $[x,x+h]\subset I$, $f\in C^2$ on a neighborhood of that interval, and $\mathrm{\mid }{f}^{\mathrm{\prime }\mathrm{\prime }}(y)\mathrm{\mid }\le M_2$ there. For the central bound assume $[x-h,x+h]\subset I$, $f\in C^3$ on a neighborhood of that interval, and $\mathrm{\mid }{f}^{\mathrm{\prime }\mathrm{\prime }\mathrm{\prime }}(y)\mathrm{\mid }\le M_3$ there. Then

$\mathrm{\mid }{D}_h^{+}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{M_{2}h}{2},\mathrm{\mid }{D}_h^{0}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{M_3{h}^2}{6}\mathrm{.}$

If each function evaluation has absolute error at most $\delta \ge 0$ as in Forward and central difference stencils with bounded evaluation errors, then

$\mathrm{\mid }{\widehat{D}}_h^{+}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{M_{2}h}{2}+\frac{2\delta }{h},\mathrm{\mid }{\widehat{D}}_h^{0}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{M_3{h}^2}{6}+\frac{2\delta }{h}\mathrm{.}$

For the central stencil the evaluation contribution is in fact at most $\delta \mathrm{/}h$; the displayed $2\delta \mathrm{/}h$ is a uniform conservative bound for both stencils.

Let $m,n\ge 1$, $U\subseteq {\mathbb{R}}^n$ be open, $x\in U$, $v\in {\mathbb{R}}^n$, and $h>0$ with the required line segments contained in $U$. Use the estimators of Forward and central finite-difference Jacobian-vector and Hessian-vector products. For a vector-valued map, $F\in C^k(U)$ means that each component $F_i$ is in $C^k(U)$. For $F\in C^2$, suppose each component $g_i(t):=F_i(x+tv)$ satisfies $\mathrm{\mid }{g}_i^{\mathrm{\prime }\mathrm{\prime }}(t)\mathrm{\mid }\le b_{2,i}$ for $0\le t\le h$, with $b_{2,i}\ge 0$; put $M_{2,F}:=({\sum }_{i<m}{b}_{2,i}^2{)}^{1\mathrm{/}2}$. For $F\in C^3$, suppose $\mathrm{\mid }{g}_i^{\mathrm{\prime }\mathrm{\prime }\mathrm{\prime }}(t)\mathrm{\mid }\le b_{3,i}$ for $-h\le t\le h$, with $b_{3,i}\ge 0$; put $M_{3,F}:=({\sum }_{i<m}{b}_{3,i}^2{)}^{1\mathrm{/}2}$. Then

$\mathrm{\parallel }{J}_h^{+}(F,x,v)-JF(x)v{\mathrm{\parallel }}_2\le M_{2,F}h\mathrm{/}2,\mathrm{\parallel }{J}_h^0(F,x,v)-JF(x)v{\mathrm{\parallel }}_2\le M_{3,F}{h}^2\mathrm{/}6.$

For $f\in C^3$, set $q_j(t):={\mathrm{\partial }}_{j}f(x+tv)$, suppose $\mathrm{\mid }{q}_j^{\mathrm{\prime }\mathrm{\prime }}(t)\mathrm{\mid }\le a_{2,j}$ on $[0,h]$, with $a_{2,j}\ge 0$, and put $M_{2,H}:=({\sum }_{j<n}{a}_{2,j}^2{)}^{1\mathrm{/}2}$. For $f\in C^4$, suppose $\mathrm{\mid }{q}_j^{\mathrm{\prime }\mathrm{\prime }\mathrm{\prime }}(t)\mathrm{\mid }\le a_{3,j}$ on $[-h,h]$, with $a_{3,j}\ge 0$, and put $M_{3,H}:=({\sum }_{j<n}{a}_{3,j}^2{)}^{1\mathrm{/}2}$. Then

$\mathrm{\parallel }{H}_h^{+}(f,x,v)-H_f(x)v{\mathrm{\parallel }}_2\le M_{2,H}h\mathrm{/}2,\mathrm{\parallel }{H}_h^0(f,x,v)-H_f(x)v{\mathrm{\parallel }}_2\le M_{3,H}{h}^2\mathrm{/}6.$

If each vector evaluation has Euclidean error at most $\delta$, the forward evaluation contribution is at most $2\delta \mathrm{/}h$ and the central contribution at most $\delta \mathrm{/}h$ (hence also at most $2\delta \mathrm{/}h$). For a scalar output coordinate on a common probability space, suppose its two square-integrable evaluation noises are centered, with variances ${\sigma }_1^2,{\sigma }_0^2$ and covariance $c$. Their forward-difference noise variance is $({\sigma }_1^2+{\sigma }_0^2-2c)\mathrm{/}{h}^2$; the central formula has the same numerator divided by $4h^2$. If these centered noises are independent and both have variance ${\sigma }^2$, then the formulas reduce to $2{\sigma }^2\mathrm{/}{h}^2$ forward and ${\sigma }^2\mathrm{/}(2h^2)$ central. No variance sum or root-sum-square total-error claim is made without the stated centering and covariance information.

If every evaluation call used in a given estimate has Euclidean error at most a common $\delta \ge 0$, then the inexact estimates of Forward and central finite-difference Jacobian-vector and Hessian-vector products satisfy

$\mathrm{\parallel }{\widehat{J}}_h^{+}(F,x,v)-JF(x)v{\mathrm{\parallel }}_2\le M_{2,F}h\mathrm{/}2+2\delta \mathrm{/}h,\mathrm{\parallel }{\widehat{J}}_h^0(F,x,v)-JF(x)v{\mathrm{\parallel }}_2\le M_{3,F}{h}^2\mathrm{/}6+\delta \mathrm{/}h,$

$\mathrm{\parallel }{\widehat{H}}_h^{+}(f,x,v)-H_f(x)v{\mathrm{\parallel }}_2\le M_{2,H}h\mathrm{/}2+2\delta \mathrm{/}h,\mathrm{\parallel }{\widehat{H}}_h^0(f,x,v)-H_f(x)v{\mathrm{\parallel }}_2\le M_{3,H}{h}^2\mathrm{/}6+\delta \mathrm{/}h\mathrm{.}$

The central evaluation contribution is at most $\delta \mathrm{/}h$, so the conservative $2\delta \mathrm{/}h$ contribution is also valid for the central estimates.

Let the objects and notation be those of The column-interference graph and collision-free coloring for a sparse Jacobian, and suppose $\chi$ is a proper coloring of the graph of a safe structural row pattern. For every row $0\le i<m$, every color $c$, and every column $j\in S_i$ with $\chi (j)=c$,

$J_{ij}=(J{s}^c{)}_i\mathrm{.}$

Entries with $j\mathrm{\notin }{S}_i$ are zero by the safe-pattern hypothesis. Thus all entries are recovered from one compressed product per color. A color collision in a row does not give a general recovery guarantee: two distinct Jacobians supported on that colliding row can have the same compressed product.

For a function $f$ with the values required by the central stencils at $x\pm h$ and $x\pm h\mathrm{/}2$, the Richardson extrapolated central difference is

$R_{h}f(x):=\frac{4D_{h\mathrm{/}2}^{0}f(x)-D_h^{0}f(x)}{3}\mathrm{.}$

Its fourth-order error guarantee is conditional on a uniform expansion with a remainder bounded by a constant times $h^4$; that condition is stated and proved from a bounded fifth derivative in A bounded fifth derivative gives the central-difference expansion required by Richardson extrapolation.

Let $H>0$. Consider the conservative error bounds from Forward and central differences have first- and second-order error bounds on a permitted step interval $(0,H]$, where every stencil with $h\le H$ stays inside the domain. If $M_2>0$ and $\delta >0$, the forward bound

$E_{+}(h)=M_{2}h\mathrm{/}2+2\delta \mathrm{/}h$

has the unique unconstrained minimizer $h_{+}=\sqrt{4\delta \mathrm{/}{M}_2}$ and minimum $2\sqrt{M_2\delta }$. On $(0,H]$ the minimizer is $\min \{h_{+},H\}$. If $M_3>0$ and $\delta >0$, the central conservative bound

$E_0(h)=M_3{h}^2\mathrm{/}6+2\delta \mathrm{/}h$

has the unique unconstrained minimizer $h_0=\sqrt[3]{6\delta \mathrm{/}{M}_3}$ and minimum $M_3{h}_0^2\mathrm{/}2$; on $(0,H]$ the minimizer is $\min \{h_0,H\}$. These optimize the displayed error bounds, not the unknown actual error.

If $\delta =0$ and the corresponding truncation constant is positive, the bound has infimum $0$ as $h\downarrow 0$ and no minimizer for $h>0$. If the corresponding truncation constant is $0$ and $\delta >0$, the bound decreases with $h$, so its capped minimum is at $H$. If $\delta$ and the corresponding truncation constant are both zero, every permitted step is optimal.

Let $I\subseteq \mathbb{R}$ be open, $x\in I$, and $h>0$ with $[x-h,x+h]\subset I$. Suppose that for every $0<t\le h$ the central difference has the actual uniform expansion

$D_t^{0}f(x)=f^{\mathrm{\prime }}(x)+c_2{t}^2+r_t,\mathrm{\mid }{r}_t\mathrm{\mid }\le M{t}^4,$

where $c_2$ and $M\ge 0$ do not depend on $t$. Then the Richardson estimate of Richardson extrapolation for a central difference with a fourth-order remainder satisfies

$\mathrm{\mid }{R}_{h}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{5M{h}^4}{12}\mathrm{.}$

In particular, under A bounded fifth derivative gives the central-difference expansion required by Richardson extrapolation the bound holds with $c_2=f^{\mathrm{\prime }\mathrm{\prime }\mathrm{\prime }}(x)\mathrm{/}6$ and $M=M_5\mathrm{/}120$, giving $\mathrm{\mid }{R}_{h}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le M_5{h}^4\mathrm{/}288$. If each of the four function values at $x\pm h\mathrm{/}2$ and $x\pm h$ is evaluated as $\widehat{f}(y)=f(y)+e_y$ with $\mathrm{\mid }{e}_y\mathrm{\mid }\le \delta$, define ${\widehat{R}}_{h}f(x)$ by replacing those four values in $R_{h}f(x)$ by their evaluations. Then

$\mathrm{\mid }{\widehat{R}}_{h}f(x)-f^{\mathrm{\prime }}(x)\mathrm{\mid }\le \frac{5M{h}^4}{12}+\frac{3\delta }{h}\mathrm{.}$