[sympy] fixed bug to make sure empty sums return zero. (#1696)
Pieter
fixed bug to make sure empty sums return zero. Issue 3175
You can merge this Pull Request by running:
git pull https://github.com/uctpphd/sympy summations
Or view, comment on, or merge it at:
https://github.com/sympy/sympy/pull/1696
Commit Summary
- fixed bug to make sure empty sums return zero.
File Changes
- M sympy/concrete/summations.py (2)
Patch Links
—
Reply to this email directly or view it on GitHub.
Aaron Meurer
Looks good. Just needs a test.
Anders Kaseorg
The right definition is ∑[n=a … a-1] = 0 and ∑[n=a … b] = -∑n=b+1 … a-1. This is required so that the identity ∑[n=a … b-1] + ∑[n = b … c-1] = ∑[n=a … c-1] holds in all cases, not just a < b < c.
raoulb
once we take into account that he summs over m <= i < n, upper bound not
included. (When he says n = m we should read b = a-1, no?).
BTW: Where did you copy this from? I think I have read is written
exactly that way, but I can not find it right now.
Anders Kaseorg
I wrote it from memory, but yes, it’s equivalent to Karr’s definition, which reads:
1.4 Conventions and Notation. Consider ∑[m≤i<n] f(i). If m < n, this has an obvious meaning. If m ≥ n, this is summation over the null set, which is customarily defined to be zero. We shall follow this convention when m = n, but when m > n, we shall say that
Definition 3. ∑[m≤i<n] f(i) ≜ −∑[n≤i<m] f(i), where m > n.
N.B. This abuse of notation means that “m ≤ i < n”, when written under a ∑, does not imply that m < n.
Aaron Meurer
Karr's definition concerns formal summations. This convention might be useful for an internal representation, but for the user facing Sum, this represents the more standard definition. So I think it should be 0. See http://en.wikipedia.org/wiki/Empty_sum.
Anders Kaseorg
Both definitions agree that, say, ∑[n=10…9] n = 0. But what reference (Wikipedia or otherwise) suggests to you that ∑[n=10…5] n = 0 is a more standard definition than ∑[n=10…5] n = 5·6/2 − 9·10/2 = −30?
Aaron Meurer
I didn't notice that the Wikipedia article limits itself to b = a - 1. I'm not sure why it does. http://en.wikipedia.org/wiki/Empty_product seems a little better to me. The basic idea is that ∑[n=a..b] f(n) means sum({f(n)|n=a..b}), and we define the sum of an empty set to be the additive identity (we also define the sum of a singleton set to be the single element of that set). That article talks about products, but obviously sums should work exactly the same way.
Which reminds me, is Prodoct in SymPy doing the right thing with respect to this?
Aaron Meurer
What do other computer algebra systems do for this?
Anders Kaseorg
I don’t see anything in the empty product article about the notation ∏[n=a…b], much less the case where b < a − 1.
Other CASs give mixed results.
Maple 16 agrees with Karr:
> sum(n, n=10..5);
-30
Mathematica 9.0 may give either answer:
In[1]:= Sum[n,{n,a,b}]
Out[1]= -((-1 + a - b)*(a + b))/2
In[2]:= Sum[n,{n,a,b}] /. {a -> 10, b -> 5}
Out[2]= -30
In[3]:= Sum[n,{n,10,5}]
Out[3]= 0
Matlab R2012b (8.0.0.783) may give either answer:
>> syms n a b
>> symsum(n, n, a, b)
ans =
(b*(b + 1))/2 - (a*(a - 1))/2
>> subs(symsum(n, n, a, b), [a, b], [10, 5])
ans =
-30
>> symsum(n, n, 10, 5)
ans =
0
Sage 5.4 may give either answer (I’m not sure if this counts as “other”):
>>> var('a', 'b', 'n')
(a, b, n)
>>> sum(n, n, a, b)
-1/2*a^2 + 1/2*b^2 + 1/2*a + 1/2*b
>>> sum(n, n, a, b).substitute(a=10, b=5)
-30
>>> sum(n, n, 10, 5)
0
Maxima 5.28.0 gives zero (it’s not smart enough to evaluate the general form):
(%i1) display2d:false;
(%o1) false
(%i2) sum(n, n, a, b);
(%o2) 'sum(n,n,a,b)
(%i3) sum(n, n, 10, 5);
(%o3) 0
Anders Kaseorg
(I’m rather surprised by the self-inconsistency of Mathematica, Matlab, and Sage, by the way—regardless of the outcome of our notational quibbling, I would think it’s a pretty fundamental property that substituting a=10, b=5 into ∑[n=a…b] n gives ∑[n=10…5] n.)
Aaron Meurer
That's because when it evaluates the general form, it computes it assuming that a < b.
Aaron Meurer
In sympy/concrete/summations.py:
> @@ -141,7 +141,7 @@ def doit(self, **hints): > i, a, b = limit > dif = b - a > if dif.is_Integer and dif < 0: > - a, b = b, a > + return S.zero
Should be S.Zero, with a capital Z.
Anders Kaseorg
I found reference to this question in Graham, Knuth, Patashnik, Concrete Mathematics, exercise 2.1 (p. 62):
What does the notation ∑[k=4 to 0] q_k mean?
The given solution (p. 487) is:
2.1. There’s no agreement about this; three answers are defensible:
(1) We can say that ∑[k=m to n] q_k is always equivalent to ∑[m≤k≤n] q_k; then the stated sum is zero.
(2) A person might say that the given sum is q₄ + q₃ + q₂ + q₁ + q₀, by summing over decreasing values of k. But this conflicts with the generally accepted convention that ∑[k=1 to n] q_k = 0 when n = 0.
(3) We can say that ∑[k=m to n] q_k = ∑[k≤n] q_k − ∑[k<m] q_k; then the stated sum is −q₁ − q₂ − q₃. This convention may appear strange, but it obeys the useful law ∑[k=a to b] + ∑[k=b+1 to c] = ∑[k=a to c] for all a, b, c.
It’s best to use the notation ∑[k=m to n] only when n − m ≥ −1; then both conventions (1) and (3) agree.
(Karr’s definition is (3).)
Aaron Meurer
Maybe I'm missing something, but how is it consistant with n - m = -1? If it's n - m ≥ 0, I can understand that. That basically corresponds to a fourth option, which is to make summation(n, (n, 1, 0)) and summation(n, (n, 4, 0)) raise ValueError.
Whatever we end up choosing, we ought to document it in the docstring, with a little discussion of the options, and why we chose the one we did.
@uctpphd sorry that you chose such a controversial bug to fix.
Anders Kaseorg
Maybe I'm missing something, but how is it consistant with
n - m = -1?
In Karr’s paper, ∑[k=m to m−1] = ∑[m≤k<m] was defined to be zero. In CM’s definition (3), ∑[k=m to m−1] = ∑[k≤m−1] − ∑[k<m] = 0, since k ≤ m − 1 and k < m are the same range.
Aaron Meurer
This has gotten me thinking about why this rule holds in integration, whether it is also just a convention there, or if there is a more fundamental reason for it. Since summations can be viewed as a special case of integration, the answer to that question would help clear up this one.
I asked this question at http://math.stackexchange.com/q/261244/781, which should hopefully get an answer from someone who knows a thing or two about the mathematics.
Aaron Meurer
In Karr’s paper, ∑[k=m to m−1] = ∑[m≤k<m] was defined to be zero. In CM’s definition (3), ∑[k=m to m−1] = ∑[k≤m−1] − ∑[k<m] = 0, since k ≤ m − 1 and k < m are the same range.
OK, I thought it might just be an issue of whether the last inequality was < or ≤.
Anders Kaseorg
Your StackExchange question is slightly wrong (and perhaps this is what you were confused about above): you claim that “some authors define ∑[n=4 to 1] n as −∑[n=1 to 4] n”, while the definition that Karr, CM (3), and I are advocating is ∑[n=4 to 1] n = −∑[n=2 to 3] n.
Aaron Meurer
Oh, that's just my bad (I confused the way integrals work, but with sums you have to have start with b + 1).
Anders Kaseorg
I didn’t realize, by the way, that you were looking for mathematical evidence, not just published papers. The mathematical evidence for Karr’s definition is overwhelming, and there’s no need to appeal to an analogy with Lebesgue integration. Basically every identity you ever learned about summation also holds for Karr summation, even when n is zero or negative; for example:
- ∑[k=1 to n] 1 = n
- ∑[k=1 to n] k = n(n + 1)/2
- ∑[k=1 to n] k² = n(n + 1)(2n + 1)/6
- ∑[k=1 to n] k³ = n²(n + 1)²/4, etc.
- ∑[k=0 to n-1] r^k = (1 − r^n)/(1 − r) [r ≠ 1]
- ∑[k=1 to n] f(k) = ∑[k=1 to n−1] f(k) + f(n) = f(1) + ∑[k=2 to n] f(k)
- Every identity listed here (except the binomial coefficient identities, since those are secretely infinite sums to which this discussion doesn’t apply).
And the Fundamental theorem of discrete calculus is directly equivalent to Karr’s definition.
Aaron Meurer
Right. So one question we need to decide is, do we want Sum to represent Karr summation, or should we have a separate object for that? I remember this was discussed before, and we at least decided that we wouldn't use Sum to represent indefinite summation. I guess there are several different definitions of a summation, which are all related, but technically different:
- A sum over a finite subset of a ring (we define the sum of a singleton to be the element in the singleton and the sum over the empty set to be 0, or just define
sum_A f(n)as0 + f(a_1) + f(a_2) + ... + f(a_n), whereA = {a_1, a_2, ..., a_n}). - An infinite series.
- A formal power series. The difference here from the last point is that there is no notion of convergence.
- Summation as the antidifference operator, which is the inverse of Δ.
In SymPy, Sum/summation already try to represent at least the first three of these (finite and infinite summations are computed, and you are allowed to construct and manipulate series which converge nowhere). We don't really have a notion of the fourth, so it doesn't represent it. It's not clear if it should, especially if it would mean conflicting results from the first three definitions, which is what most users would expect.
Anders Kaseorg
Do you see a conflict between the four definitions? (Since you keep bringing up the empty sum, I want to be clear again that all proposed definitions agree that ∑[n=10 to 9] f(n) is the empty sum and equals 0.)
Ronan Lamy
@asmeurer: I'm convinced now that we should use Karr's convention. The other possibilities are just too painful to consider for symbolic limits. Also, writing explicit limits with the wrong ordering on purpose is excessively rare, so I wouldn't consider that there are user expectations to break.
Aaron Meurer
Yep. The current master doit behavior is ok then, right? We just need to fix evalf.
Ronan Lamy
No, it's wrong:
In [1]: Sum(n, (n, 5, 3)).doit()
Out[1]: 12
Aaron Meurer
(Just realized I never actually posted this yesterday)
Let's forget about sums where the limits differ by one. That's just a special case as far as I'm concerned. It's clear to see from the code change that this applies to any b - a < 0.
The conflict is obvious against Sum(x, (x, 4, 0)). For the first definition, this should be 0, because it represents the sum over an empty set. For the fourth case, it is -Sum(x, (x, 1, 3)). I added points two and three merely to point out that there are multiple definitions of a sum. For case 2, I guess whatever we do for finite sums should also apply to Sum(1/x**2, (x, -oo, 1)). Case 3 is actually an interesting one, since in that case, summations do not actually represent numbers per se. Rather, they are formal objects. In that case, the negation convention is probably the most convenient. On the other hand, everything I've seen for formal sums involves formal power series, which always have the same limits. So this issue has never come up.
By the way, there's another interesting angle to all this. Consider summations where the limits are not integers. How should we define those. Should summation(n, (n, 0, pi)) be the same as summation(n, (n, 0, 3)), or should it be pi*(pi + 1)/2 (or should it be an error)?
Ultimately, we need to determine what convention will be the most useful to users, and stick to it. Raising an error seems like a copout, like we just couldn't agree on a convention so we just picked to have none. It also is not possible to tell if limits are integral or how to compare them in general if they are symbolic. Using the Karr definition makes things consistant when substituting before and after summation computation. The disadvantage is that it means that the summation is not actually what the user might expect, namely, summation(f(n), (n, a, b)) does not mean sum({f(n)|a ≤ n ≤ b, n integer}). But it already doesn't mean that in the formal power series sense (i.e., if you enter a nonsensical summation and evaluate it). I'm actually starting to lean toward the Karr definition, both for non-integral limits and for reversed limits.
Regardless, there is still something to do here, because as I noted on the issue page, doit and evalf should be consistant. And we also need to document why we chose whatever we chose.
Anders Kaseorg
I don’t see how the failure of (x, 4, 0) to represent a finite subset of a ring poses any more of a conflict than the failure of (x, -oo, 1) to represent a finite subset of a ring, or even the failure of sqrt(-1) to represent a real number.
(n, 0, pi) is an interesting question. Maple gives this fascinating result:
> sum(n, n=0..Pi);
1/2*(Pi+1)^2-1/2*Pi-1/2
> evalf(sum(n, n=0..Pi));
6.505598528
for which one can in fact consruct a reasonable mathematical justification by generalizing the antidifference operator with Ramanujan’s formula. But it also seems reasonable to defer that question and raise an error for now.
Aaron Meurer
That's the same result that SymPy gives. It's again a consistency between substituting values before or after computing the summation. You can also think of it as an analytic continuation if you are so inclined. I guess you could also compute summation(n, (n, 0, I)).
raoulb
For non-integral limits maybe we should look at this text:
"How to Add a Non-Integer Number of Terms, and How to Produce Unusual Infinite Summations" by Markus Müller.
http://www.math.tu-berlin.de/~mueller/HowToAdd.pdf
I found it in the references of http://en.wikipedia.org/wiki/Indefinite_sum
right now and had no time to read even parts of it. However it looks
very interesting!
Aaron Meurer
And also "How to Add a Noninteger Number of Terms: From Axioms to New Identities" by the same authors, which seems to be a more formal followup. That's very interesting stuff.
Pieter
@asmeurer Hi, I wasn't looking into this too much yesterday, because I was preoccupied with something. This actually looks pretty interesting and I would like to improve it/collaborate once we decide what to settle for. The Karr definition seems manageable to implement. I will look into the evalf function today to see how that can be consistent.
Aaron Meurer
@uctpphd we are definitely going with the Karr definition. Both doit and evalf need to be fixed, and a note added to the docstring (I can write the note if you're not sure what to put).
Julien Rioux
SymPy Bot Summary: 
Failed after merging uctpphd/summations (d956cf9) into master (d503614).
@uctpphd: Please fix the test failures.
PyPy 2.0.0-beta-1; 2.7.3-final-42: fail
Python 2.7.2-final-0: fail
Python 3.2.1-final-0: fail
Sphinx 1.1.3: pass
Docs build command: make clean && make html-errors && make latex && cd _build/latex && xelatex sympy-*.tex
Aaron Meurer
SymPy Bot Summary: Failed after merging uctpphd/summations (d956cf9) into master (d503614).
@uctpphd: Please fix the test failures.
Python 2.5.0-final-0: fail
Python 2.6.6-final-0: fail
Python 2.7.2-final-0: fail
Python 2.6.8-final-0: fail
Python 2.7.3-final-0: fail
PyPy 2.0.0-beta-1; 2.7.3-final-42: fail
Python 3.2.2-final-0: fail
Python 3.3.0-final-0: fail
Python 3.2.3-final-0: fail
Python 3.3.0-final-0: fail
Python 3.3.0-final-0: fail
**Sphinx 1.1.3:** pass