Stable special function recurrences and certified evaluation

For $y(z)={\sum }_{k\ge 0}{c}_k{z}^k$, set $c_2=0$ and $c_{k+3}=c_k\mathrm{/}((k+3)(k+2))$ for $k\ge 0$. Define $\mathrm{Ai}$ by this series with

$c_0=\frac{1}{3^{2\mathrm{/}3}\mathrm{\Gamma }(2\mathrm{/}3)},c_1=-\frac{1}{3^{1\mathrm{/}3}\mathrm{\Gamma }(1\mathrm{/}3)},$

and define $\mathrm{Bi}$ by the same recurrence with

$c_0=\frac{1}{3^{1\mathrm{/}6}\mathrm{\Gamma }(2\mathrm{/}3)},c_1=\frac{3^{1\mathrm{/}6}}{\mathrm{\Gamma }(1\mathrm{/}3)}\mathrm{.}$

Here $\mathrm{\Gamma }$ is Euler's integral on the right half-plane from Euler's Gamma function on the right half-plane. The standard Airy initial values are those in DLMF §9.2(ii); in particular, the power of $3$ in $\mathrm{Bi}(0)$ is in the denominator.

A certified enclosure for a real value is a nonempty closed interval $B=[a,b]$ with finite dyadic endpoints; for a complex value it is a nonempty closed rectangle $B=[a_1,b_1]+i[a_2,b_2]$ with finite dyadic endpoints. Its radius is half its length in the real case and half its Euclidean diagonal in the complex case; $d(B,0)$ is the Euclidean distance from $B$ to $0$. An evaluator item declares whether its input is an exact dyadic point or a dyadic interval/rectangle enclosure; a successful range result contains the exact function image of every admissible point in an input enclosure. Constants not exactly dyadic must be supplied with certified enclosing intervals or rectangles.

The interval backend uses the endpoint operations in Outward-rounded interval, rectangle and ball operations; the soundness of finite arithmetic expression DAGs follows from Outward-rounded operations are inclusion-isotone and enclose DAG evaluations. Real elementary function ranges use Certified real elementary-function evaluations enclose their ranges, and complex exponential and principal-logarithm ranges use Certified complex exponential and principal logarithm enclose their ranges. An operation that returns a non-finite status supplies no finite enclosure claim.

A finite-budget evaluator receives a branch input enclosure wholly within the stated domain, absolute tolerance ${\epsilon }_a>0$, optional relative tolerance ${\epsilon }_r>0$, finite integer work and precision caps, and the fixed format and truncation parameters required by the backend. Its result has one of four statuses:

• success: the returned enclosure contains the exact requested value and its radius is at most ${\epsilon }_a$. If $0\mathrm{\notin }B$ and relative tolerance was requested, the evaluator may also report relative_certified: true only when $\mathrm{rad}(B)\mathrm{/}d(B,0)\le {\epsilon }_r$.
• inconclusive: the procedure exhausted a finite cap before its requested enclosure test passed.
• domain: the input is outside the stated branch domain, violates a stated parameter condition, or lies on a branch cut excluded by that branch.
• range: a certified finite enclosure cannot be represented in the requested output format; the procedure returns no unsound finite interval.

If a success enclosure contains $0$, it still supplies the stated absolute bound, but it carries no relative-accuracy claim. The work and precision caps do not imply that success occurs for every input.

For $z\in \mathbb{C}$, define the Fresnel functions and error function by

$C(z)={\sum }_{k=0}^{\mathrm{\infty }}\frac{(-1{)}^k(\pi \mathrm{/}2{)}^{2k}{z}^{4k+1}}{(2k)!(4k+1)},S(z)={\sum }_{k=0}^{\mathrm{\infty }}\frac{(-1{)}^k(\pi \mathrm{/}2{)}^{2k+1}{z}^{4k+3}}{(2k+1)!(4k+3)},\mathrm{erf}(z)=\frac{2}{\sqrt{\pi }}{\sum }_{k=0}^{\mathrm{\infty }}\frac{(-1{)}^k{z}^{2k+1}}{k!(2k+1)},\mathrm{erfc}(z)=1-\mathrm{erf}(z)\mathrm{.}$

The following local lemmas prove convergence and show that $C,S,\mathrm{erf}$ equal the standard entire integral functions with their usual normalizations.

