HDU2823 The widest road

The widest road

给出两个点集，求出两个凸包后，求两个凸包的最近距离。需要特判两个凸包内含和相交各种乱七八糟的情况。

(1 ≤ m,n ≤ 1000)

题解

https://blog.csdn.net/clover_hxy/article/details/54022026

旋转卡壳求两凸包的最近距离

考虑如下的算法， 算法的输入是两个分别有(m)和(n)个顺时针给定顶点的凸多边形P和Q。

1. 计算P上y坐标值最小的顶点（称为 yminP ）和Q上y坐标值最大的顶点（称为 ymaxQ）。

2. 为多边形在 yminP 和 ymaxQ 处构造两条切线 LP 和 LQ 使得他们对应的多边形位于他们的右侧。

   此时 LP 和 LQ 拥有不同的方向， 并且 yminP 和 ymaxQ 成为了多边形间的一个对踵点对。

3. 计算距离(yminP,ymaxQ) 并且将其维护为当前最小值。

4. 顺时针同时旋转平行线直到其中一个与其所在的多边形的边重合。

5. 如果只有一条线与边重合， 那么只需要计算“顶点-边”对踵点对和“顶点-顶点”对踵点对距离。 都将他们与当前最小值比较， 如果小于当前最小值则进行替换更新。如果两条切线都与边重合，那么情况就更加复杂了。如果边“交叠”，也就是可以构造一条与两条边都相交的公垂线（但不是在顶点处相交）， 那么就计算“边-边”距离。 否则计算三个新的“顶点-顶点”对踵点对距离。 所有的这些距离都与当前最小值进行比较， 若小于当前最小值则更新替换。

6. 重复执行步骤4和步骤5， 直到新的点对为(yminP,ymaxQ)。

7. 输出最小距离。

时间复杂度(O(nlog n))。

struct point {float128 x,y;};

IN bool operator==(CO point&a,CO point&b){
	return a.x==b.x and a.y==b.y;
}
IN bool operator<(CO point&a,CO point&b){
	return a.x!=b.x?a.x<b.x:a.y<b.y;
}
IN point operator+(CO point&a,CO point&b){
	return (point){a.x+b.x,a.y+b.y};
}
IN point operator-(CO point&a,CO point&b){
	return (point){a.x-b.x,a.y-b.y};
}
IN point operator*(CO point&a,float128 b){
	return (point){a.x*b,a.y*b};
}
IN point operator/(CO point&a,float128 b){
	return (point){a.x/b,a.y/b};
}
IN float128 cross(CO point&a,CO point&b){
	return a.x*b.y-a.y*b.x;
}
IN float128 dot(CO point&a,CO point&b){
	return a.x*b.x+a.y*b.y;
}
IN float128 len(CO point&a){
	return sqrt(dot(a,a));
}

IN bool on_seg(CO point&a,CO point&b,CO point&c){
	if(a==b or a==c) return 1;
	return cross(a-b,c-b)==0 and dot(b-a,c-a)<0;
}
IN bool seg_inter(CO point&a,CO point&b,CO point&c,CO point&d){
	if(on_seg(a,c,d) or on_seg(b,c,d) or on_seg(c,a,b) or on_seg(d,a,b)) return 1;
	return cross(c-a,b-a)*cross(d-a,b-a)<0 and cross(a-c,d-c)*cross(b-c,d-c)<0;
}
IN float128 dist_seg(CO point&a,CO point&b,CO point&c){
	if(b==c) return len(a-b);
	if(dot(a-b,c-b)<0) return len(a-b);
	if(dot(a-c,b-c)<0) return len(a-c);
	return abs(cross(a-b,c-b))/len(c-b);
}

struct line {point x,y;};

IN bool on_left(CO point&a,CO line&b){
	return cross(b.y,a-b.x)>0;
}
IN int quad(CO point&a){
	if(a.x>0 and a.y>=0) return 1;
	else if(a.x<=0 and a.y>0) return 2;
	else if(a.x<0 and a.y<=0) return 3;
	else return 4;
}
IN bool operator<(CO line&a,CO line&b){
	if(quad(a.y)!=quad(b.y)) return quad(a.y)<quad(b.y); // edit 1
	return cross(a.y,b.y)==0?on_left(a.x,b):cross(a.y,b.y)>0;
}
IN point intersect(CO line&a,CO line&b){
	return b.x+b.y*cross(b.x-a.x,a.y)/cross(a.y,b.y);
}

CO int N=2e3+10;
CO float128 inf=1e9;
point p[N],q[N],ch[N],it[N];
line ln[N];

