Legendre Polynomials are defined for arguments lying in -1 ⩽ x ⩽ 1

[sloisel]

I got this error:

julia> Pl(1.2,3)
DomainError with 1.2:
Legendre Polynomials are defined for arguments lying in -1 ⩽ x ⩽ 1

I am somewhat perplexed. A polynomial defined on -1 ⩽ x ⩽ 1 is of course defined for all complex x, and in particular for x=1.2. Evaluating a Legendre polynomial at x=1.2 is particularly straightforward when using recurrence relations. I checked the source code of this library and it looks to me like it uses the recurrence relation. I'm not sure why it refuses to work at x=1.2.

Is this a bug?

[jishnub]

What result are you expecting? Legendre polynomials are solutions to the differential equation within [-1,1]. Is there a well-defined extension to other domains? I mean the polynomial may be evaluated elsewhere for sure, but what properties are we expecting at x=1.2?

[sloisel]

Please see the Wikipedia article. The definition of the polynomials by way of differential equations is one of the possible definitions. The differential equation you speak of has singularities at x=-1 and x=+1, but the solution is nevertheless defined for the whole complex plane.

The Legendre polynomial of degree 3 is 1/2(5x^3-3x), so

julia> legendreP(1.2,3)
(2.5199999999999996, 9.299999999999997)

This is from my own implementation, which returns both P and its derivative P'.

Another implementation can be found in FastGaussQuadrature as the unexported function innerRec:

julia> FastGaussQuadrature.innerRec(3,1.2)
([2.5199999999999996], [9.299999999999999])

You can see in the Wikipedia article, a table of the first 11 Legendre polynomials, and note that they are bona fide polynomials defined for arbitrary complex values of the parameter x. (Actually, the polynomials are defined for x in any ring, e.g. x could also be a square matrix.)

If you implemented Legendre polynomials by numerically solving the ODE, then presumably you'd have to take some care crossing the singular points x=±1, but the definition/implementation of the polynomials by recurrence formula has no such difficulties and behaves as desired for arbitrary complex x.

[jishnub]

A PR to remove the restriction is welcome

[sloisel]

I don't actually know how to do that, and the incentive to learn is diminished by the fact that I've now implemented my own? Also, I wouldn't want to mess up your code. Hopefully you will take my information as a courtesy?

[jishnub]

Yes, thanks anyway

[sloisel]

I just realized I had failed to answer one of your questions, which I do now for posterity. The Gauss-Legendre quadrature nodes are the roots of the Legendre polynomials. These roots can be located by polynomial rootfinding algorithms (e.g. Newton iteration, Durand-Kerner). These algorithms occasionally drift outside of [-1,1] in search for the roots.