Numeric approximation of infinite sum

How should I get a numeric approximation of the sum of the legendre polynomials expansion symbolically (The following formula is part of a bigger symbolic formula in which the only part that does not get calculated is the legendre polynomials expansion).

syms n 
symsum(subs(legendreP(n,1)),n,0,Inf)

ans = 

As in the following example

syms k 
symsum(1^k/factorial(k),k,1,Inf)

ans = 

[Torsten]

How do you arrive at the above formula ?

Answer 1

[John D'Errico]

Um, first, legendreP will return either 1 or -1, depending on the value of n. So effectively, we can replace that call to LegendreP using this identity

legendreP(n,-1) == (-1)^n

The reason is, that is the endpoint of the interval a Legendre polynomial typically lives on, and those Legendre polynomials are normalized. You may not believe me, so a simple experiment might help to convince you.

for n = 0:10
    legendreP(n,-1)
end

ans = 1

ans = -1

ans = 1

ans = -1

ans = 1

ans = -1

ans = 1

ans = -1

ans = 1

ans = -1

ans = 1

As I said, it can be directly replaced with (-1)^n.

Next, as I look at the terms in that sum, I see it reduces trivially. Almost all of it goes away. For example, sigma1 is in there, but do you not see that you can factor sigma1 out? In fact, then you should see that if this is sigma1

sigma1 = 1/4^(2*n+1)

then

4^(2*n+1) * sigma1 == 1

So sigma1 goes away completely, as does that extra factor of 4^(2*n+1).

Next, what is 2^n/(2^(2*n+1)?

2^n/(2^(2*n+1) == 1/2^(n+1)

But then you have another factor of 4^n. Can you write 4^n in a different way? Yes.

4^n = (2^n)*(2^n)

But since we have a term of the form 1/2^(n+1), then again, we have a massive cancellation occur. Seriously. This is all just basic algebra. In the end, the expression inside that infinite sum reduces to this product:

(-1)^n*2^(n-1)*(1+n)/n

Which you should realize, even without resorting to the symbolic toolbox, that infinite sum diverges. Look carefully at the terms. Note that you have a factor of 2^(n-1) in the denominator. Even though this is a series with alternating signs, that won't help. If you need to do so, look up the rules for when an infinite series will converge.

But then even worse, the term with n==0? That first term results in a divide by zero. In case you don't believe me though, we can ask MATLAB:

syms n
symsum((-1)^n*2^(n-1)*(1+n)/n,[0,inf])

Error using symengine
First summand has a zero denominator.

Error in sym/symsum (line 70)
rSym = mupadmex('symobj::map',fsym.s,'symobj::symsum',x.s,a.s,b.s);

I'm sorry. But that sum is meaningless as you have written it. Perhaps you have made multiple mistakes in writing it down. Perhaps when happens outside of that expression is not as you say. Perhaps there are some terms you have left out, or you have miswritten it. But what you have written here has no finite value we can attribute to it.