[ SHOI 2014 ] 概率充电器

(\)

(Description)

---

一个含(N)个元器件的树形结构充电器，第(i)个元器件有(P_i)的概率直接从外部被充电，连接(i,j)的边有(P_{i,j})的概率导电，元器件只有外部充电和从已充电元器件导电两种方式被充电，求最后被充电的元器件个数的期望。

• (Nin [0,5 imes10^5])，(P_i,P_{ij}in [0,1])

(\)

(Solution (1))

---

• 期望(=)贡献( imes)概率，而个数的贡献只能为(1)，所以被点亮元器件个数的期望(=)每个元器件被点亮概率之和。

• 因为是无根树，所以要换根进行两次(DP)，假设树的形态定为以(1)为根。设(f[i])表示只考虑当前节点和所有当前节点子树的节点，当前节点被点亮的概率，设(g[i])表示只考虑除掉当前节点和所有当前节点子树的节点剩下的节点，当前节点被点亮的概率。两遍(DFS)分别处理(f[i])和(g[i])。

• 考虑(f[i])的求法，应当是自底向上的。首先当前节点自己点亮的概率显然是(P_i)，剩下需要考虑的就是子树对根节点的影响。需要注意的是，不同子树对根的影响虽然是"或"的关系，但却不是互斥的。也就是说，如果暴力的累加每一个子树的贡献，多个子树同时做出贡献的部分会被计算多次，因为直接累加一棵子树的贡献正确的前提是其他子树都不能做出贡献，即都不能将根节点点亮。一种做法是容斥计算，还有一种更简便的方法是正难则反。注意到当前点被点亮的概率(=1-)不被点亮的概率，而不被点亮的条件是"与"的关系，也就是说不被点亮的部分是可以连乘得到的，于是有：

  [f[i]=1-(1-P_i) imes prod_{vin son[i]}(1-f[v] imes P_{i,v}) ]

• 对于(g[i])的求法就是换根那一套了。但是这次只有(1)号节点的(g)是确定为(0)，所以自上向下更新。考虑对于一个子节点，其(g)值可以表示为两部分，被父节点的(g)部分点亮，或是父节点被其他子树点亮，再点亮当前节点，同样的，它们之间是"或"的关系，容斥起来非常麻烦。转化之后只需考虑如何从父节点的(f)值中去掉当前节点的部分，发现在累乘时关于当前节点的一项是独立的，可以直接去掉，于是有：

  [g[i]=ig(1-(1-g[fa]) imes frac{1-f[fa]}{1-f[i] imes P_{i,fa}}ig) imes P_{i,fa} ]

  注意去掉当前点贡献时，如果当前点原来就没有贡献，算得的分母是(0)，所以应当忽略贡献，分母变为(1)。

• 统计答案的时候是一样的道理，两者是"或"的关系，转化为"与"的关系后再取补即可：

  [ans=sum_{i=1}^N 1-(1-f[i]) imes (1-g[i]) ]

(\)

(Code (1))

---

#include<cmath>
#include<cstdio>
#include<cctype>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define N 500010
#define R register
#define gc getchar
using namespace std;

int n,tot,hd[N],p[N];
double ans,tmp,f[N],g[N];
struct edge{int to,nxt,p;}e[N<<1];

inline void add(int u,int v,int p){
  e[++tot].to=v; e[tot].p=p;
  e[tot].nxt=hd[u]; hd[u]=tot;
}

inline int rd(){
  int x=0; bool f=0; char c=gc();
  while(!isdigit(c)){if(c=='-')f=1;c=gc();}
  while(isdigit(c)){x=(x<<1)+(x<<3)+(c^48);c=gc();}
  return f?-x:x;
}

inline void dfs1(int u,int fa){
  f[u]=0.01*(100-p[u]);
  for(R int i=hd[u],v;i;i=e[i].nxt)
    if((v=e[i].to)!=fa){
      dfs1(v,u);
      f[u]*=(1.0-f[v]*(0.01*e[i].p));
    }
  f[u]=1.0-f[u];
}

inline void dfs2(int u,int fa){
  for(R int i=hd[u],v;i;i=e[i].nxt)
    if((v=e[i].to)!=fa){
      double pf=(1.0-f[v]*(0.01*e[i].p));
      if(pf<1e-8) pf=1.0;
      g[v]=(1.0-(1.0-g[u])*(1.0-f[u])/pf)*(0.01*e[i].p); 
      dfs2(v,u);
    }
}