For $a,b,c,z\in \mathbb{C}$ with $c\mathrm{\notin }\{0,-1,-2,\dots \}$ and $\mathrm{\mid }z\mathrm{\mid }<1$, define

$(a{)}_0=1,(a{)}_k={\prod }_{j=0}^{k-1}(a+j),{}_2{F}_1(a,b;c;z)={\sum }_{k=0}^{\mathrm{\infty }}\frac{(a{)}_k(b{)}_k}{(c{)}_{k}k!}{z}^k\mathrm{.}$

Only this Gauss-series branch of ${}_2{F}_1$ is defined on this page; no value outside $\mathrm{\mid }z\mathrm{\mid }<1$ or across a continuation cut is implicit.

For an integer $n\ge 0$ and $z\in \mathbb{C}$, define

$J_n(z)={\sum }_{k=0}^{\mathrm{\infty }}\frac{(-1{)}^k(z\mathrm{/}2{)}^{2k+n}}{k!(k+n)!},I_n(z)={\sum }_{k=0}^{\mathrm{\infty }}\frac{(z\mathrm{/}2{)}^{2k+n}}{k!(k+n)!}\mathrm{.}$

These are the integer-order Bessel and modified Bessel functions used on this page. For real $x\ge 0$, the scaled modified Bessel value means $I_n^{\mathrm{s}\mathrm{c}}(x)=e^{-x}{I}_n(x)$. This notation does not assert a finite evaluation method or a truncation bound; those are proved separately.

On $\mathrm{Re}z>0$, define

${\mathrm{\ell }}_{+}(z)=-\mathrm{Log}z-{\gamma }_{E}z+{\sum }_{n=1}^{\mathrm{\infty }}(\frac{z}{n}-\mathrm{Log}(1+\frac{z}{n}))\mathrm{.}$

Here $\mathrm{Log}$ is the principal logarithm from Complex logarithms, the principal logarithm, and principal and multivalued complex powers, and ${\gamma }_E$ denotes Euler's constant $\gamma$ from The Euler-Mascheroni constant. Both $z$ and $1+z\mathrm{/}n$ lie in the principal-logarithm domain. The local uniform convergence of this series and the identity $e^{{\mathrm{\ell }}_{+}(z)}=\mathrm{\Gamma }(z)$ are established by Weierstrass log-Gamma series converges to a Gamma branch ↗.

For each integer $M\ge 1$, put $H_M={\sum }_{k=1}^{M}1\mathrm{/}k$. Euler's constant $\gamma$ satisfies $0\le H_M-\log M-\gamma \le \frac{1}{M}\mathrm{.}$ Thus $[H_M-\log M-1\mathrm{/}M,H_M-\log M]$ is a finite enclosure for $\gamma$.

For every $z\in \mathbb{C}\setminus \{0,-1,-2,\dots \}$, $\mathrm{\Gamma }(z+1)=z\mathrm{\Gamma }(z),$ where $\mathrm{\Gamma }$ is the meromorphic continuation of Euler's Gamma function.

For real $z>0$, let $D_z:=z^{-1}\frac{d}{dz}$ and define, for integers $\mathrm{\ell }\ge 0$, $j_{\mathrm{\ell }}(z):=(-z{)}^{\mathrm{\ell }}{D}_z^{\mathrm{\ell }}\left(\frac{\sin z}{z}\right),h_{\mathrm{\ell }}^{(1)}(z):=(-z{)}^{\mathrm{\ell }}{D}_z^{\mathrm{\ell }}(-i\frac{e^{iz}}{z})\mathrm{.}$ The definition fixes the outgoing $e^{iz}$ convention. The functions solve the radial Helmholtz equation and have the recurrence, origin, Wronskian, and infinity behavior by Spherical Bessel and outgoing Hankel recurrences solve the radial Helmholtz equation ↗.

These symbolic definitions assert no stable or certified numerical evaluation algorithm for large $\mathrm{\ell }$ or $z$.

Fix real $x>0$. Let $D_0(x)=-J_1(x)$ and $D_1(x)=J_0(x)$, and for every integer $j\ge 1$ define

$D_{j+1}(x)=\frac{2j}{x}{D}_j(x)-D_{j-1}(x)\mathrm{.}$

Thus $D$ is the second sequence obtained from the same order recurrence as $J_j(x)$. It is a local auxiliary solution, not the Bessel function of the second kind.

