bzoj 4037: [HAOI2015]Str 矩阵乘法&动态规划

       这道题目确实是一道好题啊。。（据ydc说还是抄ZJ的题？）

       首先根据定义推出f()的递推式：f(i)=Σ(j=1,m)f(i-j)，考虑到m<=5，因此这个可以用矩阵乘法加速运算。（补充一下以前一直没有写：矩阵乘法有一种比较简单的写法，就是用p[1][i]（第一位不一定是1但需要相同）表示状态i的答案，然后用a[i][j]=1表示上一步的状态i可以转移到下一步，这样用q[1][i]表示下一步的状态，那么q=p*a，这个应该是最好理解的写法了。可以自己推一下）

       然后来关注g()数组，如果用dp[i]表示g()的前i位构成的数字的答案（注意是一个矩阵），那么可以得到递推是dp[i]=Σdp[j]*w(g(j+1,i))，其中g(i,j)表示g()的第i-j位构成的数组，w(x)表示数字x的转移矩阵（因此dp[i]=B*w(g(1,i))，矩阵B仅有B[1][m]=1），注意递推式那里用了乘法，原因是w()和dp()都是矩阵，而一个矩阵向后转移x步就相当于*w(x)。

       实际实现中用c[i][j]表示w(i*10^j)，方便操作。

AC代码如下：

#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
#define mod 998244353
#define N 1005
using namespace std;

struct matrix{ int p[7][7]; }c[15][N],dp[N],a;
int s[N],n,m; char ch[N];
matrix operator *(matrix x,matrix y){
	int i,j,k; matrix z;
	for (i=1; i<=m; i++)
		for (j=1; j<=m; j++){
			ll tmp=0;
			for (k=1; k<=m; k++) tmp+=(ll)x.p[i][k]*y.p[k][j];
			z.p[i][j]=tmp%mod;
		}
	return z;
}
matrix operator +(matrix x,matrix y){
	int i,j;
	for (i=1; i<=m; i++)
		for (j=1; j<=m; j++){
			x.p[i][j]+=y.p[i][j]; if (x.p[i][j]>=mod) x.p[i][j]-=mod;
		}
	return x;
}
matrix ksm(matrix x,int y){
	matrix z=x;
	for (y--; y; y>>=1,x=x*x) if (y&1) z=z*x; return z;
}
int main(){
	scanf("%s",ch+1); n=strlen(ch+1); int i,j;
	scanf("%d",&m);
	for (i=1; i<=n; i++) s[i]=ch[i]-'0';
	for (i=1; i<=m; i++){
		a.p[i][m]=c[0][1].p[i][i]=1;
		if (i>1) a.p[i][i-1]=1;
	}
	for (i=1; i<=9; i++){
		c[i][1]=c[i-1][1]*a;
		for (j=2; j<=n; j++) c[i][j]=ksm(c[i][j-1],10);
	}
	dp[0].p[1][m]=1; matrix tmp;
	for (i=1; i<=n; i++){
		tmp=c[s[i]][1];
		for (j=i-1; j>=0; j--){
			dp[i]=dp[i]+dp[j]*tmp;
			if (j && s[j]) tmp=tmp*c[s[j]][i-j+1];
		}
	}
	printf("%d\n",dp[n].p[1][m]);
	return 0;
}

by lych

2016.3.16