bzoj1778

高斯消元+矩阵的逆

来自popoqqq大神

求矩阵的逆：把I-T放在左边，P/Q*S放在右边，这样就形成了一个n*2n的矩阵，然后把左边高斯消元，右边就是求完逆的矩阵，其实就是ans，矩阵的逆跟乘法逆元是一样的，只不过是矩阵的逆元

然后输出a[i][n+1]，事实上矩阵只有n*(n+1)

构造转移概率矩阵是a[u][v]=1.0/d[v]*(1-p/q),就是v->u的概率乘上在v不爆炸的概率，我们想一想，假设我们从1->n,1->2,有1/d[1]的概率转移，并且不能爆炸才能走过去，要乘上(1-p/q),然后2->3，要乘上1/d[2]的概率走过去，再乘上(1-p/q)，不爆炸才能走过去，这就是转移的概率，每次矩阵自乘，就是b[i][j] += a[i][k]*a[k][j],求出又走一步的概率，最先开始S便是行向量，表示从1开始，还没走就爆炸的概率就是自己在这里爆炸，就是乘T^0,第一次转移就是乘上T,

#include<bits/stdc++.h>

using namespace std;

const double eps = 1e-;

const int N = ;

int n, m;

double a[N][N];

double p, q, t;

vector<int> G[N];

void gauss_jordan()

{

    a[][n + ] = t;

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

    {

        a[i][i] += 1.0;

        for(int j = ; j < G[i].size(); ++j)

        {

            int u = G[i][j];

            a[i][u] -= (1.0 - t) / (double)(G[u].size());

        }

    }

    for(int now = ; now <= n; ++now)

    {

        int x = now;

        for(int i = now; i <= n; ++i) if(fabs(a[i][now]) > fabs(a[x][now])) x = i;

        for(int i = ; i <= n + ; ++i) swap(a[now][i], a[x][i]);

        double t = a[now][now];

        for(int i = ; i <= n + ; ++i) a[now][i] /= t;

        for(int i = ; i <= n; ++i) if(fabs(a[i][now]) > eps && now != i)

        {

            t = a[i][now];

            for(int j = ; j <= n + ; ++j) a[i][j] -= t * a[now][j];

        }

    }

    for(int i = ; i <= n; ++i) printf("%.9f\n", a[i][n + ]);

}

int main()

{

    scanf("%d%d%lf%lf", &n, &m, &p, &q);

    t = p / q;

    for(int i = ; i <= m; ++i)

    {

        int u, v;

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

        G[u].push_back(v);

        G[v].push_back(u);

    }

    gauss_jordan();

    return ;

}

巴特西

bzoj1778

最新文章

热门文章