今天的专题是矩阵快速幂，其实就是把快速幂算法应用到矩阵中去，把普通的乘法重载成矩阵相乘，从而解决一些实际问题。基本操作如下：
快速幂取模：应用的是取模运算对乘法的可结合性，以及二进制的原理——一个数可以被拆分成若干个2的幂相加之和。具体代码如下：
int qpow(int a, int b, int p) { int res = 1; while(b) {//循环到b为0时，所有的因子都已相乘，算法结束 if(b & 1) {//相当于b%2 == 1 res = (res * a) % c;//b为奇数，需补乘一个a } b >>= 1;//相当于b /= 2 a = a*a % c;//相当于int k = a*a } return res; } 应用到矩阵中，则定义一个矩阵乘法运算，并将power运算稍作变化。代码如下： struct mat { ll a[maxn][maxn]; mat() {//构造方法 memset(a, 0, sizeof(a)); } }; mat init() {//初始化 mat tt; for(int i = 0;i < n;i++) { for(int j = 0;j < n;j++) { if(i == j) tt.a[i][j] = 1; else tt.a[i][j] = 0; } } return tt; } mat mul(mat x, mat y) {//乘法 mat tt; for(int i = 0;i < n;i++) { for(int j = 0;j < n;j++) { tt.a[i][j] = 0; for(int k = 0;k < n;k++) { tt.a[i][j] += x.a[i][k] * y.a[k][j]; } tt.a[i][j] %= mod; } } return tt; } mat mul(mat x, mat y) {//乘法 mat tt; for(int i = 0;i < n;i++) { for(int j = 0;j < n;j++) { tt.a[i][j] = 0; for(int k = 0;k < n;k++) { tt.a[i][j] += x.a[i][k] * y.a[k][j]; } tt.a[i][j] %= mod; } } return tt; }而本题就是一个应用的实例。题目要求这样一个以递归定义的表达式第n项的值（将结果mod2147493647）：
F(n) = F(n-1) + 2*F(n-2) + n^4.
由于数据范围很大，直接用循环算或用递归算的话必然会爆long long，用大数运算的话时空效率均很低。故利用线性代数中基的思想，采取矩阵相乘的形式将其拆开运算。
具体过程......这个辣鸡博客不让贴图片好气啊...!
因此仅通过对系数矩阵A做n-2次乘幂运算，再与一个已知的向量相乘，得到的向量第一个元素即为F（n）。具体ac代码如下：
#include<cstdio> #include<cstring> #include<iostream> #include<algorithm> using namespace std; typedef long long ll; const ll mod = 2147493647; const ll maxn = 2147483647; int cas; ll n, a, b; struct mat { ll a[7][7]; }; void print(mat ans) {//打印一个矩阵，一直wa，调试时方便 for(int i = 0;i < 7;i++) { for(int j = 0;j < 7;j++) printf("%d ", ans.a[i][j]); printf("\n"); } return; } mat init0() {//初始化为0 mat tt; for(int i = 0;i < 7;i++) for(int j = 0;j < 7;j++) tt.a[i][j] = 0; return tt; } mat init1() {//初始化为单位矩阵 mat tt; for(int i = 0;i < 7;i++) for(int j = 0;j < 7;j++) tt.a[i][j] = (i==j); return tt; } mat mul(mat x, mat y) {//矩阵 乘法 mat tt; for(int i = 0;i < 7;i++) { for(int j = 0;j < 7;j++) { tt.a[i][j] = 0; for(int k = 0;k < 7;k++) { tt.a[i][j] += x.a[i][k] * y.a[k][j]; } tt.a[i][j] %= mod; } } return tt; } mat power(mat tt, ll b) {//快速幂取模 mat res = init1(); while(b) { if(b & 1) { res = mul(res, tt); } b >>= 1; tt = mul(tt, tt); } return res; } int main() { while(~scanf("%d", &cas)) { mat coef = init0(); coef.a[0][0] = 1; coef.a[0][1] = 2; coef.a[0][2] = 1; coef.a[0][3] = 4; coef.a[0][4] = 6; coef.a[0][5] = 4; coef.a[0][6] = 1; coef.a[1][0] = 1; coef.a[2][2] = 1; coef.a[2][3] = 4; coef.a[2][4] = 6; coef.a[2][5] = 4; coef.a[2][6] = 1; coef.a[3][3] = 1; coef.a[3][4] = 3; coef.a[3][5] = 3; coef.a[3][6] = 1; coef.a[4][4] = 1; coef.a[4][5] = 2; coef.a[4][6] = 1; coef.a[5][5] = 1; coef.a[5][6] = 1; coef.a[6][6] = 1; //print(coef); while(cas--) { scanf("%I64d %I64d %I64d", &n, &a, &b); if(n == 1) printf("%I64d\n", a); else if(n == 2) printf("%I64d\n", b); else { mat aa = power(coef, n-2); mat bb = {b, 0, 0, 0, 0, 0, 0, a, 0, 0, 0, 0, 0, 0, 16,0, 0, 0, 0, 0, 0, 8, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0}; /*ll ans = 0; ans = (ans + aa.a[0][0]*b) % mod; ans = (ans + aa.a[0][1]*a) % mod; ans = (ans + aa.a[0][2]*16)% mod; ans = (ans + aa.a[0][3]*8) % mod; ans = (ans + aa.a[0][4]*4) % mod; ans = (ans + aa.a[0][5]*2) % mod; ans = (ans + aa.a[0][6]*1) % mod;*/ mat cc = mul(aa, bb); //printf("d\n", ans); printf("%lld\n", cc.a[0][0]); } } } return 0; }
