快速幂矩阵小结
【目录】
1.快速幂（取模）
2.矩阵快速幂
一、快速幂（取模）
1.快速幂取模算法
const int mod = 1e9+7; long long qPow(long long a ,long long n){ long long ans = 1; //如果n的二进制最后一位是1 结果参与运算 //因为如果是0，那么幂运算之后该项是1,1乘任何数还是那个数，对结果不影响 while(n > 0){ if(n & 1) ans = (ans* a) % mod; a = (a*a) % mod;//底数加倍 n >>= 1;//移位 } return ans; }
二、快速幂矩阵
1.适用场景
a.清楚状态转移但无法求出通项公式或者递推公式。
（1）初始状态可设为单位阵
（2）通过简单试验可以得到状态转移矩阵
（3）后继状态可以通过当前状态转移之，其实就是在当前矩阵基础上乘上一个转移矩阵
（4）最终终态就是转移矩阵的n-k次方（k根据初态到终态需要转移多少次确定），考虑所有非零项之和或者是特定位置的数值
b.已知递推公式求第N项但是数值范围过大
用数组递推不仅会会超空间（比如109），而且会超时（朴素的一步一步递推速度太慢规模一大必然超时）
此时用类似于快速幂的思想加速转移，可以类比求21000时速度是越来越快的，求矩阵A的K次方速度也是越来越快的。
2.算法实现
int n = 4;//矩阵维数 //矩阵定义 struct Matrix { ll ma[5][5]; }; int L,M; //矩阵乘法运算 取模 Matrix multi(Matrix a, Matrix b) { Matrix ans;//存放结果矩阵 memset(ans.ma, 0, sizeof(ans.ma)); ll i, j, k; for(i=1; i<=n; ++i) for(j=1; j<=n; ++j) for(k=1; k<=n; ++k) ans.ma[i][j] = (ans.ma[i][j] + a.ma[i][k]*b.ma[k][j]%M) % M; return ans; } //矩阵快速幂 取模 Matrix pow(Matrix a, ll x) { Matrix ans; memset(ans.ma, 0, sizeof(ans.ma)); ll i, j; //单位阵初始化 for(i=1; i<=n; ++i) ans.ma[i][i] = 1; while(x) { if(x & 1) ans = multi(ans, a); a = multi(a, a); x >>= 1; } return ans; }
3.递推式的转移矩阵的构造
引例：斐波那契 递推公式 f(n) = f(n-1) + f(n-2)
这里的转移矩阵是A矩阵，注意如果横着写f(n),f(n-1)那么A就放在后边，如果竖着写A就放在前面，一定要注意，竖前横后（书签很厚）。
一般性构造原理。先确定矩阵维数，对应矩阵相乘原则m*n矩阵要和n*k矩阵相乘才可以，把转移矩阵A的所有项全部设成未知数，然后根据转移方程式进行矩阵相乘，得到若干方程，令每个方程左右两边系数分别相等，就能求出所有未知数进而确定转移方程A。
快速构造矩阵: 从矩阵Sn,到矩阵Sn-1。逐行构造转移矩阵T，该行Ti与矩阵Sn-1的向量积为Sni。
【核心思想】答案矩阵里含有所有生成数据的基础，或者称为基元素，转移矩阵根据两个相邻的答案矩阵的转移构造。
4.例题
根据前两种状态构造状态转移矩阵即可。
#include<cstdio> #include<cstring> #include<algorithm> #include<iostream> #define maxn using namespace std; typedef long long ll; int n = 4; struct Matrix { ll ma[5][5]; }; int L,M; Matrix multi(Matrix a, Matrix b) { Matrix ans; memset(ans.ma, 0, sizeof(ans.ma)); ll i, j, k; for(i=1; i<=n; ++i) for(j=1; j<=n; ++j) for(k=1; k<=n; ++k) ans.ma[i][j] = (ans.ma[i][j] + a.ma[i][k]*b.ma[k][j]%M) % M; return ans; } Matrix pow(Matrix a, ll x) { Matrix ans; memset(ans.ma, 0, sizeof(ans.ma)); ll i, j; for(i=1; i<=n; ++i) ans.ma[i][i] = 1; while(x) { if(x & 1) ans = multi(ans, a); a = multi(a, a); x >>= 1; } return ans; } int main() { Matrix ans; while(cin>>L>>M) { if(L == 0){ cout<<0<<endl; continue; } if(L == 1){ cout<<2%M<<endl; continue; } if(L == 2){ cout<<4%M<<endl; continue; } memset(ans.ma , 0, sizeof ans.ma); ans.ma[1][2] = 1; ans.ma[2][4] = 1; ans.ma[3][1] = 1; ans.ma[3][2] = 1; ans.ma[4][3] = 1; ans.ma[4][4] = 1; ans = pow(ans, L-2); ll sum = 0; for(int i=1; i<=n; i++){ for(int j=1; j<=n; j++){ //cout<< ans.ma[i][j]<<" "; sum += ans.ma[i][j]; sum %= M; } } cout<<sum<<endl; } return 0; }
之前的题目都是点到点的状态转移，这道题是二维的，可以将其看成是列到列转移，构造转移矩阵时运用快速构造法。
#include<iostream> #include<cstring> #include<string> #include<vector> #include<cstdio> using namespace std; typedef long long ll; int n = 14;//矩阵维数 const ll M = 10000007; //矩阵定义 struct Matrix { ll ma[15][15]; }; //矩阵乘法运算 取模 Matrix multi(Matrix a, Matrix b) { Matrix ans;//存放结果矩阵 memset(ans.ma, 0, sizeof(ans.ma)); ll i, j, k; for(i=1; i<=n; ++i) for(j=1; j<=n; ++j) for(k=1; k<=n; ++k) ans.ma[i][j] = (ans.ma[i][j] + a.ma[i][k]*b.ma[k][j]%M) % M; return ans; } //矩阵快速幂 取模 Matrix qpow(Matrix a, ll x) { Matrix ans; memset(ans.ma, 0, sizeof(ans.ma)); ll i, j; //单位阵初始化 for(i=1; i<=n; ++i) ans.ma[i][i] = 1; while(x) { if(x & 1) ans = multi(ans, a); a = multi(a, a); x >>= 1; } return ans; } int N,m; int main(){ while(cin>>N>>m){ ll cnt = 1; ll num[15]; num[1] = 23; ll s; for(int i=1; i<=N; i++){ cin>>s; num[++cnt] = s; } num[cnt+1] = 3; if(N == 0 && m == 0){ cout<<0<<endl; continue; } // cout<<"打印初始矩阵"<<endl; // for(int i=1; i<=N+2; i++){ // cout<<num[i]<<" "; // } // cout<<"\n"<<endl; Matrix ans; memset(ans.ma, 0, sizeof(ans.ma)); for(int i=1; i<=N+1; i++){ ans.ma[i][1] = 10; ans.ma[i][N+2] = 1; } ans.ma[N+2][N+2] = 1; for(int i=1; i<=N+1; i++){ for(int j=2; j<=i; j++){ ans.ma[i][j] = 1; } } // cout<<"打印转移矩阵......"<<endl; // for(int i=1; i<=N+2; i++){ // for(int j=1; j<=N+2; j++){ // cout<<ans.ma[i][j]<<" "; // } // cout<<endl; // } // cout<<"打印完毕。"<<endl; ans = qpow(ans , m); ll res = 0; for(int i=1; i<=N+2; i++){ res += ans.ma[N+1][i] * num[i]; } res = (res % M + M) % M; cout<<res<<endl; } return 0; }
