矩阵快速幂

矩阵快速幂在动态规划，快速计算矩阵有着非常重要的作用，可以很大程度上减少计算消耗的时间。它的实质其实就是：矩阵计算＋快速幂。

1.矩阵乘法

对于矩阵A＊B，要遵循的是A的列数要等于B的行数，其他要求没有：

#include <iostream>#include <algorithm>#include <vector>#include <cstdio>#include <ctime>using namespace std;const int M = 10000;typedef long long LL;typedef vector<int> vec;typedef vector<vec> Matrix;//计算矩阵相乘Matrix mul(Matrix& A, Matrix& B) {if (A[0].size() != B.size()) {cerr << "Wrong !" << endl;exit(0);}Matrix C(A.size(), vec(B[0].size()));for (int i = 0; i < A.size(); i++) {for (int j = 0; j < B.size(); j++) {for (int k = 0; k < B[0].size(); k++) {C[i][k] = (C[i][k] + A[i][j]*B[j][k]) % M;}}}return C;}

2.快速幂方法

这在求阶乘中用得最多，其实借用的就是计算机内二进制计数的原则，因为任何一个整数都可以写成Sigma(2^i)，所以在实际计算中，对结果会有影响（或者说对结果有贡献）的就是二进制数位上为1的权值，所以，对于一个数的N次方求解可以写成：

//求A的n次方int res = 1;while (n > 0) {if (n & 1) res *= A;A *= A;n >> =1;}

这里举个例子：

要求斐波那契数列的第N项（1<= 10^16)，直接暴力肯定会超时（按ACM的标准）。

斐波那契数列的规律是：F(n+2) = F(n+1) + F(n)。怎么和矩阵联系在一起？我们可以发现：

|F(n+1)|   |1   1| | F(n) |      |1   1|^n    |F(1)||      | = |     |*|      | ==>  |     |   *  |    || F(n) |   |1   0| |F(n-1)|      |1   0|      |F(0)|

问题解决，而且10^16以内耗时不超过2.7MS，一下是C++代码：

#include <iostream>#include <algorithm>#include <vector>#include <cstdio>#include <ctime>using namespace std;const int M = 10000;typedef long long LL;typedef vector<int> vec;typedef vector<vec> Matrix;//计算矩阵相乘Matrix mul(Matrix& A, Matrix& B) {if (A.size() != B[0].size()) {cerr << "Wrong !" << endl;exit(0);}Matrix C(A.size(), vec(B[0].size()));for (int i = 0; i < A.size(); i++) {for (int j = 0; j < B.size(); j++) {for (int k = 0; k < B[0].size(); k++) {C[i][k] = (C[i][k] + A[i][j]*B[j][k]) % M;}}}return C;}//计算矩阵阶乘Matrix Pow(Matrix A, LL n) {Matrix temp(A.size(), vec(A.size()));//单位化处理for (int i = 0; i < A.size(); i++) temp[i][i] = 1;//快速幂处理while (n > 0) {if (n & 1) temp = mul(temp, A);A = mul(A, A);n >>= 1;}return temp;}int main() {LL n;double totalTime;clock_t start, finish;puts("Please enter a num:");scanf("%lld", &n);start = clock();Matrix A(2, vec(2));A[0][0] = 1; A[0][1] = 1;A[1][0] = 1; A[1][1] = 0;A = Pow(A, n);finish = clock();totalTime = (double)(finish - start) / CLOCKS_PER_SEC;printf("Proceeding time is: %.5lf, The answer is: %d\n", totalTime, A[1][0]);return 0;}