c++ - 通过平方进行模幂运算的溢出可能性

Question

我希望实现 fermat's little theorem用于主要测试。这是我编写的代码:

lld expo(lld n, lld p) //2^p mod n
{
    if(p==0)
        return 1;
    lld exp=expo(n,p/2);
    if(p%2==0)
        return (exp*exp)%n;
    else
        return (((exp*exp)%n)*2)%n;
}

bool ifPseudoPrime(lld n)
{
    if(expo(n,n)==2)
        return true;
    else
        return false;
}

注意:我取了a(<=n-1)的值作为2 .

现在，数字 n 可以大到 10^18 .这意味着变量 exp可以达到 10^18 附近的值.这进一步意味着表达式 (exp*exp)最高可达10^36从而导致溢出。我该如何避免这种情况。

我对此进行了测试，它运行良好直到 10^9 .我正在使用 C++

Answer 1

如果模数接近您可以使用的最大整数类型的限制，事情就会变得有些复杂。如果您不能使用实现双整数运算的库，您可以通过将因子拆分为低阶和高阶部分来自行进行模乘法。

如果模数m太大了 2*(m-1)溢出，事情变得非常挑剔，但如果 2*(m-1)不会溢出，可以忍受。

让我们假设您拥有并使用 64 位无符号整数类型。

您可以通过将因子拆分为低 32 位和高 32 位来计算模积，然后乘积拆分为

a = a1 + (a2 << 32)    // 0 <= a1, a2 < (1 << 32)
b = b1 + (b2 << 32)    // 0 <= b1, b2 < (1 << 32)
a*b = a1*b1 + (a1*b2 << 32) + (a2*b1 << 32) + (a2*b2 << 64)

计算a*b (mod m)与 m <= (1 << 63) , 减少四个产品中的每一个模 m ,

p1 = (a1*b1) % m;
p2 = (a1*b2) % m;
p3 = (a2*b1) % m;
p4 = (a2*b2) % m;

合并这些转变的最简单方法是

for(i = 0; i < 32; ++i) {
    p2 *= 2;
    if (p2 >= m) p2 -= m;
}

p3 相同p4 的 64 次迭代.然后

s = p1+p2;
if (s >= m) s -= m;
s += p3;
if (s >= m) s -= m;
s += p4;
if (s >= m) s -= m;
return s;

这种方式不是很快，但是对于这里需要的少数乘法，它可能足够快了。应该通过减少轮类次数来获得小的加速；先计算(p4 << 32) % m ,

for(i = 0; i < 32; ++i) {
    p4 *= 2;
    if (p4 >= m) p4 -= m;
}

然后所有p2 , p3和 p4 的当前值需要乘以 2³² 模 m ,

p4 += p3;
if (p4 >= m) p4 -= m;
p4 += p2;
if (p4 >= m) p4 -= m;
for(i = 0; i < 32; ++i) {
    p4 *= 2;
    if (p4 >= m) p4 -= m;
}
s = p4+p1;
if (s >= m) s -= m;
return s;