CODEWITHSUNDEEP
pythonnumpy
Duck
Duck

Reputation: 11

Associated Legendre Function

Hi I am writing Python code which returns the associated Legendre function. Using numpy poly1d function on this part,

firstTerm = (np.poly1d([-1,0,1]))**(m/2.0) # HELP!

It yields an error since it can only be raised to integer.

Is there any other alternative where I can raise the desired function to power 1/2 and etc.?

Upvotes: 1

Views: 1318

Answers (1)

Andras Deak -- Слава Україні
Andras Deak -- Слава Україні

Reputation: 35109

The reason you can't raise your poly1d to half-integer power is that that would not be a polynomial, since it would contain square roots.

While in principle you could orthogonalize the functions yourself, or construct the functions from something like sympy.special.legendre, but your safest bet is symbolic math. And hey, we already have sympy.functions.special.polynomials.assoc_legendre! Since symbolic math is slow, you should probably use sympy.lambdify to turn each function into a numerical one:

import sympy as sym

x = sym.symbols('x')
n = 3
m = 1
legfun_sym = sym.functions.special.polynomials.assoc_legendre(n,m,x)
legfun_num = sym.lambdify(x,legfun_sym)

print(legfun_sym)
print(legfun_num)
x0 = 0.25
print(legfun_sym.evalf(subs={x:x0}) - legfun_num(x0))

This prints

-sqrt(-x**2 + 1)*(15*x**2/2 - 3/2)
<function <lambda> at 0x7f0a091976e0>
-1.11022302462516e-16

which seems to make sense (the first is the symbolic function at x, the second shows that lambdify indeed creates a lambda from the function, and the last one is the numerical difference of the two functions at the pseudorandom point x0 = 0.25, and is clearly zero within machine precision).

Upvotes: 3

Related Questions