平面最近点对问题
概述
给定 \(n\) 个二维平面上的点，求一组欧几里得距离最近的点对。
其中 \(n\leq 10^5\)，对于每个点对 \((x,y)\) 满足 \(0\leq x,y\leq 10^9\)。
算法分析
显然暴力枚举的时间复杂度为 \(O(n^2)\)。
对于这种看上去很不可做的题目，或许我们可以考虑分治。
首先每个点按 \(x\) 的大小升序排序。
然后我们可以每次将需要处理的节点分为两部分 \((l,mid),(mid+1,r)\)。
考虑合并区间，显然 \(ans(l,r)=min\{ans(l,mid),ans(mid+1,r),ans(跨越 mid)\}\)。
首先我们可以分治下去，直接得到 \(min\{ans(l,mid),ans(mid+1,r)\}\)。
我们将这个值作为半径，以中点为圆心画圆，显然只有圆内的点有可能刷新答案。
但是圆不太好处理，或许我们可以将它放大为圆的外接正方形，显然对答案不会有影响。
我们假设集合 \(A\) 表示点的横坐标 \(x\) 满足 \(mid_x-r\leq x\leq mid_x+r\) 的点的集合。
对于每个点 \(p(x,y)\in A\)，我们在集合 \(B\) 中寻找答案，其中集合 \(B=\{q(x,y)|q\in A,|y_p-y_q|\leq r\}\)。
尝试在分治内部将点按照 \(y\) 升序排列，那么我们就有了合并时的算法：设当前区间为 \((l,r),n=l-r+1\)。
- 对区间内的点按 \(y\) 升序排序。快排只能到 \(O(n\log n)\)，归并排序可以达到 \(O(n)\)。
- 构建集合 \(A\)。\(O(n)\)
- 对每个 \(p(x,y)\in A\) 构建相应的集合 \(B\) 并尝试刷新答案。\(O(?)\)
时间复杂度
前两步的时间复杂度已经简单得出了。
但是看上去似乎第三步的复杂度与集合 \(B\) 的大小线性相关，实际不然。
根据简单的鸽巢原理/抽屉原理可以得知，集合 \(B\) 的大小将无法超过一个很小的数。（好像是13）
所以时间复杂度仍然是很好的 \(O(n)\)。（常数略大）
所以分治背景下的总时间复杂度是 \(O(n\log n)\)。
代码实现
值得注意的是集合 \(B\) 可能不需要真正被建立，只需要思想上存在就好。
这样会大大降低代码实现复杂度。
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define N 100010
using namespace std;
const double INF=1e15;
int T,n;
struct node{
double x,y;
bool operator<(const node &a)const{return x<a.x;}
}a[N],tmp[N];
int read(){
int x=0,f=1;char c=getchar();
while(c<'0' || c>'9') f=(c=='-')?-1:1,c=getchar();
while(c>='0' && c<='9') x=x*10+c-48,c=getchar();
return x*f;
}
double get_dis(double x1,double y1,double x2,double y2){
double f1=x1-x2,f2=y1-y2;
return sqrt(f1*f1+f2*f2);
}
double dfs(int l,int r){
if(l==r) return INF;
int mid=(l+r)>>1;
double mid_x=a[mid].x;
double ans=min(dfs(l,mid),dfs(mid+1,r));
//以下是归并排序
int i=l,j=mid+1,cnt=0;
while(i<=mid && j<=r){
if(a[i].y<a[j].y) tmp[++cnt]=a[i++];
else tmp[++cnt]=a[j++];
}
while(i<=mid) tmp[++cnt]=a[i++];
while(j<=r) tmp[++cnt]=a[j++];
for(int i=l;i<=r;i++) a[i]=tmp[i-l+1];
//以下是构建集合 A
cnt=0;
for(int i=l;i<=r;i++)
if(a[i].x>=mid_x-ans && a[i].x<=mid_x+ans) tmp[++cnt]=a[i];
//以下是模拟集合 B
for(int i=1;i<=cnt;i++){
for(int j=i-1;j>=1;j--){
if(tmp[i].y-tmp[j].y>ans) break;
ans=min(ans,get_dis(tmp[i].x,tmp[i].y,tmp[j].x,tmp[j].y));
}
}
return ans;
}
void work(){
n=read();
for(int i=1;i<=n;i++)
scanf("%lf %lf",&a[i].x,&a[i].y);
sort(a+1,a+n+1);
printf("%.3lf\n",dfs(1,n));
}
int main(){
T=read();
while(T--) work();
return 0;
}
算法拓展
我们还可以加上一个限制条件：
给定 \(n\) 个属于集合 \(N\) 的点和 \(n\) 个属于集合 \(M\) 的点，求 属于不同集合的点 的平面最近点对。
看上去复杂了很多，但其实只需要加上几句话就好。
主要思想是将统一集合内的点打上同一标记，如果刷新答案时发现两点属于同一集合就不刷新。
时间复杂度仍然是 \(O(n\log n)\)。
struct node{
double x,y;bool typ;
bool operator<(const node &a)const{return x<a.x;}
}a[N<<1],tmp[N<<1];
//以下是模拟集合 B
for(int i=1;i<=cnt;i++){
for(int j=i-1;j>=1;j--){
if(tmp[i].y-tmp[j].y>ans) break;
if(tmp[i].typ != tmp[j].typ) ans=min(ans,get_dis(tmp[i].x,tmp[i].y,tmp[j].x,tmp[j].y));
//加了个 if
}
}
void work(){
n=read();
for(int i=1;i<=n;i++)
scanf("%lf %lf",&a[i].x,&a[i].y),a[i].typ=false;//打标记
for(int i=n+1;i<=2*n;i++)
scanf("%lf %lf",&a[i].x,&a[i].y),a[i].typ=true;//打标记
sort(a+1,a+2*n+1);
printf("%.3lf\n",dfs(1,n*2));
}
总结
完结撒花。
