Computing the Gamma function using contour integrals and rational approximations

Some of the best methods for computing the gamma function are based on numerical evaluation of Hankel's contour integral. For example, Temme evaluates this integral based on steepest-decent contours by the trapezoid rule. Here we investigate a different approach to the integral: the application of the trapezoid rule on Talbot-type contours using optimal parameters recently derived by Weideman for computing inverse Laplace transforms. Relatedly, we also investigate quadrature formulas derived from best approximations to exp(z) on the negative real axis, following Cody, Meinardus and Varga. The two methods are closely related and both converge geometrically. We find that the new methods are competitive with existing ones, even though they are based on generic tools rather than on specific analysis of the gamma function

The exponentially convergent trapezoidal rule

It is well known that the trapezoidal rule converges geometrically when applied to analytic functions on periodic intervals or the real line. The mathematics and history of this phenomenon are reviewed and it is shown that far from being a curiosity, it is linked with computational methods all across scientific computing, including algorithms related to inverse Laplace transforms, special functions, complex analysis, rational approximation, integral equations, and the computation of functions and eigenvalues of matrices and operators

Evaluating matrix functions for exponential integrators via Carathéodory-Fejér approximation and contour integrals

Among the fastest methods for solving stiff PDE are exponential integrators, which require the evaluation of $f(A)$, where $A$ is a negative definite matrix and $f$ is the exponential function or one of the related ``$\varphi$ functions'' such as $\varphi_1(z) = (e^z-1)/z$. Building on previous work by Trefethen and Gutknecht, Gonchar and Rakhmanov, and Lu, we propose two methods for the fast evaluation of $f(A)$ that are especially useful when shifted systems $(A+zI)x=b$ can be solved efficiently, e.g. by a sparse direct solver. The first method method is based on best rational approximations to $f$ on the negative real axis computed via the Carathéodory-Fejér procedure, and we conjecture that the accuracy scales as $(9.28903\dots)^{-2n}$, where $n$ is the number of complex matrix solves. In particular, three matrix solves suffice to evaluate $f(A)$ to approximately six digits of accuracy. The second method is an application of the trapezoid rule on a Talbot-type contour

On the Power Series Expansion of the Reciprocal Gamma Function

Using the reflection formula of the Gamma function, we derive a new formula for the Taylor coefficients of the reciprocal Gamma function. The new formula provides effective asymptotic values for the coefficients even for very small values of the indices. Both the sign oscillations and the leading order of growth are given.Comment: Corrected a sign in equation (3.21) due to a minor error in (3.19) where the fraction was inadvertently inverted. Now the rough approximation provides an elementary proof that the order of the reciprocal gamma function is 1 and that its type is maxima

Multi-point Taylor Expansions of Analytic Functions

Taylor expansions of analytic functions are considered with respect to several points, allowing confluence of any of them. Cauchy-type formulas are given for coefficients and remainders in the expansions, and the regions of convergence are indicated. It is explained how these expansions can be used in deriving uniform asymptotic expansions of integrals. The method is also used for obtaining Laurent expansions in several points as well as Taylor-Laurent expansions.Comment: 20 pages, 7 figures. Keywords: multi-point Taylor expansions, Cauchy's theorem, analytic functions, multi-point Laurent expansions, uniform asymptotic expansions of integral

Basic Methods for Computing Special Functions

This paper gives an overview of methods for the numerical evaluation of special functions, that is, the functions that arise in many problems from mathematical physics, engineering, probability theory, and other applied sciences. We consider in detail a selection of basic methods which are frequently used in the numerical evaluation of special functions: converging and asymptotic series, including Chebyshev expansions, linear recurrence relations, and numerical quadrature. Several other methods are available and some of these will be discussed in less detail. We give examples of recent software for special functions where these methods are used. We mention a list of new publications on computational aspects of special functions available on our website

Some properties of WKB series

We investigate some properties of the WKB series for arbitrary analytic potentials and then specifically for potentials $x^N$ ($N$ even), where more explicit formulae for the WKB terms are derived. Our main new results are: (i) We find the explicit functional form for the general WKB terms $\sigma_k'$, where one has only to solve a general recursion relation for the rational coefficients. (ii) We give a systematic algorithm for a dramatic simplification of the integrated WKB terms $\oint \sigma_k'dx$ that enter the energy eigenvalue equation. (iii) We derive almost explicit formulae for the WKB terms for the energy eigenvalues of the homogeneous power law potentials $V(x) = x^N$, where $N$ is even. In particular, we obtain effective algorithms to compute and reduce the terms of these series.Comment: 18 pages, submitted to Journal of Physics A: Mathematical and Genera

More Special Functions Trapped

We extend the technique of using the Trapezoidal Rule for efficient evaluation of the Special Functions of Mathematical Physics given by integral representations. This technique was recently used for Bessel functions, and here we treat Incomplete Gamma functions and the general Confluent Hypergeometric Function.Comment: 6 page