【CF446D】DZY Loves Games 高斯消元+矩阵乘法

【CF446D】DZY Loves Games

题意：一张n个点m条边的无向图，其中某些点是黑点，1号点一定不是黑点，n号点一定是黑点。问从1开始走，每次随机选择一个相邻的点走过去，经过恰好k个黑点到达n的概率。

$nle 500,mle 500000,kle 10^9$，黑点个数不超过100.

题解：一眼就知道是高斯消元和矩乘什么的。我们先预处理出f[i][j]表示从第i个黑点开始走到的下一个黑点是j的概率。这个用高斯消元容易搞定。然后上矩乘即可。但是问题在于如果这样做的话我们要做n遍高斯消元。不过我们发现每次消元时左边的系数矩阵都是不变的，所以我们可以将n个方程组放到一起消元，复杂度就变成$O(n^3+10^6log k)$了。

#include <cstring>
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
using namespace std;
typedef double db;
int n,m,K,tot;
db f[510][610];
int dan[510],d[510],p[510],pa[100010],pb[100010];
struct M
{
	db v[105][105];
	M () {memset(v,0,sizeof(v));}
	db * operator [] (const int &a) {return v[a];}
	inline M operator * (const M &a) const
	{
		M b;
		int i,j,k;
		for(i=1;i<=tot;i++)	for(j=1;j<=tot;j++)	for(k=1;k<=tot;k++)	b.v[i][j]+=v[i][k]*a.v[k][j];
		return b;
	}
}S,T;
inline void pm(int y)
{
	while(y)
	{
		if(y&1)	S=S*T;
		T=T*T,y>>=1;
	}
}
inline int rd()
{
	int ret=0,f=1;	char gc=getchar();
	while(gc<'0'||gc>'9')	{if(gc=='-')	f=-f;	gc=getchar();}
	while(gc>='0'&&gc<='9')	ret=ret*10+(gc^'0'),gc=getchar();
	return ret*f;
}
int main()
{
	n=rd(),m=rd(),K=rd();
	int i,j,k,a,b;
	for(i=1;i<=n;i++)
	{
		dan[i]=rd();
		if(dan[i])	dan[i]=++tot,p[tot]=i;
	}
	for(i=1;i<=m;i++)	pa[i]=rd(),pb[i]=rd(),d[pa[i]]++,d[pb[i]]++;
	for(i=1;i<=m;i++)
	{
		a=pa[i],b=pb[i];
		if(!dan[a])	f[b][a]-=1.0/d[a];
		else	f[b][dan[a]+n]+=1.0/d[a];
		if(!dan[b])	f[a][b]-=1.0/d[b];
		else	f[a][dan[b]+n]+=1.0/d[b];
	}
	f[1][n+tot+1]=1;
	for(i=1;i<=n;i++)	f[i][i]+=1;
	for(i=1;i<=n;i++)
	{
		for(j=k=i;j<=n;j++)	if(fabs(f[j][i])>fabs(f[k][i]))	k=j;
		if(i!=k)	for(j=i;j<=n+tot+1;j++)	swap(f[i][j],f[k][j]);
		db tmp=f[i][i];
		for(j=i;j<=n+tot+1;j++)	f[i][j]/=tmp;
		for(j=1;j<=n;j++)	if(j!=i)
		{
			tmp=f[j][i];
			for(k=1;k<=n+tot+1;k++)	f[j][k]-=f[i][k]*tmp;
		}
	}
	for(i=1;i<=tot;i++)	for(j=1;j<=tot;j++)	T[i][j]=f[p[j]][i+n];
	for(i=1;i<=tot;i++)	S[1][i]=f[p[i]][n+tot+1];
	pm(K-2);
	printf("%.10lf",S[1][tot]);
	return 0;
}