题意
Sol
设(f[i][j])表示Petya在(i),(Vasya)在(j)的概率,我们要求的是(f[i][i])
直接列方程高斯消元即可,由于每个状态有两维,因此时间复杂度为(O(n^6))
注意不能从终止节点转移而来
#include<bits/stdc++.h>
using namespace std;
const int MAXN = 2333;
inline int read() {
char c = getchar(); int x = 0, f = 1;
while(c < '0' || c > '9') {if(c == '-') f = -1; c = getchar();}
while(c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
return x * f;
}
int N, M, a, b, Lim;
double p[MAXN], f[MAXN][MAXN], E[MAXN][MAXN], deg[MAXN];
vector<int> v[MAXN];
int id[MAXN][MAXN], tot;
void Pre() {
f[id[a][b]][Lim + 1] = -1;
for(int i = 1; i <= N; i++) {
for(int j = 1; j <= N; j++) {
int now = id[i][j];
--f[now][now];//tag
if(i != j) f[now][now] += p[i] * p[j];
for(auto &x : v[i]) {
for(auto &y : v[j]) {
if(x == y) continue;
int nxt = id[x][y];
f[now][nxt] += (1.0 - p[x]) * (1.0 - p[y]) / deg[x] / deg[y];
}
}
for(auto &x : v[i]) {
int nxt = id[x][j];
if(x == j) continue;
f[now][nxt] += (1.0 - p[x]) * p[j] / deg[x];
}
for(auto &y : v[j]) {
int nxt = id[i][y];
if(i == y) continue;
f[now][nxt] += p[i] * (1.0 - p[y]) / deg[y];
}
}
}
}
void Gauss() {
for(int i = 1; i <= Lim; i++) {
int mx = i;
for(int j = i + 1; j <= Lim; j++) if(f[j][i] > f[mx][i] && f[j][i] != 0) swap(j, mx);
if(mx != i) swap(f[i], f[mx]);
// assert(fabs(f[i][i] < 1e-13));
for(int j = 1; j <= Lim; j++) {
if(i == j) continue;
double p = f[j][i] / f[i][i];
for(int k = i; k <= Lim + 1; k++) f[j][k] -= f[i][k] * p;
}
}
for(int i = 1; i <= Lim; i++) f[i][i] = f[i][Lim + 1] / f[i][i];
}
int main() {
N = read(); M = read(); a = read(); b = read(); Lim = N * N;
for(int i = 1; i <= M; i++) {
int x = read(), y = read();
v[x].push_back(y); v[y].push_back(x);
deg[x]++; deg[y]++;
}
for(int i = 1; i <= N; i++) scanf("%lf", &p[i]);
for(int i = 1; i <= N; i++)
for(int j = 1; j <= N; j++)
id[i][j] = ++tot;
Pre();
Gauss();
for(int i = 1; i <= N; i++) printf("%.10lf ", f[id[i][i]][id[i][i]]);
return 0;
}