[Scipy-tickets] [SciPy] #1296: associated legendre polynomials in scipy.special fail for calculations with degree higher than approx. 22

SciPy Trac scipy-tickets@scipy....
Sat Oct 9 16:18:33 CDT 2010


#1296: associated legendre polynomials in scipy.special fail for calculations with
degree higher than approx. 22
------------------------------+---------------------------------------------
 Reporter:  josip.mihaljevic  |       Owner:  pv   
     Type:  defect            |      Status:  new  
 Priority:  high              |   Milestone:  0.9.0
Component:  scipy.special     |     Version:  0.7.0
 Keywords:  legendre          |  
------------------------------+---------------------------------------------

Comment(by warren.weckesser):

 (It looks like you have 'order' and 'degree' reversed in your description.
 The first argument of lpmv is the order, and the second is the degree.
 Your title of the bug report matches this usage.)

 It does appear that lpmv is not a very good implementation.  I have
 attached a file that contains two different implementations of the
 associated Legendre polynomial.  One is based on a recursion formula (see
 Numerical Recipes for a similar implementation), and one is based on an
 integral representation.  Both are much more accurate than lpmv.  (The
 first two arguments of assoc_legendre_p and assoc_legendre_p2 are the
 degree then the order.)

 Order 0, degree 16:
 {{{
 In [325]: lpmv(0, 16, 0.001)
 Out[325]: 0.19635390721305157

 In [326]: assoc_legendre_p(16, 0, 0.001)
 Out[326]: 0.19635390806272032

 In [327]: assoc_legendre_p2(16, 0, 0.001)
 Out[327]: (0.19635390806272018, 2.3468911764622238e-15)
 }}}
 Maxima's assoc_legendre_p:
 {{{
 assoc_legendre_p(16, 0, 0.001);
 interval(.1963539080627201,5.223692285805479*10^-15)
 }}}

 Order 0, degree 22:
 {{{
 In [328]: lpmv(0, 22, 0.001)
 Out[328]: -0.16814844246255234

 In [329]: assoc_legendre_p(22, 0, 0.001)
 Out[329]: -0.16814554527766912

 In [330]: assoc_legendre_p2(22, 0, 0.001)
 Out[330]: (-0.16814554527766906, 4.7065059058218157e-14)
 }}}
 Maxima:
 {{{
 assoc_legendre_p(22, 0, 0.001);
 interval(-.1681455452776688,5.982697311505185*10^-15)
 }}}

 Order 0, degree 40 (lpmv is garbage):
 {{{
 In [331]: lpmv(0, 40, 0.001)
 Out[331]: -21449.235366821289

 In [332]: assoc_legendre_p(40, 0, 0.001)
 Out[332]: 0.1252678976534484

 In [333]: assoc_legendre_p2(40, 0, 0.001)
 Out[333]: (0.12526789765344834, 1.1965288185813496e-13)
 }}}
 Maxima:
 {{{
 assoc_legendre_p(40, 0, 0.001);
 interval(.1252678976534483,7.666725305397109*10^-15)
 }}}

-- 
Ticket URL: <http://projects.scipy.org/scipy/ticket/1296#comment:1>
SciPy <http://www.scipy.org>
SciPy is open-source software for mathematics, science, and engineering.


More information about the Scipy-tickets mailing list