Lucas-Lehmer 素性测试的更快的按位模数

Question

Lucas-Lehmer 素性测试 测试素数以确定它们是否也是梅森素数。瓶颈之一是计算 (s**2 − 2) % (2**p - 1) 时的模运算。

使用按位运算可以大大加快速度（请参阅 LL 链接），迄今为止我所拥有的最好的结果是：

def mod(n,p):
    """ Returns the value of (s**2 - 2) % (2**p -1)"""
    Mp = (1<<p) - 1
    while n.bit_length() > p: # For Python < 2.7 use len(bin(n)) - 2 > p
        n = (n & Mp) + (n >> p)
    if n == Mp:
        return 0
    else:
        return n

一个简单的测试用例是 p 有 5-9 位数字和 s 有 10,000 多个数字（或更多；它们是什么并不重要）。解决方案可以通过 mod((s**2 - 2), p) == (s**2 - 2) % (2**p -1) 进行测试。请记住，LL 测试中需要 p - 2 次模运算迭代，每次迭代都会以指数方式增加 s，因此需要优化。

有没有办法使用纯Python（包括Python 3）进一步加快速度？有更好的办法吗？

Answer 1

我能找到的最好的改进是从模函数中完全删除 Mp = (1< ，并在开始 LL 测试迭代之前在 LL 函数中预先计算它。使用 while n > Mp: 而不是 while n.bit_length() > p: 也节省了一些时间。

Answer 2

在 n 比 2^p 长得多的情况下，您可以通过执行以下操作来避免一些二次时间的痛苦：

def mod1(n,p):
    while n.bit_length() > 3*p:
        k = n.bit_length() // p
        k1 = k>>1
        k1p = k1*p
        M = (1<<k1p)-1
        n = (n & M) + (n >> k1p)
    Mp = (1<<p)-1
    while n.bit_length() > p:
        n = (n&Mp) + (n>>p)
    if n==Mp: return 0
    return n

[已编辑，因为我之前搞砸了格式；感谢本杰明指出了这一点。寓意：不要从空闲窗口复制粘贴到 SO 中。抱歉！]

（注意：将 n 的长度减半而不是去掉 p 的标准，以及 k1 的确切选择，都有点错误，但这并不重要，所以我没有费心去修复它们。）

如果我取 p=12345 和 n=9**200000 （是的，我知道 p 不是素数，但这在这里并不重要），那么速度大约快 13 倍。

不幸的是，这对您没有帮助，因为在 LL 测试中，n 永远不会大于 (2^p)^2。对不起。