void convex(point p[],int&n){
	if(n==1) return; // edit 1
	sort(p+1,p+n+1);
	int m=0;
	for(int i=1;i<=n;++i){
		for(;m>=2 and cross(p[i]-ch[m-1],ch[m]-ch[m-1])>=0;--m);
		ch[++m]=p[i];
	}
	int o=--m;
	for(int i=n;i>=1;--i){
		for(;m>=o+2 and cross(p[i]-ch[m-1],ch[m]-ch[m-1])>=0;--m);
		ch[++m]=p[i];
	}
	n=--m,copy(ch+1,ch+n+1,p+1),p[n+1]=p[1];
//	cerr<<"p="<<endl;
//	for(int i=1;i<=n;++i)
//		cerr<<" ("<<p[i].x<<","<<p[i].y<<")"<<endl;
}
bool in_poly(point a,point p[],int n){
	for(int i=1;i<=n;++i)
		if(cross(a-p[i],p[i+1]-p[i])>0) return 0;
	return 1;
}
bool seg_poly_inter(point a,point b,point p[],int n){
	for(int i=1;i<=n;++i)
		if(seg_inter(a,b,p[i],p[i+1])) return 1;
	return 0;
}
bool halfplane(point p[],int n,point q[],int m){
	for(int i=1;i<=n;++i)
		ln[i].x=p[i],ln[i].y=p[i+1]-p[i];
	for(int i=1;i<=m;++i)
		ln[n+i].x=q[i],ln[n+i].y=q[i+1]-q[i];
	sort(ln+1,ln+n+m+1);
	int num=1;
	for(int i=2;i<=n+m;++i){
		if(cross(ln[i].y,ln[num].y)==0) continue;
		ln[++num]=ln[i];
	}
	int l=1,r=1;
	for(int i=2;i<=num;++i){
		for(;l<r and !on_left(it[r],ln[i]);--r);
		for(;l<r and !on_left(it[l+1],ln[i]);++l);
		ln[++r]=ln[i];
		if(l<r) it[r]=intersect(ln[r-1],ln[r]);
	}
	for(;l<r and !on_left(it[r],ln[l]);--r);
	return r-l+1>=3;
}
float128 rotate(point p[],int n,point q[],int m){
	int x=1,y=1;
	for(int i=2;i<=n;++i)if(p[i].y<p[x].y) x=i;
	for(int i=2;i<=m;++i)if(q[i].y>q[y].y) y=i;
	float128 ans=inf;
	for(int i=1;i<=n;++i,x=x%n+1){
		for(;cross(p[x+1]-p[x],q[y+1]-q[y])>=0;y=y%m+1); // edit 2
//		cerr<<"now ("<<p[x].x<<","<<p[x].y<<") ("<<p[x+1].x<<","<<p[x+1].y<<")"<<endl;
//		cerr<<" to ("<<q[y].x<<","<<q[y].y<<") ("<<q[y+1].x<<","<<q[y+1].y<<")"<<endl;
		if(cross(p[x+1]-p[x],q[y+1]-q[y])==0){
			ans=min(ans,dist_seg(p[x],q[y],q[y+1]));
			ans=min(ans,dist_seg(p[x+1],q[y],q[y+1]));
			ans=min(ans,dist_seg(q[y],p[x],p[x+1]));
			ans=min(ans,dist_seg(q[y+1],p[x],p[x+1]));
		}
		else ans=min(ans,dist_seg(q[y],p[x],p[x+1]));
	}
//	cerr<<"ans="<<ans<<endl;
	return ans;
}
void real_main(int n,int m){
	for(int i=1;i<=n;++i) scanf("%lf%lf",&p[i].x,&p[i].y);
	for(int i=1;i<=m;++i) scanf("%lf%lf",&q[i].x,&q[i].y);
	convex(p,n);
	convex(q,m);
	if(n>m) swap(n,m),swap(p,q);
	float128 ans=inf;
	if(n==1){
		if(m==1) ans=len(p[1]-q[1]);
		else if(m==2) ans=dist_seg(p[1],q[1],q[2]);
		else{
			if(in_poly(p[1],q,m)) ans=0;
			else{
				for(int i=1;i<=m;++i)
					ans=min(ans,dist_seg(p[1],q[i],q[i+1]));
			}
		}
		printf("%.4lf
",ans);
		return;
	}
	if(n==2){
		if(m==2){
			if(seg_inter(p[1],p[2],q[1],q[2])) ans=0;
			else{
				ans=min(ans,dist_seg(p[1],q[1],q[2]));
				ans=min(ans,dist_seg(p[2],q[1],q[2]));
				ans=min(ans,dist_seg(q[1],p[1],p[2]));
				ans=min(ans,dist_seg(q[2],p[1],p[2]));
			}
		}
		else{
			if(in_poly(p[1],q,m) or in_poly(p[2],q,m)) ans=0;
			else if(seg_poly_inter(p[1],p[2],q,m)) ans=0;
			else{
				for(int i=1;i<=m;++i){
					ans=min(ans,dist_seg(p[1],q[i],q[i+1]));
					ans=min(ans,dist_seg(p[2],q[i],q[i+1]));
					ans=min(ans,dist_seg(q[i],p[1],p[2]));
					ans=min(ans,dist_seg(q[i+1],p[1],p[2]));
				}
			}
		}
		printf("%.4lf
",ans);
		return;
	}
	if(halfplane(p,n,q,m)) ans=0;
	else{
		ans=min(ans,rotate(p,n,q,m));
		ans=min(ans,rotate(q,m,p,n));
	}
	printf("%.4lf
",ans);
}
int main(){
	for(int n,m;scanf("%d%d",&n,&m)!=EOF;) real_main(n,m);
	return 0;
}