数学期望/高斯消元/马尔可夫过程
刘汝佳老师白书上的例题- -b
本体不满足拓扑关系，但马尔可夫过程是可以高斯消元解的……
用「高斯·约当消元」更方便！
1 //UVA 10828 2 #include<cmath> 3 #include<vector> 4 #include<cstdio> 5 #include<cstring> 6 #include<cstdlib> 7 #include<iostream> 8 #include<algorithm> 9 #define rep(i,n) for(int i=0;i<n;++i) 10 #define F(i,j,n) for(int i=j;i<=n;++i) 11 #define D(i,j,n) for(int i=j;i>=n;--i) 12 #define pb push_back 13 using namespace std; 14 int getint(){ 15 int v=0,sign=1; char ch=getchar(); 16 while(ch<'0'||ch>'9'){ if (ch=='-') sign=-1; ch=getchar();} 17 while(ch>='0'&&ch<='9'){ v=v*10+ch-'0'; ch=getchar();} 18 return v*=sign; 19 } 20 const int N=110; 21 const double eps=1e-8; 22 typedef double Matrix[N][N]; 23 /******************tamplate*********************/ 24 25 void gauss_jordan(Matrix A,int n){ 26 int i,j,k,r; 27 rep(i,n){ 28 r=i; 29 for(j=i+1;j<n;++j) 30 if (fabs(A[j][i]) > fabs(A[r][i])) r=j; 31 if (fabs(A[r][i]) < eps)continue; 32 if (r!=i) F(j,0,n) swap(A[r][j],A[i][j]); 33 34 rep(k,n) if (k!=i) 35 D(j,n,i) A[k][j]-=A[k][i]/A[i][i]*A[i][j]; 36 } 37 } 38 Matrix A; 39 int n,d[N]; 40 vector<int> pre[N]; 41 bool inf[N]; 42 43 int main(){ 44 #ifndef ONLINE_JUDGE 45 freopen("10828.in","r",stdin); 46 freopen("10828.out","w",stdout); 47 #endif 48 int cs=0; 49 while(scanf("%d",&n)==1 && n){ 50 memset(d,0,sizeof d); 51 rep(i,n) pre[i].clear(); 52 int a,b; 53 while(scanf("%d%d",&a,&b)==2 && a && b){ 54 a--; b--; 55 d[a]++; 56 pre[b].pb(a); 57 } 58 memset(A,0,sizeof (A)); 59 rep(i,n){ 60 A[i][i]=1; 61 rep(j,pre[i].size()) 62 A[i][pre[i][j]]-=1.0/d[pre[i][j]]; 63 if (i==0) A[i][n]=1; 64 } 65 66 gauss_jordan(A,n); 67 memset(inf,0,sizeof inf); 68 for(int i=n-1;i>=0;i--){ 69 if ( fabs(A[i][i])<eps && fabs(A[i][n])>eps ) inf[i]=1; 70 for(int j=i+1;j<n;j++) 71 if (fabs(A[i][j])>eps && inf[j]) inf[i]=1; 72 } 73 74 int q,u; 75 scanf("%d",&q); 76 printf("Case #%d: ",++cs); 77 while(q--){ 78 scanf("%d",&u); u--; 79 if (inf[u]) printf("infinity "); 80 else printf("%.3lf ",fabs(A[u][u]) < eps ? 0.0 : A[u][n]/A[u][u]); 81 } 82 } 83 return 0; 84 }