Ask Your Question
1

Why does list(primes(2,10)) behave different from [2,3,5,7]?

asked 2012-12-13 10:57:27 +0100

anstei gravatar image

I'm trying to run the following code: http://pastebin.com/UrBHRxr7

I experience a very strange behaviour. With this code I'll hit a memory limit even for the first integral that is being calculated. If I comment out line 17 and use line 16 instead (so primes = [2,3,5,7] instead of primes = list(primes(2,10))) the computation goes through quickly and without any problems!

What could be the source for this? How can I rewrite the code so that I can use parameters instead of a hard-coded list?

I'm using Sage Version 4.8, Release Date: 2012-01-20.

edit retag flag offensive close merge delete

Comments

Of course making the list of primes is 28 microseconds versus less than a microsecond for just the list, but this can't be the problem. Notice that `sage: list(primes(2,10)) == [2,3,5,7]` yields `True`! Huh.

kcrisman gravatar imagekcrisman ( 2012-12-13 11:21:23 +0100 )edit

1 Answer

Sort by ยป oldest newest most voted
2

answered 2012-12-13 12:56:20 +0100

DSM gravatar image

updated 2012-12-13 12:56:53 +0100

One difference between [2,3,5,7] and list(primes(2, 10)) when you're running Sage in Python mode is that the types will be different, Python ints vs. Sage Integers:

from sage.all import *

i = [2,3,5,7]
j = list(primes(2, 10))

print i, map(type, i)
print j, map(type, j)

print i == j

gives

[2, 3, 5, 7] [<type 'int'>, <type 'int'>, <type 'int'>, <type 'int'>]
[2, 3, 5, 7] [<type 'sage.rings.integer.Integer'>, <type 'sage.rings.integer.Integer'>, <type 'sage.rings.integer.Integer'>, <type 'sage.rings.integer.Integer'>]
True

This will affect gamma_sq, for example:

sage: gamma_sq(2, theta1)                  
1/4*(e^(2*I*pi*theta1) - 2)*(e^(-2*I*pi*theta1) - 2)/((e^(2*I*pi*theta1) - 1)*(e^(-2*I*pi*theta1) - 1))
sage: gamma_sq(2, theta1).subs(theta1=1/2)
9/16
sage: gamma_sq(int(2), theta1)            
1/((e^(2*I*pi*theta1) - 1)*(e^(-2*I*pi*theta1) - 1))
sage: gamma_sq(int(2), theta1).subs(theta1=1/2)
1/4

due to truncating division: 1/2 gives 1/2 in QQ in Sage, but 0 in Python. (Note that I had to use the Python-loaded version of the function so that Sage's preparser didn't wrap the 2; preparser(False) would have worked too.)

I can't reproduce a MemoryError but I only have 5.5 at hand.

Using Integers, I get:

Representation Polynomial: x0 + x1 + x2
Prime number 2: 
First integral: -16*(e^(2*I*pi*theta2) - 1)*(e^(-2*I*pi*theta2) - 1)*(2*e^(6*I*pi*theta2) + 7*e^(4*I*pi*theta2) - 15*e^(2*I*pi*theta2) + 1)*pi^2*e^(-2*I*pi*theta2)/((e^(2*I*pi*theta2) - 2)*(e^(-2*I*pi*theta2) - 2)) 
Second Integral: 116*pi^2 
Prime number 3: 
First integral: -36*(e^(2*I*pi*theta2) - 1)*(e^(-2*I*pi*theta2) - 1)*(12*e^(6*I*pi*theta2) + 52*e^(4*I*pi*theta2) - 20*e^(2*I*pi*theta2) + 31)*pi^2*e^(-2*I*pi*theta2)/((e^(2*I*pi*theta2) - 3)*(e^(-2*I*pi*theta2) - 3)) 
Second Integral: -256*pi^2 
Prime number 5: 
First integral: -100*(e^(2*I*pi*theta2) - 1)*(e^(-2*I*pi*theta2) - 1)*(80*e^(6*I*pi*theta2) + 496*e^(4*I*pi*theta2) + 156*e^(2*I*pi*theta2) + 391)*pi^2*e^(-2*I*pi*theta2)/((e^(2*I*pi*theta2) - 5)*(e^(-2*I*pi*theta2) - 5)) 
Second Integral: -28144*pi^2 
Prime number 7: 
First integral: -196*(e^(2*I*pi*theta2) - 1)*(e^(-2*I*pi*theta2) - 1)*(252*e^(6*I*pi*theta2) + 2052*e^(4*I*pi*theta2) + 1100*e^(2*I*pi*theta2) + 1751)*pi^2*e^(-2*I*pi*theta2)/((e^(2*I*pi*theta2) - 7)*(e^(-2*I*pi*theta2) - 7)) 
Second Integral: -282944*pi^2
edit flag offensive delete link more

Comments

Thanks, that's it. It also explains why the result for different p was the same whereas the math says it can't be equal.

anstei gravatar imageanstei ( 2012-12-14 06:38:44 +0100 )edit

Oh, of course - I even checked the types, but in `sage:` so completely useless... good work, DSM.

kcrisman gravatar imagekcrisman ( 2012-12-14 09:02:52 +0100 )edit

Your Answer

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

Add Answer

Question Tools

Stats

Asked: 2012-12-13 10:57:27 +0100

Seen: 476 times

Last updated: Dec 13 '12