洛谷 P3706 - [SDOI2017]硬币游戏(高斯消元)

发布时间 2023-05-22 17:56:33作者: tzc_wk

听说是 PGF 板板题,但是不会 PGF,怎么办捏(

暴力做法显然是建出 AC 自动机但是高斯消元,但是状态数高达 \(nm\),有没有优化的余地呢?

注意到终止状态只有 \(n\) 个,AC 自动机上其他节点表示的状态其实都可以归结为“非终止状态”,因此我们考虑设 \(n\) 个变量 \(x_1\sim x_n\) 表示最先出现字符串 \(s_i\) 的概率,和一个辅助变量 \(x_{n+1}\) 表示当前位于不合法状态的概率(当然如果将这个解释为“概率”其实有点牵强,因为我们设立这个变量并不是为了解出 \(x_{n+1}\) 这个没啥实际意义的东西,而是计算出 \(x_1\sim x_n\) 之间的比例关系,然后根据 \(\sum x_i=1\) 解出所有 \(x_i\),但是其实我也找不到更好的解释方式,只好将这个东西用”概率“来解释)

考虑如何建立这 \(n+1\) 个变量之间的方程,对于每个字符串 \(i\),我们考虑现在我们在不合法状态,然后以 \(2^{-m}\) 的概率恰好在后面添加了字符串 \(s_i\),那么这个概率应该用什么来衡量呢?你可能说是 \(x_i\),但是问题在于,不一定非得添加完全部 \(m\) 个字符才能第一次得到 \(s_i\),比方说 \(s_i=101\),然后现在字符串末尾刚好有个 \(10\),这样添加一个 \(1\) 以后就到了 \(x_i\) 的状态。换句话说,添加完 \(m\) 个字符以后恰好第一次出现 \(s_i\) 的概率 \(x_i\) 只是 \(x_{n+1}·2^{-m}\) 的一部分。为了覆盖所有的情况,我们考虑所有字符串 \(s_j\)(当然,\(j\) 也可以等于 \(i\)),如果 \(s_i\) 长度为 \(k\) 的前缀与 \(s_j\) 长度为 \(k\) 的后缀完全相等,说明有可能是从 \(x_k\) 的状态以 \(2^{-(m-k)}\) 的概率添加了 \(s_i\) 的后 \(m-k\) 的字符得到当前状态的,因此我们有

\[2^{-m}x_{n+1}=\sum\limits_{j=1}^n\sum\limits_{k=1}^m[s_i[1...k]=s_j[m-k+1...m]]x_j·2^{-(m-k)} \]

当然最后还有个 \(\sum x_i=1\)。直接高斯消元即可。

有点卡精度,很不明白为什么不取模,有点复古。下头/ruo

const int MAXN=300;
const u64 BS=131;
const ld EPS=1e-10;
int n,m;char s[MAXN+5][MAXN+5];
u64 pw_hsh[MAXN+5],pre[MAXN+5][MAXN+5],suf[MAXN+5][MAXN+5];
ld a[MAXN+5][MAXN+5],pw2[MAXN+5],res[MAXN+5];
int main(){
	scanf("%d%d",&n,&m);
	for(int i=(pw_hsh[0]=pw2[0]=1);i<=m;i++)pw_hsh[i]=pw_hsh[i-1]*BS,pw2[i]=pw2[i-1]*0.5;
	for(int i=1;i<=n;i++)scanf("%s",s[i]+1);
	for(int i=1;i<=n;i++){
		for(int j=1;j<=m;j++)pre[i][j]=pre[i][j-1]*BS+s[i][j];
		for(int j=m;j;j--)suf[i][j]=suf[i][j+1]+pw_hsh[m-j]*s[i][j];
	}
	for(int i=1;i<=n;i++){
		a[i][n+1]=-pw2[m];
		for(int j=1;j<=n;j++)for(int k=1;k<=m;k++)
			if(suf[j][m-k+1]==pre[i][k])a[i][j]+=pw2[m-k];
	}
	for(int i=1;i<=n;i++)a[n+1][i]=1;a[n+1][n+2]=1;
	for(int i=1;i<=n+1;i++){
		int t=i;
		for(int j=i+1;j<=n+1;j++)if(fabs(a[j][i]-a[t][i])>EPS)t=j;
		for(int j=i;j<=n+2;j++)swap(a[i][j],a[t][j]);
		for(int j=i+1;j<=n+2;j++)a[i][j]/=a[i][i];a[i][i]=1;
		for(int j=i+1;j<=n+1;j++){
			for(int k=i+1;k<=n+2;k++)a[j][k]-=a[j][i]*a[i][k];
			a[j][i]=0;
		}
	}
	for(int i=n+1;i;i--){
		res[i]=a[i][n+2];
		for(int j=i+1;j<=n+1;j++)res[i]-=res[j]*a[i][j];
	}
	for(int i=1;i<=n;i++)printf("%.10Lf\n",res[i]);
	return 0;
}