The generalized hypergeometric function is formally defined by
It can be interpreted using analytic continuation or regularization when the sum does not converge. In a looser sense, we understand “hypergeometric functions” to be linear combinations of generalized hypergeometric functions with prefactors that are products of exponentials, powers, and gamma functions.
In this section, we define
and
For the conventional generalized hypergeometric function \({}_pF_{q}\), compute \({}_pH_{q+1}\) with the explicit parameter \(b_q = 1\), or remove a 1 from the \(a_i\) parameters if there is one.
Computes a factor C such that
We check that \(\operatorname{Re}(b+n) > 0\) for all lower parameters b. If this does not hold, C is set to infinity. Otherwise, we cancel out pairs of parameters \(a\) and \(b\) against each other. We have
and
for all \(k \ge n\). This gives us a constant D such that \(T(k+1) \le D T(k)\) for all \(k \ge n\). If \(D \ge 1\), we set C to infinity. Otherwise, we take \(C = \sum_{k=0}^{\infty} D^k = (1-D)^{-1}\).
As currently implemented, the bound becomes infinite when \(n\) is too small, even if the series converges.
Heuristically attempts to choose a number of terms n to sum of a hypergeometric series at a working precision of prec bits.
Uses double precision arithmetic internally. As currently implemented, it can fail to produce a good result if the parameters are extremely large or extremely close to nonpositive integers.
Numerical cancellation is assumed to be significant, so truncation is done when the current term is prec bits smaller than the largest encountered term.
This function will also attempt to pick a reasonable truncation point for divergent series.
Computes \(s = \sum_{k=0}^{n-1} T(k)\) and \(t = T(n)\). Does not allow aliasing between input and output variables. We require \(n \ge 0\).
The forward version computes the sum using forward recurrence.
The bs version computes the sum using binary splitting.
The rs version computes the sum in reverse order using rectangular splitting. It only computes a magnitude bound for the value of t.
The fme version uses fast multipoint evaluation.
The default version automatically chooses an algorithm depending on the inputs.
Like acb_hypgeom_pfq_sum(), but taking advantage of \(w = 1/z\) possibly having few bits.
Computes
directly from the defining series, including a rigorous bound for the truncation error \(\varepsilon\) in the output.
If \(n < 0\), this function chooses a number of terms automatically using acb_hypgeom_pfq_choose_n().
Computes \({}_pH_{q}(z)\) directly using the defining series, given parameters and argument that are power series. The result is a power series of length len.
An error bound is computed automatically as a function of the number of terms n. If \(n < 0\), the number of terms is chosen automatically.
If regularized is set, the regularized hypergeometric function is computed instead.
Let \(U(a,b,z)\) denote the confluent hypergeometric function of the second kind with the principal branch cut, and let \(U^{*} = z^a U(a,b,z)\). For all \(z \ne 0\) and \(b \notin \mathbb{Z}\) (but valid for all \(b\) as a limit), we have (DLMF 13.2.42)
Moreover, for all \(z \ne 0\) we have
which is equivalent to DLMF 13.2.41 (but simpler in form).
We have the asymptotic expansion
where \({}_2F_0(a,b,z)\) denotes a formal hypergeometric series, i.e.
The error term \(\varepsilon_n(z)\) is bounded according to DLMF 13.7. A case distinction is made depending on whether \(z\) lies in one of three regions which we index by \(R\). Our formula for the error bound increases with the value of \(R\), so we can always choose the larger out of two indices if \(z\) lies in the union of two regions.
Let \(r = |b-2a|\). If \(\operatorname{Re}(z) \ge r\), set \(R = 1\). Otherwise, if \(\operatorname{Im}(z) \ge r\) or \(\operatorname{Re}(z) \ge 0 \land |z| \ge r\), set \(R = 2\). Otherwise, if \(|z| \ge 2r\), set \(R = 3\). Otherwise, the bound is infinite. If the bound is finite, we have
in terms of the following auxiliary quantities
Computes \(U(a,b,z)\) as a power series truncated to length len, given \(a, b, z \in \mathbb{C}[[x]]\). If \(b[0] \in \mathbb{Z}\), it computes one extra derivative and removes the singularity (it is then assumed that \(b[1] \ne 0\)). As currently implemented, the output is indeterminate if \(b\) is nonexact and contains an integer.
Computes \(U(a,b,z)\) as a sum of two convergent hypergeometric series. If \(b \in \mathbb{Z}\), it computes the limit value via acb_hypgeom_u_1f1_series(). As currently implemented, the output is indeterminate if \(b\) is nonexact and contains an integer.
Computes \(U(a,b,z)\) using an automatic algorithm choice. The function acb_hypgeom_u_asymp() is used if \(a\) or \(a-b+1\) is a nonpositive integer (in which case the asymptotic series terminates), or if z is sufficiently large. Otherwise acb_hypgeom_u_1f1() is used.
Computes the error function respectively using
and an automatic algorithm choice. The asymp version takes a second precision to use for the U term.
Computes the complementary error function \(\operatorname{erfc}(z) = 1 - \operatorname{erf}(z)\). This function avoids catastrophic cancellation for large positive z.
Computes the imaginary error function \(\operatorname{erfi}(z) = -i\operatorname{erf}(iz)\). This is a trivial wrapper of acb_hypgeom_erf().
Computes the Bessel function of the first kind via acb_hypgeom_u_asymp(). For all complex \(\nu, z\), we have
where
Nicer representations of the factors \(A_{\pm}\) can be given depending conditionally on the parameters. If \(\nu + \tfrac{1}{2} = n \in \mathbb{Z}\), we have \(A_{\pm} = (\pm i)^{n} (2 \pi z)^{-1/2}\). And if \(\operatorname{Re}(z) > 0\), we have \(A_{\pm} = \exp(\mp i [(2\nu+1)/4] \pi) (2 \pi z)^{-1/2}\).
Computes the Bessel function of the first kind from
Computes the Bessel function of the first kind \(J_{\nu}(z)\) using an automatic algorithm choice.
Computes the Bessel function of the second kind \(Y_{\nu}(z)\) from the formula
unless \(\nu = n\) is an integer in which case the limit value
is computed. As currently implemented, the output is indeterminate if \(\nu\) is nonexact and contains an integer.
Computes the modified Bessel function of the first kind \(I_{\nu}(z) = z^{\nu} (iz)^{-\nu} J_{\nu}(iz)\) respectively using asymptotic series (see acb_hypgeom_bessel_j_asymp()), the convergent series
or an automatic algorithm choice.
Computes the modified Bessel function of the second kind via via acb_hypgeom_u_asymp(). For all \(\nu\) and all \(z \ne 0\), we have
Computes the modified Bessel function of the second kind \(K_{\nu}(z)\) as a power series truncated to length len, given \(\nu, z \in \mathbb{C}[[x]]\). Uses the formula
If \(\nu[0] \in \mathbb{Z}\), it computes one extra derivative and removes the singularity (it is then assumed that \(\nu[1] \ne 0\)). As currently implemented, the output is indeterminate if \(\nu[0]\) is nonexact and contains an integer.
Computes the modified Bessel function of the second kind from
if \(\nu \notin \mathbb{Z}\). If \(\nu \in \mathbb{Z}\), it computes the limit value via acb_hypgeom_bessel_k_0f1_series(). As currently implemented, the output is indeterminate if \(\nu\) is nonexact and contains an integer.
Computes the upper incomplete gamma function respectively using
and an automatic algorithm choice. The automatic version also handles other special input such as \(z = 0\) and \(s = 1, 2, 3\). The singular version evaluates the finite sum directly and therefore assumes that s is not too large. If modified is set, computes the exponential integral \(z^{-s} \Gamma(s,z) = E_{1-s}(z)\) instead.
The branch cut conventions of the following functions match Mathematica.
Computes the generalized exponential integral \(E_s(z)\). This is a trivial wrapper of acb_hypgeom_gamma_upper().
Computes the exponential integral \(\operatorname{Ei}(z)\), respectively using
and an automatic algorithm choice.
Computes the sine integral \(\operatorname{Si}(z)\), respectively using
and an automatic algorithm choice.
Computes the cosine integral \(\operatorname{Ci}(z)\), respectively using
and an automatic algorithm choice.
Computes the hyperbolic sine integral \(\operatorname{Shi}(z) = -i \operatorname{Si}(iz)\). This is a trivial wrapper of acb_hypgeom_si().
Computes the hyperbolic cosine integral \(\operatorname{Chi}(z)\), respectively using
and an automatic algorithm choice.
If offset is zero, computes the logarithmic integral \(\operatorname{li}(z) = \operatorname{Ei}(\log(z))\).
If offset is nonzero, computes the offset logarithmic integral \(\operatorname{Li}(z) = \operatorname{li}(z) - \operatorname{li}(2)\).