矩阵快速幂
矩阵是种神奇的东西,无奈我这种数学渣对于他的理解只能是表面。
快速幂应该都很熟悉,实际上把快速幂里面的乘法换成矩阵乘法就可以变为矩阵快速幂
上份模板:
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 typedef vector<int>vec; 2 typedef vector<vec>mat; 3 typedef long long ll; 4 mat mul(mat &A,mat &B){ 5 mat C (A.size(),vec(B[0].size())); 6 for(int i=0;i<A.size();++i){ 7 for(int k=0;k<B.size();++k){ 8 for(int j=0;j<B[0].size();++j){ 9 C[i][j]=(C[i][j]+(A[i][k]*B[k][j])%mod)%mod; 10 } 11 } 12 } 13 return C; 14 } 15 mat pow(mat A,ll n){ 16 mat B (A.size(),vec(A.size())); 17 for(int i=0;i<A.size();++i)B[i][i]=1; 18 while(n>0){ 19 if(n&1)B = mul(B,A); 20 A = mul(A,A); 21 n>>=1; 22 } 23 return B; 24 }
矩阵快速幂有什么作用呢?优化递推关系!!!!
先看一下求第n项费布那切数列f[n]=f[n-1]+f[n-2]
这样我们就可以快速的求得第n项
题意:用红蓝绿黄去涂n个块,问红绿都为偶数块的方案数
假设当前前i个已经涂完其中红绿都为偶数有ai种,一奇一偶为bi,都为奇为ci ,那么对于i+1 有
a(i+1)=2*a(i)+b(i)
b(i+1)=2*a(i)+2*b(i)+2*c(i)
c(i+1)=b(i+1)+2*c(i)
这样就转化为求a(n)了
利用矩阵快速幂可以很好的解决。
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include <iostream> 2 #include <cstdio> 3 #include <cstring> 4 #include <set> 5 #include <algorithm> 6 #include <map> 7 #include <queue> 8 #include<cmath> 9 #include<vector> 10 #define maxn 50010 11 #define maxm 100010 12 #define mod 10007 13 #define INF 0x3f3f3f3f 14 using namespace std; 15 typedef vector<int>vec; 16 typedef vector<vec>mat; 17 typedef long long ll; 18 mat mul(mat &A,mat &B){ 19 mat C (A.size(),vec(B[0].size())); 20 for(int i=0;i<A.size();++i){ 21 for(int k=0;k<B.size();++k){ 22 for(int j=0;j<B[0].size();++j){ 23 C[i][j]=(C[i][j]+(A[i][k]*B[k][j])%mod)%mod; 24 } 25 } 26 } 27 return C; 28 } 29 mat pow(mat A,ll n){ 30 mat B (A.size(),vec(A.size())); 31 for(int i=0;i<A.size();++i)B[i][i]=1; 32 while(n>0){ 33 if(n&1)B = mul(B,A); 34 A = mul(A,A); 35 n>>=1; 36 } 37 return B; 38 } 39 int n; 40 void solve(){ 41 mat A(3,vec(3)); 42 A[0][0]=2;A[0][1]=1;A[0][2]=0; 43 A[1][0]=2;A[1][1]=2;A[1][2]=2; 44 A[2][0]=0;A[2][1]=1;A[2][2]=2; 45 A = pow(A,n); 46 printf("%d\n",A[0][0]); 47 } 48 int main (){ 49 int t; 50 scanf("%d",&t); 51 while(t--){ 52 scanf("%d",&n); 53 solve(); 54 } 55 }
题意:给你一个4*n的方格,问用1*2的去填有多少种方案。
这里给出Matrix67大神的一篇博客点我
里面的问题九就是这样一个问题,主要还是转移矩阵的构造,对于行数比较小的情况我们可以自己手动推,也更有助于理解。
但是实际上可以通过计算机来计算这个转移的矩阵,有点类似DFA。
我们只需要给他转移的规则就可以构造出这个矩阵(记忆化搜索)具体的看代码吧
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include <iostream> 2 #include <cstdio> 3 #include <cstring> 4 #include <set> 5 #include <algorithm> 6 #include <map> 7 #include <queue> 8 #include<cmath> 9 #include<vector> 10 #define maxn 50010 11 #define maxm 100010 12 #define INF 0x3f3f3f3f 13 using namespace std; 14 typedef vector<int>vec; 15 typedef vector<vec>mat; 16 typedef long long ll; 17 int mod; 18 mat mul(mat &A,mat &B){ 19 mat C (A.size(),vec(B[0].size())); 20 for(int i=0;i<A.size();++i){ 21 for(int k=0;k<B.size();++k){ 22 for(int j=0;j<B[0].size();++j){ 23 C[i][j]=(C[i][j]+(A[i][k]*B[k][j])%mod)%mod; 24 } 25 } 26 } 27 return C; 28 } 29 mat pow(mat A,ll n){ 30 mat B (A.size(),vec(A.size())); 31 for(int i=0;i<A.size();++i)B[i][i]=1; 32 while(n>0){ 33 if(n&1)B = mul(B,A); 34 A = mul(A,A); 35 n>>=1; 36 } 37 return B; 38 } 39 int n; 40 int dp[17][17]; 41 void dfs(int i,int s){ 42 if(dp[i][s]==1)return; 43 dp[i][s]=1; 44 for(int k=0;k<4;++k){ 45 if(k<3&&!(s&(1<<k))&&!(s&(1<<(k+1)))){ 46 dfs(i,s+(1<<k)+(1<<(k+1))); 47 } 48 } 49 } 50 void solve(){ 51 memset(dp,0,sizeof(dp)); 52 for(int i=0;i<(1<<4);++i)dfs(i,(1<<4)-1-i); 53 mat A(16,vec(16)); 54 mat I(16,vec(16)); 55 for(int i=0;i<16;++i){ 56 for(int j=0;j<16;++j){ 57 if(i==j)I[i][j]=1; 58 else I[i][j]=0; 59 A[i][j]=dp[i][j]; 60 // cout<<A[i][j]<<" "; 61 } 62 //cout<<endl; 63 } 64 A = pow(A,n+1); 65 A = mul(I,A); 66 printf("%d\n",A[15][0]); 67 } 68 int main (){ 69 int t; 70 while(scanf("%d%d",&n,&mod)!=EOF){ 71 if(mod==0&&n==0)break; 72 solve(); 73 } 74 }
题意给你矩阵A 求sum(A^0+A^1+...+A^k)%m
有了这个就可以直接搞了。实际上有更快的算法,基本原理就是Matrix里面写的
A^1+A^2+...+A^6 = (A^1+A^2+A^3)+A^3* (A^1+A^2+A^3)
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include <iostream> 2 #include <cstdio> 3 #include <cstring> 4 #include <set> 5 #include <algorithm> 6 #include <map> 7 #include <queue> 8 #include<cmath> 9 #include<vector> 10 #define maxn 50010 11 #define maxm 100010 12 #define INF 0x3f3f3f3f 13 using namespace std; 14 typedef vector<int>vec; 15 typedef vector<vec>mat; 16 typedef long long ll; 17 int mod; 18 mat mul(mat &A,mat &B){ 19 mat C (A.size(),vec(B[0].size())); 20 for(int i=0;i<A.size();++i){ 21 for(int k=0;k<B.size();++k){ 22 for(int j=0;j<B[0].size();++j){ 23 C[i][j]=(C[i][j]+(A[i][k]*B[k][j])%mod)%mod; 24 } 25 } 26 } 27 return C; 28 } 29 mat pow(mat A,ll n){ 30 mat B (A.size(),vec(A.size())); 31 for(int i=0;i<A.size();++i)B[i][i]=1; 32 while(n>0){ 33 if(n&1)B = mul(B,A); 34 A = mul(A,A); 35 n>>=1; 36 } 37 return B; 38 } 39 int n; 40 mat M; 41 void solve(int k){ 42 mat A(2*n,vec(2*n)); 43 for(int i=0;i<n;++i){ 44 for(int j=0;j<n;++j){ 45 A[i][j]=M[i][j]; 46 } 47 A[i+n][i]=A[i+n][i+n]=1; 48 } 49 A = pow(A,k+1); 50 for(int i=0;i<n;++i){ 51 for(int j=0;j<n;++j){ 52 int a = A[n+i][j]%mod; 53 if(i==j)a = (a+mod-1)%mod; 54 printf("%d%c",a,j+1==n?'\n':' '); 55 } 56 57 } 58 } 59 int main (){ 60 int k; 61 while(scanf("%d%d%d",&n,&k,&mod)!=EOF){ 62 M = mat(n,vec(n)); 63 for(int i=0;i<n;++i){ 64 for(int j=0;j<n;++j){ 65 scanf("%d",&M[i][j]); 66 } 67 } 68 solve(k); 69 } 70 }
矩阵快速幂的题目主要还是构造矩阵,如何把题目转化为矩阵这一步才是关键!!!!