【面试相关】python实现快速幂取余算法详解

假设我们要计算 2102^{10}210 对1000取模的结果，可以很简单的得到24。但是如果要求 210002^{1000}21000 对1000取模的结果，常规方法就行不通了，因为常规的变量无法容纳这么大的数值。为此，需要借助数学技巧求解。

循环求余法

首先引入一个好理解的取模运算公式：
(a×b)%p=(a%p×b%p)%p(a \times b) \% p = (a \% p \times b \% p) \% p (a×b)%p=(a%p×b%p)%p
举个例子：
a=9,b=9(9×9)%6=(9%6×9%6)%6=（3×3)%6=3a=9, b=9\\ (9 \times 9) \% 6 = (9\%6 \times 9\%6)\%6 = （3 \times 3) \% 6 = 3 a=9,b=9(9×9)%6=(9%6×9%6)%6=（3×3)%6=3
受上面计算过程的启发，不难得到下面的取模运算公式：
xa%p=((x%p)a)%px^a \% p = \left( (x \% p)^a \right) \% p xa%p=((x%p)a)%p
即 a 个 x 的乘积可以分解为 a 个 x 对 p 余数的乘积再对 p 取余，即 ((x%p)a)%p\left( (x \% p)^a \right) \% p((x%p)a)%p 。

举个例子：
x=9,a=2,p=6(92)%6=(9%6)2%6=32%6=3x=9, a=2, p=6\\ (9^2) \% 6 =(9\%6 )^2\%6 = 3 ^2 \% 6 = 3 x=9,a=2,p=6(92)%6=(9%6)2%6=32%6=3
对于上面的例子，如果 a = 1000，难道最后还得计算 31000%63^{1000} \% 631000%6 吗？其实可以这么做：

在进行 a 次迭代计算 x 对 p 余数的过程中，先将第 iii 次迭代得到的余数用变量 rem 保存起来，在进行第 i+1i + 1i+1 次迭代的时候，由 (rem * (x%p)) % p 所得到的结果就是 xi+1%px^{i+1} \%pxi+1%p 的结果。

举个例子，求 x3%px^3 \% px3%p 的过程如下：

1，rem = x % p
2，rem = x^2 % p = ((x%p) * (x%p)) % p = (rem * (x%p)) % p
3，rem = x^3 % p = ((x^2%p) * (x%p)) % p = (rem * (x%p)) % p

可以看到这个过程的时间复杂度是 O(N)，其中 N = a，代码实现如下：

# 求 (x^a) % p —— 循环求余法
def reminder1(x, a, p):rem = 1for _ in range(a):rem = (rem * (x%p)) % preturn remx, a, p = 2,100000000, 3
print(reminder1(x, a, p)) # 用时：9.06s

观察(rem * (x%p)) % p，发现既然括号的外面要用 p 取模，那么内部的(x%p)就显得没有必要了，因此可以简化为：(rem * x) % p ：

# 求 (x^a) % p —— 循环求余法
def remainder2(x, a, p):rem = 1for _ in range(a):rem = (rem * x) % preturn remx, a, p = 2,100000000, 3
print(reminder2(x, a, p)) # 用时：6.76s

快速幂求余

上面的过程中，我们每一次迭代都是计算 x%px\% px%p，所以需要迭代 a 次。其实这个过程还有优化的空间：

如果我们每次都计算 x2%px^2 \% px2%p，那只需要迭代 a/2 次。也就是说，可以把指数 a 级的问题转换为 a/2 级的问题，如果继续分解这个计算过程，又可以得到：
xa%p=（x2%p）a/2%p=((x2%p)2%p)a/4%p=...x^a \% p = （x^2 \% p）^{a/2} \% p = \left( (x^2 \% p)^2\%p \right)^{a/4} \% p = ... xa%p=（x2%p）a/2%p=((x2%p)2%p)a/4%p=...
最终计算 xa%px^a \%pxa%p 的时间复杂度仅为 O(logN)O(logN)O(logN) 。

实际上，上面的式子忽略了一个问题， a 可能是奇数也可能是偶数，因此对 a/2 的处理上需要遵循下面的公式：
xa%p={(x2%p)a//2%p,a为偶数[(x%p)(xa−1%p)]%p=[x(x2%p)a//2]%p,a为奇数x^{a} \% p=\left\{\begin{array}{l} \left(x^{2} \% p\right)^{a / / 2} \% p \qquad \qquad \qquad \qquad \qquad \qquad \qquad, a 为偶数\\ {\left[(x \% p)\left(x^{a-1} \% p\right)\right] \% p=\left[x\left(x^{2} \% p\right)^{a / / 2}\right] \% p} \quad\; ,a为奇数 \end{array}\right. xa%p={(x2%p)a//2%p,a为偶数[(x%p)(xa−1%p)]%p=[x(x2%p)a//2]%p,a为奇数
式子中的 // 表示除法运算的结果向下取整。此计算过程的代码实现如下：

# 求 (x^a) % p —— 快速幂求余
def remainder3(x, a, p):rem = 1while a > 0:# 无论 a 是奇数或偶数，a//2的结果是一样的# 因此，对于奇数要进行特殊处理if a % 2: rem = (rem * x) % p# 底数变大为原来的平方：x = x ** 2 % p# 指数变为原来的一半a //= 2return rem
x, a, p = 2,100000000, 3
print(remainder3(x, a, p)) # 几乎等于0

为什么 rem 只在 a 为奇数的时候才计算呢？

假设 a = 4，那么

第一次迭代后 a=2， x = x^2 % p；
第二次迭代后 a=1，x = (x^2 % p)^2 %p；
第三次迭代时，可以求得 rem = (1 * x) % p，且 a = 0;
结束循环，返回了 rem。

假设 a = 5，那么：

第一次迭代后 rem = x%p，a=2， x = x^2 % p；
第二次迭代后 a=1，x = (x^2 % p)^2 %p；
第三次迭代时，可以求得 rem = (rem * x) % p，且 a = 0;
结束循环，返回了 rem。

也就是说无论 a 是奇数还是偶数， a//2的计算结果都是交替出现的奇数和偶数，并且最终都会得到 a = 1 计算出 rem，然后结束循环，返回 rem。因此只需要在 a 为奇数的时候计算 rem 即可保证最后得到的 rem 就是我们想要的值。

参考文章：

面试题14- II. 剪绳子 II（数学推导 / 贪心思想 + 快速幂求余，清晰图解）

快速幂算法（全网最详细地带你从零开始一步一步优化）