最近点对问题

http://blog.csdn.net/lonelycatcher/article/details/7973046/

在二维平面上的n个点中，如何快速的找出最近的一对点，就是最近点对问题。

    一种简单的想法是暴力枚举每两个点，记录最小距离，显然，时间复杂度为O(n^2)。

    在这里介绍一种时间复杂度为O(nlognlogn)的算法。其实，这里用到了分治的思想。将所给平面上n个点的集合S分成两个子集S1和S2，每个子集中约有n/2个点。然后在每个子集中递归地求最接近的点对。在这里，一个关键的问题是如何实现分治法中的合并步骤，即由S1和S2的最接近点对，如何求得原集合S中的最接近点对。如果这两个点分别在S1和S2中，问题就变得复杂了。

    为了使问题变得简单，首先考虑一维的情形。此时，S中的n个点退化为x轴上的n个实数x1，x2，...，xn。最接近点对即为这n个实数中相差最小的两个实数。显然可以先将点排好序，然后线性扫描就可以了。但我们为了便于推广到二维的情形，尝试用分治法解决这个问题。

    假设我们用m点将S分为S1和S2两个集合，这样一来，对于所有的p(S1中的点)和q(S2中的点)，有p<q。

    递归地在S1和S2上找出其最接近点对{p1,p2}和{q1,q2}，并设

d = min{ |p1-p2| , |q1-q2| }

    由此易知，S中最接近点对或者是{p1,p2}，或者是{q1,q2}，或者是某个{q3,p3}，如下图所示。

 

    如果最接近点对是{q3,p3}，即|p3-q3|<d，则p3和q3两者与m的距离都不超过d，且在区间(m-d,d]和(d,m+d]各有且仅有一个点。这样，就可以在线性时间内实现合并。

    此时，一维情形下的最近点对时间复杂度为O(nlogn)。

    在二维情形下，类似的，利用分治法，但是难点在于如何实现线性的合并？

 

    由上图可见，形成的宽为2d的带状区间，最多可能有n个点，合并时间最坏情况下为n^2,。但是，P1和P2中的点具有以下稀疏的性质，对于P1中的任意一点，P2中的点必定落在一个d X 2d的矩形中，且最多只需检查六个点（鸽巢原理）。

    这样，先将带状区间的点按y坐标排序，然后线性扫描，这样合并的时间复杂度为O(nlogn)，几乎为线性了。

 

    光说不练也不行，经过自己的思考和参考网上的程序，完成了最近点对的程序，并在各OJ上成功AC了。

    POJ3714 ZOJ2107 HDU1007

[cpp] view plain copy

1. /** 
2. 最近点对问题，时间复杂度为O(n*logn*logn) 
3. */  
4. #include <iostream>  
5. #include <cstdio>  
6. #include <cstring>  
7. #include <cmath>  
8. #include <algorithm>  
9. using namespace std;  
10. const double INF = 1e20;  
11. const int N = 100005;  
12.   
13. struct Point  
14. {  
15.     double x;  
16.     double y;  
17. }point[N];  
18. int n;  
19. int tmpt[N];  
20.   
21. bool cmpxy(const Point& a, const Point& b)  
22. {  
23.     if(a.x != b.x)  
24.         return a.x < b.x;  
25.     return a.y < b.y;  
26. }  
27.   
28. bool cmpy(const int& a, const int& b)  
29. {  
30.     return point[a].y < point[b].y;  
31. }  
32.   
33. double min(double a, double b)  
34. {  
35.     return a < b ? a : b;  
36. }  
37.   
38. double dis(int i, int j)  
39. {  
40.     return sqrt((point[i].x-point[j].x)*(point[i].x-point[j].x)  
41.                 + (point[i].y-point[j].y)*(point[i].y-point[j].y));  
42. }  
43.   
44. double Closest_Pair(int left, int right)  
45. {  
46.     double d = INF;  
47.     if(left==right)  
48.         return d;  
49.     if(left + 1 == right)  
50.         return dis(left, right);  
51.     int mid = (left+right)>>1;  
52.     double d1 = Closest_Pair(left,mid);  
53.     double d2 = Closest_Pair(mid+1,right);  
54.     d = min(d1,d2);  
55.     int i,j,k=0;  
56.     //分离出宽度为d的区间  
57.     for(i = left; i <= right; i++)  
58.     {  
59.         if(fabs(point[mid].x-point[i].x) <= d)  
60.             tmpt[k++] = i;  
61.     }  
62.     sort(tmpt,tmpt+k,cmpy);  
63.     //线性扫描  
64.     for(i = 0; i < k; i++)  
65.     {  
66.         for(j = i+1; j < k && point[tmpt[j]].y-point[tmpt[i]].y<d; j++)  
67.         {  
68.             double d3 = dis(tmpt[i],tmpt[j]);  
69.             if(d > d3)  
70.                 d = d3;  
71.         }  
72.     }  
73.     return d;  
74. }  
75.   
76.   
77. int main()  
78. {  
79.     while(true)  
80.     {  
81.         scanf("%d",&n);  
82.         if(n==0)  
83.             break;  
84.         for(int i = 0; i < n; i++)  
85.             scanf("%lf %lf",&point[i].x,&point[i].y);  
86.         sort(point,point+n,cmpxy);  
87.         printf("%.2lf\n",Closest_Pair(0,n-1)/2);  
88.     }  
89.     return 0;  
90. }