HDU 1007 Quoit Design【計算幾何/分治/最近點對】
Quoit Design
Time Limit: 10000/5000 MS (Java/Others) Memory Limit: 65536/32768 K (Java/Others)
Total Submission(s): 58566 Accepted Submission(s): 15511
Assume that all the toys are points on a plane. A point is encircled by the ring if the distance between the point and the center of the ring is strictly less than the radius of the ring. If two toys are placed at the same point, the radius of the ring is considered to be 0.
Input The input consists of several test cases. For each case, the first line contains an integer N (2 <= N <= 100,000), the total number of toys in the field. Then N lines follow, each contains a pair of (x, y) which are the coordinates of a toy. The input is terminated by N = 0.
Sample Input 2 0 0 1 1 2 1 1 1 1 3 -1.5 0 0 0 0 1.5 0
Sample Output 0.71 0.00 0.75
Author CHEN, Yue
Source ZJCPC2004 【分析】：
平面最近點對，即平面中距離最近的兩點
分治算法：
int SOLVE(int left,int right)//求解點集中區間[left,right]中的最近點對
{
double ans; //answer
0) 調用前的預處理:對所有點排序,以x為第一關鍵詞y為第二關鍵字 , 從小到大;
1) 將所有點按x坐標分成左右兩部分;
/* 分析當前集合[left,right]中的最近點對，有兩種可能：
1. 當前集合中的最近點對，點對的兩點同屬於集合[left,mid]或同屬於集合[mid,right]
則ans = min(集合1中所有點的最近距離, 集合2中所有點的最近距離)
2. 當前集合最近點對中的兩點分屬於不同集合：[left,mid]和[mid,right]
則需要對兩個集合進行合並，找出是否存在p∈[left,mid],q∈[mid,right],使得distance(p,q)小於當前ans(即步驟1中求得的ans);
*/
2) Mid = (left+right)/2;
ans = min( SOLVE(left,mid), SOLVE(mid,right) );
即：遞歸求解左右兩部分中的最近距離，並取最小值；
//此步驟實現上文分析中的第一種情況
/*
再次進行分析
我們將集合[left,right]用x = mid這條直線分割成兩部分
則如果畫出直線l1:x=mid-ans 和 l2:x=mid+ans，顯然如果有p∈[left,mid], q∈[mid,right]且distance(p,q) < ans則p,q一定在直線l1和直線l2之間，否則distance(p,q)必定大於ans。
於是掃描出在l1和l2之間的點
*/
3) 建立緩存數組temp[];
for i = left TO right
{
如果 abs(Point[i].x - Point[mid].x) <= ans
則向temp中加入點Point[i];
}
/*
對於temp中的點，枚舉求所有點中距離最近兩點的距離，然後與ans比較即可。
枚舉的時候不必兩兩枚舉。觀察下圖中的點p
不難發現，若有q∈[mid,mid+ans]使得distance(p,q) <
ans，則q點的位置一定在圖中畫出的一個2ans×ansd的矩形中。可以證明點集[mid,mid+ans]中的、矩形外的點與p點的距離一定大於ans。於是我們可以對temp以y為唯一關鍵字從小到大排序,進行枚舉, 更新ans,然後在枚舉時判斷：一旦枚舉到的點與p點y值之差大於ans，停止枚舉。最後就能得到該區間的最近點對。
*/
4) sort(temp);
for i = 0 TO k-1
{
for j = i+1 TO k-1
如果 temp[j].y - temp[i].y >= ans break;
ans = min( ans, distance(temp[i], temp[j]) );
}
5) return ans;
}
算法的時間復雜度
由鴿巢原理，代碼中第四步的枚舉實際上最多只會枚舉6個點，效率極高（一種蒟蒻的證明請看下方的評論）
#include <bits/stdc++.h> using namespace std; typedef long long LL; const int N = 200000+5; const int mod = 2009; const double eps=1e-8; const int INF=0x7fffffff; int n; struct point{ double x,y; point(double x=0,double y=0):x(x),y(y) {} bool operator < (const point& p) const { if(x!=p.x) return x<p.x; else return y<p.y; } }p[N],tmp[N]; bool cmpy(point a,point b){ return a.y<b.y; } double dis(point a,point b){ return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)); } double closest_pair(int left,int right){ double d=INF; if(left==right)//如果此時區間查到了單點，直接返回最大值 return d; if(left+1==right) //如果左區間與右區間只差一，就直接返回他們的距離 return dis(p[left],p[right]); int mid=(left+right)>>1; //計算中點，位運算加速 double d1=closest_pair(left,mid); //分別計算出兩個區間的值 double d2=closest_pair(mid,right); d=min(d1,d2); //取兩個區間中的最小值 int k=0; //計算數量 for(int i=left;i<=right;i++){ if(fabs(p[mid].x-p[i].x)<=d){ //如果x軸的坐標相減滿足<=d載入數組 tmp[k++]=p[i]; } } //按照y軸的值排序，默認按照升序 sort(tmp,tmp+k,cmpy); for(int i=0;i<k;i++) for(int j=i+1;j<k&&tmp[j].y-tmp[i].y<d;j++){ //在已經篩出的數中計算最小值 double d3=dis(tmp[i],tmp[j]); d=min(d,d3); //如果有最小值更新 } return d; //直接返回最小值 } int main() { while(scanf("%d",&n),n){ for(int i=0;i<n;i++){ double a,b; scanf("%lf%lf",&a,&b); p[i]=point(a,b); } sort(p,p+n); //化成有序數列 printf("%.2lf\n",closest_pair(0,n-1)/2); //求的是半徑，若求直徑不必/2 } return 0; }平面點對
HDU 1007 Quoit Design【計算幾何/分治/最近點對】