平面最近点对

cdq分治可以很好地处理平面点对间具有某种性质的值或数量,最近点对也不例外。
参考OI Wiki
先对x排序,cdq返回点集内部最近点对的距离记为mindis。
考虑如何合并左右区间,对于处于点集A的点a和B的点b,显然
(|a.x-b.x|<=mindis)
设a,b属于C,找到C后,对于点p,k属于C,显然也要满足
(|p.y-k.y|<=mindis)
可以证明,对于点p,C中满足上述性质的点不会超过7个,因为
(|p.x-k.x|<=mindis)(|p.y-k.y|<=mindis),把区间分为2个正方形,正方形边长为mindis,对于每个正方形,再分为4个小正方形,而且每个小正方形内部最多有1个点,因为如果有大于两个点,那么,存在点x,y,,使得
(dis(x,y)<sqrt(2(frac{mindis}{2})^2)<mindis)
也就是说,我们每次合并两个点集,计算点x时,与点y的dist中,y与x属于同一个集合的个数最多为7个,即计算时无效的计算最多为7次。
所以,为了方便计算,我们将cdq后的点按y的大小归并。
Raid
给定两个点集,计算点集之间的最小距离。
给每个点加个flag,flag相同的点为同一个点集,距离为inf

#include<bits/stdc++.h>
using namespace std;
char buf[1<<20],*P1=buf,*P2=buf;
#define gc() (P1==P2&&(P2=(P1=buf)+fread(buf,1,1<<20,stdin),P1==P2)?EOF:*P1++)
#define TT template<class T>inline
TT bool read(T &x){
    x=0;char c=gc();bool f=0;
    while(c<48||c>57){if(c==EOF)return 0;f^=(c=='-'),c=gc();}
    while(47<c&&c<58)x=(x<<3)+(x<<1)+(c^48),c=gc();
    if(f)x=-x;return 1;
}
TT void print(T x){
    if(x<0)putchar('-'),x=-x;
    if(x>9)print(x/10);
    putchar('0'+x%10);
}
TT bool read(T&a,T&b){return read(a)&&read(b);}
TT bool read(T&a,T&b,T&c){return read(a)&&read(b)&&read(c);}
typedef long long ll;
const ll MAXN=1e5+8;
#define lowbit(x) (x&(-x))
struct Point{
    int x,y,flag;
}point[MAXN<<1],tmp[MAXN<<1];
double dis(Point&a,Point&b){
    if(a.flag==b.flag)return 1e30;//同一点集
    ll dx=a.x-b.x,dy=a.y-b.y;
    return sqrt(dx*dx+dy*dy);
}
bool cmp(Point&a,Point&b){
    if(a.x^b.x)return a.x<b.x;
    return a.y<b.y;
}
int n;
double cdq(int l,int r){
    if(l==r)return 1e30;//单个点
    int mid=(l+r)>>1;
    int midx=point[mid].x;
    //中点的x值,如果abs(p.x-midx)<mindis,就是符合点集C的条件
    double mindis=cdq(l,mid);//好像不能直接min,调试的时候好像先递归的you右区间。
    mindis=min(mindis,cdq(mid+1,r));
    //按y的值归并
    for(int i=l,j=mid+1,p=l;i<=mid||j<=r;++p)
        if(i<=mid&&(j>r||point[i].y<=point[j].y))tmp[p]=point[i++];
        else tmp[p]=point[j++];
    int cnt=0;
    for(int i=l;i<=r;++i){
        point[i]=tmp[i];//copy
        if(abs(point[i].x-midx)<mindis){//即点x属于点集C
        //点集C中满足abs(y)小于mindis的点,遍历一遍。
            for(int j=cnt;j&&point[i].y-tmp[j].y<mindis;--j)
                mindis=min(mindis,dis(point[i],tmp[j]));
            tmp[++cnt]=point[i];//点集C存放在tmp中,但是不会修改有用的值。
        }
    }return mindis;
}
int main() {
    int t;read(t);
    while(t--){
        read(n);
        for(int i=1;i<=n;++i)read(point[i].x,point[i].y);
        for(int i=n+1;i<=n<<1;++i){
            read(point[i].x,point[i].y);
            point[i].flag=1;
        }
        sort(point+1,point+1+n*2,cmp);
        printf("%.3f
",cdq(1,n<<1));
    }
    return 0;
}
原文地址:https://www.cnblogs.com/foursmonth/p/14161869.html