Ask Your Question
1

A simple sum causes division by zero exception

asked 2021-02-19 12:53:02 +0200

7sDream gravatar image

I'm very new to SageMath and did some small experiments yesterday. This is my Code:

x, n, r, i = var('x, n, r, i')
f(x, n, r) = x * sum((1 + r) ^ i, i, 1, n)
f(1, 30, 0.0)

It works well if r != 0, but will raise a exception of "division by zero" if equal.

My ipynb file: nbviewer.jupyter.org/gist/7sDream/5db3cfd153269fd1a1cacaf0b60f69bb

It seems some optimizations of sum did not consider the range of variable r. But there is also no change if I added assume(r >= 0)before define f.

So how can I disable this sum optimizations? Or what is wrong with the way I use it?

edit retag flag offensive close merge delete

Comments

Curiously :

sage: def ef(x, n, r): return sum([(1+r)^u for u in (1..n)])
sage: ef(1, 30, 0)
30
sage: ef(1, 30, 0.0)
30.0000000000000
Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2021-02-19 16:30:40 +0200 )edit

Interestingly, if we define fc = fast_callable( x * sum((1 + r) ^ i, i, 1, n), vars=[x,n,r] ), then fc(1, 30, 0.0) evaluates to NaN.

Max Alekseyev gravatar imageMax Alekseyev ( 2021-02-20 04:02:15 +0200 )edit
1

Note that the "optimization" takes place at the creation of f :

sage: var("x, n, r, i")
(x, n, r, i)
sage: f(x,n,r,i)=x*sum((1+r)^i,i,1,n)
sage: f
(x, n, r, i) |--> ((r + 1)^(n + 1) - r - 1)*x/r

FWIW, sympy has an interesting translation :

>>> x, n, r, i = symbols("x, n, r, i")
>>> f = Lambda((x, n, r), x*summation((1+r)**i, (i, 1, n)))
>>> f
Lambda((x, n, r), x*Piecewise((n, Eq(r + 1, 1)), (-(r - (r + 1)**(n + 1) + 1)/r, True)))
Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2021-02-20 10:54:24 +0200 )edit

1 Answer

Sort by ยป oldest newest most voted
2

answered 2021-02-20 11:32:37 +0200

Emmanuel Charpentier gravatar image

updated 2021-02-20 11:47:36 +0200

Sage's implementation of sum is incorrect in this case :

sage: sum((1+r)^i, i, 1, n)
((r + 1)^(n + 1) - r - 1)/r

which is incorrect if r==0.

Compare this with sympy's implementation of summation :

>>> x, n, r, i = symbols("x, n, r, i")
>>> x*summation((1+r)**i, (i, 1, n))
x*Piecewise((n, Eq(r + 1, 1)), (-(r - (r + 1)**(n + 1) + 1)/r, True))

Interestingly, Mathematica falls in the same trap as sage :

In[2]:= InputForm[Sum[(1+r)^i, {i, 1, n}]]    
Out[2]//InputForm= ((1 + r)*(-1 + (1 + r)^n))/r

The workaround is to use the algorithm="sympy" optional argument of sum :

sage: sum((1+r)^i, i, 1, n, algorithm="sympy")
cases(((r + 1 == 1, n), (1, ((r + 1)^(n + 1) - r - 1)/r)))
sage: f(x,n,r,i)=x*sum((1+r)^i,i,1,n, algorithm="sympy")
sage: f
(x, n, r, i) |--> x*cases(((r + 1 == 1, n), (1, ((r + 1)^(n + 1) - r - 1)/r)))
sage: f(1, 30, 0)
30
sage: f(1, 30, 0.0)
30

HTH,

EDIT : This seems to be previously unreported, and is now Trac#31418.

edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2021-02-19 12:53:02 +0200

Seen: 61 times

Last updated: Feb 20