高斯消元解线性期望方程的妙用

BZOJ3143

http://www.lydsy.com/JudgeOnline/problem.php?id=3143

E_i表示经i点的期望
E_（u，v）表示经（u，v）的期望

特别地有

构造线性期望方程，高斯消元即可.

根据排序不等式，贪心即可

#include<cstdio>
#include<algorithm>
using namespace std;
double a[511][511],x[511],w[250011],ans;
int deg[511];
int u[250011],v[250011];
int n,m;
inline void init(){
    for(register int i=1;i<=m;++i){
        a[u[i]][v[i]]+=1.00/deg[v[i]];
        a[v[i]][u[i]]+=1.00/deg[u[i]];
    }
    for(register int i=1;i<n;++i)
        a[i][i]=-1.00;
    a[1][n+1]=-1.00;
    for(register int i=1;i<=n+1;++i)
        a[n][i]=0.00;
    a[n][n]=1.00;
     
}
inline void gauss(){
    double tmp;
    for(register int i=1,p;i<=n;++i)
        for(register int j=i+1;j<=n;++j){
            tmp=1.00*a[j][i]/a[i][i];
            for(register int k=i;k<=n+1;++k)
                a[j][k]-=1.00*a[i][k]*tmp;
        }
    for(register int i=n;i;--i){
        for(register int j=i+1;j<=n;++j)a[i][n+1]-=x[j]*a[i][j];
        x[i]=a[i][n+1]/a[i][i];
    }
}
inline bool cmp(double a,double b){
    return a>b;
}
int main(){
    scanf("%d%d",&n,&m);
    for(register int i=1;i<=m;++i){
        scanf("%d%d",v+i,u+i);
        ++deg[v[i]];++deg[u[i]];
    }
    init();
    gauss();
    for(register int i=1;i<=m;++i)w[i]=1.00*x[u[i]]/(1.00*deg[u[i]])+1.00*x[v[i]]/(1.00*deg[v[i]]);
    sort(w+1,w+m+1,cmp);
    for(register int i=1;i<=m;++i)
        ans+=1.00*i*w[i];
    printf("%.3lf",ans);
    return 0;
}

---

BZOJ2337

http://www.lydsy.com/JudgeOnline/problem.php?id=2337

按二进制枚举每位，然后做法类似与上面

E_i表示i到n的二进制第k位为1的期望

对于每次消元后对E₁*2^k求个和

#include<cstdio>
#include<algorithm>
using namespace std;
const int maxn=1e6+7;
double a[501][501],_x[501],ans;
int nxt[maxn],las[maxn],to[maxn],w[maxn],deg[maxn];
int n,m,x,y,z,tot;
inline void add(int x,int y,int z){
    nxt[++tot]=las[x];
    las[x]=tot;
    to[tot]=y;
    w[tot]=z;
    ++deg[y];
}
inline void gauss(){
    double tmp;
    for(register int i=1,p;i<=n;++i){
        p=i;
        while(!a[p][i])++p;
        if(p!=i)swap(a[p],a[i]);
        for(register int j=i+1;j<=n;++j){
            tmp=a[j][i]/a[i][i];
            for(register int k=i;k<=n+1;++k)
                a[j][k]-=1.00*a[i][k]*tmp;
        }
    }
    for(register int i=n;i;--i){
        for(register int j=i+1;j<=n;++j)
            a[i][n+1]-=1.00*_x[j]*a[i][j];
        _x[i]=a[i][n+1]/a[i][i];
    }
}
int main(){
    scanf("%d%d",&n,&m);
    for(register int i=1;i<=m;++i){
        scanf("%d%d%d",&x,&y,&z);
        if(x!=y)add(x,y,z),add(y,x,z);
        else add(x,y,z);
    }
    for(register int i=0;i<=30;++i){
        for(register int i=1;i<=n;++i)
            for(register int j=1;j<=n+1;++j)
                a[i][j]=0.00;
        for(register int j=1;j<n;++j){
            a[j][j]=1.00;
            for(register int e=las[j];e;e=nxt[e]){
                x=to[e];
                if((w[e]>>i)&1){
                    a[j][x]+=1.00/deg[j];
                    a[j][n+1]+=1.00/deg[j];
                }
                else
                    a[j][x]-=1.00/deg[j];
            }
        }
        a[n][n]=1.00;
        gauss();
        ans+=1.00*_x[1]*(1<<i);
    }
    printf("%.3lf
",ans);
    return 0;
} 

　　

---

BZOJ3534

http://www.lydsy.com/JudgeOnline/problem.php?id=3534

就是求这一坨，新边权为

结合矩阵树定理，构造行列式，做高斯消元，最后乘上后面那一坨就好了

#include<cstdio>
#include<algorithm>
using namespace std;
const int mod=998244353;
typedef double db;
db a[5001][5001],ans;
int x,y,n,m;
inline db fabs(db x){
	return x>=0?x:-x;
}
inline db gauss(){
	db tmp;
	for(register int i=1,t;i<n;++i){
		t=i;
		for(register int j=i+1;j<n;++j)
			if(fabs(a[t][i])<fabs(a[j][i]))
				t=j;
		if(!a[t][i])return 0.00;
		swap(a[t],a[i]);
		for(register int j=i+1;j<n;++j){
			tmp=1.00*a[j][i]/a[i][i];
			for(register int k=i;k<n;++k)
				a[j][k]-=1.00*a[i][k]*tmp;
		}
	}
	tmp=1.00;
	for(register int i=1;i<n;++i)
		tmp=tmp*a[i][i];
	return fabs(tmp);
}
int main(){
	scanf("%d",&n);
	ans=1;
	for(register int i=1;i<=n;++i)
		for(register int j=1;j<=n;++j){
			scanf("%lf",&a[i][j]);
			if(i==j)continue;
			if(i<j)ans=ans*(1.00-a[i][j]);
			a[i][j]=a[i][j]/(1.00-a[i][j]);
		}
	for(register int i=1;i<=n;++i)
		for(register int j=1;j<=n;++j)
			if(i!=j)
				a[i][i]-=a[i][j];
	ans=ans*gauss();
	printf("%.8lf",ans);
	return 0;
}