企业官网建设流程全解析

快速幂这个名字，我每次跟人解释的时候都习惯说一句话：它不是一个奇技淫巧，而是把“连乘”这件事从线性思维里解放出来的一种底层套路。整数快速幂、矩阵快速幂、Strassen矩阵乘法，本质上都在回答同一个问题——能不能少算几次乘法？今天这篇就把整数快速幂和Strassen放一起写，因为它们的思路是连通的。

如果你写过求 a^n 的朴素循环，或者被 n=10^18 的模幂运算卡到怀疑人生；如果你在某个神经网络或者数值计算的Demo里发现矩阵乘法成了性能瓶颈，那么这篇文章值得你看完。我尽量把每一步推导都写透，顺带把我在实际编码里踩过的坑也抖出来。

1. 整数快速幂：从朴素循环到倍增压缩

1.1 朴素幂运算的瓶颈到底在哪

先看最直观的写法。求 a^n，循环乘 n 次：

long long power_naive(long long a, long long n) { long long result = 1; for (long long i = 0; i < n; ++i) result *= a; return result; }

这段代码的逻辑没有任何问题，问题在复杂度。假设底数 a 和指数 n 都是很大的数，尤其当题目要求“结果对某个质数取模”，n 可以轻易到 10^18。这个量级下，循环 10^18 次乘法，普通机器上是要跑到天荒地老的。

有人可能会说：那我把循环展开、用位运算加速不就行了？没用。O(n) 的时间复杂度不改变，循环再优化也只是把常数缩小，指数一大照样崩。这里的核心矛盾是：我们不需要把 a 连乘 n 次，因为连乘过程里绝大多数中间结果是后面用不到的。

举个很俗的例子。让你算 2 的 100 次方，你不会真的在纸上连续乘 100 次。你会先算出 2^10=1024，然后用 1024^2 推出 2^20，再用 2^20 平方得到 2^40，这样指数翻着倍往上涨。快速幂的灵感就来自这个朴素的观察。

1.2 二进制拆分与倍增原理

快速幂的核心操作叫“二进制拆分”。把指数 n 写成二进制形式，比如 n=13，二进制是 1101，即：

13 = 8 + 4 + 1 = 2^3 + 2^2 + 2^0

那么 a^13 可以拆成：

a^13 = a^8 × a^4 × a^1

你注意到没有，原本需要 13 次乘法，现在只要把 a^1、a^2、a^4、a^8 逐个算出来，然后把二进制位里是 1 的那些挑出来乘在一起。a^1 已知，a^2 是 a 的平方，a^4 是 a^2 的平方，a^8 是 a^4 的平方——每一步只用一次乘法，指数翻倍。这就是“倍增”的含义。

整个过程需要的乘法次数，大约等于 n 的二进制位数加上 1 的个数。13 的二进制有 4 位，其中 3 个 1，所以大约 4 次平方操作加上 2 次结果累乘，总共不超过 6 次乘法。相比原始的 13 次，已经不是一个量级。

关键公式写出来就是：

result = 1 base = a while n > 0: if n & 1 == 1: result *= base base *= base n >>= 1

每次循环里做两件事：检查当前二进制最低位是不是 1，是就把当前 base 乘进结果；然后把 base 平方，相当于指数一倍一倍地翻。最后 n 右移一位，处理下一个二进制位。

1.3 带模运算的快速幂实现

实际比赛和工程里，很少直接算很大的幂，因为数字早就爆了。绝大多数场景是算模幂，即 (a^n) % m。取模运算对乘法有很好的性质：乘积的模等于模的乘积再取模，所以计算过程中每一步都能取模，保证中间数值不会溢出。

我给出一版我常用的 C++ 实现：

typedef long long ll; ll mod_pow(ll a, ll n, ll m) { ll result = 1; a %= m; while (n > 0) { if (n & 1) result = (result * a) % m; a = (a * a) % m; n >>= 1; } return result; }

这里有个特别容易踩的坑：如果 m 本身接近 LLONG_MAX，那么 (result * a) 这一步可能在取模之前就溢出了。C++ 乘法溢出是未定义行为，经常导致结果莫名其妙。解决方式有几种：

• 如果编译器支持__int128，直接用__int128做中间乘法。
• 如果只有 64 位环境，可以用“快速乘”的加法拆乘法，把 a*b 拆成类似快速幂的形式，用加法+取模完成。
• 必要时引入高精度库。

实际工作中我建议直接上__int128，性能损失可接受，但代码安全系数高很多。

再给一个 Python 版本，Python 自带大整数，写起来最直接：

