题目描述
小L有一座环形花园，沿花园的顺时针方向，他把各个花圃编号为1~N(2<=N<=10^15)。他的环形花园每天都会换一个新花样，但他的花园都不外乎一个规则，任意相邻M(2<=M<=5,M<=N)个花圃中有不超过K(1<=K<M)个C形的花圃，其余花圃均为P形的花圃。
例如，N=10,M=5,K=3。则
CCPCPPPPCC 是一种不符合规则的花圃；
CCPPPPCPCP 是一种符合规则的花圃。
请帮小L求出符合规则的花园种数Mod 1000000007
由于请您编写一个程序解决此题。
输入输出格式
输入格式：
一行，三个数N,M,K。
输出格式：
花园种数Mod 1000000007
输入输出样例
【样例输入1】 10 5 3 【样例输入2】 6 2 1
【样例输出1】 458 【样例输出2】 18
说明
【数据规模】
40%的数据中，N<=20；
60%的数据中，M=2；
80%的数据中，N<=10^5。
100%的数据中，N<=10^15。
题解：
先讲一讲70分的做法，首先一看数据范围，M<=5，想到从m下手。然后由于它是一个环，想到破环成链，又因为只跟后m位有关系，所以只需将序列拓宽m个单位即可。综合以上几点，可以想到比较暴力的做法，对后m位进行状压，f[i][j]表示当前处理到第i位，后m位的状态位j的方案数。然后预处理出一个v[i][j]表示状态i转移到状态j知否可行，ok[i]表示状态i是否合法。于是状压DP的解法就出来了。70分代码如下：
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<stack> 9 #include<ctime> 10 #include<vector> 11 #define mod (1000000007) 12 using namespace std; 13 typedef long long lol; 14 lol n,m,l,ans,f[100010][64],to; 15 bool v[64][64],ok[64]; 16 lol gi() 17 { 18 lol ans=0,f=1; 19 char i=getchar(); 20 while(i<'0'||i>'9'){if(i=='-')f=-1;i=getchar();} 21 while(i>='0'&&i<='9'){ans=ans*10+i-'0';i=getchar();} 22 return ans*f; 23 } 24 bool judge(lol a,lol b) 25 { 26 lol x[10],y[10],cnt1=0,cnt2=0,c=a,d=b; 27 for(lol i=1;i<=m;i++){x[i]=a&1;cnt1+=a&1;a>>=1;y[i]=b&1;cnt2+=b&1;b>>=1;} 28 for(lol i=1;i<m;i++)if(x[i]!=y[i+1])return false; 29 return true; 30 } 31 bool pd(lol x) 32 { 33 lol ans=0; 34 while(x){ans+=x&1;x>>=1;} 35 if(ans<=l)return true; 36 else return false; 37 } 38 void dp(int x) 39 { 40 memset(f,0,sizeof(f)); 41 f[m][x]=1; 42 for(lol i=m+1;i<=n+m;i++) 43 for(lol j=0;j<=to;j++) 44 for(lol k=0;k<=to;k++) 45 if(v[j][k])f[i][k]=(f[i][k]+f[i-1][j])%mod; 46 ans=(ans+f[n+m][x])%mod; 47 } 48 int main() 49 { 50 n=gi();m=gi();l=gi(); 51 to=(1<<m)-1; 52 for(lol i=0;i<=to;i++)if(pd(i))ok[i]=1; 53 for(lol i=0;i<=to;i++) 54 for(lol j=0;j<=to;j++) 55 if(ok[i]&&ok[j]&&judge(i,j))v[i][j]=1; 56 for(lol i=0;i<=to;i++)if(ok[i])dp(i); 57 printf("%lld ",ans); 58 return 0; 59 }
很显然的是这份代码的空间复杂度是O(n*2^m)，根据数据范围会爆内存。时间复杂度是O(n*(2^m)^2)很显然也会超时。
那么该怎么办呢？考虑矩阵快速幂优化。（实际上根本不是优化，因为用了矩阵快速幂之后就根本不是DP了）
想想，我们预处理的v数组表示的是不是每一种状态之间的连通性，而根据Floyed矩阵的性质，把这个矩阵n次方后，表示的就是从i状态走n次到j的路径数。那么就可以用快速幂解题了。直接把v数组n次方，最后统计一下ok[i]==1（i状态合法）的情况下v[i][i]的和就可以了。（它是一个环，所最前的状态应该和最后的状态相同）
代码如下：
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<stack> 9 #include<ctime> 10 #include<vector> 11 #define mod (1000000007) 12 using namespace std; 13 typedef long long lol; 14 lol n,m,l,ans,f[100010][64],to; 15 bool ok[64]; 16 struct matrix 17 { 18 lol a[64][64]; 19 matrix(){for(int i=0;i<=to;i++)for(int j=0;j<=to;j++)a[i][j]=0;} 20 matrix(lol b[3][3]){for(int i=0;i<=to;i++)for(int j=0;j<=to;j++)a[i][j]=b[i][j];} 21 friend matrix operator * (const matrix a,const matrix b) 22 { 23 matrix ans; 24 for(int i=0;i<=to;i++) 25 for(int j=0;j<=to;j++) 26 for(int k=0;k<=to;k++) 27 ans.a[i][j]=(ans.a[i][j]+(a.a[i][k]%mod)*(b.a[k][j]%mod)%mod)%mod; 28 return ans; 29 } 30 }v,a; 31 lol gi() 32 { 33 lol ans=0,f=1; 34 char i=getchar(); 35 while(i<'0'||i>'9'){if(i=='-')f=-1;i=getchar();} 36 while(i>='0'&&i<='9'){ans=ans*10+i-'0';i=getchar();} 37 return ans*f; 38 } 39 bool judge(lol a,lol b) 40 { 41 lol x[10],y[10],cnt1=0,cnt2=0,c=a,d=b; 42 for(lol i=1;i<=m;i++){x[i]=a&1;cnt1+=a&1;a>>=1;y[i]=b&1;cnt2+=b&1;b>>=1;} 43 for(lol i=1;i<m;i++)if(x[i]!=y[i+1])return false; 44 return true; 45 } 46 bool pd(lol x) 47 { 48 lol ans=0; 49 while(x){ans+=x&1;x>>=1;} 50 if(ans<=l)return true; 51 else return false; 52 } 53 matrix make(matrix s,lol x) 54 { 55 matrix ans=s; 56 x--; 57 while(x) 58 { 59 if(x&1)ans=ans*s; 60 s=s*s; 61 x>>=1; 62 } 63 return ans; 64 } 65 int main() 66 { 67 n=gi();m=gi();l=gi(); 68 to=(1<<m)-1; 69 for(lol i=0;i<=to;i++)if(pd(i))ok[i]=1; 70 for(lol i=0;i<=to;i++) 71 for(lol j=0;j<=to;j++) 72 if(ok[i]&&ok[j]&&judge(i,j))v.a[i][j]=1; 73 a=make(v,n); 74 for(lol i=0;i<=to;i++)if(ok[i])ans=(ans+a.a[i][i])%mod; 75 printf("%lld ",ans); 76 return 0; 77 }