For real $x>0$, define $\mathrm{erfcx}(x)=e^{x^2}\mathrm{erfc}(x)$. For each integer $k\ge 0$, define the finite positive moment

$M_k(x)={\int }_0^{\mathrm{\infty }}{t}^k{e}^{-t^2-2xt}dt\mathrm{.}$

The moments are auxiliary quantities for a finite-depth interval enclosure of $\mathrm{erfcx}$. They are not the modified Bessel functions $I_k$.

Let $\mathrm{\Omega }=\mathbb{C}\setminus (-\mathrm{\infty },0]$. For every $z\in \mathrm{\Omega }$, an admissible integer $m\ge 0$ with $\mathrm{Re}(z+m)>0$ exists by the Archimedean property Dedekind completeness implies the Archimedean property. Choose any such $m$ and set

$\mathrm{\ell }(z)={\mathrm{\ell }}_{+}(z+m)-{\sum }_{j=0}^{m-1}\mathrm{Log}(z+j),$

where the sum is empty if $m=0$ and ${\mathrm{\ell }}_{+}$ and $\mathrm{Log}$ are defined in Log-Gamma branch on the right half-plane. The next lemma proves that this value is independent of the admissible $m$ and satisfies $e^{\mathrm{\ell }(z)}=\mathrm{\Gamma }(z)$. This is a continued branch of log-Gamma, not a separate principal $\mathrm{Log}(\mathrm{\Gamma }(z))$.

The coefficient recurrence in Airy initial-value series makes both series solve $y^{\mathrm{\prime }\mathrm{\prime }}=zy$. In a residue class $k\equiv r(\mathrm{m}\mathrm{o}\mathrm{d}3)$, successive nonzero terms have absolute ratio

$\frac{\mathrm{\mid }z{\mathrm{\mid }}^3}{(k+3)(k+2)}\mathrm{.}$

If the ratio of each term after a chosen first omitted term is bounded by $q<1$, that class's tail is at most the first omitted term's absolute value divided by $1-q$. The series and their first two derivative series converge locally uniformly on $\mathbb{C}$.

For the $k$th term $e_k$ of the $\mathrm{erf}$ series in Fresnel, error, and Gauss hypergeometric series,

$\frac{\mathrm{\mid }{e}_{k+1}\mathrm{\mid }}{\mathrm{\mid }{e}_k\mathrm{\mid }}=\frac{\mathrm{\mid }z{\mathrm{\mid }}^2}{k+1}\frac{2k+1}{2k+3}\le \frac{\mathrm{\mid }z{\mathrm{\mid }}^2}{k+1}\mathrm{.}$

The displayed quotient is asserted when $e_k\mathrm{\ne }0$. At $z=0$, every $e_k=0$, so the series and every tail are zero; the corresponding division-free inequality $\mathrm{\mid }{e}_{k+1}\mathrm{\mid }\le \mathrm{\mid }z{\mathrm{\mid }}^2\mathrm{\mid }{e}_k\mathrm{\mid }\mathrm{/}(k+1)$ holds there as well.

Thus, if $\mathrm{\mid }z{\mathrm{\mid }}^2\mathrm{/}(K+2)=q<1$, the tail after terms $0,\dots ,K$ is bounded by $\mathrm{\mid }{e}_{K+1}\mathrm{\mid }\mathrm{/}(1-q)$. The series converges locally uniformly. For every piecewise-$C^1$ path $\gamma :[0,1]\to \mathbb{C}$ with $\gamma (0)=0$ and $\gamma (1)=z$, it equals $\frac{2}{\sqrt{\pi }}{\int }_{\gamma }{e}^{-t^2}dt,$ where the complex line integral is defined by The complex line integral over a rectifiable path as a componentwise Riemann–Stieltjes integral. In particular, the straight-segment integral has this value, and the integral is path-independent.

For the terms $c_k$ of $C$ and $s_k$ of $S$ in Fresnel, error, and Gauss hypergeometric series,

$\frac{\mathrm{\mid }{c}_{k+1}\mathrm{\mid }}{\mathrm{\mid }{c}_k\mathrm{\mid }}\le \frac{{\pi }^2\mathrm{\mid }z{\mathrm{\mid }}^4\mathrm{/}4}{(2k+1)(2k+2)},\frac{\mathrm{\mid }{s}_{k+1}\mathrm{\mid }}{\mathrm{\mid }{s}_k\mathrm{\mid }}\le \frac{{\pi }^2\mathrm{\mid }z{\mathrm{\mid }}^4\mathrm{/}4}{(2k+2)(2k+3)}\mathrm{.}$

