hdu1588 矩阵快速幂

//看了很多的博客 后来队友指点才懂
//
sum=f(g(0))+f(g(1))+.... //sum=A^(b-1)*|...|.... //要将b-1换,防止出现b=0时有负一,用A^b代替,取下面的即可 //这样问题成了 sum=A^b(A+A^(2k)+A^(3k)+...+A^(k(n-1))); //令B=A^k次,就简单了。 /* 主要要求1+A+A^2+A^3+...+A^(n-1)次方 | A A | | A^2 A^2+A | | A^(k-1) A^(k-1)+A^(k-2)+A^(k-3)... | 令B= | | B^2=| | 这样可以得到 B^(n-1)=| | | 0 1 | | 0 1 | | 0 1 | */ #include<stdio.h> #include<string.h> #define maxn 30 #define ll __int64 ll n,mod; struct Mat { ll mat[maxn][maxn]; }; Mat cal1(Mat a,Mat b,int nn)//矩阵乘法 { Mat c; memset(c.mat,0,sizeof(c.mat)); int i,j,k; for(i=0;i<nn;i++) for(j=0;j<nn;j++) for(k=0;k<nn;k++) { c.mat[i][j]+=((a.mat[i][k]*b.mat[k][j])%mod); c.mat[i][j]%=mod; } return c; } Mat cal2(Mat a,ll k,int nn)//矩阵幂 { int i,j; Mat c; for(i=0;i<nn;i++) for(j=0;j<nn;j++) if(i==j)c.mat[i][j]=1; else c.mat[i][j]=0; while(k) { if(k&1) c=cal1(c,a,nn); k=k>>1; a=cal1(a,a,nn); } return c; } Mat A,B,S; void initA() { A.mat[0][0]=1; A.mat[0][1]=1; A.mat[1][0]=1; A.mat[1][1]=0; } void initB() { int i,j; for(i=0;i<2;i++) for(j=0;j<2;j++) { B.mat[i][j]=A.mat[i][j]; B.mat[i][j+2]=A.mat[i][j]; } for(i=2;i<4;i++) for(j=0;j<2;j++) B.mat[i][j]=0; for(i=2;i<4;i++) for(j=2;j<4;j++) if(i==j) B.mat[i][j]=1; else B.mat[i][j]=0; } Mat getmat() { int i,j; Mat c; for(i=0;i<2;i++) for(j=2;j<4;j++) c.mat[i][j-2]=B.mat[i][j]; for(i=0;i<2;i++) for(j=0;j<2;j++) if(i==j) c.mat[i][j]+=1; return c; } int main() { ll i,j,k,b; while(scanf("%I64d %I64d %I64d %I64d",&k,&b,&n,&mod)!=EOF) { initA(); S=cal2(A,b,2); A=cal2(A,k,2); initB(); B=cal2(B,n-1,4); Mat temp=getmat(); S=cal1(S,temp,2); /*for(i=0;i<2;i++) { for(j=0;j<2;j++) printf("%I64d ",S.mat[i][j]); printf(" "); }*/ printf("%d ",S.mat[1][0]); } }
原文地址:https://www.cnblogs.com/sweat123/p/4666254.html