在 C 中求幂和除数非常大的模数

Question

我需要 C 中的一个函数来计算 a^n mod q，其中除数 q 被确定为非常大 (15,383,399,235,709,406,497) 并且指数 n 可能也为这么大。

基于模乘属性，即(a * b) mod n = ((a mod n) * (b mod n)) mod n，我的尝试如下：

#include <stdio.h>

unsigned long long int modExp(unsigned long long int base, unsigned long long int expo, unsigned long long int mod)
{
    unsigned long long int out = 1;
    unsigned long long int cnt;
    for (cnt=expo; cnt>0; cnt--)
    {
        out = (out * base) % mod;
    }
    return out;
}

int main()
{
    printf("%llu", modExp(3, (189 + 50 * 223), 15383399235709406497));
    return 0;
}

如上面的 main 函数所示，我用基数 3、指数 (189 + 50 * 223) 和除数 15383399235709406497 测试了我的函数 modExp。它给出了输出 3915400295876975163 和一些警告

In function ‘findxA’:
findxA.c:17:32: warning: integer constant is so large that it is unsigned [enabled by default]
     unsigned long long int h = 12036625823877237123;
                                ^
findxA.c:17:5: warning: this decimal constant is unsigned only in ISO C90 [enabled by default]
     unsigned long long int h = 12036625823877237123;
     ^
findxA.c:18:32: warning: integer constant is so large that it is unsigned [enabled by default]
     unsigned long long int q = 15383399235709406497;
                                ^
findxA.c:18:5: warning: this decimal constant is unsigned only in ISO C90 [enabled by default]
     unsigned long long int q = 15383399235709406497;
     ^
findxA.c: In function ‘main’:
findxA.c:34:48: warning: integer constant is so large that it is unsigned [enabled by default]
     printf("%llu", modExp(3, (189 + 50 * 223), 15383399235709406497));
                                                ^
findxA.c:34:5: warning: this decimal constant is unsigned only in ISO C90 [enabled by default]
     printf("%llu", modExp(3, (189 + 50 * 223), 15383399235709406497));
     ^

为了验证这个答案，我将结果（由 C 函数给出）与通过计算表达式 3^(189 + 50 * 223) `mod` 15383399235709406497 给出的输出进行比较，写在 Haskell 中。此 Haskell 表达式计算为不同的小数 12349118475990906085。我认为是我的 C 函数出错了，因为我相信 Haskell 在处理如此大的小数方面做得更好。

如何改进我的 C 函数 modExp？

编辑：我尝试了的选项但是，当我试图处理unsigned long long int的十进制时，我已经切换了每个输入return 类型从 int 到 unsigned long long int。这导致了分段错误。

Edit2: 我错误地使用了上面link描述的功能。它不会给出分段错误；但它仍然没有输出正确的小数。

Answer 1

这种方式永远行不通！

这里单独乘法的这个结果

out = (out * base) % mod;

可能需要超过 64 位的基础数据类型。如果发生整数溢出（即最高有效位被截断），mod 操作就会出错！

使用更大的数据类型，或使用不同的方法:-)

顺便说一句，如果您愿意，请使用此测试代码来查看您的输入是否确实发生了两次整数溢出：

#include <stdio.h>

unsigned long long int modExp(unsigned long long int base, unsigned long long int expo, unsigned long long int mod)
{
    unsigned long long int out = 1;
    while(expo>0)
    {
        if(expo % 2 == 1)
            out = (out * base) % mod;
        expo = expo >> 1;
        if (base * base < base) printf("WARNING! INTEGER OVERFLOW!!!!\n");
        base = (base * base) % mod;
    }
    return out;
}

int main()
{
    printf("%llu", modExp(3, (189 + 50 * 223), 15383399235709406497));
    return 0;
}

Answer 2

减少溢出的几率，可以靠(x*y)%z == ((x%z) * (y%z)) % z.

例如（未经测试）：

unsigned long long int modExp(unsigned long long int base, unsigned long long int expo, unsigned long long int mod)
{
    unsigned long long int baseMod = base%mod;
    unsigned long long int out = 1;
    unsigned long long int cnt;
    for (cnt=expo; cnt>0; cnt--)
    {
        out = ( (out%mod) * baseMod) % mod;
    }
    return out;
}

另请注意，求幂可以像 "product of squares" 那样更有效地完成。例如，x ** 5 == 1*(x**4) * 0*(x**2) * 1*(x**1) 因为 5 == 1*4 + 0*2 + 1*1（或 5 == 101b）。

例如（未经测试）：

unsigned long long int modExp(unsigned long long int base, unsigned long long int expo, unsigned long long int mod)
{
    unsigned long long int out = 1;
    unsigned long long int temp = base%mod;

    while(exp0 > 0) {
        if( (exp0 & 1) != 0) {
            out = (out * temp) % mod;
        }
        temp = (temp*temp) % mod;
        exp0 >>= 1;
    }
    return out;
}

对于大指数，这应该会在性能上产生巨大差异（例如，如果指数为 123456，则循环将有 17 次迭代而不是 123456 次迭代）。

还有；如果这是为了某种密码学（并非不可能）；那么你应该清楚地说明这一点，因为你可能想要 "constant time" （以减少计时边信道的机会 - 例如，从 modExp() 执行所需的时间推断有关指数的信息）。

最后；即使进行了更改，最大 15,383,399,235,709,406,497 的数字对于 unsigned long long 来说可能太大了（您需要确保 mod * mod <= ULLONG_MAX）；这意味着您可能需要 use/implement "big integer" 操作（例如，可能 typedef struct { uint32_t digit[4]; } MY_NUMBER_TYPE 处理 128 位整数，具有乘法和模函数）。

在 C 中求幂和除数非常大的模数

Finding modulus with very large exponentiation and divisor in C

c

modulo

unsigned-long-long-int