algorithm - 找到 a^b^c^... mod m

Question

我想计算：

a ^{b ^{c ^{d ^{。^{. ^.}}}}}模数

你知道任何有效的方法，因为这个数字太大但是 a ， b ， c ， ... 和 m 适合一个简单的 32 位 int。

有任何想法吗？

警告：这个问题与寻找^b mod m 不同。

另请注意， a ^{b ^c}与 (a ^b ) ^c不同。后者等于 a ^bc。求幂是右结合的。

score 24 · Accepted Answer

a ^{b ^c} mod m = a ^{b ^c mod n} mod m，其中 n = φ(m)欧拉总函数。

如果 m 是素数，则 n = m-1。

编辑：正如 Nabb 指出的那样，这仅在 a 与 m 互质时成立。所以你必须先检查一下。

score 19 · Accepted Answer

答案不包含正确性的完整正式数学证明。我认为这里没有必要。此外，它在 SO 上会非常难以辨认（例如，没有 MathJax）。
我将使用（只是一点点）特定的素数分解算法。这不是最好的选择，但足够了。

tl;博士

我们要计算a^x mod m。我们将使用函数modpow(a,x,m)。如下面所描述的。

如果x足够小（不是指数形式或存在p^x | m）只需计算它并返回
使用函数拆分为素数并p^x mod m为每个素数分别计算modpow
1. 计算c' = gcd(p^x,m)和t' = totient(m/c')
2. 计算w = modpow(x.base, x.exponent, t') + t'
3. 保存pow(p, w - log_p c', m) * c'在A表格中
多个 A 中的所有元素并以 m 为模返回

这里pow应该看起来像 python 的 pow。

主要问题：

因为目前的最佳答案只是关于特殊情况gcd(a,m) = 1，并且 OP 没有考虑这个假设，所以我决定写这个答案。我还将使用欧拉的全部定理。引用维基百科：

Euler's totient theorem：
如果n和a是互质正整数，则其中 φ(n) 是Euler's totient 函数。

正如 Nabb在评论中所显示的那样，这个假设numbers are co-prime非常重要。所以，首先我们需要确保这些数字是互质的。（为了更清楚假设。）因为，我们可以在其中进行因式分解，分别计算，然后以简单的方式计算答案。数最多有质因数，所以我们将被迫执行最多的计算。x = b^(c^...) a^x mod m = ((p1^alpha)^x mod m)*(p2... a = p1^α * p2^... aq1 = (p1^alpha)^x mod m,q2 = (p2^beta)^x mod m...(q1 * q2 * q3 * ... mod m)o(log a)o(log a)

事实上，我们不必拆分到每个素因数a（如果不是全部m与其他指数一起出现），我们可以与相同的指数组合，但现在并不值得注意。

现在看看(p^z)^x mod m问题，p质数在哪里。注意一些重要的观察：

如果a,b是小于m和c的正整数是某个正整数和，那么 true 是句子。

使用上述观察，我们可以得到实际问题的解决方案。我们可以很容易地计算出来gcd((p^z)^x, m)。m如果 x*z 很大，它是我们可以除以多少次的数字p。让m' = m /gcd((p^z)^x, m). （注意(p^z)^x = p^(z*x)。）让c = gcd(p^(zx),m)。现在我们可以很容易地（看下面）w = p^(zx - c) mod m'使用欧拉定理进行计算，因为这个数字是互质的！之后，使用上述观察，我们可以收到p^(zx) mod m. 根据上述假设wc mod m'c = p^(zx) mod m，所以现在的答案p^(zx) mod m = wc很w,c容易计算。

因此我们可以很容易地计算出a^x mod m。

`a^x mod m`使用欧拉定理计算

现在假设a,m是互质的。如果我们要计算a^x mod m，我们可以计算t = totient(m)并注意a^x mod m = a^(x mod t) mod m。它可能会有所帮助，如果x很大并且我们只知道的特定表达x，例如x = 7^200。

看例子x = b^c。我们可以通过平方算法及时计算t = totient(m)和x' = b^c mod t使用幂Θ(log c)。之后（使用相同的算法）a^x' mod m，等于解决方案。

如果x = b^(c^(d^...)我们将递归地解决它。首先计算t1 = totient(m)，之后t2 = totient(t1)等等。例如采取x=b^(c^d). 如果t1=totient(m), a^x mod m = a^(b^(c^d) mod t1), 并且我们能够说b^(c^d) mod t1 = b^(c^d mod t2) mod t1, 在哪里t2 = totient(t1). 我们使用平方算法求幂计算的所有内容。注意：如果某个totient不是指数的互质数，则有必要使用与主问题相同的技巧（实际上，我们应该忘记它是指数并递归解决问题，就像在主问题中一样）。在上面的例子中，如果t2不是与 c 的互质，我们必须使用这个技巧。

计算`φ(n)`

注意简单的事实：

如果gcd(a,b)=1, 那么φ(ab) = φ(a)*φ(b)
如果p是素数φ(p^k)=(p-1)*p^(k-1)

因此，我们可以分解n(ak. n = p1^k1 * p2^k2 * ...) 并使用事实 2 单独计算φ(p1^k1),φ(p2^k2),...。然后使用事实 1 将其组合。φ(n)=φ(p1^k1)*φ(p2^k2)*...

值得记住的是，如果我们要重复计算，我们可能要使用埃拉托色尼筛法并将素数保存在表中。它会减少常数。

python示例：（它是正确的，原因与此分解算法相同）

def totient(n) :          # n - unsigned int
    result = 1
    p = 2                 #prime numbers - 'iterator'
    while p**2 <= n :
        if(n%p == 0) :    # * (p-1)
            result *= (p-1)
            n /= p
        while(n%p == 0) : # * p^(k-1)
            result *=  p
            n /= p
        p += 1
    if n != 1 :
        result *= (n-1)
    return result         # in O(sqrt(n))

案子：`a`^{b^c}`mod m`

因为它实际上多次做同样的事情，我相信这个案例会告诉你如何解决这个问题。
首先，我们必须分裂a成主要力量。最佳代表将是 pair <number, exponent>。
c++11示例：

std::vector<std::tuple<unsigned, unsigned>> split(unsigned n) {
  std::vector<std::tuple<unsigned, unsigned>> result;
  for(unsigned p = 2; p*p <= n; ++p) {
    unsigned current = 0;
    while(n % p == 0) {
      current += 1;
      n /= p;
     }
    if(current != 0)
     result.emplace_back(p, current);
   }
  if(n != 1)
   result.emplace_back(n, 1);
  return result;
 }

拆分后，我们必须计算(p^z)^(b^c) mod m=p^(z*(b^c)) mod m每一对。首先我们应该检查，如果p^(z*(b^c)) | m。如果，是的，答案就是 (p^z)^(b^c)，但只有在z,b,c非常小的情况下才有可能。我相信我不必向它展示代码示例。
最后，如果p^(z*b^c) > m我们必须计算答案。首先，我们必须计算c' = gcd(m, p^(z*b^c))。之后我们就可以计算了t = totient(m')。和 (z*b^c - c' mod t)。这是获得答案的简单方法。

function modpow(p, z, b, c, m : integers) # (p^z)^(b^c) mod m
    c' = 0
    m' = m
    while m' % p == 0 :
        c' += 1
        m' /= p
    # now m' = m / gcd((p^z)^(b^c), m)
    t = totient(m')
    exponent = z*(b^c)-c' mod t
    return p^c' * (p^exponent mod m')

在Python工作示例下面：

def modpow(p, z, b, c, m) : # (p^z)^(b^c) mod m
    cp = 0
    while m % p == 0 :
        cp += 1
        m /= p              # m = m' now
    t = totient(m)
    exponent = ((pow(b,c,t)*z)%t + t - (cp%t))%t
                            # exponent = z*(b^c)-cp mod t
    return pow(p, cp)*pow(p, exponent, m)

使用这个函数，我们可以很容易地计算(p^z)^(b^c) mod m，在我们只需将所有结果 ( mod m) 相乘后，我们还可以持续计算所有内容。下面的例子。（我希望我没有犯错，写作。）唯一的假设，b，c 足够大（b^c > log(m)ak。每个p^(z*b^k)不划分m），这是简单的检查，我看不出有什么意义。

def solve(a,b,c,m) : # split and solve
    result = 1
    p = 2            # primes
    while p**2 <= a :
        z = 0
        while a % p == 0 :
                     # calculate z
            a /= p
            z += 1
        if z != 0 :
            result *=  modpow(p,z,b,c,m)
            result %= m
        p += 1
    if a != 1 :      # Possible last prime
        result *= modpow(a, 1, b, c, m)
    return result % m

看起来，像它的作品。
DEMO，它是正确的！

score 1 · Accepted Answer

因为对于任何关系a=x^y，关系对于您使用的数字基数（基数 2、基数 6、基数 16 等）都是不变的。
由于 mod N 操作相当于提取基数 N 中的最低有效位 (LSD)
由于以 N 为底的结果 A 的 LSD 只能受以 N 为底的 X 的 LSD 影响，而不受高位数字的影响。（例如 34*56 = 30*50+30*6+50*4+4*5 = 10*(3+50+3*6+5*4)+4*6）

因此，从LSD(A)=LSD(X^Y)我们可以推断

LSD(A)=LSD(LSD(X)^Y)

所以

A mod N = ((X mod N) ^ Y) mod N

和

(X ^ Y) mod N = ((X mod N) ^ Y) mod N)

因此，您可以在每个功率步骤之前进行 mod，这使您的结果保持在整数范围内。

这假设 a 不是负数，并且对于任何 x^y，a^y < MAXINT

这个答案回答了错误的问题。（亚历克斯）

score 0 · Accepted Answer

模幂是解决这个问题的正确方法，这里有一点提示：

要找到 a ^{b ^{c ^d}} % m 你必须从计算 a % m 开始，然后是^{a b} % m，然后是^{a b ^c} % m，然后a ^{b ^{c ^d}} % m ...（你明白了）

要找到a ^b % m，基本上需要两个思路：[Let B=floor(b/2)]

a ^b = (a ^B ) ²如果 b 是偶数或 a ^b = (a ^B ) ² *a 如果 b 是奇数。

(X*Y)%m = ((X%m) * (Y%m)) % m
(% = mod)

因此，
如果 b 是偶数
a ^b % m = (a ^B % m) ² % m
或者如果 b 是奇数
a ^b % m = (((a ^B % m) ² ) * (a % m)) % m

因此，如果您知道^B的值，则可以计算该值。

要找到^B，请应用类似的方法，将 B 划分直到达到 1。

例如计算 16 ¹³ % 11：

16 ¹³ % 11 = (16 % 11) ¹³ % 11 = 5 ¹³ % 11 = (5 ⁶ % 11) * (5 ⁶ % 11) * (5 % 11) <---- (I)

求 5 ⁶ % 11:
5 ⁶ % 11 = ((5 ³ % 11) * (5 ³ % 11)) % 11 <----(II)
求5 ³ %11:
5 ³ % 11 = ((5 ¹ % 11) * (5 ¹ % 11) * (5 % 11)) % 11
= (((5 * 5) % 11) * 5) % 11 = ((25 % 11) * 5) % 11 = (3 * 5) % 11 = 15 % 11 = 4
将此值代入 (II) 得到
5 ⁶% 11 = (((4 * 4) % 11) * 5) % 11 = ((16 % 11) * 5) % 11 = (5 * 5) % 11 = 25 % 11 = 3
将此值代入 (I ) 给出
5 ¹³ % 11 = ((3 % 11) * (3 % 11) * 5) % 11 = ((9 % 11) * 5) % 11 = 45 % 11 = 4

这样 5 ¹³ % 11 = 4
这样你就可以计算任何形式的^{5 ¹³} % 11 等等......

score 0 · Accepted Answer

查看A^X mod MasX增加的行为。它最终必须进入一个循环。假设循环有长度P并且在N步骤之后开始。则X >= N暗示A^X = A^(X+P) = A^(X%P + (-N)%P + N) (mod M)。因此我们可以通过计算A^B^C来计算y=B^C, z = y < N ? y : y%P + (-N)%P + N, return A^z (mod m)。

请注意，我们可以在幂树上递归地应用此策略，因为导出的方程要么有一个指数 <，要么有一个指数M，它涉及一个较小的指数塔和一个较小的红利。

唯一的问题是您是否可以有效地计算N和P给出A和M。请注意，高估N是可以的。我们可以开始N，M一切都会好起来的。P有点难。如果A和M是不同的素数，那么P=M-1。如果A有的所有M质因数，那么我们会卡在 0 和P=1。我会把它作为一个练习来解决这个问题，因为我不知道怎么做。

///Returns equivalent to list.reverse().aggregate(1, acc,item => item^acc) % M
func PowerTowerMod(Link<int> list, int M, int upperB = M)
    requires M > 0, upperB >= M
    var X = list.Item
    if list.Next == null: return X
    var P = GetPeriodSomehow(base: X, mod: M)
    var e = PowerTowerMod(list.Next, P, M)
    if e^X < upperB then return e^X //todo: rewrite e^X < upperB so it doesn't blowup for large x
    return ModPow(X, M + (e-M) % P, M)

score 0 · Accepted Answer

Tacet 的回答很好，但有可能进行大量简化。

x 的幂 mod m 是前周期的。如果 x 与 m 互质，则 x 的幂是周期性的，但即使没有这个假设，周期之前的部分也不长，最多为 m 的素数分解中的指数的最大值，最多为 log_2 m . 周期的长度除以 phi(m)，实际上是 lambda(m)，其中 lambda 是Carmichael 函数，最大乘法阶模 m。这可以明显小于 phi(m)。Lambda(m) 可以从 m 的素数分解中快速计算出来，就像 phi(m) 一样。Lambda(m) 是 lambda(p_i^e_i) 在 m 的素数分解中对所有素数幂 p_i^e_i 的 GCD，对于奇数素数幂，lambda(p_i^e_i) = phi(p_i^e^i)。lambda(2)=1, lamnda(4)=2, lambda(2^n)=2^(n-2) 对于 2 的较大幂。

定义 modPos(a,n) 为 {0,1,..,n-1} 中 a 的同余类的代表。对于非负 a，这只是 a%n。对于负数，由于某种原因，a%n 被定义为负数，因此 modPos(a,n) 为 (a%n)+n。

将 modMin(a,n,min) 定义为与至少为 min 的 mod n 一致的最小正整数。如果是肯定的，您可以将其计算为 min+modPos(a-min,n)。

如果 b^c^... 小于 log_2 m （我们可以通过递归取对数来检查这个不等式是否成立），那么我们可以简单地计算 a^b^c^... 否则，a^b^c^ ... mod m = a^modMin(b^c^..., lambda(m), [log_2 m])) mod m = a^modMin(b^c^... mod lambda(m), lambda (m),[log_2 m])。

例如，假设我们要计算 2^3^4^5 mod 100。请注意，3^4^5 只有 489 位，所以这可以通过其他方法实现，但它足够大，您不想计算它直接。但是，通过我在这里给出的方法，您可以手动计算 2^3^4^5 mod 100。

由于 3^4^5 > log_2 100，

2^3^4^5 mod 100 
= 2^modMin(3^4^5,lambda(100),6) mod 100 
= 2^modMin(3^4^5 mod lambda(100), lambda(100),6) mod 100
= 2^modMin(3^4^5 mod 20, 20,6) mod 100.

让我们计算 3^4^5 mod 20。由于 4^5 > log_2 20，

3^4^5 mod 20 
= 3^modMin(4^5,lambda(20),4) mod 20 
= 3^modMin(4^5 mod lambda(20),lambda(20),4) mod 20
= 3^modMin(4^5 mod 4, 4, 4) mod 20
= 3^modMin(0,4,4) mod 20
= 3^4 mod 20
= 81 mod 20
= 1

我们可以将其代入前面的计算：

2^3^4^5 mod 100 
= 2^modMin(3^4^5 mod 20, 20,6) mod 100
= 2^modMin(1,20,6) mod 100
= 2^21 mod 100
= 2097152 mod 100
= 52.

请注意，2^(3^4^5 mod 20) mod 100 = 2^1 mod 100 = 2，这是不正确的。你不能减少到基地权力的前周期部分。

algorithm - 找到 a^b^c^... mod m

a b c d 。. . 模数

6 回答 6

tl;博士

主要问题：

a^x mod m使用欧拉定理计算

计算φ(n)

案子：abc mod m

Related

Reference

a ^{b ^{c ^{d ^{。^{. ^.}}}}}模数

`a^x mod m`使用欧拉定理计算

计算`φ(n)`

案子：`a`^{b^c}`mod m`