bzoj4520K远点对

Posted 2020-12-24 yoyoball

Portal --> bzoj4520

Description

　　给你平面内(n)个点的坐标，求欧氏距离下第(k)远的点对

Solution

　　因为kd其实。。严格来说挺不熟的用的太少了qwq

　　然后不知道为啥第一反应凸包直径取(k)次qwq然而这样有一个问题就是。。取完一次之后删点不知道要删直径中两个点中的哪一个，所以。。不太靠谱

　　正解应该是kd-tree

　　其实这题挺暴力的，时间复杂度也很玄学（不会算qwq）貌似kd的题复杂度就没有不玄学的。。

　　因为不知道答案是哪两个点，初步的想法是我们干脆维护一个大小为(k)的小根堆，对于每一个点，都在kd-tree里面查以其作为两个点之一的点对的前(k)大距离，如果说比小根堆的堆顶更优那就用当前的结果把小根堆的堆顶替换掉，这样对每一个点都操作一遍之后，堆顶就是答案了

?　　然而实际上，我们会发现一个点对会被计算两次，所以我们其实应该维护一个(k*2)的堆，查找也是(k*2)而不是(k)

　　再稍微具体一点的话在kd中查找的流程大概是这样的：

1、计算当前点与固定点的距离，如果比堆顶优就替换

2、用一个估值函数分别计算两个子树的可能最远距离(lval)和(rval)

3、如果说(lval>rval)则优先遍历左子树，否则优先遍历右子树

4、遍历一个子树的前提条件是：当前堆中不足(k*2)个元素或者该子树的估值函数返回值优于堆顶

　　接下来就是这个估值函数要怎么写了：

　　（这里提供的这种写法比较。。水。。其实如果有心去卡的话貌似是可以卡掉的qwq）

　　我们考虑分别记录该子树内的(x)的最大最小值和(y)的最大最小值，然后估值函数就返回(x)的最大最小值与固定点的(x)的最大差值的平方+(y)的最大最小值与固定点的(y)的最大差值的平方，这样我们就能得到一个最优情况下的最大值了

　　

　　代码大概长这个样子

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<queue>
#define ll long long
using namespace std;
const int N=100010,inf=2147483647;
int which;
ll sqr(ll x){return x*x;}
struct Dot{
    int data[2];
    Dot(){}
    Dot(int x,int y){data[0]=x; data[1]=y;}
    friend ll dis(Dot a,Dot b){return sqr(a.data[0]-b.data[0])+sqr(a.data[1]-b.data[1]);}
    friend bool operator < (Dot a,Dot b){return a.data[which]<b.data[which];}
}a[N];
int n,k,tot;
bool cmp(ll x,ll y){return x>y;}
priority_queue<ll,vector<ll>,greater<ll> > q;
namespace Kd{/*{{{*/
    int ch[N][2],mnx[N],mxx[N],mny[N],mxy[N];
    int rt;
    void pushup(int x){
        mxx[x]=mnx[x]=a[x].data[0]; mxy[x]=mny[x]=a[x].data[1];
        if (ch[x][0]){
            mxx[x]=max(mxx[x],mxx[ch[x][0]]);
            mxy[x]=max(mxy[x],mxy[ch[x][0]]);
            mnx[x]=min(mnx[x],mnx[ch[x][0]]);
            mny[x]=min(mny[x],mny[ch[x][0]]);
        }
        if (ch[x][1]){
            mxx[x]=max(mxx[x],mxx[ch[x][1]]);
            mxy[x]=max(mxy[x],mxy[ch[x][1]]);
            mnx[x]=min(mnx[x],mnx[ch[x][1]]);
            mny[x]=min(mny[x],mny[ch[x][1]]);
        }
    }
    int _build(int l,int r,int now){
        if (l>r) return 0;
        int mid=l+r>>1;
        which=now;
        nth_element(a+l,a+mid,a+r+1);
        ch[mid][0]=_build(l,mid-1,now^1);
        ch[mid][1]=_build(mid+1,r,now^1);
        pushup(mid);
        return mid;
    }
    ll val(int x,Dot &delta){
        if (!x) return -1;
        ll ret=max(sqr(delta.data[0]-mnx[x]),sqr(delta.data[0]-mxx[x]))+
               max(sqr(delta.data[1]-mny[x]),sqr(delta.data[1]-mxy[x]));
        return ret;
    }
    void build(int n){rt=_build(1,n,0);}
    void _query(int x,int k,Dot &delta){
        if (!x) return;
        ll d=dis(delta,a[x]),lval=val(ch[x][0],delta),rval=val(ch[x][1],delta);
        if (q.size()<k) q.push(d);
        else{
            if (d>q.top())
                q.pop(),q.push(d);
        }
        if (lval>rval){
            if (lval>q.top()||q.size()<k) _query(ch[x][0],k,delta);
            if (rval>q.top()||q.size()<k) _query(ch[x][1],k,delta);
        }
        else{
            if (rval>q.top()||q.size()<k) _query(ch[x][1],k,delta);
            if (lval>q.top()||q.size()<k) _query(ch[x][0],k,delta);
        }
    }
    void query(Dot &delta,int k){_query(rt,k,delta);}
}/*}}}*/

int main(){
#ifndef ONLINE_JUDGE
    freopen("2.in","r",stdin);
#endif
    scanf("%d%d",&n,&k);
    for (int i=1;i<=n;++i)
        scanf("%d%d",&a[i].data[0],&a[i].data[1]);
    k*=2;
    Kd::build(n);
    tot=0;
    for (int i=1;i<=n;++i) Kd::query(a[i],k);
    printf("%lld
",q.top());
}