Hello all, Most current implementations of associated Legendre functions (ALFs), including the current GSL version fail due to underflow above degree and order (L,M) 1900. A 2002 paper by Holmes and Featherstone shows that with a very simple scaling modification this limit can be extended to 2700 in double precision. I have recently needed to compute high degree ALFs so I have implemented their algorithm and made a new GSL extension if others are interested.
There are several other advantages in this extension over the current GSL implementation: 1) I have included support for several normalization conventions: - Schmidt semi-normalized ALFs - Spherical harmonic normalized ALFs (GSL already has this) - Fully normalized ALFs - Unnormalized ALFs (GSL has this) 2) The user can choose whether to compute the Condon-Shortley phase of (-1)^m 3) The GSL array versions currently only compute ALFs for a fixed order m and multiple degrees l. This extension now efficiently computes all ALFs for a given maximum degree lmax for all (l,m) up to lmax. 4) The routines in this extension use a workspace which precomputes various factors in the recurrence relations which should lead to a (minor) speedup in the calculations. This extension includes a complete test suite and documentation. The docs contain a plot which shows how the current GSL implementation fails for L,M = 2700. I have put a link on the main GSL webpage for the extension. Patrick Alken
