HDU 1588 矩阵快速幂 嵌套矩阵
这个题目搞了我差不多一个下午,之前自己推出一个公式,即 f[n+k]=k*f[n]+f[n-1]结果发现根本不能用,无法降低复杂度。
后来又个博客的做法相当叼,就按他的做法来了
即 最终求得是 S(n)=f[b]+f[b+k]+f[b+2*k]....f[b+n*k] (原题的意思好像是不用加到第n项,但实测确实要加到该项)
然后我们令 A={1,1}(标准的斐波那契矩阵)
{1,0}发现 f[b]=A^b,f[b+k]=A^(b+k),....f[b+nk]=A^(b+nk);
提取公共因子 A^b.S(n)=A^b*(E+A^K+A^K^2....A^K^n)
再令K=A^K (K和E都是矩阵,E为单位矩阵)。于是S(n)=A^b*(E+K+k^2+K^3+...K^n);
这样我们的目的大概出来了,就是要尽可能短的算出 上式括号中的内容(A^b可以用经典斐波那契矩阵很快算出来)
于是构造一个嵌套矩阵
令 SK={K,E},(K,E,0均为矩阵),SK^N={K^N , E+K+K^2+K^3+...K^N}
{0,E} {0 , E }
这样只要求出SK^N,取其第一行的第二项 再与 A^b相乘,即可得出最终结果
嵌套矩阵跟普通矩阵差不多,调用写好的矩阵乘法即可实现嵌套矩阵的乘法和乘方
#include <iostream> #include <cstdio> #include <cstring> #define ll __int64 using namespace std; ll k,b,n,M; struct Mat { ll mat[3][3]; } Fb,E,zero,A; struct Smart { Mat d[3][3]; } SK,SE; void test(Mat t) { for (int i=0;i<2;i++) { for (int j=0;j<2;j++) cout<<t.mat[i][j]<<" "; cout<<endl; } } void test2(Smart tmp) { for (int i=0;i<2;i++) { for (int j=0;j<2;j++) test(tmp.d[i][j]); } } Mat operator +(Mat a,Mat b) { Mat c; for (int i=0;i<2;i++) { for (int j=0;j<2;j++) { c.mat[i][j]=a.mat[i][j]+b.mat[i][j]; c.mat[i][j]%=M; } } return c; } Mat operator *(Mat a,Mat b) { Mat c; memset(c.mat,0,sizeof c.mat); for (int i=0;i<2;i++) { for (int j=0;j<2;j++) { for (int k=0;k<2;k++) { c.mat[i][j]+=a.mat[i][k]*b.mat[k][j]; if (c.mat[i][j]>=M) c.mat[i][j]%=M; } } } return c; } Mat operator ^(Mat a,ll x) { //if (x==0) return a; // cout<<x<<" the x"<<endl; Mat c=E; for (;x;x>>=1) { if (x&1) c=c*a; a=a*a; } return c; } Smart operator *(Smart a,Smart b) //大矩阵的乘法 调用普通矩阵的乘法和加法公式即可 { Smart c; Mat tmp; for (int i=0;i<2;i++) { for (int j=0;j<2;j++) { c.d[i][j]=zero; for (int k=0;k<2;k++) { tmp=a.d[i][k]*b.d[k][j]; c.d[i][j]=c.d[i][j]+tmp; } } } return c; } Smart operator ^(Smart a,ll x) //大矩阵的乘法 调用普通矩阵的乘法公式即可 { Smart c=SE; // if (x==0) return a; //cout<<x<<" the s x "<<endl; for (;x;x>>=1) { if (x&1) c=c*a; a=a*a; // test2(c); } // test2(c); return c; } void init() { memset(zero.mat,0,sizeof zero); memset(E.mat,0,sizeof E.mat); for (int i=0;i<2;i++) E.mat[i][i]=1; Mat A; memset(A.mat,0,sizeof A.mat); SE.d[0][0]=SE.d[1][1]=E; SE.d[0][1]=SE.d[1][0]=zero; A.mat[0][0]=A.mat[0][1]=A.mat[1][0]=1; Fb=A^b; //test(Fb); SK.d[0][0]=A^k; SK.d[0][1]=SK.d[1][1]=E; SK.d[1][0]=zero; } int main() { while (scanf("%I64d %I64d %I64d %I64d",&k,&b,&n,&M)!=EOF) { init(); Smart s=SK^n; Mat tmp=s.d[0][1]; //test(tmp); Mat ans=Fb*tmp; printf("%I64d\n",ans.mat[1][0]);//因为Fb求出来 真正的 A^b在[1][0]处,但是最终结果为何也是1 0处 我觉得还是有待考究 } return 0; }