高斯-约当消元法（随机程序，UVA 10828）

就是有一类题目：

就是给你一个图，和一个起始点（可以很抽象，比如UVA 11762需要你自己建模。也可以很具体，比如就是一个网格）

然后每个节点都有一定概率转移到另一些节点。

具体解法就是：

如果是有向无环图就用动态规划或者记忆化搜索。

否则就是高斯消元。

本题就是高斯消元。

然后如果使用高斯-约当消元法就可以省去回代的过程。

大白书上的代码是用double类型的，虽然有一些技巧可以提高数值稳定性，但是多少会有精度损失的。想要完全准确得重写分数运算。

double类型有一个技巧技巧就是定义一个const double eps=1e-8;如果小于eps就当做0，否则就当做常数。

还有求解的过程，比如ax=b

如果a==0，b！=0，那就是无解，此时可以当做x=inf。只需开一个数组记录一下即可。

如果a==b==0，那解就是0。

否则答案就是b/a。

然后所有跟inf关联的量都是inf（跟无解的关联就是无解）。

输出时特判输出就好了。

代码

#include<bits/stdc++.h>
#define maxn 110
using namespace std;
const double eps=1e-8;

typedef double MAT[maxn][maxn];

void gsyd(MAT A,int n)
{
    for(int i=0;i<n;i++)
    {
        int r=i;
        for(int j=i+1;j<n;j++)
            if(fabs(A[j][i])>fabs(A[r][i]))
                r=j;
        if(fabs(A[r][i])<eps) continue;
        if(r!=i) for(int j=0;j<=n;j++) swap(A[r][j],A[i][j]);
        for(int j=0;j<n;j++) if(j!=i)
            for(int k=n;k>=i;k--) A[j][k]-=A[j][i]/A[i][i]*A[i][k];
    }
}

int n;
vector<int>pre[maxn];
int od[maxn];
int inf[maxn];
int kase;

int main()
{
    while(scanf("%d",&n)==1&&n)
    {
        for(int i=0;i<n;i++) pre[i].clear();
        memset(od,0,sizeof(od));
        int u,v;
        while(scanf("%d %d",&u,&v)==2&&(u+v))
        {
            u--;v--;
            pre[v].push_back(u);
            od[u]++;
        }
        MAT A;
        memset(A,0,sizeof(A));
        memset(inf,0,sizeof(inf));
        A[0][n]=1;
        for(int i=0;i<n;i++)
        {
            A[i][i]=1;
            for(unsigned int j=0;j<pre[i].size();j++)
                A[i][pre[i][j]]-=1.0/od[pre[i][j]];
        }
        gsyd(A,n);
        for(int i=n-1;~i;i--)
        {
            if(fabs(A[i][i])<eps&&fabs(A[i][n])>eps) inf[i]=1;
            for(int j=i+1;j<n;j++) if(fabs(A[i][j])>eps&&inf[j]) inf[i]=1;
        }
        printf("Case #%d:\n",++kase);
        int q;
        scanf("%d",&q);
        while(q--)
        {
            int x;
            scanf("%d",&x);
            x--;
            if(inf[x]) puts("infinity");
            else printf("%.3lf\n",abs(A[x][x])<eps?0.0:A[x][n]/A[x][x]);
        }
    }
    return 0;
}