Spherical harmonics and Laplace operator

The spherical harmonic function $$Y_l^m(\theta,\phi)$$ is defined to be an eigenfunction of the angular part of the Laplace operator with eigenvalue $$-l(l+1)$$. In other words, it solves the PDE: $$\Delta_{\theta \phi}f := \frac{\partial^2 f}{\partial \theta^2} + \frac{1}{\tan \theta} \frac{\partial f}{\partial \theta} + \frac{1}{\sin^2 \theta} \frac{\partial f}{\partial \phi} = -l(l+1) f.$$ I tried to check this using the Mathematica implementation of the spherical harmonics SphericalHarmonicY. If I plug this into the left-hand side of the above equation,

D[SphericalHarmonicY[l, m, \[Theta], \[Phi]] , {\[Theta], 2}] +
1/Tan[\[Theta]] D[SphericalHarmonicY[l, m, \[Theta], \[Phi]] , \[Theta]]  + 1/Sin[\[Theta]]^2 D[
SphericalHarmonicY[l, m, \[Theta], \[Phi]] , {\[Phi], 2}]

I get: $$-m(m+1)Y_l^m + \cdots.$$ In the dots, there are complicated expressions involving gamma functions and spherical harmonics. Now, it turns out that this actually does equal $$-l(l+1) Y_l^m$$ simply by trying out different values of $$l$$ and $$m$$, keeping $$\theta$$ and $$\phi$$ arbitrary. However, this is not at first obvious and I was wondering if there was a way to make Mathematica recognize this simplification?

I tried do use Refine and Element to tell Mathematica that $$l$$ and $$m$$ and integers, but this did not help.

• Please provide Mathematica code that you've used. Oct 1 '21 at 15:47

To simplify the resulting expression, for general arguments, essentially one needs to apply the following identity:

https://functions.wolfram.com/HypergeometricFunctions/LegendreP2General/17/01/01/0003/

• Could you elaborate your direction, giving us details and a Mathematica code? TIA. Oct 10 '21 at 14:14

Too long for a comment. Indeed,

FullSimplify[D[SphericalHarmonicY[l, m, \[Theta], \[Phi]], {\[Theta], 2}] +
1/Tan[\[Theta]] D[SphericalHarmonicY[l, m, \[Theta], \[Phi]], \[Theta]] +
1/Sin[\[Theta]]^2 D[SphericalHarmonicY[l, m, \[Theta], \[Phi]], {\[Phi], 2}],
Assumptions -> {l, m} \[Element] Integers]

-m (1 + m) SphericalHarmonicY[l, m, \[Theta], \[Phi]] + (E^(-2 I \[Phi]) Sqrt[ Gamma[1 + l - m]] ((1/Sqrt[Gamma[l - m]]) 2 E^(I \[Phi]) (1 + m) Cot[\[Theta]] Sqrt[Gamma[2 + l + m]] SphericalHarmonicY[l, 1 + m, \[Theta], \[Phi]] + ( Sqrt[Gamma[3 + l + m]] SphericalHarmonicY[l, 2 + m, \[Theta], \[Phi]])/Sqrt[ Gamma[-1 + l - m]]))/(Sqrt[Gamma[1 + l + m]])

, but

N[(E^(-2 I \[Phi]) Sqrt[
Gamma[1 + l -
m]] ((1/Sqrt[Gamma[l - m]]) 2 E^(I \[Phi]) (1 +
m) Cot[\[Theta]] Sqrt[Gamma[2 + l + m]] SphericalHarmonicY[
l, 1 + m, \[Theta], \[Phi]] + (Sqrt[
Gamma[3 + l + m]] SphericalHarmonicY[l,
2 + m, \[Theta], \[Phi]])/Sqrt[Gamma[-1 + l - m]]))/(Sqrt[
Gamma[1 + l + m]]) /. {l -> 3, m -> 2, \[Phi] -> Pi/4, \[Theta] -> Pi/3}]

0. - 2.29947 I

Maybe, a bug in FullSimplify.