洛谷P3343 [ZJOI2015]地震后的幻想乡 [DP，概率期望]
思路
题目给了一个提示：对于\(n\)个\([0,1]\)的随机变量，其中第\(k\)小的期望大小是\(\frac{k}{n+1}\)。
这引导我们枚举边的相对大小的全排列，然后求最小生成树
设\(P(x)\)表示最小生成树中最大一条边的排名是\(x\)的概率，那么有
\[ans=\frac 1 {m+1}\sum_x xP(x) \Leftrightarrow (m+1)ans=\sum_x xP(x)
\]
恰好是\(x\)比较麻烦，再设\(f_x\)表示最大排名大于\(x\)的概率，那么\(P(x)=f_{x-1}-f_x\)。
于是有
\[(m+1)ans=\sum_{x=1}^m x(f_{x-1}-f_x)=\sum_{x=0}^{m-1} f_x
\]
然而，由于概率容易被卡精度，将\(f_x\)的意义改为最大排名大于\(x\)的方案数，即用排名小于等于\(x\)的边无法使图连通的方案数。
再由于排名随机，这也就等价于选\(x\)条边使得图不连通的方案数。
然后就可以开心地DP了。
记\(f_{S,i}\)表示点集为\(S\)，选了\(i\)条边时，这些点不连通的方案数。
而\(g_{S,i}\)表示连通的方案数。
根据套路，有
\[f_{S,i}+g_{S,i}={cnt_S\choose i}\\
f_{S,i}=\sum_{k\in T \subset S} \sum_{j} g_{T,j}{cnt_{S-T}\choose i-j}
\]
其中\(cnt_S\)表示点集\(S\)中边数，\(k\)为\(S\)中编号最小的点。
最后计算答案：
\[ans=\frac{1}{m+1} \sum_{i=0}^{m-1} \frac {f_{U,i}} {m \choose i}
\]
代码
#include<bits/stdc++.h>
clock_t t=clock();
namespace my_std{
using namespace std;
#define pii pair<int,int>
#define fir first
#define sec second
#define MP make_pair
#define rep(i,x,y) for (int i=(x);i<=(y);i++)
#define drep(i,x,y) for (int i=(x);i>=(y);i--)
#define go(x) for (int i=head[x];i;i=edge[i].nxt)
#define templ template<typename T>
#define sz 150
#define S 40000
typedef long long ll;
typedef double db;
mt19937 rng(chrono::steady_clock::now().time_since_epoch().count());
templ inline T rnd(T l,T r) {return uniform_int_distribution<T>(l,r)(rng);}
templ inline bool chkmax(T &x,T y){return x<y?x=y,1:0;}
templ inline bool chkmin(T &x,T y){return x>y?x=y,1:0;}
templ inline void read(T& t)
{
t=0;char f=0,ch=getchar();double d=0.1;
while(ch>'9'||ch<'0') f|=(ch=='-'),ch=getchar();
while(ch<='9'&&ch>='0') t=t*10+ch-48,ch=getchar();
if(ch=='.'){ch=getchar();while(ch<='9'&&ch>='0') t+=d*(ch^48),d*=0.1,ch=getchar();}
t=(f?-t:t);
}
template<typename T,typename... Args>inline void read(T& t,Args&... args){read(t); read(args...);}
void file()
{
#ifndef ONLINE_JUDGE
freopen("a.in","r",stdin);
#endif
}
inline void chktime()
{
#ifndef ONLINE_JUDGE
cout<<(clock()-t)/1000.0<<'\n';
#endif
}
#ifdef mod
ll ksm(ll x,int y){ll ret=1;for (;y;y>>=1,x=x*x%mod) if (y&1) ret=ret*x%mod;return ret;}
ll inv(ll x){return ksm(x,mod-2);}
#else
ll ksm(ll x,int y){ll ret=1;for (;y;y>>=1,x=x*x) if (y&1) ret=ret*x;return ret;}
#endif
// inline ll mul(ll a,ll b){ll d=(ll)(a*(double)b/mod+0.5);ll ret=a*b-d*mod;if (ret<0) ret+=mod;return ret;}
}
using namespace my_std;
int n,m;
pii edge[sz*sz];
int cnt[S];
db f[S][sz],g[S][sz];
db C[sz][sz];
void init(){rep(i,0,50) C[i][0]=1;rep(i,1,50) rep(j,1,50) C[i][j]=C[i-1][j-1]+C[i-1][j];}
int main()
{
file();
init();
read(n,m);
int x,y;
rep(i,1,m) read(x,y),edge[i]=MP(x,y);
rep(id,0,(1<<n)-1) rep(i,1,m) if (((1<<(edge[i].fir-1))&id) && ((1<<(edge[i].sec-1))&id)) ++cnt[id];
rep(i,1,n) f[1<<(i-1)][0]=0,g[1<<(i-1)][0]=1;
#define lowbit(x) (x&(-x))
rep(id,1,(1<<n)-1) if (id!=lowbit(id)) rep(i,0,cnt[id])
{
x=lowbit(id);
for (int k=(id-1)&id;k;k=(k-1)&id)
if (k&x)
rep(j,0,min(cnt[k],i))
f[id][i]+=g[k][j]*C[cnt[id^k]][i-j];
g[id][i]=C[cnt[id]][i]-f[id][i];
}
db ans=0;
rep(i,0,m) ans+=f[(1<<n)-1][i]/C[m][i];
ans/=(m+1);
printf("%.6lf",ans);
return 0;
}
