Tough Legendre Integral I am currently fighting with the following integral. I have simplified it to this one here:
$\int_{-1}^{\cos(\alpha)} P_l(t)P_{l'}(t) dt$, where $P_l$ is the l-th Legendre polynomial.
unfortunately you cannot use orthogonality, so this is somewhat hard to do, but maybe somebody here has an idea.
 A: Assume $n \ne m$, we have
$$
\frac{d}{dx}\left[(1-x^2)P'_n(x))\right] + n(n+1) P_n(x) = 0\\
\frac{d}{dx}\left[(1-x^2)P'_m(x)\right] + m(m+1) P_m(x) = 0
$$
This implies
$$\begin{align}
 & ( n(n+1) - m(m+1) )P_n(x)P_m(x)\\
=& -\left( \frac{d}{dx}\left[(1-x^2)P_n'(x)\right] \right) P_m(x)
   +P_n(x) \left( \frac{d}{dx}\left[(1-x^2)P_m'(x)\right] \right)\\
=& \frac{d}{dx}\left[(1-x^2)\left(P_n(x)P'_m(x) - P'_n(x)P_m(x)\right)\right]\\
\end{align}$$
As a result,
$$
\int_{-1}^{t}P_n(x)P_m(x) dx = \frac{(1-t^2)\left(P_n(t)P'_m(t) - P'_n(t)P_m(t)\right)}{(n+m+1)(n-m)}
\tag{*1}$$
Using the recurrence relation:
$$ (x^2-1)P'_n(x) = nxP_n(x) - nP_{n-1}(x)$$
We get
$$(*1) = \frac{(n-m)tP_n(t)P_m(t)-nP_{n-1}(t)P_m(t) + mP_n(t)P_{m-1}(t)}{(n+m+1)(n-m)}$$
For $n = m$ and general $t$, no idea how to express the integral in compact form. However,
for $t = 0$, we have:
$$\int_{-1}^0 P_n(x)^2 dx =  \frac12 \int_{-1}^1 P_n(x)^2 dx = \frac{1}{2n+1}$$

Note For those with sharp eyes, one will notice in immediate steps of above
  derivations, there are terms invoking the function $P_{-1}$. There are
  there to simplify the expression. Please treat all $P_{-1}$ as
  constant $-1$.

Update
For general $t \in (-1,1)$, let $\theta_t(x)$ be the step function defined by:
$$\theta_t(x) = \begin{cases}1,&\text{ if } x < t\\1/2,&\text{ if } x = t\\0,&\text{ if } x > t.\end{cases}$$
Expand $\theta_t(x)$ in terms of $P_l(x)$ over $[-1,1]$, we have:
$$\begin{align}\theta_t(x) 
= & \sum_{l=0}^{\infty}\frac{2l+1}{2} \left[\int_{-1}^{1} \theta_t(y) P_l(y) dy \right] P_l(x)\\
= & \sum_{l=0}^{\infty}\frac{2l+1}{2} \left[\int_{-1}^{t} P_l(y) dy \right] P_l(x)\\
= & \frac12 \sum_{l=0}^{\infty} \left[\int_{-1}^{t} (P'_{l+1}(y) - P'_{l-1}(y)) dy \right] P_l(x)\\
= & \frac12 \sum_{l=0}^{\infty} \left( P_{l+1}(t) - P_{l-1}(t) \right) P_l(x)
\end{align}$$
Multiply both sides by $P_n(x)^2$ and integrate, we get:
$$\int_{-1}^{t} P_{n}(x)^2 dx =\int_{-1}^{1} \theta_t(x) P_{n}(x)^2 dx
=\frac12 \sum_{l=0}^{\infty} \left[ \int_{-1}^{1} P_n(x)^2 P_l(x) dx \right]
\left( P_{l+1}(t) - P_{l-1}(t) \right)
$$
The integral appeared in above expansion can be expressed in terms of the Wigner 3-j symbols:
$$\int_{-1}^{1} P_{j_1}(x) P_{j_2}(x) P_{j_3}(x) dx 
= 2 \begin{pmatrix}j_1 & j_2 & j_3\\0 & 0 & 0\end{pmatrix}^2$$
In general, the Wigner 3-j symbols $\begin{pmatrix}j_1 & j_2 & j_3\\m_1 & m_2 & m_3\end{pmatrix}$ vanishes unless


*

*$j_1, j_2, j_3$ satisfy the triangular inequalities.

*$m_1 + m_2 + m_3 = 0$

*if $m_1 = m_2 = m_3 = 0$, then $j_1+j_2+j_3$ need to be even.


For the special case we need where $j_1 + j_2 + j_3 = 2g$ is even, $\begin{pmatrix}j_1 & j_2 & j_3\\0 & 0 & 0\end{pmatrix}$ is equal to:
$$(-1)^g \sqrt{ \frac{(2g - 2j_1)!(2g-2j_2)!(2g-2j_3)!}{(2g+1)!}}
        \frac{g!}{(g-j_1)!(g-j_2)!(g-j_3)!}
$$
Substitute this back into $(*2)$, we get:
$$\begin{align}
\int_{-1}^{t} P_{n}(x)^2 dx 
= & \sum_{k=0}^{n}
\begin{pmatrix}n & n & 2k\\0 & 0 & 0\end{pmatrix}^{\!2} \left( P_{2k+1}(t) - P_{2k-1}(t) \right)\\
= & \sum_{k=0}^{n}\frac{\binom{2k}{k}^2\binom{2n-2k}{n-k}}{(2n+2k+1)\binom{2n+2k}{n+k}}\left( P_{2k+1}(t) - P_{2k-1}(t) \right)
\end{align}
$$
