首页 > 代码库 > poj3233—Matrix Power Series
poj3233—Matrix Power Series
题目链接:http://poj.org/problem?id=3233
题目意思:给一个矩阵n*n的矩阵A和一个k,求一个式子 S = A + A2 + A3 + … + Ak。
这个需要用到等比数列和的二分加速。
当n为奇数的时候,Sn=Sn-1+A^k;
当n为偶数的时候,Sn=(S[n/2]+E)*A^(k/2)
自己xjb推一下就知道等比数列和的二分加速是咋回事了。我举个例子,我们假设求等比数列2,4,8,16,32,64的和s=(8+1)*(2+4+8),而2+4+8=(2+4)+8,而(2+4)=(2+1)*2,其他的可以一次类推。
这个过程我们可以直接用一个递归过程算出来,其余我们套用矩阵快速幂的模板就好了。
代码:
1 //Author: xiaowuga 2 #include<iostream> 3 #include<cstring> 4 #define maxx INT_MAX 5 #define minn INT_MIN 6 #define inf 0x3f3f3f3f 7 #define size 35 8 using namespace std; 9 typedef long long ll; 10 int n,k,mod; 11 struct Mat{ 12 int mat[size][size]; 13 void clear(){ 14 memset(mat,0,sizeof(mat)); 15 } 16 17 Mat operator *(const Mat &e) const{ 18 Mat tmp; 19 tmp.clear(); 20 for(int k=0;k<n;k++) 21 for(int i=0;i<n;i++){ 22 if(mat[i][k]==0) continue; 23 for(int j=0;j<n;j++){ 24 if(e.mat[k][j]==0) continue; 25 tmp.mat[i][j]+=mat[i][k]*e.mat[k][j]%mod; 26 tmp.mat[i][j]%=mod; 27 } 28 } 29 return tmp; 30 } 31 Mat operator +(const Mat &e) const{ 32 Mat tmp; 33 tmp.clear(); 34 for(int i=0;i<n;i++){ 35 for(int j=0;j<n;j++){ 36 tmp.mat[i][j]=(mat[i][j]%mod+e.mat[i][j]%mod)%mod; 37 } 38 } 39 return tmp; 40 } 41 }; 42 Mat m,E; 43 Mat pow(Mat ma,ll num){ 44 Mat ans; 45 ans.clear(); 46 for(int i=0;i<n;i++) ans.mat[i][i]=1; 47 while(num){ 48 if(num&1) ans=ans*ma; 49 num/=2; 50 ma=ma*ma; 51 } 52 return ans; 53 } 54 Mat fun(int x){ 55 if(x==1) return m; 56 if(x&1) return fun(x-1)+pow(m,x); 57 else return (pow(m,x/2)+E)*fun(x/2); 58 } 59 int main() { 60 ios::sync_with_stdio(false);cin.tie(0); 61 cin>>n>>k>>mod; 62 E.clear(); 63 for(int i=0;i<n;i++) E.mat[i][i]=1; 64 for(int i=0;i<n;i++) 65 for(int j=0;j<n;j++) cin>>m.mat[i][j]; 66 Mat ans=fun(k); 67 for(int i=0;i<n;i++){ 68 for(int j=0;j<n;j++) 69 cout<<ans.mat[i][j]<<" "; 70 cout<<endl; 71 } 72 return 0; 73 }
poj3233—Matrix Power Series
声明:以上内容来自用户投稿及互联网公开渠道收集整理发布,本网站不拥有所有权,未作人工编辑处理,也不承担相关法律责任,若内容有误或涉及侵权可进行投诉: 投诉/举报 工作人员会在5个工作日内联系你,一经查实,本站将立刻删除涉嫌侵权内容。