快速幂算法
目录
- 简介
- 矩阵幂
- 线性递推
- 一个ctf题目
一、简介
快速幂，也称二进制取幂，可以做到在$O(\log n)$ 时间内计算$a^n$ 。
如果用连乘的方法，将$a^n$ 写成$\underbrace{ a\cdot a\cdot\cdot\cdot a }_{ n个 }$ 的形式，则需要$O(n)$的时间来计算，在n很大时无法计算
二进制取幂的想法是，将取幂的任务按照指数的二进制表示来分割成更小的任务．
[!example] 假设要计算 $3^{13}$。如果将它展开为连乘式，需要 $13-1=12$ 次乘法。但是，因为
\[3^{13}=3^{(1101)_2}=3^8 \times 3^4 \times 3^1,\]所以，只要能快速计算出 $3^1,3^2,3^4,3^8$，就能通过 2 次乘法计算出 $3^{13}$ 的值。于是，只需要知道一个快速的方法来计算上述 3 的 $2^k$ 次幂的序列。这是容易的，因为序列中（除第一个）任意一个元素都是其前一个元素的平方。
根据这些分析，可以得到 $3^{13}$ 的计算过程如下：
\[\begin{align*} 3^1 &= 3, \\ 3^2 &= (3^1)^2 = 3^2 = 9, \\ 3^4 &= (3^2)^2 = 9^2 = 81, \\ 3^8 &= (3^4)^2 = 81^2 = 6561, \\ 3^{13} &= 6561 \times 81 \times 3 = 1594323. \end{align*}\]过程中，只进行了 5 次乘法运算。
具体的实现，有递归和非递归两种形式
迭代实现
设 $n$ 的二进制表示为 $(n_tn_{t-1}\cdots n_1n_0)_2$，也就是说，有
\[n = n_t 2^t + n_{t-1} 2^{t-1} + \cdots + n_1 2^1 + n_0 2^0,\]其中，$n_i \in {0,1}$。那么，就有
\[\begin{align*} a^n &= a^{n_t 2^t + n_{t-1} 2^{t-1} + \cdots + n_1 2^1 + n_0 2^0} \\ &= a^{n_0 2^0} \times a^{n_1 2^1} \times \cdots \times a^{n_{t-1} 2^{t-1}} \times a^{n_t 2^t}. \end{align*}\]注意，只有 $n_i = 1$ 的项才会真正出现在乘积的计算中。
根据这一表达式，可以首先在 $O(\log n)$ 时间内计算出 $a$ 的 $O(\log n)$ 个 $2^k$ 次幂的取值，然后花费 $O(\log n)$ 的时间选择等于 1 的二进制位对应的幂次乘到最终结果中。
1
2
3
4
5
6
7
8
def binpow(a, b, p):
res = 1
while b > 0:
if b & 1:
res = res * a % p
a = a * a % p
b >>= 1
return res
这就是快速幂的迭代版本实现。
递归实现
这一过程同样可以通过递归形式实现。注意到，指数 $n$ 的二进制展开可以递归地写作
\[(n_tn_{t-1}\cdots n_1n_0)_2 = 2 \times (n_tn_{t-1}\cdots n_1)_2 + n_0.\]因此，幂次 $a^n$ 可以递归地计算为
\[a^n = \begin{cases} 1, & n = 0, \\ \left(a^{\lfloor n/2 \rfloor}\right)^2, & n > 0 \text{ and } n \text{ is even}, \\ \left(a^{\lfloor n/2 \rfloor}\right)^2 \cdot a, & n > 0 \text{ and } n \text{ is odd}. \end{cases}\]1
2
3
4
5
6
7
8
def binpow(a, b, p):
if b == 0:
return 1
res = binpow(a, b // 2, p)
if (b % 2) == 1:
return res * res * a % p
else:
return res * res % p
这就是快速幂的递归版本实现。
[!note] 在模意义下取幂，只要在快速幂过程中取模数就行了
二、矩阵幂
实际上，任何满足乘法结合律的幂运算都可以使用快速幂进行高效计算，比如矩阵的幂。
同理，利用二进制拆分指数，每次平方底数矩阵并在指数为奇数时累乘结果矩阵。代码实现如下：
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
def mat_mul(A, B, mod=None):
n = len(A)
C = [[0] * n for _ in range(n)]
for i in range(n):
for k in range(n):
if A[i][k] == 0:
continue
aik = A[i][k]
for j in range(n):
C[i][j] += aik * B[k][j]
if mod is not None:
C[i][j] %= mod
return C
def mat_pow(mat, power, mod=None):
n = len(mat)
res = [[0] * n for _ in range(n)]
for i in range(n):
res[i][i] = 1
base = mat
while power > 0:
if power & 1:
res = mat_mul(res, base, mod)
base = mat_mul(base, base, mod)
power >>= 1
return res
时间复杂度为$O(m^3\log n)$ ,m为维数
三、应用：解决线性递推方程
我们先看最简单的线性递推：Fibonacci
\[F_{1}=F_{2}=1,F_{n}=F_{n-1}+F_{n-2},n\geq3\]显然递推公式能写成矩阵形式
\[(F_{n}\ \ F_{n-1})=(F_{n-1} \ \ F_{n-2})\begin{pmatrix} 1\ \ 1 \\ 1\ \ 0 \end{pmatrix}\] \[(F_{n}\ \ F_{n-1})=(1,1) \begin{pmatrix} 1\ \ 1 \\ 1\ \ 0 \end{pmatrix} ^{n-2}\]此处的矩阵幂就可以使用快速幂了。
更一般的，对于递推
\[f_{n}=a_{1}f_{n-1}+a_{2}f_{n-2}+\dots+a_{k}f_{n-k}\]可以构造矩阵
\[A=\begin{pmatrix} 0 & 1 & 0 & \dots & 0 \\ 0 & 0 & 1 & \dots & 0 \\ \dots \\ 0 & 0 & 0 & \dots & 1 \\ a_{k} & a_{k-1} & a_{k-2} & \dots & a_{1} \\ \end{pmatrix}\]做一次乘法，得到
\[A \begin{pmatrix} f_{i-k} \\ f_{i-k+1} \\ f_{i-k+2} \\ \dots \\ f_{i-1} \end{pmatrix}= \begin{pmatrix} f_{i-k+1} \\ f_{i-k+2} \\ f_{i-k+3} \\ \dots \\ f_{i} \end{pmatrix}\]累乘得
\[\begin{pmatrix} f_{n-k+1} \\ f_{n-k+2} \\ f_{n-k+3} \\ \dots \\ f_{n} \end{pmatrix}=A^{n-k} \begin{pmatrix} f_{1} \\ f_{2} \\ f_{3} \\ \dots \\ f_{k} \end{pmatrix}\]四、TJCTF 2026
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
from pathlib import Path
from Crypto.Util.number import long_to_bytes
M = 0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f
def encrypt(val, key):
amp = val
p = M
prev = 1
current = amp
for _ in range(key - 1):
new_val = (2 * amp * current - prev) % p
prev = current
current = new_val
return current
flag_bytes = Path("flag.enc").read_bytes()
base_val = 0x1337C0DE
frequency_key = 10**25
secret = encrypt(base_val, frequency_key)
encrypted = bytes([a ^ b for a, b in zip(flag_bytes, long_to_bytes(secret))])
print(encrypted)
可以看到，解密脚本就是一个线性递推
\[T_0(x) = 1,\quad T_1(x) = x,\quad T_{n+1}(x) = 2x\,T_n(x) - T_{n-1}(x)\]这个实际上是著名的切比雪夫多项式。
题目中，解密脚本没能得到结果。主要的问题在于，使用了复杂度$O(n)$的纯递推算法。
我们构造矩阵：
\[A=\begin{pmatrix} 0& 1 \\ -1& 2x \end{pmatrix}\] \[\begin{pmatrix} T_{n-1}(x) \\ T_{n}(x) \\ \end{pmatrix}=A^{n-2} \begin{pmatrix} T_{1}(x) \\ T_{2}(x) \\ \end{pmatrix}\]转换成快速幂就算的很顺利了。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
from pathlib import Path
from Crypto.Util.number import long_to_bytes
M = 0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f
def encrypt(n, x, mod):
if n == 0:
return 1
x %= mod
a, b = 1, x
bits = bin(n)[2:]
for bit in bits:
a2 = (2 * a * a - 1) % mod
b2 = (2 * a * b - x) % mod
if bit == '0':
a, b = a2, b2
else:
a = b2
b = (2 * x * b2 - a2) % mod
return a
flag_bytes = Path("flag.enc").read_bytes()
base_val = 0x1337C0DE
frequency_key = 10**25
secret = encrypt(frequency_key,base_val, M)
encrypted = bytes([a ^ b for a, b in zip(flag_bytes, long_to_bytes(secret))])
print(encrypted)
参考：