计算几何-最近点对

在二维平面上，有n个点，求距离最近的一对点

方法：

1.暴力：将n个点先按横坐标排序，再从第一个点开始，依次比较所有比它横坐标大的点的距离，找到最小距离d，再从第二个点开始，依次比较所有比它横坐标大的点的距离，一直比到倒数第二个点。时间复杂度O(n^2)

2.分治法：将n个点排序后，以横坐标的大小进行平分为左子集和右子集，分别找出左子集和右子集的最短距离dl，dr，取d=min(dl，dr)，再依次比较左子集中横坐标较大的点和右子集中横坐标较小的点之间的距离，如果有小于d的点对，则更新d值。时间复杂度O(nlgn)

                                    

//hdu 1007

先说下题意，很简单，给n个点的坐标，求距离最近的一对点之间距离的一半。

第一行是一个数n表示有n个点，接下来n行是n个点的x坐标和y坐标。实数。


这个题目其实就是求最近点对的距离。《算法导论》上有详细讲解，王晓东的书上也有代码。主要思想就是分治。先把n个点按x坐标排序，然后求左边n/2个和右边n/2个的最近距离，最后合并。

合并要重点说一下，比较麻烦。


首先，假设点是n个，编号为1到n。我们要分治求，则找一个中间的编号m，先求出1到m点的最近距离设为d1，还有m+1到n的最近距离设为d2。这里的点需要按x坐标的顺序排好，并且假设这些点中，没有2点在同一个位置。（若有，则直接最小距离为0了）。


然后，令d为d1, d2中较小的那个点。如果说最近点对中的两点都在1-m集合中，或者m+1到n集合中，则d就是最小距离了。但是还有可能的是最近点对中的两点分属这两个集合，所以我们必须先检测一下这种情况是否会存在，若存在，则把这个最近点对的距离记录下来，去更新d。这样我们就可以得道最小的距离d了。


关键是要去检测最近点对，理论上每个点都要和对面集合的点匹配一次，那效率还是不能满足我们的要求。所以这里要优化。怎么优化呢？考虑一下，假如以我们所选的分割点m为界，如果某一点的横坐标到点m的横坐标的绝对值超过d1并且超过d2，那么这个点到m点的距离必然超过d1和d2中的小者，所以这个点到对方集合的任意点的距离必然不是所有点中最小的。



所以我们先把在m为界左右一个范围内的点全部筛选出来，放到一个集合里。

筛选好以后，当然可以把这些点两两求距离去更新d了，不过这样还是很慢，万一满足条件的点很多呢。这里还得继续优化。首先把这些点按y坐标排序。假设排序好以后有p个点，编号为1到p。那么我们用1号去和2到p号的点求一下距离，然后2号和3到p号的点求一下距离。。。还没完，因为这样比，求的次数还是O(p^2), 所以其实和没优化没区别。


假设有一个点q,坐标是xq, yq。可以证明在以q为底边中点，长为2d，宽为d的矩形区域内不会有超过6个点。具体证明过程可以参看算法导论。


利用这个结论我们就可以来继续优化比较的过程了。刚刚我们是用用1号点去和2到p号的点求一下距离，我们知道以1号点构造图中矩形内，不会有超过6个点存在。但是我们又不能直接从1号求到6号，因为这里的p个点是按y坐标排序而不是按距离排序的，有可能在y坐标上，前10个点距离1号点都很近，但是前6个点的x坐标很远，而第10个点的x坐标和1号点的x坐标很进，这样第10个点到1号点的距离反而更近。


那么我们怎么利用这个结论呢？应该这样，假设1号点和2到p号点比较，由于y坐标排序的缘故，假设第p个点的y坐标距离1号点的y坐标大于当前能求出的最小值，那么这点以及这点后的所有点距离1号点的距离必然大于当前已经获得的最小值，所以直接不用比较后面的了。
又因为满足比较条件的点很少，不会超过6个，所以这里可以看成O(1)的效率。那么整个算法的效率大概是在O(nlogn)，非常快！

具体见代码：

//hdu 1007

#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <iostream>
#include <math.h>
using namespace std;
const double eps = 1e-6;
const int MAXN = 100010;
const double INF = 1e20;
struct Point
{
    double x, y;
};
double dist(Point a, Point b)
{
    return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
Point p[MAXN];
Point tmpt[MAXN];
bool cmpxy(Point a, Point b)
{
    if (a.x != b.x)
    {
        return a.x < b.x;
    }
    else
    {
        return a.y < b.y;
    }
}
bool cmpy(Point a, Point b)
{
    return a.y < b.y;
}
double Closest_Pair(int left, int right)
{
    double d = INF;
    if (left == right)
    {
        return d;
    }
    if (left + 1 == right)
    {
        return dist(p[left], p[right]);
    }
    int mid = (left + right) / 2;
    double d1 = Closest_Pair(left, mid);
    double d2 = Closest_Pair(mid + 1, right);
    d = min(d1, d2);
    
    int k = 0;
    for (int i = left; i <= right; i++)
    {
        if (fabs(p[mid].x - p[i].x) <= d)
        {
            tmpt[k++] = p[i];
        }
    }
    sort(tmpt, tmpt + k, cmpy);
    for (int i = 0; i < k; i++)
    {
        for (int j = i + 1; j < k && tmpt[j].y - tmpt[i].y < d; j++)
        {
            d = min(d, dist(tmpt[i], tmpt[j]));
        }
    }
    return d;
}
int main()
{
    int n;
    while (scanf("%d", &n) == 1 && n)
    {
        for (int i = 0; i < n; i++)
        {
            scanf("%lf%lf", &p[i].x, &p[i].y);
        }
        sort(p, p + n, cmpxy);
        printf("%.2lf\n", Closest_Pair(0, n - 1) / 2);
    }
    return 0;
}