Отказ от ответственности: я поддерживаю gmpy и gmpy2.
Я не смогу проверить это до вечера. Но пара комментариев и вопросов.
Вместо использования (s * s)% M, используйте pow (s, 2, M). Это должно быть быстрее.
Что произойдет, если вы используете gmpy2.mpz () вместо gmpy2.xmpz ()?
Вы используете 64-битную версию Python и gmpy2? (Полагаю, что так, но я просто хочу подтвердить.)
Что касается диапазона против xrange, в Python 3.x диапазон заменил xrange.
Редактировать с дополнительной информацией.
Причиной сбоя было переполнение внутренних структур в 32-разрядной сборке. Использование 64-битной версии Python и gmpy или gmpy2 - правильное решение.
Функция unpack (x, n) похожа на split () для строки: она делит число на серию n-битных значений. Это эквивалентно, но намного быстрее, чем:
def unpack(x,n):
r = []
m = 2**n
while x:
x, temp = divmod(x,m)
r.append(temp)
return r
Некоторая документация доступна через help(gmpy2.unpack)
, но лучшая документация есть в моем списке дел.
Причина, по которой unpack () может использоваться для исключения операции%, аналогична добавлению цифр числа base-10 для проверки делимости на 9. В этом случае unpack () создает числа p-битов и мы делим на 2 ** р - 1.
Вот тестовый код:
import gmpy2
import time
def mersenne1(p):
'''Primality test for Mersenne prime: 2**p -1.
Uses native Python longs. Does not verify that p is prime.'''
s = 4
M = 2**p - 1
for i in range(p-2):
s = ((s*s)-2) % M
return False if s else True
def mersenne2(p):
'''Primality test for Mersenne prime: 2**p -1.
Uses gmpy2.mpz. Does not verify that p is prime.'''
s = gmpy2.mpz(4)
M = gmpy2.mpz(2)**p - 1
for i in range(p-2):
s = ((s*s)-2) % M
return False if s else True
def mersenne3(p):
'''Primality test for Mersenne prime: 2**p -1.
Uses gmpy2.mpz and no mod. Does not verify that p is prime.'''
s = gmpy2.mpz(4)
M = gmpy2.mpz(2)**p - 1
for i in range(p-2):
s = (s*s)
s = sum(gmpy2.unpack(s, p))
s = sum(gmpy2.unpack(s, p))
if s < 2:
s = M - 2 + s
else:
s = s - 2
return False if s else True
if __name__ == "__main__":
p = 44497
start = time.time()
result1 = mersenne1(p)
print("Elapsed time: {:6.3f}".format(time.time() - start))
start = time.time()
result2 = mersenne2(p)
print("Elapsed time: {:6.3f}".format(time.time() - start))
start = time.time()
result3 = mersenne3(p)
print("Elapsed time: {:6.3f}".format(time.time() - start))
if result1 == result2 == result3:
print("All three tests are equal!")
else:
print("Oops, something has gone wrong.")
И немного времени работы ...
C:\x64\Python32>python.exe mersenne.py
Elapsed time: 163.683
Elapsed time: 12.782
Elapsed time: 3.630
All three tests are equal!