抗拒黄泉

题目大意

一个n*m的网格上，有一些格子可操作。
每一轮等概率选择一个可操作的格子操作，若选择了(i,j)，那么第i行与第j行会被标记为操作过。
求期望轮数使得所有行和所有列均操作过。
n,m<=20，n*m<=200

容斥

我们用P(i)表示已经过了i轮，还没有使得所有行列均操作过的概率。
那么答案即为
$\sum_{i=0}^\infty(P(i)-P(i+1))*(i+1)$
其中P(i)-P(i+1)表示第i+1轮所有行列均操作过的概率
展开这个式子得
$\sum_{i=0}^\infty P(i)$
$\sum_{i=0}^\infty (\sum_{|S|是奇数}(1-p(S))^i-\sum_{|S|是偶数}(1-p(S))^i)$
S表示二进制状态代表哪些行列一定不被操作过。
这显然就是一个容斥，我限定了一些行列一定不被操作过，但是其实不只这些行列不被操作过，所以是至少|S|个行列没被操作过，那么就是个简单容斥。
p(S)代表什么？代表每一轮选择的可操作格在S中不被操作的行列中的概率。因为我们不希望这些行列被操作，所以每一轮都不能选择这些可操作格，共进行了i轮，那么就是上面写的那样。
我们交换主体，就是先枚举S再枚举i，就变成了等比数列求和！
设q=1-p(S)
那我们要求$\sum_{i=0}^\infty q^i$
设$t=\sum_{i=0}^\infty q^i$
那么$qt=\sum_{i=1}^\infty q^i$
相减移项得$t=\frac{q^\infty -1}{q-1}=\frac{-1}{q-1}=\frac{1}{1-q}$
把q用1-p(S)代入得答案是$\frac{1}{p(S)}$
那么想一个暴力算法，可以2^(n+m)枚举集合，然后就可以统计了。
这显然过不了，我们思考优化。

背包DP

我们可以只枚举行的二进制状态，然后同样是搜索，但是我们记忆化，这里就直接说dp。
我们的贡献和什么有关？不可选的可操作格的奇偶性以及p(S)，而$p(S)=\frac{k}{tot}$，tot表示总共的可操作格数量，k表示S中的可操作格数量。
因此我们需要知道可操作格的数量。
对列做背包dp，设dpi,j,k表示做完了前i列，奇偶性为j，可操作格的数量为k的方案数。
那么两种转移，第i+1列是否一定不能被操作。
最后算一发即可。具体见代码。

#include<cstdio>
#include<algorithm>
#define fo(i,a,b) for(i=a;i<=b;i++)
using namespace std;
typedef double db;
const int maxn=20+10,maxtot=200+10;
int dp[maxn][2][maxtot],num[maxn];
bool bz[maxn][maxn],pd[maxn][maxn],fl[maxn];
int i,j,k,l,t,n,m,tot,top,cnt;
db ans;
void calc(){
    int i,j,k;
    fo(i,0,m)
        fo(j,0,1)
            fo(k,0,tot)
                dp[i][j][k]=0;
    dp[0][cnt%2][top]=1;
    fo(i,0,m-1)
        fo(j,0,1)
            fo(k,top,tot)
                if (dp[i][j][k]){
                    dp[i+1][1-j][k+num[i+1]]+=dp[i][j][k];
                    dp[i+1][j][k]+=dp[i][j][k];
                }
    fo(j,0,1)
        fo(k,top,tot)
            if (k&&dp[m][j][k]){
                if (j%2) ans+=((db)tot/(db)k)*dp[m][j][k];else ans-=((db)tot/(db)k)*dp[m][j][k];
            }
}
void dfs(int x){
    if (x==n+1){
        calc();
        return;
    }
    int i;
    fl[x]=1;
    cnt++;
    fo(i,1,m)
        if (bz[x][i]){
            num[i]--;
            top++;
        }
    dfs(x+1);
    cnt--;
    fo(i,1,m)
        if (bz[x][i]){
            num[i]++;
            top--;
        }
    fl[x]=0;
    dfs(x+1);
}
int main(){
    freopen("refuse.in","r",stdin);freopen("refuse.out","w",stdout);
    scanf("%d%d",&n,&m);
    fo(i,1,n)
        fo(j,1,m){
            scanf("%d",&t);
            bz[i][j]=t;
            tot+=t;
        }
    if (n>m){
        fo(i,1,n)
            fo(j,1,m)
                pd[j][i]=bz[i][j];
        swap(n,m);
        fo(i,1,n)
            fo(j,1,m)
                bz[i][j]=pd[i][j];
    }
    fo(i,1,n)
        fo(j,1,m)
            if (bz[i][j]) num[j]++;
    dfs(1);
    printf("%.5lf\n",ans);
}