def mod_pow(a, n, m): result = 1 a %= m while n > 0: if n & 1: result = (result * a) % m a = (a * a) % m n >>= 1 return result

1.4 二进制位遍历的直观演示

为了让你彻底吃透，我手写一遍 n=13 的完整过程。设 a=3，m 先不管（直接算）：

初始：result=1，base=3，n=13（二进制 1101）

• 第一步：n 的最低位是 1，result = 1×3 = 3；base = 3×3 = 9；n 右移变成 6。
• 第二步：n=6（二进制 110），最低位是 0，result 不变；base = 9×9 = 81；n 右移变成 3。
• 第三步：n=3（二进制 11），最低位是 1，result = 3×81 = 243；base = 81×81 = 6561；n 右移变成 1。
• 第四步：n=1，最低位是 1，result = 243×6561；base = 6561×6561；n 右移变成 0，循环结束。

而 3^13 = 1594323，上面分步取到 result 的时候正好是 3×81×6561 = 1594323，验证通过。

你发现规律没有：平方操作始终在"无条件"地执行，它不管当前二进制位是否为 1；只有二进制位为 1 时才把 base 乘进 result。这就是为什么算法能稳定压缩复杂度——无论指数多大，都只用 O(log n) 次循环，每次循环只有平方和条件累乘。

1.5 从单个数值到递推关系

整数快速幂最经典的延伸是“矩阵快速幂”。比如计算斐波那契数列第 n 项，如果用递推 f(n) = f(n-1) + f(n-2)，线性时间在 n 到 10^18 时依然不可用。但可以把递推写成矩阵形式：

| f(n+1) | | 1 1 |^n | f(1) | | f(n) | = | 1 0 | × | f(0) |

求矩阵的 n 次幂，直接用“整数快速幂”的思路，只不过底数从整数换成矩阵，乘法换成矩阵乘法。这块后面讲矩阵快速幂的时候再细说。

2. Strassen矩阵乘法：把 O(n³) 拉进 2.81 次方的世界

2.1 矩阵乘法的痛点和朴素分治思路

矩阵乘法是很多数值计算的核心操作。三层循环直写的复杂度是 O(n³)，如果有两个 1024×1024 的矩阵相乘，光内层运算就到 10 亿次乘加，再跑个图像处理或科学模拟，性能立刻成为瓶颈。

朴素的优化思路是分治。把 A、B 各分成四个 n/2 大小的子块：

A = | A11 A12 | | A21 A22 | B = | B11 B12 | | B21 B22 |

那么 C = A×B 可以表示为：

C11 = A11×B11 + A12×B21 C12 = A11×B12 + A12×B22 C21 = A21×B11 + A22×B21 C22 = A21×B12 + A22×B22

这样需要 8 次 n/2 规模的矩阵乘法。递推式 T(n) = 8T(n/2) + O(n²)，算出来的复杂度还是 O(n³)，并没有本质改善。纯粹的大块拆小块，只是把乘法和加法换了个顺序，工作量一点没少。

那么问题来了：能不能减少子矩阵乘法的次数？Strassen 的回答是：能，用加法次数换乘法次数。

2.2 Strassen 的 7 次乘法推导

Strassen 的核心是构造 7 个矩阵 M1 到 M7，每个 M 只需要一次子矩阵乘法，然后通过加减法组合得到 C 的四个子块。我先给出定义：

M1 = (A11 + A22)(B11 + B22) M2 = (A21 + A22)B11 M3 = A11(B12 - B22) M4 = A22(B21 - B11) M5 = (A11 + A12)B22 M6 = (A21 - A11)(B11 + B12) M7 = (A12 - A22)(B21 + B22)

然后 C 的四个子块用这些 M 组合出来：

C11 = M1 + M4 - M5 + M7 C12 = M3 + M5 C21 = M2 + M4 C22 = M1 - M2 + M3 + M6

乍一看这组式子像魔术，但逐项展开就能验证。以 C12 为例，按定义是 A11×B12 + A12×B22。代人 M3 + M5：

M3 + M5 = A11(B12 - B22) + (A11 + A12)B22 = A11B12 - A11B22 + A11B22 + A12B22 = A11B12 + A12B22

完美吻合。C11 的展开稍微长一点，但每一项加减最后正好约到 A11B11 + A12B21。我第一次推这个推导的时候，最大的感受是：这不是什么高深的魔法，它就是很精细的“乘法次数压缩计划”，把原本 8 次矩阵乘法合并成 7 次，代价是增加了大量矩阵加减法。

2.2.1 C11 的详细展开

