[最小方差生成树 Kruskal] SRM 611 500pts

小R要给小Z修建庄园公路,目标是最小化道路长度的标准差。通过Kruskal算法,根据边的平方和平均值的关系进行贪心选择,通过预计算所有可能的平均值来找到最优解。关键在于理解边的平方和平均值的关系,并暴力尝试所有可能的‘关键点’作为平均值进行Prim算法构建最小生成树。
摘要由CSDN通过智能技术生成

题目大意:

小R计划给小Z买一个很大的庄园,庄园里的建筑可以抽象为二
维平面坐标系上的N<=20个点。对于任意两个点,他都可以修建一条笔
直的公路使它们相连(该公路的长度即为两点的欧式距离) 。
为了能腾出更多的资金给小Z买礼物,小R要恰好修建N-1条公
路,使得任意两个建筑都能直接或者间接的到达。
但小R有强迫症。他要使他修建的道路长度的这个数列的标准差
尽量的小。设他修建的N-1条道路的长度的平均数为ave。
其方差S=∑(Ai-ave)^2 标准差为:sqrt(S/(N-1))

题解:

假设我们已经知道了 ave,我们就可以把
方差的式子拆开,变成: ai^2+ave^2+2*ai*ave 。我们把所有边按
ai^2+2*ai*ave 排序,贪心地用 kruskal 连成树即可。
那 ave 究竟是多少呢?没事,我们可以多随机几次试试。根据标
算的做法,随机多次后搜到正确结果的期望还是挺高的。

我们考虑之前的做法在快排的时候的 cmp 函数。如果边 i 优于边
j,说明 ai^2+2*ai*ave<aj^2+2*aj*ave,化简后得到(ai+aj)/2<ave。也
就是说,我们不必确切知道 ave 的值,只需知道所有边组的(ai+aj)/2
与 ave 的关系。假设我们把 m^2 的(ai+aj)/2 全部提取出来标在数
轴上,如果一个 ave 是落在两个点之间的,它在哪里都会造成一样的
结果。
所以我们只需暴力提取 n^4 个“关键点”当做 ave 去做一遍 MST
即可。做的时候一定要用 prim 算法,否则加一个 log 可能会跪。


#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;

inline char nc()
{
	static char buf[100000],*p1=buf,*p2=buf;
	if (p1==p2) { p2=(p1=buf)+fread(buf,1,100000,stdin); if (p1==p2) return EOF; }
	return *p1++;
}

inline void read(int &x)
{
	char c=nc(),b=1;
	for (;!(c>='0' && c<='9');c=nc()) if (c=='-') b=-1;
	for (x=0;c>='0' && c<='9';x=x*10+c-'0',c=nc()); x*=b;
}

int Unique(double *a,int l,int r){
	int pnt=l;
	for (int i=l+1;i<=r;i++) if (fabs(a[i]-a[i-1])>1e-7) a[++pnt]=a[i];
	return pnt-l+1;
}

const int N=25;

int n;
int x[N],y[N];

inline double sqr(int x){ return (double)x*x; }
inline double Dist(int a,int b){ return sqrt(sqr(x[a]-x[b])+sqr(y[a]-y[b])); }

double ans=1e130,AVE;

double a[N*N];

double sx[N*N*N*N]; int icnt;

inline bool cmp(int x,int y){
	if (fabs(a[x]-a[y])<1e-7) return x<y;
	if (a[x]>a[y]) return a[x]+a[y]<2*AVE;
	return a[x]+a[y]>2*AVE;
}

struct edge{
	int u,v,idx; double w;
	edge() { }
	edge(int u,int v,int idx,double w):u(u),v(v),idx(idx),w(w) { }
	bool operator < (const edge &B) const{
		return cmp(idx,B.idx);
	}
}E[N*N];
int m;

namespace TSet{
	int fat[N],rank[N];
	inline void init(int n) {
		for (int i=1;i<=n;i++) fat[i]=i,rank[i]=0;
	}
	inline int Fat(int u){
		return u==fat[u]?u:fat[u]=Fat(fat[u]);
	}
	inline bool Union(int x,int y){
		x=Fat(x); y=Fat(y); if (x==y) return 0;
		if (rank[x]>rank[y]) swap(x,y);
		if (rank[x]==rank[y]) rank[y]++;
		fat[x]=y; return 1;
	}
}

inline void Krus(){
	using namespace TSet;
	double ave=0,ret=0;
	sort(E+1,E+m+1); init(n);
	for (int i=1;i<=m;i++)
		if (Union(E[i].u,E[i].v))
			ave+=E[i].w,ret+=E[i].w*E[i].w;
	ave/=(n-1);
	ret=ret/(n-1)-ave*ave;
	ans=min(ans,ret);
}

int main()
{
	freopen("mst.in","r",stdin);
	freopen("mst.out","w",stdout);
	read(n);
	for (int i=1;i<=n;i++) read(x[i]);
	for (int i=1;i<=n;i++) read(y[i]);
	for (int i=1;i<=n;i++)
		for (int j=i+1;j<=n;j++)
			a[++m]=Dist(i,j),E[m]=edge(i,j,m,a[m]);
	for (int i=1;i<=m;i++)
		for (int j=i+1;j<=m;j++)
			sx[++icnt]=(a[i]+a[j])/2.0;
	sort(sx+1,sx+icnt+1); icnt=Unique(sx,1,icnt);
	for (int i=1;i<=icnt;i++)
	{
		AVE=sx[i];
		Krus();
	}
	printf("%.11lf\n",sqrt(ans));	
	return 0;
}


评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值