矩阵快速幂

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

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;
}