Each displayed quotient is asserted only when its denominator is nonzero. At $z=0$, every $c_k=s_k=0$, so both series and every tail are zero. The corresponding division-free inequalities obtained by multiplying by $\mathrm{\mid }{c}_k\mathrm{\mid }$ or $\mathrm{\mid }{s}_k\mathrm{\mid }$ hold for every $z$.

Once the applicable bound is at most $q<1$ for all terms after a first omitted term, the absolute tail is bounded by that term divided by $1-q$. Both series converge locally uniformly on $\mathbb{C}$. For every piecewise-$C^1$ path $\gamma :[0,1]\to \mathbb{C}$ with $\gamma (0)=0$ and $\gamma (1)=z$, they equal $C(z)={\int }_{\gamma }\cos (\pi t^2\mathrm{/}2)dt,S(z)={\int }_{\gamma }\sin (\pi t^2\mathrm{/}2)dt,$ where the complex line integral is defined by The complex line integral over a rectifiable path as a componentwise Riemann–Stieltjes integral. In particular, these integrals are path-independent and agree with the straight-segment integrals.

Let $a,b,c,z\in \mathbb{C}$, $\mathrm{\mid }z\mathrm{\mid }<1$, and $c\mathrm{\notin }\{0,-1,-2,\dots \}$. If $z=0$, then $h_0=1$ and $h_k=0$ for $k\ge 1$, so the series equals $1$ and every tail after $h_0$ is zero. Suppose instead that $z\mathrm{\ne }0$ and neither $a$ nor $b$ is a nonpositive integer. Then every term $h_k=(a{)}_k(b{)}_k{z}^k\mathrm{/}((c{)}_{k}k!)$ is nonzero, and for any integer $N>\mathrm{\mid }c\mathrm{\mid }$,

$\frac{\mathrm{\mid }{h}_{k+1}\mathrm{\mid }}{\mathrm{\mid }{h}_k\mathrm{\mid }}=\mathrm{\mid }z\mathrm{\mid }\frac{\mathrm{\mid }a+k\mathrm{\mid }\mathrm{\mid }b+k\mathrm{\mid }}{\mathrm{\mid }c+k\mathrm{\mid }(k+1)}\le q_N:=\mathrm{\mid }z\mathrm{\mid }\frac{(1+\mathrm{\mid }a\mathrm{\mid }\mathrm{/}N)(1+\mathrm{\mid }b\mathrm{\mid }\mathrm{/}N)}{1-\mathrm{\mid }c\mathrm{\mid }\mathrm{/}N}(k\ge N)\mathrm{.}$

For sufficiently large $N$, $q_N<1$; then the tail starting at $h_N$ is at most $\mathrm{\mid }{h}_N\mathrm{\mid }\mathrm{/}(1-q_N)$. If $a$ or $b$ is a nonpositive integer, the series terminates and its tail after the last nonzero term is zero.

For $n\ge 0$, $z\mathrm{\ne }0$, and either series in Integer-order Bessel series, the absolute term ratio is

$r_k=\frac{\mathrm{\mid }z{\mathrm{\mid }}^2}{4(k+1)(k+n+1)}\mathrm{.}$

At $z=0$, the order-zero series have only the constant term $a_0=1$ and positive orders have only zero terms; in every case the tail after terms $k=0,\dots ,K$ is zero. The series for $J_n$ and $I_n$ converge absolutely and locally uniformly in $z\in \mathbb{C}$. If $K\ge 0$ and $q=\mathrm{\mid }z{\mathrm{\mid }}^2\mathrm{/}[4(K+2)(K+n+2)]<1$, then the tail after terms $k=0,\dots ,K$ has absolute value at most $\mathrm{\mid }{a}_{K+1}\mathrm{\mid }\mathrm{/}(1-q)$, where $a_k$ is the corresponding series term. The first two formal termwise derivative series also converge locally uniformly.

For every $z$ with $\mathrm{Re}z>0$, the series defining ${\mathrm{\ell }}_{+}(z)$ in Log-Gamma branch on the right half-plane converges locally uniformly, defines a holomorphic function, and satisfies $e^{{\mathrm{\ell }}_{+}(z)}=\mathrm{\Gamma }(z)$ for Euler's Gamma function on that half-plane.