[BZOJ4820]硬币游戏 KMP+高斯消元

4820: [Sdoi2017]硬币游戏

Time Limit: 10 Sec  Memory Limit: 128 MB

Description

周末同学们非常无聊，有人提议，咱们扔硬币玩吧，谁扔的硬币正面次数多谁胜利。大家纷纷觉得这个游戏非常符

合同学们的特色，但只是扔硬币实在是太单调了。同学们觉得要加强趣味性，所以要找一个同学扔很多很多次硬币

，其他同学记录下正反面情况。用H表示正面朝上，用T表示反面朝上，扔很多次硬币后，会得到一个硬币序列。比

如HTT表示第一次正面朝上，后两次反面朝上。但扔到什么时候停止呢？大家提议，选出n个同学，每个同学猜一个

长度为m的序列，当某一个同学猜的序列在硬币序列中出现时，就不再扔硬币了，并且这个同学胜利，为了保证只

有一个同学胜利，同学们猜的n个序列两两不同。很快，n个同学猜好序列，然后进入了紧张而又刺激的扔硬币环节

。你想知道，如果硬币正反面朝上的概率相同，每个同学胜利的概率是多少。

Input

第一行两个整数n,m。

接下里n行，每行一个长度为m的字符串，表示第i个同学猜的序列。

1<=n,m<=300

Output

输出n行，第i行表示第i个同学胜利的概率。

输出与标准输出的绝对误差不超过10^-6即视为正确。

Sample Input

3 3
THT
TTH
HTT

Sample Output

0.3333333333
0.2500000000
0.4166666667

 

题解：

生无可恋.jpg

在我做过的题里继小凸玩密室又一神题啊……我可能学了假的概率DP

上来我们可以发现这本是一道字符集为2的AC自动机的DP，迅速码好读入和处理fail指针。

然后我们继续套路，发现这个获胜概率在Trie图上互相牵制不能递推，所以我们要用高斯消元。

……高斯消元？

如果我们对于每个AC自动机的节点都这样做的话，仅时间复杂度就变成了O(300⁶），直接螺旋上天了。

那么我们只能不考虑每个节点，而是考虑每个串了。

设p[i]为第i个同学获胜的概率，也就是说第一个匹配到第i个串的概率。

首先我们可以列出第一个方程:Σp[i]=1.0

那么我们考虑，在一个不是单词节点的节点，假设已经经过的字符状态为S，

如果我们向后面添加m个字符，那么显然，由于我们是完全随机添加的，所以匹配到每个串i的概率都是一样的，我们设为H。

不难看出，如果我们匹配到i，那么这种添加方式可以包含所有匹配到i的情况。

但是这里H不一定等于p[i]：我们经过每个非单词节点的概率也不一定一样；

并且同时，由于我们之前已经匹配了状态为S的一些字符，我们可能在添加不到m个字符时就匹配到了某个串j（j可以等于i）

如果在某个节点加上字符串i的前k个字符后就已经到达了字符串j的终止节点,那么j的后k个字符必然等于i的前k个字符.

在匹配上j后,(虽然继续匹配是非法的，但是我们要减去这种非法状态，所以我们还要计算这种状态发生的概率.）

我们还要继续生成字符使得接下来m-k的字符等于串i的后m-k个字符，也就是说，p[i]应该在H的基础上减去p[j]*(1/2)^m-k

这里的“前k个字符重叠”，我们可以利用KMP的失配指针来处理。

（其实这里的处理方法有很多，除了kmp，hash也可以，只要能找出两串的重叠部分即可）

那么最后，我们可以列出方程：p[i]=H-Σp[j]*(1/2)m-k，移项得p[i]+Σp[j]*(1/2)^m-k-H=0

再加上一开始的方程Σp[i]=1.0，我们就可以对n+1个变量列出n+1个方程来解方程了！

代码见下：

#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
using namespace std;
const int N=310;
int n,l,cnt,fail[N<<1];
char s[N][N],str[N<<1];
double f[N][N];
inline void kmp()
{
    register int i,j,len;
    for(i=2,j=0,len=(l<<1);i<=len;++i)
    {
        while(j&&str[j+1]!=str[i])j=fail[j];
        j=(str[j+1]==str[i])?j+1:j,fail[i]=j;
    }
}
void Swap(int a,int b)
    {for(register int i=1;i<=cnt+1;++i)swap(f[a][i],f[b][i]);}
void Execution(int a,int b,double t)
    {for(register int i=1;i<=cnt+1;++i)f[a][i]+=f[b][i]*t;}
inline void gauss()
{
    register int i,j,k;
    for(i=1;i<=cnt;++i)
    {
        for(j=i+1;j<=cnt;++j)if(fabs(f[i][i])<fabs(f[j][i]))Swap(i,j);
        for(j=1;j<=cnt;++j)if(j!=i)Execution(j,i,-f[j][i]/f[i][i]);
    }
    for(i=1;i<=cnt;++i)f[i][cnt+1]/=f[i][i];
}
int main()
{
    scanf("%d%d",&n,&l);register int i,j,k;
    for(i=1;i<=n;++i)scanf("%s",s[i]+1);
    for(i=1;i<=n;++i)
        for(j=1;j<=n;++j)
        {
            for(k=1;k<=l;++k)str[k]=s[i][k],str[l+k]=s[j][k];
            for(kmp(),k=fail[l<<1];k;k=fail[k])
                if(k<l)f[i][j]+=pow(0.5,l-k);
        }
    for(i=1;i<=n;++i)f[i][i]+=1.0,f[i][n+1]-=1.0;//n+1代表未知变量H
    for(i=1,f[n+1][n+2]=1.0;i<=n;++i)f[n+1][i]=1.0;
    cnt=n+1;gauss();
    for(i=1;i<=n;++i)printf("%.10lf\n",f[i][cnt+1]);
}