「NOIP十连赛day3」涂色游戏 - 动态规划+组合数学+矩阵快速幂

题目大意

他们找到了一个n行m列呈网格状的画板。小A拿出了p支不同颜色的画笔，开始在上面涂色。看到小A涂好的画板，小B觉得颜色太单调了，于是把画板擦干净，希望涂上使它看起来不单调的颜色（当然，每个格子里只能涂一种颜色）。小B想知道一共有多少种不单调的涂色方案。我们定义一个涂色方案是不单调的，当且仅当任意相邻两列都出现了至少q种颜色。

---

题目分析

根据数据范围，我们肯定是对$m$用$\log$的算法解决，相邻两行的处理又很像状压Dp。
因此我们得出了初步思路：状压Dp+矩阵快速幂。
但是很显然，因为数据范围很大，状压是无法实现的。我们需要用一个更小的状态空间表示一列的状态。

设$f[i,S]$表示前$i$列，第$i$列有$S$种不同颜色的方案数。
我们可以考虑从前一列$f[i-1,T]$转移得到，转移的条件取决于这两列的相同颜色个数。
设相同颜色个数为$w$，需满足条件$S+T-w\ge q$。
因此我们可以这样转移：　

$$f[i,S]=f[i-1,T]\times \frac{g[n,S]}{C_p^S}\sum_w C_T^wC_{p-T}^{S-w}$$

其中$g[i,j]$的意思是一列有$i$个格子，染上$j$种不同颜色的方案数。
转移可以这样理解：我们用$\sum_w C_T^wC_{p-T}^{S-w}$表示从$p$中选择颜色满足条件，然后把颜色放入$n$行的格子中，即乘上$g[n,S]$，但是我们实际上我们已经将$p$种颜色映射为了我们选择的颜色，而$g[n,S]$仍然是以前未映射的颜色，这样计算有误，我们通过除以$C_p^S$将所有选出$S$个颜色的算作一类就正确了。

这个转移显然可以使用矩阵快速幂优化，不细讲了。

---

代码

#include<algorithm>
#include<iostream>
#include<iomanip>
#include<cstring>
#include<cstdlib>
#include<vector>
#include<cstdio>
#include<cmath>
#include<queue>
using namespace std;
typedef long long LL;
inline const LL Get_Int() {
	LL num=0,bj=1;
	char x=getchar();
	while(x<'0'||x>'9') {
		if(x=='-')bj=-1;
		x=getchar();
	}
	while(x>='0'&&x<='9') {
		num=num*10+x-'0';
		x=getchar();
	}
	return num*bj;
}
const int maxn=105;
const LL mod=998244353;
struct Matrix {
	LL n,m,a[maxn][maxn];
	Matrix(LL n,LL m) {
		init(n,m);
	}
	Matrix(LL n,LL m,char E) { //单位矩阵
		init(n,m);
		for(int i=1; i<=n; i++)a[i][i]=1;
	}
	void init(LL n,LL m) {
		this->n=n;
		this->m=m;
		memset(a,0,sizeof(a));
	}
	LL* operator [] (const LL x) {
		return a[x];
	}
	Matrix operator * (Matrix& b) {
		Matrix c(n,b.m);
		for(int i=1; i<=n; i++)
			for(int j=1; j<=b.m; j++)
				for(int k=1; k<=m; k++)
					c[i][j]=(c[i][j]+a[i][k]*b[k][j])%mod;
		return c;
	}
	void operator *= (Matrix& b) {
		*this=*this*b;
	}
	Matrix operator ^ (LL b) {
		Matrix ans(n,m,'e'),a=*this;
		while(b>0) {
			if(b&1)ans=ans*a;
			a*=a;
			b>>=1;
		}
		return ans;
	}
};
LL Quick_Pow(LL a,LL b) {
	LL ans=1;
	while(b>0) {
		if(b&1)ans=ans*a%mod;
		a=a*a%mod;
		b>>=1;
	}
	return ans;
}
LL n,m,p,q,C[105][105],g[105][105];
int main() {
	n=Get_Int();
	m=Get_Int();
	p=Get_Int();
	q=Get_Int();
	C[0][0]=1;
	for(int i=1; i<=100; i++) {
		C[i][0]=C[i][i]=1;
		for(int j=1; j<i; j++)
			C[i][j]=(C[i-1][j-1]+C[i-1][j])%mod;
	}
	g[0][0]=1;
	for(int i=1; i<=n; i++)
		for(int k=1; k<=p; k++)
			g[i][k]=(k*g[i-1][k]%mod+(p-(k-1))*g[i-1][k-1]%mod)%mod;
	Matrix B(p,p);
	for(LL S=1; S<=p; S++)
		for(LL T=1; T<=p; T++) {
			LL tmp=0;
			for(int w=0; w<=min(min(S,T),min(p,S+T-q)); w++)tmp=(tmp+C[T][w]*C[p-T][S-w]%mod)%mod;
			tmp=tmp*g[n][S]%mod*Quick_Pow(C[p][S],mod-2)%mod;
			B[S][T]=tmp;
		}
	B=B^(m-1);
	LL ans=0;
	for(int i=1; i<=p; i++)
		for(int j=1; j<=p; j++)
			ans=(ans+B[j][i]*g[n][i]%mod)%mod;
	printf("%lld\n",ans);
	return 0;
}