I have encountered the following sum involving Legendre polynomials, which I hope to reduce to something involving a $\delta$-function: $$ \frac{d^2}{dx^2} \sum_{\ell=0}^{\infty} \frac{2 \ell + 1}{2 \ell (\ell+1)} P_{\ell}(x) P_{\ell}(y) $$

I have played around with the recurrence relation for derivatives of Legendre polynomials without much success ($(1-x^2) P_{\ell}'(x) = -\ell x P_{\ell}(x) + \ell P_{\ell-1}(x)$). I have done some numerical investigation to try and gain intuition, and it does appear that the result vanishes for $x \neq y$, as expected for something proportional to $\delta(x-y)$. However, I've run into tricky divergences and other problems for the $x=y$ case.

If anyone has come across a sum or identity similar to this, I would be most obliged.