计算几何模板——平面最近点对（附C++源代码） - [转载经典]

   转：http://anonympine.blogbus.com/logs/33968379.html

求平面最近点对的核心思想乃是二分，用递归实现。具体操作如下：

     若点的个数很少（比如小于3或者小于5），就直接枚举。

     如点的个数很多，按现将所有点按X排序，并按X坐标平均的分成左右两个部分（假设分割线为X=nx），分别求出两边的最短距离minl与minr并令ans=min（minl，minr）。

     求出左右两边的最小值之后，剩下的工作就是合并。易见若该点集存在点对（a，b）的最近距离小于ans，则a，b一定分别在x=nx的两边，切nx-a.x与nx-b.x的绝对值肯定小于ans。

     据此我们可以将点集中所有X值在（nx-ans，nx+ans）的点都选出来，那么满足条件的（a，b）肯定都在其中。

     易见若存在（a，b）两点他们之间的距离小于ans，那么a.y-b.y的绝对值也肯定小于ans。

     综上存在（a，b）两点他们之间的距离小于ans那，（a，b）一定在一个长为2*ans宽为ans的矩形之中。而且这个矩形被X=nx平分成两个ans*ans的矩形，由于无论是在左边还是在右边，任意两点的之间的距离总是小于等于ans的，所以两个ans*ans 的矩形中最多只有4个点（分别在四个顶点上），长为2*ans宽为ans的矩形最多有8个点。

     据此我们将所有X值在（nx-ans，nx+ans）的点按他们的Y值进行排序。依次看每个点与它之后的7个点的距离是否小于ans，若小于则更新ans，最后求出来的结果就是平面最近点对的距离。保留产生该距离的两个点即可得到最近点对。

     练手题目：Pku2107，Vijos1012

附C++代码（Pku2107）：

#include <iostream>
#include <cmath>

const long maxsize = 100000;

typedef struct 
{ 
    double x, y; 
} PointType;

long list[maxsize], listlen,n;
PointType point[maxsize];

int sortcmp(const void *,const void *); 
double dis(PointType,PointType);
double getmin(double,double);
int listcmp(const void *,const void *); 
double shortest(long,long);
int init(void);

int main() 
{ 
    while (init())
       printf("%.2lf\n",shortest(0, n - 1)/2);    
    return 0;
}

int sortcmp(const void *a, const void *b) 
{ 
    if (((PointType*)a)->x < ((PointType*)b)->x)    
       return -1;   
    else 
       return 1; 
}

double dis(PointType a, PointType b) 
{ 
    return sqrt((a.x - b.x)*(a.x - b.x) + (a.y - b.y)*(a.y - b.y)); 
}

double getmin(double a, double b) 
{ 
    return a<b?a:b;
}

int listcmp(const void *a, const void *b) 
{ 
    if (point[*(int*)a].y < point[*(int*)b].y)    
       return -1;   
    else 
       return 1; 
}

double shortest(long l, long r) 
{ 
    if (r - l == 1) 
        return dis(point[l], point[r]);   
    if (r - l == 2)    
       return getmin(getmin(dis(point[l], point[l+1]), dis(point[l], point[r])), dis(point[l+1], point[r]));   
    long i, j, mid = (l + r) >> 1;   
    double curmin = getmin(shortest(l, mid), shortest(mid + 1, r));   
    listlen = 0;   
    for (i = mid; i >= l && point[mid+1].x - point[i].x <= curmin; i --)    
       list[listlen++] = i;   
    for (i = mid + 1; i <= r && point[i].x - point[mid].x <= curmin; i ++)    
       list[listlen++] = i;   
    qsort(list, listlen, sizeof(list[0]), listcmp);   
    for (i = 0; i < listlen; i ++)    
       for (j = i + 1; j < listlen && point[list[j]].y - point[list[i]].y <= curmin; j ++)     
          curmin = getmin(curmin, dis(point[list[i]], point[list[j]]));   
    return curmin; 
}

int init(void)
{
    int i;
    scanf("%d", &n);      
    for (i=0;i<n;i++) 
       scanf("%lf%lf",&point[i].x,&point[i].y);      
    qsort(point,n,sizeof(point[0]),sortcmp);
    return n;
}