最近点对
设\(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我们有了一个算法:
- 构建\(B\)
- 将\(B\)按照纵坐标排序
- 对于每一个\(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;
}