Legendre expansion of the Dirac delta function

There is a known expansion for the Dirac delta function in the interval $$(-1, 1)$$ in terms of the Legendre polynomials as

$$\delta(x) = \sum_{k = 0}^{\infty} (-1)^k \frac{(4k + 1) (2k)!}{2^{2k + 1} (k!)^2} P_{2k}(x) .$$

I would like to verify this identity numerically at the endpoint, i.e., $$x = 1$$. We know that $$\delta(1) = 0$$, so inserting $$x = 0$$ in the right-hand side of the above may result in zero. If we add the first two hundreds of the terms in the sum using Mathematica, results in $$7.994$$. As we add more terms the value of the sum increases, for example, for the first three hundred terms, we obtain $$9.784$$; and it becomes larger as we add more terms.

Thus, it seems that the identity is not valid at $$x = 1$$. However, if we set the number of terms to infinity as

N[Sum[(-1)^k (4 k + 1) (2 k)! LegendreP[2 k, 1]/(2^(2 k + 1) k! k!), {k, 0, Infinity}]]


now, it returns $$0$$. How can I understand the Mathematica's behavior?

• The so-called $\delta$- function is not a usual function, but a distribution. The convergence is valid in the weak topology only (see encyclopediaofmath.org/wiki/Generalized_function). Jun 19 '20 at 14:46
• @user64494 actually, it makes perfect sense to evaluate a generalized function away from its singularities. What doesn't make sense is $\delta(0)$, but this is not what the OP is trying to do. Jun 19 '20 at 22:16
• @user64494 that same link doesn't hesitate to talk about Heaviside function (which is also a distribution!) having some values at $x\ge0$ and $x<0$ (see the Examples section after Differentiation). And the modulus function, being a distribution, is (even better) defined everywhere on $\mathbb R$. It's straightforward to define how you'd evaluate a generalized function $f$ at (non-singular) $x_0$: just consider a sequence of ever-narrowing around $x_0$ test functions, similar to the sequence of nascent delta functions. The actions of $f$ on them will converge to the value $f(x_0)$. Jun 20 '20 at 9:55
• @user64494 please don't re-post the same comment just for the sake of editing it to add several words. I get pinged for each of its copies. Jun 20 '20 at 12:33
• @user64494 as for what defines what, maybe in some cases it's useful to distinguish the distribution and the associated function defined at that distribution's non-singular support. But the very term "generalized function" already hints that it's often a needless distinction. It may be useful when you're trying to prove something related to the properties of distributions as the functionals that can't be defined as integrals of a smooth function times a test function, but that isn't needed too often in applied maths contexts. Jun 20 '20 at 12:40

As I understood you start from the completeness relation $$\sum_{\ell=0}^\infty \frac{2\ell + 1}{2} P_\ell(x)P_\ell(y) = \delta(x-y)$$ and use that $$P_n(0) = \begin{cases} \frac{(-1)^{m}}{4^m} \tbinom{2m}{m} = \frac{(-1)^{m}}{2^{2m}} \frac{(2m)!}{\left(m!\right)^2} & \text{for} \quad n = 2m \\ 0 & \text{for} \quad n = 2m+1 \,. \end{cases}$$

The convergence is, however, not uniform

d[x_,n_]:=Sum[(-1)^k (4 k + 1) (2 k)! LegendreP[2 k, x]/(2^(2 k + 1) k! k!), {k, 0,n}]
Plot[d[x,#],{x,-1,1},PlotRange->All]&/@Range[5,25,5]


• Thanks for your answer! What does Range[5, 25, 5] do in your code?
– user72595
Jun 19 '20 at 13:06
• It just creates a list {5,10,15,20,25} serving as upper limits for summations Jun 19 '20 at 13:08
• The so-called $\delta$- function is not a usual function, but a distribution. The convergence is valid in the weak topology only (see encyclopediaofmath.org/wiki/Generalized_function). Jun 19 '20 at 14:45
• @user64494 This is an important clarification, thank you for providing the link. Jun 19 '20 at 14:55

I think you're misinterpreting Mathematica output. First of all, notice that the first output you get from Mathematica is a message that precludes the subsequent numerical results being sensible (or at the very least must make you very cautious of them):

Sum::div: Sum does not converge.

If you omit the N[] wrapper, you'll get the actual sum simplified to take into account that $$P_n(1)\equiv1$$:

$$\sum\limits_{k=0}^\infty \frac{(-1)^k 2^{-1-2k}(1+4k)(2k)!}{(k!)^2}.\tag1$$

Plotting the terms of this sum, you'll see that they increase indefinitely in magnitude and alternate in sign. This should already make it obvious that it can't converge.

Now, with your N[] wrapper you also get another message:

NumericalMathNSequenceLimit::seqlim: The general form of the sequence could not be determined, and the result may be incorrect.

What you get after this (on my system it's about $$1.82\times10^{-14}$$) is a result of an attempt at numerical summation of the divergent series.

While trying to naively sum this series, hoping that after a large enough value of $$k$$ changes will become negligible, is futile (since the sum is divergent), one may take a different approach, and use a resummation, e.g. the following scheme (known as Abel summation): define

$$S(a)=\sum\limits_{k=0}^\infty \frac{(-1)^k 2^{-1-2k}(1+4k)(2k)!}{(k!)^2}e^{-ak}=\frac{e^a-1}{2\sqrt{e^{-a}+1}(e^a+1)}.\tag2$$

Then eliminate the $$a$$ by an operation that would, if taken on the terms of the series, yield $$(1)$$:

$$\lim\limits_{a\to0} S(a)=0.\tag3$$

This is one of the possible values that could be "assigned" to the series $$(1)$$. I suppose the N[]` function was also trying to evaluate the series using some kind of regulator to speedup convergence. This would indeed work for a convergent sum, since resummation schemes normally preserve the values of convergent sums, only supplying additional values for divergent ones. But in this case, since the series is divergent, you got a resummation, which may or may not be useful to you.

• Many thanks for your detailed analysis of the question.
– user72595
Jun 20 '20 at 12:18
• A question: Shouldn't we have $\sqrt{1 + e^{-a}}$ instead of $\sqrt{1 - e^{-a}}$ in Eq. (2)?
– user72595
Jun 21 '20 at 13:45
• @Novice indeed, that was a typo, fixed now. Jun 21 '20 at 14:05