【UVA】【10828】随机程序

时间：2023-03-08 20:05:22

数学期望/高斯消元/马尔可夫过程

　　刘汝佳老师白书上的例题- -b

　　本体不满足拓扑关系，但马尔可夫过程是可以高斯消元解的……

　　用「高斯·约当消元」更方便！

 //UVA 10828

 #include<cmath>

 #include<vector>

 #include<cstdio>

 #include<cstring>

 #include<cstdlib>

 #include<iostream>

 #include<algorithm>

 #define rep(i,n) for(int i=0;i<n;++i)

 #define F(i,j,n) for(int i=j;i<=n;++i)

 #define D(i,j,n) for(int i=j;i>=n;--i)

 #define pb push_back

 using namespace std;

 int getint(){

     int v=,sign=; char ch=getchar();

     while(ch<''||ch>''){ if (ch=='-') sign=-; ch=getchar();}

     while(ch>=''&&ch<=''){ v=v*+ch-''; ch=getchar();}

     return v*=sign;

 }

 const int N=;

 const double eps=1e-;

 typedef double Matrix[N][N];

 /******************tamplate*********************/

 void gauss_jordan(Matrix A,int n){

     int i,j,k,r;

     rep(i,n){

         r=i;

         for(j=i+;j<n;++j)

             if (fabs(A[j][i]) > fabs(A[r][i])) r=j;

         if (fabs(A[r][i]) < eps)continue;

         if (r!=i) F(j,,n) swap(A[r][j],A[i][j]);

         rep(k,n) if (k!=i)

             D(j,n,i) A[k][j]-=A[k][i]/A[i][i]*A[i][j];

     }

 }

 Matrix A;

 int n,d[N];

 vector<int> pre[N];

 bool inf[N];

 int main(){

 #ifndef ONLINE_JUDGE

     freopen("10828.in","r",stdin);

     freopen("10828.out","w",stdout);

 #endif

     int cs=;

     while(scanf("%d",&n)== && n){

         memset(d,,sizeof d);

         rep(i,n) pre[i].clear();

         int a,b;

         while(scanf("%d%d",&a,&b)== && a && b){

             a--; b--;

             d[a]++;

             pre[b].pb(a);

         }

         memset(A,,sizeof (A));

         rep(i,n){

             A[i][i]=;

             rep(j,pre[i].size())

                 A[i][pre[i][j]]-=1.0/d[pre[i][j]];

             if (i==) A[i][n]=;

         }

         gauss_jordan(A,n);

         memset(inf,,sizeof inf);

         for(int i=n-;i>=;i--){

             if ( fabs(A[i][i])<eps && fabs(A[i][n])>eps ) inf[i]=;

             for(int j=i+;j<n;j++)

                 if (fabs(A[i][j])>eps && inf[j]) inf[i]=;

         }

         int q,u;

         scanf("%d",&q);

         printf("Case #%d:\n",++cs);

         while(q--){

             scanf("%d",&u); u--;

             if (inf[u]) printf("infinity\n");

             else printf("%.3lf\n",fabs(A[u][u]) < eps ? 0.0 : A[u][n]/A[u][u]);

         }

     }

     return ;

 }

相关文章