最近点对

最近点对

设\(p_i = (x_i,y_i)\),表示平面上的一个点。
对于给定的点集\(S\),求最近点对。

很容易想到\(O(n^2)\)的算法。
计算每一对点的距离,然后取最小值。
但今天看分治的时候看到一种\(O(nlogn)\)的算法。

我们将点集合\(S\)按照\(x\)为第一关键字,\(y\)为第二关键字的大小升序排列,然后在拆分成集合大小相同的\(S_1\)和\(S_2\)。
根据分治的想法,最近点对只会有三种情况:

  • 两个点都在\(S_1\)集合中
  • 两个点都在\(S_2\)集合中
  • 一个在\(S_1\)一个在\(S_2\)。

前面两种可以递归的寻找,所以问题在于处理第三种情况。
假设情况1的最近点对距离为\(h_1\),情况2的最近点对距离为\(h_2\),\(h = min(h_1,h_2)\)。
用\(OIwiki\)上的这张图来观察一下:
最近点对
对于情况3的情况,所选的点对最好是位于中间中间轴左右\(h\)的区域内比较好,这个区域就是上图的\(B\)。

上图是在\(x\)方向考虑,现在考虑\(y\)方向。
对于每一个在\(B\)中的\(p_i\),只需要找到另外一个点\(p_j\),且纵坐标绝对值之差小于\(h\),考虑到重复的问题,我们让\(p_j\)的纵坐标大于\(p_i\)的纵坐标。有一个神奇的证明表示\(p_j\)的个数最大为7。

于是对于情况3我们有了一个算法:

  1. 构建\(B\)
  2. 将\(B\)按照纵坐标排序
  3. 对于每一个\(p_i\),找到最近距离更新答案

接下来就是使用分治法解决问题了,时间复杂度为\(O(nlogn)\)【虽然不太会证明这个问题】

#include <iostream>
#include <algorithm>
#include <vector>
#include <cstring>
#include <map>
#include <queue>
#include <cmath>
using namespace std;
#define ll long long
#define pii pair<ll,ll>
const ll N = 1e6+50;
const ll inf = 0x3f3f3f3f;
struct Pt{
    double x,y;
}pt[N];
ll n,f[N];
bool cmpxy(const Pt&a,const Pt& b){
    return a.x == b.x?a.y<b.y:a.x<b.x;
}
bool cmpy(const ll& a,const ll& b){
    return pt[a].y < pt[b].y;
}
double dis(ll i,ll j){
    ll dx = pt[i].x-pt[j].x,dy = pt[i].y-pt[j].y;
    return sqrt(dx*dx+dy*dy);
}
double close_pair(ll l,ll r){
    double h = 1e9+50;
    if(l == r) return h;
    if(l+1 == r) return dis(l,r);
    ll mid = l+r>>1;
    double h1 = close_pair(l,mid);
    double h2 = close_pair(mid+1,r);
    //cout<<h1<<' '<<h2<<endl;
    h = min(h1,h2);
    ll k = 0;
    for(ll i = l;i <= r;i++){
        if(fabs(pt[mid].x-pt[i].x) <= h)f[++k] = i;
    }
    sort(f+1,f+k+1,cmpy);
    for(ll i = 1;i <= k;i++){
        for(ll j = i+1;j <= k;j++){
            if(pt[f[j]].y-pt[f[i]].y >= h) break;
            //printf("%lld %lld %.4lf\n",f[i],f[j],dis(f[i],f[j]));
            h = min(h,dis(f[i],f[j]));
        }
    }
    //printf("%lld %lld %.4lf\n",l,r,h);
    return h;
}
int main() {
    //ios::sync_with_stdio(false);
    scanf("%lld",&n);
    for(ll i = 1;i <= n;i++) scanf("%lf%lf",&pt[i].x,&pt[i].y);
    sort(pt+1,pt+n+1,cmpxy);
//    for(ll i = 1;i <= n;i++){
//        printf("%lf %lf\n",pt[i].x,pt[i].y);
//    }
    printf("%.4lf\n",close_pair(1,n));
    return 0;
}
上一篇:noip77


下一篇:根据进程名称杀掉进程