为了让你信服，我把 C11 展开。C11 = M1 + M4 - M5 + M7：

M1 = (A11 + A22)(B11 + B22) = A11B11 + A11B22 + A22B11 + A22B22 M4 = A22(B21 - B11) = A22B21 - A22B11 -M5 = -(A11 + A12)B22 = -A11B22 - A12B22 M7 = (A12 - A22)(B21 + B22) = A12B21 + A12B22 - A22B21 - A22B22

把四式相加，逐项约掉：A11B22 与 -A11B22 抵消，A22B11 与 -A22B11 抵消，A22B22 与 -A22B22 抵消，A22B21 与 -A22B21 抵消，-A12B22 与 +A12B22 抵消，最后只剩下 A11B11 + A12B21。刚好就是 C11 的定义。

C22 的展开也类似，我这里不全部抄完，建议你拿笔自己推一遍。推过的印象深度远高于只看文章。

2.3 复杂度分析：O(n^log2(7)) 是怎么来的

递归式变成了：

T(n) = 7T(n/2) + O(n²)

用主定理，a=7, b=2，f(n)=n²。比较 n^(log2(7)) 与 n²：log2(7)≈2.807，大于 2，所以主定理第三种情况，复杂度为：

T(n) = O(n^log2(7)) ≈ O(n^2.807)

从 3 降到 2.807，看起来幅度不算大，但 n 增大时差距是指数级的。n=1024 时，朴素算法约 10 亿次乘法，Strassen 约 2.807 次方量级，大约是 2.4 亿次乘法，理论省下 70% 以上。

但我要泼一盆冷水：这只是渐进复杂度。Strassen 的常数因子非常大，因为它额外引入 18 次矩阵加减法（每个子矩阵，加法和减法一次也是一层 O(n²) 的遍历）。实际工程里，小矩阵上 Strassen 往往还没朴素算法快。所以正规实现都会设一个阈值，比如 n < 64 或 n < 128 时直接走朴素三层循环，超过阈值才递归调用 Strassen。

2.4 Strassen 的工程实现细节

真正把 Strassen 落地到代码，有几个绕不开的问题。

2.4.1 矩阵维度不是 2 的幂

Strassen 分治要求每次都能对半分块，维度最好是 2 的幂。如果矩阵不是 2 的幂，最常见的做法是补零扩展到最近的 2 的幂。扩展出来的区域全是 0，乘积不受影响，算完后取回左上角的部分。注意补零之后 n 变大，且每层都要传子矩阵的起始坐标和尺寸，否则反复拷贝矩阵会把性能吃光。

用一个结构体加起始偏移量最方便：

struct Matrix { vector<vector<double>> data; int row_offset, col_offset, size; double& at(int i, int j) { return data[row_offset + i][col_offset + j]; } };

递归时创建四个子视图，不拷贝数据，只是调整偏移量。这一步做好了，内存开销和性能才有保证。

2.4.2 浮点数误差问题

矩阵加减法会传递误差，Strassen 由于中间结果包含大量减法，尤其是 A12 - A22、B12 - B22 这类可能出现接近 0 的情况，数值稳定性比朴素算法差。如果你在做敏感数值计算，建议用高精度浮点（比如long double或__float128）或者只用 Strassen 做加速、把误差控制作为验收指标。

我实际见过一个案例：某图像处理 Demo 里，Strassen 加速后 PS 指标（峰值信噪比）略微下降，根源就是浮点舍入误差被放大。后来把阈值调高，让中小矩阵继续用朴素算法，只在超大矩阵上走 Strassen，误差才回到可接受范围。

2.4.3 缓存友好性

朴素三层循环如果顺序不当，内存访问会非常不友好。i-j-k 循环顺序比 i-k-j 快的原因在于缓存局部性。Strassen 的递归分块天然有一定局部性，但加减法操作需要频繁遍步子矩阵，如果数据存储是行主序且不连续，性能依然受影响。

工程上更细致的选择是按块存储，比如每块内部用行主序，块的元信息存开头。这些属于性能调优的深水区，初学阶段不用强求。

3. 统一视角：快速幂、Strassen 与矩阵快速幂

3.1 从整数到矩阵：同一个倍增内核

如果只把整数快速幂和 Strassen 当两个孤立算法背下来，收获很小。它们的共同点是：乘法成本远高于加减法成本时，就用加法次数去换乘法次数。

整数快速幂的“平方复用”：每次平方把指数扩大一倍，后续结果的乘法次数对数化。Strassen 的“加法替换”：用一堆加减法把 8 次子矩阵乘法压缩成 7 次。二者都是“乘法少一点，加法多一点，整体成本下降”的思路。

