An improved algorithm and a Fortran 90 module for computing the conical function $P^m_{-1/2+i\tau}(x)$

TL;DR

This paper introduces an improved algorithm and a Fortran 90 module for accurately computing conical functions $P^m_{-1/2+i\tau}(x)$, crucial in physics and engineering, with high precision across specified parameter ranges.

Contribution

It presents a new algorithm and a Fortran 90 module that extend the computational range and accuracy for conical functions used in scientific applications.

Findings

Achieves near $10^{-12}$ relative accuracy in standard double precision.

Extends parameter ranges for stable computation.

Provides higher accuracy ($10^{-13}$ to $10^{-14}$) in monotonic regions.

Abstract

In this paper we describe an algorithm and a Fortran 90 module ({\bf Conical}) for the computation of the conical function $P^m_{-\tfrac12+i\tau}(x)$ for $x>-1$, $m \ge 0$, $\tau >0$. These functions appear in the solution of Dirichlet problems for domains bounded by cones; because of this, they are involved in a large number of applications in Engineering and Physics. In the Fortran 90 module, the admissible parameter ranges for computing the conical functions in standard IEEE double precision arithmetic are restricted to $(x,m,\tau) \in (-1,1) \times [0,\,40] \times [0,\,100]$ and $(x,m,\tau) \in (1,100) \times [0,\,100] \times [0,\,100]$. Based on tests of the three-term recurrence relation satisfied by these functions and direct comparison with Maple, we claim a relative accuracy close to $10^{-12}$ in the full parameter range, although a mild loss of accuracy can be found at some points of the oscillatory region of the conical functions. The relative accuracy increases to $10^{-13}\,-\,10^{-14}$ in the region of the monotonic regime of the functions where integral representations are computed ($-1<x<0$).