BZOJ2707 : [SDOI2012]走迷宫

#include<cstdio>
#include<cmath>
const int N=10010,M=1000010;
int n,m,x,y,i,j,S,T;
int g[3][N],nxt[3][M],v[3][M],ed,G[N],NXT[N],V[N],d[N],q[N],f[N],h,t,cnt,id[N],tot;bool vis[N],inf[N];
double e[N],a[110][110],ans[110],tmp;
inline void read(int&a){char c;while(!(((c=getchar())>='0')&&(c<='9')));a=c-'0';while(((c=getchar())>='0')&&(c<='9'))(a*=10)+=c-'0';}
inline void add(int x,int y){
  v[0][++ed]=y;nxt[0][ed]=g[0][x];g[0][x]=ed;
  v[1][ed]=x;nxt[1][ed]=g[1][y];g[1][y]=ed;
}
inline void ADD(int x,int y){d[y]++,v[2][++ed]=y;nxt[2][ed]=g[2][x];g[2][x]=ed;}
inline void addnode(int x,int y){V[++ed]=y;NXT[ed]=G[x];G[x]=ed;}
void dfs1(int x){
  vis[x]=1;
  for(int i=g[0][x];i;i=nxt[0][i])if(!vis[v[0][i]])dfs1(v[0][i]);
  q[++t]=x;
}
void dfs2(int x,int y){
  vis[x]=0,addnode(f[x]=y,x);
  for(int i=g[1][x];i;i=nxt[1][i])if(vis[v[1][i]])dfs2(v[1][i],y);
}
inline void cal(int x){
  int i,j,k,u,y,cnt;double t;
  for(tot=0,i=G[x];i;i=NXT[i])id[V[i]]=++tot;
  for(cnt=0,i=G[x];i;i=NXT[i]){
    for(u=V[i],++cnt,j=1;j<=tot+1;j++)a[cnt][j]=0;
    if(u==T){a[cnt][cnt]=1;continue;}
    for(j=g[0][u];j;j=nxt[0][j]){
      if(inf[f[y=v[0][j]]]){inf[x]=1;return;}
      a[cnt][cnt]+=1,a[cnt][tot+1]+=1;
      if(f[y]==x)a[cnt][id[y]]-=1;else a[cnt][tot+1]+=e[y];
    }
  }
  for(i=1;i<=tot;i++){
    for(k=i,j=i+1;j<=tot;j++)if(std::fabs(a[j][i])>std::fabs(a[k][i]))k=j;
    if(k!=i)for(j=i;j<=tot+1;j++)tmp=a[i][j],a[i][j]=a[k][j],a[k][j]=tmp;
    for(j=i+1;j<=tot;j++)for(t=a[j][i]/a[i][i],k=i;k<=tot+1;k++)a[j][k]-=a[i][k]*t;
  }
  for(ans[tot]=a[tot][tot+1]/a[tot][tot],i=tot-1;i;i--){
    for(ans[i]=a[i][tot+1],j=tot;j>i;j--)ans[i]-=ans[j]*a[i][j];
    ans[i]/=a[i][i];
  }
  for(i=G[x];i;i=NXT[i])e[V[i]]=ans[id[V[i]]];
}
int main(){
  read(n),read(m),read(S),read(T);
  while(m--)read(x),read(y),add(x,y);
  for(i=1;i<=n;i++)if(!vis[i])dfs1(i);
  for(ed=0,i=n;i;i--)if(vis[q[i]])dfs2(q[i],++cnt);
  for(ed=0,i=1;i<=n;i++)for(j=g[0][i];j;j=nxt[0][j])if(f[i]!=f[v[0][j]])ADD(f[v[0][j]],f[i]);
  for(i=1;i<=cnt;i++)if(!d[i])if(f[T]!=i)inf[i]=1;
  for(h=i=1,t=0;i<=cnt;i++)if(!d[i]){
    q[++t]=i;
    if(f[T]==i)cal(i);
  }
  while(h<=t)for(i=g[2][x=q[h++]];i;i=nxt[2][i])if(!(--d[v[2][i]]))cal(q[++t]=v[2][i]);
  if(inf[f[S]])puts("INF");else printf("%.3f",e[S]);
  return 0;
}