int main(){
  n=rd();
  for(R int i=1,u,v,w;i<n;++i){
      u=rd(); v=rd(); w=rd();
      add(u,v,w); add(v,u,w);
  }
  for(R int i=1;i<=n;++i) p[i]=rd();
  dfs1(1,-1);
  dfs2(1,-1);
  for(R int i=1;i<=n;++i) ans+=1.0-(1.0-f[i])*(1.0-g[i]);
  printf("%.6lf
",ans);
  return 0;
}

(\)

(Solution (2))

---

• 同样的转化问题二次扫描，但状态都从被点亮变为不被点亮，换个角度再推一遍，设(f[i])表示只考虑当前节点和所有当前节点子树的节点，当前节点不被点亮的概率，设(g[i])表示只考虑除掉当前节点和所有当前节点子树的节点剩下的节点，当前节点不被点亮的概率。

• 对于(f[i])，需要注意的是不能直接借鉴上面的写法，因为转移自子树的(f)的定义有变化。一个点不被点亮，当且仅当自己不会被外界直接点亮，也不被子树内的点点亮。不被子树内的点点亮需要分两种情况，子树直接不被点亮，或子树被点亮而边不连通，于是有：

  [f[i]=(1-P_i) imes prod_{vin son[i]}ig(f[v]+（1-f[v]) imes(1- P_{i,v})ig) ]

• 对于(g[i])，同样需要考虑父节点是否点亮，设父节点不被点亮的概率为(Q_{fa})，则有：

  [g[i]=Q_{fa}+(1-Q_{fa}) imes P_{i,fa} ]

  对与(Q_{fa})的求法跟上一种方案是一样的，讨论分成(g)和(f)两部分，有：

  [Q_{fa}=g[fa] imesfrac{f[fa]}{f[i]+(1-f[i]) imes(1-P_{i,fa})} ]

• 可以发现，反过来设计状态所需要处理的式子就容易了很多，对答案有：

  [ans=sum_{i=1}^N 1-f[i] imes g[i] ]

(\)

(Code (2))

---

#include<cmath>
#include<cstdio>
#include<cctype>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#define N 500010
#define R register
#define gc getchar
using namespace std;

double ans,tmp,f[N],g[N];
int n,tot,hd[N],p[N],fp[N];
struct edge{int to,nxt,p;}e[N<<1];

inline void add(int u,int v,int p){
  e[++tot].to=v; e[tot].p=p;
  e[tot].nxt=hd[u]; hd[u]=tot;
}

inline int rd(){
  int x=0; bool f=0; char c=gc();
  while(!isdigit(c)){if(c=='-')f=1;c=gc();}
  while(isdigit(c)){x=(x<<1)+(x<<3)+(c^48);c=gc();}
  return f?-x:x;
}

inline void dfs1(int u,int fa){
  f[u]=0.01*(100-p[u]);
  for(R int i=hd[u],v;i;i=e[i].nxt)
    if((v=e[i].to)!=fa){
      dfs1(v,u); fp[v]=e[i].p;
      f[u]*=(f[v]+(1.0-f[v])*(0.01*(100-e[i].p)));
    }
}

inline void dfs2(int u,int fa){
  if(fa==-1) g[u]=1.0;
  else{
    double pf=g[fa]*(f[fa]/(f[u]+((double)1-f[u])*(0.01*(100-fp[u]))));
    g[u]=pf+((double)1-pf)*(0.01*(100-fp[u]));
  }
  for(R int i=hd[u];i;i=e[i].nxt) if(e[i].to!=fa) dfs2(e[i].to,u);
}

int main(){
  n=rd();
  for(R int i=1,u,v,w;i<n;++i){
      u=rd(); v=rd(); w=rd();
      add(u,v,w); add(v,u,w);
  }
  for(R int i=1;i<=n;++i) p[i]=rd();
  dfs1(1,-1); dfs2(1,-1);
  for(R int i=1;i<=n;++i) ans+=1.0-f[i]*g[i];
  printf("%.6lf
",ans);
  return 0;
}