快速幂详解（幂运算与矩阵）
基础知识
矩阵乘法
一张图说明足矣：
代码实现（C）
const int N = 100; int c[N][N]; void multi(int a[][N], int b[][N], int n) { memset(c, 0, sizeof c); for(int i = 1; i <= n; i++) { for(int j = 1; j <= n; j++) { for(int k = 1; k <= n; k++) c[i][j] += a[i][k] * b[k][j]; } } return; }
另一种实现，将二层循环与三层循环对调，从而可以在第二层for循环时先判断 if(a[i][k] == 0) continue; ，对于矩阵有较多0的有一定效果。
int c[N][N]; void multi(int a[][N], int b[][N], int n) { memset(c, 0, sizeof c); for(int i = 1;i <= n; i++) { for(int k = 1; k <= n; k++) { if(a[i][k] == 0) continue; for(int j = 1; j <= n; j++) { c[i][j] += a[i][k] * b[k][j]; } } } return; }
显然矩阵乘法的复杂度是$O(n^3)$。
矩阵乘法需要A的行数与B的列数相等的（这是A*B的前提条件，可见矩阵乘法是不满足交换律的）。
矩阵快速幂只会用到方阵（除非题目是裸的矩阵乘法），矩阵快速幂都是方阵也就避免了检验相乘的前提条件，可以放心用。
快速幂
经典例题：数A的幂
我们利用分治的算法思想可以考虑如下求解一个数A的幂：
例如：X的19次方。
19的二进制为1 0 0 1 1。由$X^m * X^n = X^{m+n}$，则$X^{19} = X^{16} * X^2 * X^1$。
代码实现（C++）
double QuickPow(double x, int N) { if (N == 0) return 1; if (x == 0) return 0; if (N < 0) { x = 1 / x; N = -N; } double ans = 1.0; while(N) { if(N & 1) { ans *= x; } x *= x; N >>= 1; } return ans; }
过程解析
对于$X^19$来说，初始： ans = 1; res = x; 则
10011最后一位是1，所以是奇数：
ans = res * ans = x; res = res * res = x^2;
然后右移一位，1 0 0 1，则1001最后一位是1，所以是奇数：
ans = res * ans = x * (x^2) = x^3; res = res * res = x^2 * x^2 = x^4;
然后右移一位，1 0 0，则最后一位是0，所以当前的数为偶数：
res = res * res = x^4 * x^4 = x^8;
然后右移一位，1 0，最后一位是0，当前数是偶数：
res = res * res = x^8 * x^8 = x^16;
然后右移一位，1，最后一位是1，当前数是奇数：
ans = ans * res = (x^3) *（x^16） = x^19; res = res * res = x^32;
可以看出res = X^m，m始终是与二进制位置上的权值是相对应的。当二进制位为0时，我们只让res幂指数为2，对应下一个二进制位的权值；当二进制位为1时， ans = ans * res; ，则乘上了该乘的X幂次。
矩阵快速幂
上面介绍了整型快速幂的实现思路，而矩阵乘法可以用来优化很多线性递推式，那么将其中的乘法换成矩阵乘法，就变成了矩阵快速幂。
代码实现（C）
struct Matrix { //定义一个矩阵类型的结构体 int m[maxN][maxN]; } ans, res; //实现矩阵乘法A * B， n为方阵的阶 Matrix Mul(Matrix A, Matrix B, int n) { Matrix tmp; memset(tmp.m, 0, sizeof(tmp.m)); for (int i = 0; i < n; i++) { for (int j = 0; j < n; j++) { for (int k = 0; k < n; k++) tmp.m[i][j] = A.m[i][k] + B.m[k][j]; } } return tmp; } //快速幂算法，求矩阵res的N次幂 void QuickPower(int N, int n) { //整数快速幂ans初始化为1，对于矩阵乘法来说，则ans初始化为单位矩阵 //有A * E = A memset(ans.m, 0, sizeof(ans.m)); for (int i = 0; i < n; i++) ans.m[i][i] = 1; while (N) { if (N & 1) ans = Mul(ans, res); res = Mul(res, res); N = N >> 1; } }
经典例题：斐波拉契数列
众所周知：斐波那契数列递推公式为：$F[n] = F[n-1] + F[n-2]$， 由$F[0] = 0$，$F[1] = 1$，可以递推后面的所有数。
使用for循环是最简单直接的解法，然而对于像POJ 3070这样的题目，求斐波那契数列，其n可高达10亿，直接递推有其局限性：
（1）当n过大，递推次数过多，而测试时间有限时，可能会导致超时。
（2）想要打表实现随机访问根本不可能，先把斐波那契数列求到10亿，然后想去进行随机访问。题目未给出那么多内存，数组也开不到10亿。
因此另一个更快的解法是矩阵快速幂。
建立矩阵递推式，找到转移矩阵
显然，A就是转移矩阵。
构建矩阵递推的思路
$$A * F(n-1) = F(n)$$
一般$F(n)$与$F(n-1)$都是按照原始递推式来构建的，当然可以先猜一个$F(n)$，主要是利用矩阵乘法凑出矩阵$A$，第一行一般就是递推式，后面的行就是不需要的项就让与其相乘的系数为0。矩阵$A$就叫做转移矩阵（一定是常数矩阵），它能把$F(n-1)$转移到$F(n)$，然后这就是个等比数列，直接写出通项：$F(n) = A^{n-1} * F(1)$，此处$F1$叫初始矩阵。所以用一下矩阵快速幂然后乘上初始矩阵就能得到$F(n)$。例如，这里的$F(n)$就两个元素（两个位置），根据自己设置的$F(n)$对应位置就是对应的值。
示例
1. $f(n) = a * f(n-1) + b * f(n-2) + c$（a, b, c是常数）
$$\left( \begin{array}{ccc} a & b & 1 \\ 1 & 0 & 0 \\ 1 & 1 & 1\\ \end{array} \right) * \left( \begin{array}{c} f_{n-1} \\ f_{n-2} \\ c \\ \end{array} \right) = \left( \begin{array}{c} f_n \\ f_{n-1} \\ c \\ \end{array} \right)$$
2. $f(n) = c^n - f(n-1) $（c是常数）
$$\left( \begin{array}{cc} -1 & c \\ 0 & c \\ \end{array} \right) * \left( \begin{array}{c} f_{n-1} \\ c^{n-1} \\ \end{array} \right) = \left( \begin{array}{c} f_n \\ c^n \\ \end{array} \right)$$
（整理自网络）
参考资料：
https://www.cnblogs.com/cmmdc/p/6936196.html
https://blog.csdn.net/wust_zzwh/article/details/52058209
class Solution {public: double myPow(double x, int n) { if (n == 0) return 1; if(x==0) return 0; long N = n; if (n < 0) { x = 1/x; N = -N; } double res = 1; //快速幂 while (N > 0) { if ((N & 1) == 1) { res = res * x; } x *= x; N >>= 1; } return res; }};
