快速幂及快速幂取模

数学基础

快速幂

当我们求取公式 $a^p$ 时，传统的做法是计算 $a\times a\times \dots \times a$，这种方法的复杂度为 $O(p)$，在 $p$ 特别大时性能很差（例如RSA算法等）。

我们不妨转换一种思路，首先将 $p$ 转换为二进制，若 $p$ 共有 $n$ 位，$k_i$ 表示 $p$ 第 $i$ 位的值（$i=0,1,\dots ,n-1$），则可得出如下公式：

$$\begin{array}{rl}p& =\sum _{i=0,1,\dots ,n-1}{k}_i\cdot 2^i\\ & =k_0\cdot 2^0+k_1\cdot 2^1+k_2\cdot 2^2+\dots +k_{n-1}\cdot 2^{n-1}\end{array}$$

将其代入 $a^p$ 中可得：

$$\begin{array}{rl}{a}^p& =a^{k_0\cdot 2^0+k_1\cdot 2^1+k_2\cdot 2^2+\dots +k_{n-1}\cdot 2^{n-1}}\\ & =a^{k_0\cdot 2^0}\cdot a^{k_1\cdot 2^1}\cdot \dots \cdot a^{k_{n-1}\cdot 2^{n-1}}\end{array}$$

定义第 $i$ 项为 $a_i=b_i^{k_i}$，则 $b_i=b_{i-1}^2$ 其中 $i=1,2,\dots ,n-1$。

由于 $b_0=a$，而之后的每一项都可以由上一项推导出来，所以只需进行一个线性的计算，即可求出 $a^p$，其复杂度为 $O(\log n)$。

模运算

虽然快速幂算法能够在幂运算的指数很大时非常快地求出结果，但通常情况下，我们求如此大的指数幂，并不是为了得到目标结果，而是为了得到目标结果对某个数的余数，这在现代密码学中十分常见。

在实际的应用中，由于目标结果通常非常大，远远超出了计算机的存储能力，而我们需要的只是结果对于某个数的余数，，所以我们需要应用一些模运算的技巧来简化运算。

以下是模运算中最基本的一些公式，其中除法并不能简单的进行运算，详见乘法逆元：

$$\begin{array}{rl}(a+b)modp& =(amodp+bmodp)modp\\ (a-b)modp& =(amodp-bmodp+p)modp\\ (a\cdot b)modp& =(amodp)\cdot (bmodp)modp\end{array}$$

我们发现，快速幂算法中实际上只用到了乘法，所以可以利用模运算的特点将取模的运算直接添加到所有会影响最终结果的乘法计算中，则最终的计算结果正是我们需要的结果。

代码

幂运算的结果可能很大，应当注意数据类型以避免溢出！

快速幂

Typescript

/**
 * number类型支持的最大安全整数仅有2^53-1，大约为9e15
 * 在幂运算中非常容易溢出，建议使用更安全的bigint
 */
export function quickPower(x: bigint, p: bigint): bigint {
  let res: bigint = 1n;
  while (p) {
    if (p & 1n) res *= x;
    x *= x;
    p >>= 1n;
  }
  return res;
}

Python

def quick_power(x: int, p: int) -> int:
    res: int = 1
    while p:
        if p & 1:
            res *= x
        x *= x
        p >>= 1
    return res

快速幂取模

Typescript

/**
 * number类型支持的最大安全整数仅有2^53-1，大约为9e15
 * 在幂运算中非常容易溢出，建议使用更安全的bigint
 */
export function quickMod(x: bigint, p: bigint, mod: bigint): bigint {
  let res: bigint = 1n;
  while (p) {
    if (p & 1n) res = (res * x) % mod;
    x = (x * x) % mod;
    p >>= 1n;
  }
  return res;
}

Python

def quick_mod(x: int, p: int, mod: int) -> int:
    res: int = 1
    while p:
        if p & 1:
            res = (res * x) % mod
        x = (x * x) % mod
        p >>= 1
    return res