[期望 状压DP 子集DP 多项式] BZOJ 3925 [Zjoi2015]地震后的幻想乡

Po姐说的肯定比我清楚多

http://blog.csdn.net/popoqqq/article/details/44858691

#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<vector>
using namespace std;
typedef long long ll;
typedef __float128 ld;

inline char nc(){
  static char buf[100000],*p1=buf,*p2=buf;
  if (p1==p2) { p2=(p1=buf)+fread(buf,1,100000,stdin); if (p1==p2) return EOF; }
  return *p1++;
}

inline void read(int &x){
  char c=nc(),b=1;
  for (;!(c>='0' && c<='9');c=nc()) if (c=='-') b=-1;
  for (x=0;c>='0' && c<='9';x=x*10+c-'0',c=nc()); x*=b;
}

struct Poly{
  int n; vector<ll> a;
  Poly(){ }
  Poly(int in){
    if (in>=0) { n=in; a.resize(n+1); for (int i=0;i<=n;i++) a[i]=0; }
    if (in==-1) n=0,a.resize(1),a[0]=1;
  }
  ll &operator [](int x){
    return a[x];
  }
  friend Poly operator + (const Poly &A,const Poly &B){
    Poly ret=Poly(max(A.n,B.n));
    for (int i=0;i<=A.n;i++) ret[i]+=A.a[i];
    for (int i=0;i<=B.n;i++) ret[i]+=B.a[i];
    return ret;
  }
  friend Poly operator - (const Poly &A,const Poly &B){
    Poly ret=Poly(max(A.n,B.n));
    for (int i=0;i<=A.n;i++) ret[i]+=A.a[i];
    for (int i=0;i<=B.n;i++) ret[i]-=B.a[i];
    return ret;
  }
  friend Poly operator * (const Poly &A,const Poly &B){
    Poly ret=Poly(A.n+B.n);
    for (int i=0;i<=A.n;i++) for (int j=0;j<=B.n;j++) ret[i+j]+=A.a[i]*B.a[j];
    return ret;
  }
};

const int N=15;

int n,m,a[N];
Poly pow[N*N],f[1<<10],Ans;
int cnt[1<<10],bin[1<<10];

#define low(x) (bin[(x)&-(x)])

int main(){
  int iu,iv;
  freopen("mst.in","r",stdin);
  freopen("mst.out","w",stdout);
  read(n); read(m);
  for (int i=1;i<(1<<n);i++) cnt[i]=cnt[i^(i&-i)]+1;
  for (int i=1;i<=n;i++) bin[1<<(i-1)]=i;
  for (int i=1;i<=m;i++)
    read(iu),read(iv),a[iu]+=1<<(iv-1),a[iv]+=1<<(iu-1);
  pow[0]=Poly(-1); pow[1]=Poly(1),pow[1][0]=1,pow[1][1]=-1;
  for (int i=2;i<=m;i++) pow[i]=pow[i-1]*pow[1];
  for (int i=1;i<(1<<n);i++){
    f[i]=Poly(-1);
    int rest=i^(1<<(low(i)-1));
    for (int j=rest;j;(--j)&=rest){
      int cot=0;
      for (int t=i^j;t;t-=t&-t)
	cot+=cnt[a[low(t)]&j];
      f[i]=f[i]-f[i^j]*pow[cot];
    }
  }
  Ans=Poly(-1)-f[(1<<n)-1];
  ld ans=0;
  for (int i=0;i<=Ans.n;i++) ans+=(ld)Ans[i]/(i+1);
  printf("%.6lf\n",(double)ans);
  return 0;
}