kmjp's blog

競技プログラミング参加記です

Codeforces #231 Div2. E. Lightbulb for Minister

パッと見結構めんどくさい問題。
http://codeforces.com/contest/394/problem/E

問題

2次元座標中でN個の点が与えられている。
また、それとは別に凸なM角形が与えられる。

M角形の内部または辺上で、これらN個の点からのユークリッド距離の2乗和が最小化となる場所の2乗和を求めよ。

解法

N個の点をX[i],Y[i]とする。
ある点SX,SYとN個の点の距離の2乗和は \sum_i ((X_i-SX)^2+(Y_i-SY)^2) = SX^2 + SY^2 -2*SX*\sum_i X_i - 2*SY*\sum_i Y_i + \sum_i (X_i^2 + Y_i^2)と変形できる。
よって、A=-2*sum(X[i])、B=-2*sum(Y[i])、C=sum(X[i]^2+Y[i]^2)とすると、2乗和は(SX^2+SY^2+A*SX+B*SY+C)と変形でき、O(1)で計算できるようになる。

また、上記式を平方完成すると、SX=-A/2、SY=-B/2の時が2乗和が最小となる。
座標(SX,SY)がM角形内にあるなら、(SX,SY)が解。

(SX,SY)がM角形内にない場合は、M角形の辺上(含む点上)に解がある。
よって、M角形上の各辺を媒介変数tを用いて表現したとき、(SX^2+SY^2+A*SX+B*SY+C)が(0<=t<=1)の範囲で最小となる位置を求めればよい。

最小値を求めるだけなので二次方程式を解く必要はなく、平方完成で最小となるtを求めるだけ。

int N,M;
ll CX1,CY1,C;
int NX[100001],NY[100001];
double X[100001],Y[100001];

bool inside(double x,double y) {
	int lr=0,i;
	FOR(i,M) {
		double cr=(x-X[i])*(y-Y[(i+1)%M])-(y-Y[i])*(x-X[(i+1)%M]);
		if(cr<0) lr |= 1;
		if(cr>0) lr |= 2;
	}
	return (lr==1 || lr==2);
}

void solve() {
	int f,i,j,k,l,x,y;
	
	cin>>N;
	FOR(i,N) {
		cin>>x>>y;
		CX1 += -2*x;
		CY1 += -2*y;
		C += x*(ll)x + y*(ll)y;
	}
	cin>>M;
	FOR(i,M) cin>>X[i]>>Y[i];
	
	double mix=CX1/(-2.0)/N,miy=CY1/(-2.0)/N;
	
	if(inside(mix,miy)) {
		double r = N*(mix*mix+miy*miy) + CX1*mix + CY1*miy + C;
		_P("%.12lf\n",r);
		return;
	}
	double mi=1e20;
	
	FOR(i,M) {
		mi=min(mi, N*(X[i]*X[i]+Y[i]*Y[i]) + CX1*X[i] + CY1*Y[i] + C);
		double x1=X[i],y1=Y[i],x2=X[(i+1)%M],y2=Y[(i+1)%M];
		double aa = N*((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
		double bb = 2*N*((x2-x1)*x1+(y2-y1)*y1) + CX1*(x2-x1) + CY1*(y2-y1);
		double t = bb/(-2.0)/aa;
		if(t<=0 || t>=1) continue;
		double nx=t*(x2-x1)+x1;
		double ny=t*(y2-y1)+y1;
		mi=min(mi, N*(nx*nx+ny*ny) + CX1*nx + CY1*ny + C);
	}
	_P("%.12lf\n",mi);
}

まとめ

本番(SX,SY)に対する二乗和をO(1)で求める方法が思いつかなかった…。
これは対して難しくないし覚えておかねば。