矩阵快速幂则直接把两种思想结合——底数是矩阵，幂次是整数，矩阵乘法内部可用 Strassen 加速。理论上是双料优化，但实际工程里我要提醒你：矩阵快速幂的底数通常是 2×2 或 3×3 的小矩阵，比如斐波那契递推，这时候用 Strassen 是画蛇添足，朴素矩阵乘法已经是最优了。Strassen 的价值场景是单次矩阵乘法里 n 本身很大，而不是矩阵数量很多。

3.2 矩阵快速幂的代码骨架

给一个 2×2 矩阵快速幂求斐波那契的例子：

struct Mat2 { long long a, b, c, d; }; Mat2 mul(const Mat2& x, const Mat2& y, long long mod) { return { (x.a * y.a + x.b * y.c) % mod, (x.a * y.b + x.b * y.d) % mod, (x.c * y.a + x.d * y.c) % mod, (x.c * y.b + x.d * y.d) % mod }; } Mat2 mat_pow(Mat2 m, long long n, long long mod) { Mat2 res = {1, 0, 0, 1}; // 单位矩阵 while (n > 0) { if (n & 1) res = mul(res, m, mod); m = mul(m, m, mod); n >>= 1; } return res; }

注意单位矩阵在矩阵快速幂里对应整数快速幂里的 result=1。每次平方针对的是矩阵自身，累乘结果针对的是结果矩阵。这个代码框架可以推广到任意 k×k 矩阵，只需要把mul函数换成通用矩阵乘法。

3.3 什么时候用、什么时候不用

快速幂和 Strassen 都有明显适用边界。快速幂几乎无脑值得用，因为它的实现简单，常数因子小，不管 n 是 13 还是 10^18，都是 O(log n) 量级。Strassen 则适合 n 非常大且乘法代价占主导的场景，比如某些支持大矩阵乘法的科研计算、图像处理或者机器学习相关实验。

如果矩阵太小，Strassen 的加法开销反而拖后腿。业界通常把阈值设在 64~128 左右。我自己习惯在 128 以下用朴素三层循环加适量循环展开，128 以上才走 Strassen 分支。阈值可以通过对目标机器做基准测试来确定：分别在 32、64、128、256 的矩阵规模下各测 100 轮，看哪条分支更占优，然后设定切换点。

4. 常见问题与排查技巧实录

这部分整理我在写快速幂和 Strassen 时真实踩过的坑，按问题类型列成速查表。

快速幂最常见的“阴间 bug”是位运算优先级问题。比如把if (n & 1 == 1)写在 C++ 里，==优先级高于&，会先比较再位运算，导致判断永远为假。我建议显式加括号：if ((n & 1) == 1)，避免大脑短路。

还有一个经验之谈：测试快速幂的时候，最好准备几个固定的小数据手算对拍。比如算 2^13、3^5、5^0（结果应该是 1，很多朴素实现会忘记处理指数为 0 的情况）。模幂测试时，随便选两个质数手动算一组小样本比对，别用大样本直接看“是否非零”来判断对错。

Strassen 的排查技巧，我个人最推荐“断点分块验证”。在递归入口把参数 n、四个子矩阵的行列范围输出来，和手写的分块对照。等确认分块正确了，再把 M1~M7 分别打印出来，先验证每个 M 的计算公式对不对，再验证 C 的组合对不对。这一步代码虽然麻烦，但远比整个矩阵算完再对拍容易定位问题。

最后说一个优化小技巧：矩阵加减法在 Strassen 里占据了大量时间，可以考虑用循环展开 + 并行化（OpenMP 或 SIMD）加速。如果你公司或实验室的多核机器比较富余，加减法可以直接并行，反正子块之间互相独立。这个优化空间往往比乘法本身更可观。

写在最后：我自己的一点体会

整数快速幂和 Strassen，我第一次学的时候都背过代码，但真正“通了”是在分别手动推了一遍之后——13 的二进制展开推了一遍，M1 到 M7 的展开推了一遍。二进制展开让我理解为什么倍增能压缩循环；M 的展开让我理解为什么加法能换乘法。很多算法不是缺记忆，是缺一次亲手验证。

如果你今天只带走一个东西，记住这句：乘法贵，加法便宜，能用加法换乘法就是赢。快速幂在整数上赢了一次，Strassen 在矩阵上赢了一次。下次遇到一个性能瓶颈，先问问自己——哪里在重复乘法？能不能用加法消掉一些？这种思考方式比背下任何单一算法都值钱。