HDU 1007 kd Tree

经典的平面最近点问题。分治或kdtree都能解决。
分治的话可以去看算法导论,里面有详细的证明。
我来说说kdtree:
这个题普通做法是会TLE的,需要剪枝,可以在建树的时候维护,距离当前轴的最大的偏移。这样在query时对另一个子树起到剪枝作用,反正按理论来说,这样应该要快一些。

按下面这样写能过。。。。。这样写岂不是增枝。。

    void querymin(int rt,int D,const node &x, double &ans) {
        if(!sz[rt]) return;ID = D;
        int w = cmp(p[rt],x);	//w表示x在当前节点的左子树/右子树
        querymin(w^ls,D^1,x,ans);	//这里采用的是父亲与儿子二倍的建树关系
        if(x.id!=p[rt].id)ans = min(ans,x.dis(p[rt]));
        if((ans-sq(x.x[D]-p[rt].mc[w]))>eps) querymin(ls^w^1,D^1,x,ans);
        // p[rt].mc[0] 左子树中距离当前轴最近的的长度
    }

理论应该这样。。。却TLE了。无语了。

if((ans-sq(x.x[D]-p[rt].mc[!w]))>eps) querymin(ls^w^1,D^1,x,ans);

然后优化了一万年。。从最大方差入手优化就能过,跑的贼快。。。
每次选取当前轴的时候可以以最大方差的那个轴作为当前划分轴。

#include <bits/stdc++.h>
using namespace std;
const int maxn = (int)1e5+100;
#define sq(x) ((x)*(x))
int ID,mk;
struct node {
    double x[3];
    int id;
    inline double dis(const node&a) const{//return distance of this,a ^2
        double res=0;
        for(int i=0;i<mk;++i) res+=sq(a.x[i]-x[i]);
        return res;
    }
};
const double eps = 1e-9;
inline bool cmp(const node&a,const node&b) { return (b.x[ID]-a.x[ID])>eps;}
struct KDtree {
#define ls (rt<<1)
#define rs (rt<<1|1)
    char sz[maxn<<2];//当前轴
    node p[maxn<<2];
    void init(int n)
    {
        memset(sz,-1,4*(n+1)*sizeof(sz[0]));
    }
    int var(int l,int r,node *a) {//返回方差最大的轴
        double ave[2]={0,0};
        for(int i = l;i<=r;++i) ave[0]+=a[i].x[0],ave[1]+=a[i].x[1];
        ave[0]/=(r-l+1);ave[1]/=(r-l+1);
        double sa[2]={0,0};
        for(int i=l;i<=r;++i) sa[0]+=sq(ave[0]-a[i].x[0]),sa[1]+=sq(ave[1]-a[i].x[1]);
        return sa[0] < sa[1];
    }
    void build(int l,int r,int rt,node *a)
    {
        if(l>r)return;
        int mid=l+r>>1;
        sz[rt]=ID=var(l,r,a);
        nth_element(a+l,a+mid,a+r+1,cmp);
        p[rt]=a[mid];
        build(l,mid-1,ls,a);
        build(mid+1,r,rs,a);
    }
    void querymin(int rt,const node &x, double &ans) {
        if(sz[rt]==-1) return;
        char &D=sz[rt];ID=D;
        int w = cmp(p[rt],x);   //x在次节点的左边(0)or右边(1)
        querymin(w^ls,x,ans);
        if(x.id!=p[rt].id)ans = min(ans,x.dis(p[rt]));
        if((ans-sq(x.x[D]-p[rt].x[D]))>eps) querymin(ls^w^1,x,ans);//距离跨越当前轴
    }
}kd;
node a[maxn];
int main()
{
    freopen("1.in","r",stdin);
    freopen("1.out","w",stdout);
    int n;
    while(1) {
        clock_t c1=clock();
        scanf("%d",&n);
        if(!n)break;
        kd.init(n);mk = 2;
        for(int i = 1;i <= n; ++i)
            scanf("%lf %lf",&a[i].x[0],&a[i].x[1]),a[i].id=i;
        kd.build(1,n,1,a);
        double ans = a[1].dis(a[2]);
        for (int i = 1; i <= n; ++i) {
            kd.querymin(1,a[i],ans);
        }
        printf("%.2f\n",sqrt(ans)/2.0);
    }
    return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值