A Fast Algorithm for the Calculation of the Roots of Special Functions

Abstract

We describe a procedure for the determination of the roots of functions satisfying secondorder ordinary differential equations, including the classical special functions. The scheme is based on a combination of the Pri.ifer transform with the classical Taylor series method for the solution of ordinary differential equations, and requires 0(1) operations for the determination of each root. When the functions in question are classical orthogonal polynomials (Legendre, Hermite, etc.), the techniques presented here also provide tools for the evaluation of the weights for the associated Gaussian quadratures. The performance of the scheme for several classical special functions (prolate spheroidal wave functions, Bessel functions, and Legendre, Hermite, and Laguerre polynomials) is illustrated with numerical examples.

Noniterative Computation of Gauss--Jacobi Quadrature

“…To handle this, we can use the following expansion (see [29,Eqn. 10 20) with z = u and λ = 1 + h/u. This gives…”

“…We refer to [40] for additional examples of interpolation using Jacobi nodes and weights (some of them requiring few hundreds of nodes for double precision accuracy). In that paper the barycentric weights where computed with a Matlab version of the algorithm [20], and it is mentioned that the inclusion of the methods of [23] should improve the performance.…”

“…In [23], the method for Gauss-Jacobi quadrature relies on simple asymptotic approximations for the nodes which are iteratively refined by Newton's method, with the Jacobi polynomial (whose roots are the nodes of the formula) and its derivative computed by two types of available asymptotic expansions; the same approach was used in [39] for Gauss-Hermite quadrature. A similar strategy was considered earlier for Gauss-Legendre quadrature in [33,20,5], but using different methods for computing the Legendre polynomials; methods for Gauss-Laguerre and Gauss-Hermite were also described in [20]. A more recent approach, based in the integration of ODEs using non-oscillatory phase representations is that of [7], which appears to be competitive for Gauss rules of very high degree.…”

Asymptotic Approximations to the Nodes and Weights of Gauss–Hermite and Gauss–Laguerre Quadratures

“…(4) Next, compute η from (17). (5) Then the better value of t again follows from (19). (6) Finally, the approximation for the requested zero is x k ∼ μt, see (6).…”

“…Most iterative methods for the computation of the Gaussian nodes (with the exception of [2]) require accurate enough first approximations to ensure the convergence of the iterative method (typically the Newton method); for two recent examples, see [3,4]. An alternative approach [5], although less efficient for high degrees than iterative methods with asymptotic first approximations [3,4], consists in guessing these first approximations by integrating a Prufer-transformed ordinary differential equation (ODE) with a Runge-Kutta method, and then refining these guesses by the Newton method (however, asymptotic approximations were also used in this reference for the particular case of Gauss-Legendre quadrature). More recently, noniterative methods based on asymptotic approximations for the computation of Gauss-Legendre nodes and weights were developed in [6], which were shown to outperform iterative approaches.…”