题目描述
分析
首先判掉 (INF) 的情况
第一种情况就是不能从 (s) 走到 (t)
第二种情况就是从 (s) 出发走到了出度为 (0) 的点，这样就再也走不到 (t)
然后我们去考虑 (60) 分的做法
我们设 (dp[u]) 为当前在点 (u) 走到点 (t) 的期望步数,(du[u]) 为 (u) 的出度
那么就有
(dp[u]=sum_{u->v}^v((dp[v]+1) imes frac{1}{du[u]}))
移项之后就变成了
(dp[u]-sum_{u->v}^v(dp[v] imes frac{1}{du[u]})=1)
初始化 (dp[t]=0)
一共有 (n) 个未知数，(n)个方程可以用高斯消元解决
这样做时间复杂度是 (n^3) 的
下面我们去考虑满分做法
因为保证强连通分量的大小不超过 (100)
所以我们可以在强联通分量内使用高斯消元
对于整个 (DAG) 使用拓扑排序去解决
注意拓扑排序要建反图
对于拓扑序小于 (t) 所属强联通分量的强联通分量，我们没有必要去算
因为在到这些点之前就已经到达了 (t)
对于其它的强联通分量
如果分量内的某个点的去边所指向的点在当前的强联通分量中，我们把其作为未知数参与消元
对于不在这个强联通分量中的点，这个点肯定在之前被求出来过
我们把它作为常数
代码
#include<cstdio>
#include<cmath>
#include<algorithm>
#include<cstring>
#include<iostream>
#include<queue>
#include<vector>
inline int read(){
int x=0,fh=1;
char ch=getchar();
while(ch<'0' || ch>'9'){
if(ch=='-') fh=-1;
ch=getchar();
}
while(ch>='0' && ch<='9'){
x=(x<<1)+(x<<3)+(ch^48);
ch=getchar();
}
return x*fh;
}
typedef double db;
const db eps=1e-8;
const int maxn=1e6+5;
const int maxk=105;
const int maxm=1e4+5;
int head[maxn],tot=1;
struct asd{
int to,next;
}b[maxn],b2[maxn];
void ad(int aa,int bb){
b[tot].to=bb;
b[tot].next=head[aa];
head[aa]=tot++;
}
int h2[maxn],t2=1;
void ad2(int aa,int bb){
b2[t2].to=bb;
b2[t2].next=h2[aa];
h2[aa]=t2++;
}
int rd[maxn];
db mp[maxk][maxk];
int jl[maxk][maxk];
int low[maxn],dfn[maxn],dfnc,js,shuyu[maxn],sta[maxn],top,siz[maxn];
void tar(int now){
dfn[now]=low[now]=++dfnc;
sta[++top]=now;
for(int i=head[now];i!=-1;i=b[i].next){
int u=b[i].to;
if(!dfn[u]){
tar(u);
low[now]=std::min(low[now],low[u]);
} else if(!shuyu[u]){
low[now]=std::min(low[now],dfn[u]);
}
}
if(dfn[now]==low[now]){
js++;
while(1){
int y=sta[top--];
shuyu[y]=js;
siz[js]++;
if(y==now) break;
}
}
}
int n,m,s,t,rk[maxm],du[maxm],xh[maxm];
std::vector<int> g[maxn];
db dp[maxn];
bool jud;
void gsxy(int id){
memset(rk,0,sizeof(rk));
memset(xh,0,sizeof(xh));
int ncnt=0;
for(int i=1;i<=n;i++){
if(shuyu[i]==id){
rk[i]=++ncnt;
xh[ncnt]=i;
g[id].push_back(i);
if(i==t) jud=1;
}
}
if(jud==0){
return;
} else {
memset(mp,0,sizeof(mp));
memset(jl,0,sizeof(jl));
for(int i=1;i<=ncnt;i++){
mp[i][ncnt+1]=1;
}
for(int i=0;i<g[id].size();i++){
int now=g[id][i];
for(int i=head[now];i!=-1;i=b[i].next){
int u=b[i].to;
if(shuyu[u]==shuyu[now]){
jl[rk[now]][rk[u]]++;
} else {
mp[rk[now]][ncnt+1]+=(double)dp[u]/du[now];
}
}
}
for(int i=0;i<g[id].size();i++){
for(int j=0;j<g[id].size();j++){
int aa=g[id][i],bb=g[id][j];
if(i==j){
mp[rk[aa]][rk[aa]]=(db)(du[aa]-jl[rk[aa]][rk[aa]])/du[aa];
} else {
mp[rk[aa]][rk[bb]]=(db)(-1.0)*jl[rk[aa]][rk[bb]]/du[aa];
}
}
}
for(int i=0;i<g[id].size();i++){
if(g[id][i]==t){
int now=rk[g[id][i]];
for(int j=1;j<=ncnt+1;j++){
mp[now][j]=0;
}
mp[now][now]=1;
}
}
int now=1;
for(int i=1;i<=ncnt;i++){
double mmax=0;
int jl=0;
for(int j=now;j<=ncnt;j++){
if(std::fabs(mp[j][i])>std::fabs(mmax)){
mmax=mp[j][i];
jl=j;
}
}
if(std::fabs(mmax)<eps) continue;
if(jl!=now) std::swap(mp[jl],mp[now]);
for(int j=i+1;j<=ncnt+1;j++){
mp[now][j]/=mp[now][i];
}
mp[now][i]=1.0;
for(int j=now+1;j<=ncnt;j++){
db cs=mp[j][i];
for(int k=i;k<=ncnt+1;k++){
mp[j][k]-=cs*mp[now][k];
}
}
now++;
}
dp[xh[ncnt]]=mp[ncnt][ncnt+1];
for(int i=ncnt-1;i>=1;i--){
dp[xh[i]]=mp[i][ncnt+1];
for(int j=i+1;j<=n;j++){
dp[xh[i]]-=mp[i][j]*dp[xh[j]];
}
}
}
}
std::queue<int> q;
int huan[maxn];
void tp(){
for(int i=1;i<=js;i++){
if(rd[i]==0) q.push(i);
}
while(!q.empty()){
int now=q.front();
gsxy(now);
q.pop();
for(int i=h2[now];i!=-1;i=b2[i].next){
int u=b2[i].to;
rd[u]--;
if(rd[u]==0) q.push(u);
}
}
}
bool vis[maxn];
void dfs(int now){
vis[now]=1;
for(int i=head[now];i!=-1;i=b[i].next){
int u=b[i].to;
if(!vis[u]) dfs(u);
}
}
int main(){
memset(h2,-1,sizeof(h2));
memset(head,-1,sizeof(head));
n=read(),m=read(),s=read(),t=read();
for(int i=1;i<=m;i++){
int aa,bb;
aa=read(),bb=read();
ad(aa,bb);
du[aa]++;
if(aa==bb) huan[aa]++;
}
dfs(s);
if(vis[t]==0){
printf("INF
");
return 0;
}
for(int i=1;i<=n;i++){
if(i==t) continue;
if(vis[i] && du[i]-huan[i]==0){
printf("INF
");
return 0;
}
}
for(int i=1;i<=n;i++){
if(!dfn[i]) tar(i);
}
for(int i=1;i<=n;i++){
for(int j=head[i];j!=-1;j=b[j].next){
int u=b[j].to;
if(shuyu[i]!=shuyu[u]){
ad2(shuyu[u],shuyu[i]);
rd[shuyu[i]]++;
}
}
}
tp();
printf("%.3f
",dp[s]);
return 0;
}
更新
补一下判 (INF) 的锅
否则遇到这组数据会挂
5 6 1 5
1 1
4 4
2 1
1 5
3 4
1 3