BZOJ 1009: [HNOI2008]GT考试（DP+KMP失配+矩阵快速幂）

1009: [HNOI2008]GT考试

Time Limit: 1 Sec  Memory Limit: 162 MB
Submit: 4181  Solved: 2552
[Submit][Status][Discuss]

Description

　　阿申准备报名参加GT考试，准考证号为N位数X1X2....Xn(0<=Xi<=9),他不希望准考证号上出现不吉利的数字。
他的不吉利数学A1A2...Am(0<=Ai<=9)有M位，不出现是指X1X2...Xn中没有恰好一段等于A1A2...Am. A1和X1可以为
0

Input

　　第一行输入N,M,K.接下来一行输入M位的数。 N<=10^9,M<=20,K<=1000

Output

　　阿申想知道不出现不吉利数字的号码有多少种，输出模K取余的结果.

Sample Input

4 3 100
111

Sample Output

81

好兴奋的做到这么一道题，还发现了其中个别也不知道是否能定义为错误的数据，先说说错误数据的问题，困然我两天，终于忍不住拿别人代码比较，发现矩阵快速幂出了错，相当自信的觉得自己的模板不会错，所以就开始跑数据输出，发现AC的代码会出现负数，AC代码用的数组是int，而longlong的结果跑大数完全不对，所以可能是数据有误，也可能是我理解的不对，不过我感觉出现负数都能过就有点问题了。

然后说说这道想了两天的题，以为模拟可以，但想想这是BZOJ怎么会简单模拟，然后就开始了漫长的演算。感觉数位DP应该是可以的，但是怎么难跟m个字符串联系起来，偷偷看了别人写的，说是KMP的失配函数，然后就去看了看失配数组，果然有收获，想想也是，不就是当前数位的后缀跟m个数的前缀有多少相同的吗，这个地方可以用KMP搞定了，新的问题又来了，DP怎么搞，一看10^9，应该有公式或者矩阵快速幂，然后就在推公式，用DP[i][j]表示第i位数匹配到第j位的种类数。发现DP[i+1][p]+=Dp[i][j]，P是匹配到的任意值，也就是说这一项跟前面所有DP[i]的项都有关，通俗点讲就是当前状态的可能种类只与上一状态有关，这正是可以利用矩阵快速幂的原因。next数组储存的就是当这一位数为i时，有多少种可能情况。

好吧，还是把公式打出来吧             dp[i][j]=a0*dp[i-1][0]+a1*dp[i-1][1]+a2*dp[i-1][2]+``````+a(m-1)dp[i-1][m-1]

然后就可以矩阵快速幂了，构造一个m*m的矩阵，系数a就是失配数组求得的种类，跑一遍矩阵快速幂，最后把dp[n][0]到dp[m-1]加起来的值就是结果。

代码中longlong定义的数组，还防止出现负数，如果想A到，longlong变int，然后去掉防止负数就可以了。

代码实现：

#include<iostream>#include<algorithm>#include<cstring>#include<cmath>#include<queue>#include<cstdio>#define ll long long#define mset(a,x) memset(a,x,sizeof(a))using namespace std;const double PI=acos(-1);const int inf=0x3f3f3f3f;const double esp=1e-6;const int maxn=100000;const int mod=1e9+7;int dir[4][2]={0,1,1,0,0,-1,-1,0};ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}ll lcm(ll a,ll b){return a/gcd(a,b)*b;}ll inv(ll b){if(b==1)return 1; return (mod-mod/b)*inv(mod%b)%mod;}ll fpow(ll n,ll k){ll r=1;for(;k;k>>=1){if(k&1)r=r*n%mod;n=n*n%mod;}return r;}int next[100],n,m,q;char map[100];struct matrix{      int r,c;    ll dp[40][40];      matrix(int R=0,int C=0){          memset(dp,0,sizeof(dp));          r=R;c=C;      }  }a,b;matrix operator*(matrix A,matrix B){      matrix C(A.r,B.c);    for(int i=0;i<A.r;i++) {for(int j=0;j<B.c;j++)              for(int k=0;k<A.c;k++){                  C.dp[i][j]=(C.dp[i][j]+A.dp[i][k]*B.dp[k][j])%q;                //防止爆表出现负数                 C.dp[i][j]%=q;      C.dp[i][j]+=q;    C.dp[i][j]%=q;            }}             return C; }//得到失配函数的next数组，用来构造矩阵 void getfail(){int i=1,j=0;      next[1]=0;      while(i<=m)      {          if (j==0||map[i]==map[j]) next[++i]=++j;          else j=next[j];      }  }//矩阵快速幂 void mat_qpow(){int i,j,k;while(n){if(n&1)a=a*b;n>>=1;b=b*b;}ll ans=0;for(i=0;i<m;i++)ans=(ans+a.dp[0][i])%q;cout<<ans<<endl;}int main(){int i,j,k,l;cin>>n>>m>>q;scanf("%s",map+1);getfail();for(i=0;i<m;i++)                          //默认匹配到自己的数量为1a.dp[i][i]=1;for(i=1;i<=m;i++)                        //利用失配数组构造矩阵 {for(j=0;j<=9;j++){int temp=i;while(temp&&map[temp]-'0'!=j)temp=next[temp];b.dp[i-1][temp]++;}}//先进行模除 for(i=0;i<m;i++){for(j=0;j<m;j++)b.dp[i][j]%=q;}a.c=a.r=b.c=b.r=m;mat_qpow();return 0;}/*1000000001110000000011111111111*/