Integration of Associated Legendre Polynomial I am interested in the following integral $$I=\int_{-1}^1P_\ell^2(x)P_n(x)\mathrm{d}x,$$
where $P_n(x)$ is Legendre Polynomial of $n$th order, and $P_\ell^2$ is Associated Legendre Polynomial. Any one has any idea on how to proceed?
 A: By definition? $P_\ell^m(x) = (-1)^m (1-x^2)^{m/2} P_\ell^{(m)}(x)$ where $f^{(m)}$ stands for the $m^{th}$ derivative of $f$. We have
$$P_\ell^2(x) = (1-x^2)\frac{d^2}{dx^2}P_\ell(x)
= \frac{d}{dx}\left[(1-x^2)\frac{d}{dx}P_\ell(x)\right] + 2x P'_\ell(x)
$$
Since $P_\ell(x)$ satisfy the Legendre DE:
$$\frac{d}{dx}\left[(1-x^2)\frac{d}{dx}P_\ell(x)\right] + \ell(\ell+1) P_\ell(x) = 0$$ We find
$$\begin{align}
\int_{-1}^1 P_\ell^2(x) P_m(x) 
&= - \ell(\ell+1) \int_{-1}^1 P_\ell(x)P_m(x) dx + \int_{-1}^1 2xP'_\ell(x)P_m(x)dx\\
&= -\frac{2\ell(\ell+1)}{2\ell+1} \delta_{\ell m} + \int_{-1}^1 2xP'_\ell(x)P_m(x)dx\tag{*1}
\end{align}
$$
Let $J_{\ell,m}$ be the integral on RHS$(*1)$. There are 3 cases
Case I : $\ell < m$
$2xP_\ell(x)$ is a polynomial of degree $< m$. By orthogonality of the Legendre polynomials, $J_{\ell,m} = 0$
Case II : $\ell > m$,
$$J_{\ell,m} = \int_{-1}^1 2xP_m(x) dP_\ell(x)
= [ 2xP_m(x)P_l(x) ]_{-1}^1 - \int_{-1}^1 (2xP_m(x))'P_\ell(x) dx
$$
   Since $(2xP_m(x))'$ is a polynomial of degree $< \ell$. By orthogonality of the Legendre polynomials again, the last integral vanishes. This leads to
$$J_{\ell,m} = 2(1 + (-1)^{\ell+m}) = \begin{cases}
4, & \ell = m + 2k, k > 0,\\
0, & \text{ otherwise }
\end{cases}
$$
Case III : $\ell = m$,
$$J_{\ell,m} = \int_{-1}^1 x d P_\ell(x)^2 =
\bigg[ x P_\ell(x)^2 \bigg]_{-1}^1 - \int_{-1}^1 P_\ell(x)^2 dx
= 2 - \frac{2}{2m+1} = \frac{4m}{2m+1}$$
Substitute these 3 expressions of $J_{\ell,m}$ back into $(*1)$, we obtain:
$$\int_{-1}^1 P_\ell^2(x)P_m(x) dx = \begin{cases}
-\frac{2m(m-1)}{2m+1}, & \ell = m\\
4, & \ell = m + 2k,\; k > 0\\
0, & \text{ otherwise }
\end{cases}$$
This is the pattern pointed out in Claude Leibovici's answer.
A: This is not an answer but it is too long for a comment.
Considering $$I_{n,l}=\int_{-1}^{+1} P_l^2(x) P_n(x)\,dx$$ this integral seems to show interesting patterns I give you below (this is just based on numerical evaluation and observation).
For positive values of $n,l$

*

*for $l<n \implies I_{n,l}=0$

*for $l=n+(2k-1)\implies I_{n,l}=0$

*for $l=n+2k\implies I_{n,l}=4$

*for $l=n\implies I_{n,l}=-\frac{2n(n-1)}{2n+1}$
