快速冪取模演算法：從 O(b) 到 O(log b) 的數學推導與 C 實現

編輯：關於C語言

快速冪取模演算法：從 O(b) 到 O(log b) 的數學推導與 C 實現

計算 a^b mod c 是密碼學、數論、競賽程式的高頻操作。朴素循環相乘 O(b) 在大指數下必超時。本文詳細推導快速冪（Binary Exponentiation）的數學基礎、迭代優化過程，並給出生產級 C 實現。

一、朴素算法與溢位問題

int ans = 1;
for (int i = 1; i <= b; ++i) ans = ans * a;
ans %= c;

缺陷：

時間 O(b)，b ~ 1e9 直接 TLE
中間乘積極易溢位（即使 long long 也擋不住 a^b 指數級增長）

二、模運算核心性質：逐步取模不改結果

定理：(x * y) mod c = ((x mod c) * (y mod c)) mod c

推論：

a^b mod c
= (a * a * ... * a) mod c          (b 個 a)
= (((a mod c) * (a mod c)) mod c * ...) mod c
= ((a mod c)^b) mod c

改進版 O(b) 但無溢位：

a %= c;
for (int i = 1; i <= b; ++i)
    ans = (ans * a) % c;

三、指數二進制分解：O(log b) 核心洞見

將指數 b 展開為二進制：

b = (d_k d_{k-1} ... d_1 d_0)_2
  = d_0*2^0 + d_1*2^1 + ... + d_k*2^k  (d_i ∈ {0,1})

則：

a^b = a^{Σ d_i 2^i} = Π (a^{2^i})^{d_i}

關鍵遞推：a^{2^{i+1}} = (a^{2^i})^2

即：平方累乘法——從低位到高位，每步平方底數，當該二進制位為 1 時累乘結果。

四、迭代演算法完整推導

ans = 1
a = a % c
while (b > 0):
    if (b & 1):          // 當前最低位為 1
        ans = (ans * a) % c
    b >>= 1              // 右移，處理下一位
    a = (a * a) % c      // 底數平方，對應 2^i → 2^{i+1}
return ans

循環不變量：ans * a^b ≡ a^b_original (mod c) 恆成立。

五、C 標準庫級實現（處理邊界、型別安全）

#include <stdint.h>

// 計算 (base^exp) % mod，支援 64 位整數，避免溢位用 __int128 中間乘積
int64_t mod_pow(int64_t base, int64_t exp, int64_t mod) {
    if (mod == 1) return 0;
    int64_t result = 1;
    base %= mod;
    while (exp > 0) {
        if (exp & 1) {
            // 使用 128 位中間乘積防止溢位（GCC/Clang 支援）
            result = (int64_t)((__int128)result * base % mod);
        }
        exp >>= 1;
        if (exp)  // 避免最後一次多餘平方
            base = (int64_t)((__int128)base * base % mod);
    }
    return result;
}

// 無 __int128 環境的安全乘模（俄羅斯農夫乘法 / 二進制拆解）
int64_t mul_mod(int64_t a, int64_t b, int64_t mod) {
    int64_t res = 0;
    a %= mod; b %= mod;
    while (b) {
        if (b & 1) res = (res + a) % mod;
        a = (a + a) % mod;
        b >>= 1;
    }
    return res;
}