[高斯消元] BZOJ3640: JC的小苹果

不难想到的题。首先设 $f[i][j]$ 表示血量为 $i$，位置在 $j$ 的期望概率。
推出方程之后发现，当一个点 $a[i]=0$ 时，不能递推出来。所以要每层进行消元。
发现每次消元的系数矩阵都相同，只是常数不同。可以想到把常数看成一个多项式，先 $n^3$ 预处理，然后每层直接 $n^2$ 带入就好了。
复杂度 $O(n^3+n^2*hp)$

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int maxn=155,maxe=10005,maxhp=20005;
int n,m,N,hp,a[maxn],fir[maxn],nxt[maxe],son[maxe],tot;
struct data{
    double a[maxn];
    data(){ memset(a,0,sizeof(a)); }
    data operator - (const data &B)const{
        data c; for(int i=0;i<=n;i++) c.a[i]=a[i]-B.a[i]; return c; 
    }
    data operator * (const double &val)const{
        data c; for(int i=0;i<=n;i++) c.a[i]=a[i]*val; return c;
    }
    data operator / (const double &val)const{
        data c; for(int i=0;i<=n;i++) c.a[i]=a[i]/val; return c;
    }
} _K[maxn],Ans[maxn];
double K[maxn][maxn];
void Gs(){
    for(int i=1;i<=n;i++){
        int where; double _max=0; for(int j=i;j<=n;j++) if(K[j][i]>_max) where=j, _max=K[j][i];
        for(int j=1;j<=n;j++) swap(K[i][j],K[where][j]); swap(_K[i],_K[where]);
        for(int j=i+1;j<=n;j++){
            double t=K[j][i]/K[i][i];
            for(int k=1;k<=n;k++) K[j][k]-=t*K[i][k]; _K[j]=_K[j]-_K[i]*t;
        }
    }
    for(int i=n;i>=1;i--){
        Ans[i]=_K[i]*(-1);
        for(int j=i+1;j<=n;j++) Ans[i]=Ans[i]-Ans[j]*K[i][j];
        Ans[i]=Ans[i]/K[i][i];
    }
}
int d[maxn];
double ans,f[maxhp][maxn];
void add(int x,int y){
    son[++tot]=y; nxt[tot]=fir[x]; fir[x]=tot;
}
int main(){
    scanf("%d%d%d",&n,&m,&hp);
    if(hp==0) return printf("0.00000000"),0;
    for(int i=1;i<=n;i++) scanf("%d",&a[i]);
    for(int i=1;i<=m;i++){
        int x,y; scanf("%d%d",&x,&y);
        add(x,y); d[x]++; if(x!=y) add(y,x), d[y]++; 
    }
    for(int i=1;i<=n;i++){
        K[i][i]+=1; if(a[i]) continue; 
        if(i==1) _K[1].a[0]-=1;
        for(int j=fir[i];j;j=nxt[j]) if(son[j]!=n) K[i][son[j]]-=1.0/d[son[j]];
    }
    Gs();
    for(int i=1;i<=n;i++) f[hp][i]=Ans[i].a[0];
    ans=f[hp][n];
    memset(K,0,sizeof(K)); memset(Ans,0,sizeof(Ans)); memset(_K,0,sizeof(_K));
    for(int i=1;i<=n;i++){
        K[i][i]+=1;
        if(a[i]){ _K[i].a[i]-=1; continue; }
        for(int j=fir[i];j;j=nxt[j]) if(son[j]!=n) K[i][son[j]]-=1.0/d[son[j]];
    }
    Gs(); 
    for(int i=hp-1;i>=1;i--){
        for(int j=1;j<=n;j++) if(a[j]){
            if(i+a[j]>hp) continue;
            for(int k=fir[j];k;k=nxt[k]) if(son[k]!=n) f[i][j]+=f[i+a[j]][son[k]]/d[son[k]];
        }
        for(int j=1;j<=n;j++) if(!a[j]){
            for(int k=1;k<=n;k++) f[i][j]+=Ans[j].a[k]*f[i][k];
        }
        ans+=f[i][n];
    }
    printf("%.8lf\n",ans);
//  for(int i=hp;i>=1;i--){
//          for(int j=1;j<=n;j++) printf("%.3lf ",f[i][j]);
//          printf("\n");
//      }   
    return 0;
}