UVa10828 Back to Kernighan-Ritchie——概率转移&&高斯消元法

题意

给出一个程序控制流图，从每个结点出发到每个后继接结点的概率均相等。当执行完一个没有后继的结点后，整个程序终止。程序总是从编号为1的结点开始。你的任务是对于若干个查询结点，求出每个结点的期望执行次数。结点个数 $n < 10$.

分析

如果是有向无环图，可以直接解出递推关系，再采用记忆化搜索求解。

当这题可能有环，只能列出方程，用高斯消元解方程组。

设结点 $i$ 的出度为 $d_I$，期望的执行次数为 $x_i$。对于一个拥有三个前驱结点 $a, b, c$ 的结点 $i$，可以列出方程 $x_i = x_a / d_a + x_b / d_b + x_c / d_c$.

如果 $x_i$为无穷大或0，通过代数方法会出错的，所以我们结合实际意义考虑，

哪些结点期望执行无数次呢？就是那些无法到达终态的结点（即一直在非终态循环）；

哪些结点期望执行次数为0次呢？就是那些从起点出发到不了的点。

我们可以先用Floyd求出传递闭包，先找到无穷大点和零点，

如果使用高斯-约当消元法可以避免这一预处理，当 $A[i][i] = A[i][n] = 0$ 时 $x_i=0$，而 $A[i][i] = 0$ 但 $A[i][n] > 0$ 时 $x_i$ 为正无穷大（这个结论看似显然

#include<cstdio>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<vector>
using namespace std;

const double eps = 1e-8;
const int maxn = 100+10;
typedef double Matrix[maxn][maxn];

//结果为A[i][n]/A[i][i]
void gauss_jordan(Matrix A, int n)
{
    int i, j, k, r;
    for(i = 0;i < n;i++)
    {
         //选绝对值一行r并与第i行交换
        r = i;
        for(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(j = 0;j <= n;j++)  swap(A[r][j], A[i][j]);

        //与除第i行外的其他行进行消元
        for(k = 0;k < n;k++)  if(k != i)
            for(j = n;j >= i;j--)  A[k][j] -= A[k][i] / A[i][i] * A[i][j];
    }
}

Matrix A;
int n, d[maxn];
vector<int>pre[maxn];
int inf[maxn];

int main()
{
    int kase = 0;
    while(scanf("%d", &n) == 1 && n)
    {
        memset(d, 0, sizeof(d));
        for(int i = 0;i < n;i++)  pre[i].clear();

        int a, b;
        while(scanf("%d%d", &a, &b) == 2 && a)
        {
            a--; b--;   //编号从0开始
            d[a]++;     //结点a的出度加1
            pre[b].push_back(a);
        }

        //构造方程组
        memset(A, 0, sizeof(A));
        for(int i = 0;i < n;i++)
        {
            A[i][i] = 1;
            for(int j = 0;j < pre[i].size();j++)
                A[i][pre[i][j]] -= 1.0 / d[pre[i][j]];
            if(i == 0)  A[i][n] = 1;
        }

        //解方程组，标记无穷变量
        gauss_jordan(A, n);
        memset(inf, 0, sizeof(inf));
        for(int i = n-1; i >= 0;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;
        }

        int q, u;
        scanf("%d", &q);
        printf("Case #%d:
", ++kase);
        while(q--)
        {
            scanf("%d", &u);
            u--;
            if(inf[u])  printf("infinity
");
            else  printf("%.3f
", fabs(A[u][u]) < eps ? 0.0 : A[u][n] / A[u][u]);
        }
    }
    return 0;
}

第一次在UVA官网交上题hh，题目链接

From:

《算法竞赛入门经典训练指南》——刘汝佳、陈锋著