「题解」「APIO2018」选圆圈 Circle selection

2020年4月30日 0 亚太地区信息学奥林匹克竞赛

APIO2018 T2，一道神奇的 $\texttt{k-D Tree}$ 搜索剪枝题，或者一道毒瘤的 平衡树 + CDQ 分治 题。

题目链接：Luogu P4631/BZOJ 5465/LibreOJ 2586/APIO2018 T2。

题目

题目描述

在平面上，有 $n$ 个圆，记为 $c_1,c_2,\cdots,c_n$。我们尝试对这些圆运行这个算法：

  1. 找到这些圆中半径最大的。如果有多个半径最大的圆，选择编号最小的。记为 $c_i$。
  2. 删除 $c_i$ 及与其有交集的所有圆。两个圆有交集当且仅当平面上存在一个点，这个点同时在这两个圆的圆周上或圆内。（原文直译：如果平面上存在一个点被这两个圆所包含，我们称这两个圆有交集。一个点被一个圆包含，当且仅当它位于圆内或圆周上。）
    https://i.loli.net/2018/05/21/5b0235e3ccd74.png
  3. 重复上面两个步骤直到所有的圆都被删除。

当 $c_i$ 被删除时，若循环中第一步选择的圆是 $c_j$ ，我们说 $c_i$ 被 $c_j$ 删除。对于每个圆，求出它是被哪一个圆删除的。

数据范围

变量数据范围
$n$$[1,3\times 10^5]$
$x,y$$[-10^9,10^9]$
$r$$[1,10^9]$

时空限制

题目编号时间限制空间限制
Luogu P4631$3\text{s}$$1000\text{MiB}$
BZOJ 5465$?\text{s}$$1024\text{MiB}$
LibreOJ 2586$3\text{s}$$1024\text{MiB}$
APIO2018 T2$3\text{s}$$1024\text{MiB}$

题解

本题目有多种玄学方法。

方法一

思路

暴力

考虑暴力，肯定要先以半径为第一关键字降序排列，编号为第二关键字升序排列。然后对于每个圆，求出它能删除那些圆，并记录答案，后面只要对没有答案的进行查找即可。最坏的时间复杂度是 $\Theta(n^2)$。

下面考虑优化这个过程。

必要条件的转化

我们不难发现，如果圆 $i$ 与圆 $j$ 相交（或者说接触），那么从它们的 $x$ 轴投影和 $y$ 轴投影可以看出：区间 $[x _ i-r _ i,x _ i+r _ i]$ 与区间 $[x _ j-r _ j,x _ j+r _ j]$ 有交，$y$ 轴同理。

所以，坐标对应的投影区间有交是两圆相交的必要条件。

所以，我们可以把圆的外部抽象出一个外接正方形，如果一个圆与这个外接正方形没有公共部分，必然也与这个圆没有公共部分，而正方形和圆可以用 $\texttt{k-D Tree}$ 来储存。

用 $\texttt{k-D Tree}$ 可以通过洛谷上的数据，但是 $\texttt{k-D Tree}$ 毕竟是暴力数据结构，LibreOJ 和 UOJ 上好像有卡 $\texttt{k-D Tree}$ 的数据，不过旋转一下点就解决了。

代码

#include<bits/stdc++.h>
using namespace std;
#define reg register
#define eps 1e-6
#define getchar()(p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++)
static char buf[100000],*p1=buf,*p2=buf;
inline int read(void){
    reg bool f=false;
    reg char ch=getchar();
    reg int res=0;
    while(ch<'0'||'9'<ch)f|=(ch=='-'),ch=getchar();
    while('0'<=ch&&ch<='9')res=10*res+ch-'0',ch=getchar();
    return f?-res:res;
}

long double sina,cosa; //旋转角度的 sin 值和 cos 值

inline long double sqr(reg long double x){ //求平方和
    return x*x;
}

const int MAXN=300000+5;

int n,ans[MAXN];

namespace kD_Tree{ //k-D Tree
    #define mid ( ( (l) + (r) ) >> 1 )
    struct Point{
        long double x[2],r;
        int id;
        inline void Read(reg int i){
            static double X,Y;
            id=i;
            X=read(),Y=read();
            x[0]=X*cosa-Y*sina; //旋转
            x[1]=X*sina+Y*cosa; //旋转
            r=read(); //读入半径
            return;
        }
        inline bool operator<(const Point& a)const{
            return r==a.r?id<a.id:r>a.r;
        }
    };
    int root,tot;
    int lson[MAXN],rson[MAXN];
    long double Min[MAXN][2],Max[MAXN][2];
    Point p[MAXN],val[MAXN];
    int WD;
    inline bool cmp(const Point& a,const Point& b){
        return a.x[WD]<b.x[WD];
    }
    inline void pushup(reg int k){ //维护信息
        reg long double x=val[k].x[0],y=val[k].x[1],R=val[k].r;
        Min[k][0]=x-R;Max[k][0]=x+R;
        Min[k][1]=y-R;Max[k][1]=y+R;
        reg int l=lson[k],r=rson[k];
        if(l){
            Min[k][0]=min(Min[k][0],Min[l][0]);
            Max[k][0]=max(Max[k][0],Max[l][0]);
            Min[k][1]=min(Min[k][1],Min[l][1]);
            Max[k][1]=max(Max[k][1],Max[l][1]);
        }
        if(r){
            Min[k][0]=min(Min[k][0],Min[r][0]);
            Max[k][0]=max(Max[k][0],Max[r][0]);
            Min[k][1]=min(Min[k][1],Min[r][1]);
            Max[k][1]=max(Max[k][1],Max[r][1]);
        }
        return;
    }
    inline int Build(reg int l,reg int r,reg int wd){ //建树
        if(l>r)
            return 0;
        reg int k=++tot;
        WD=wd;
        nth_element(p+l,p+mid,p+r+1,cmp);
        val[k]=p[mid];
        lson[k]=Build(l,mid-1,wd^1);
        rson[k]=Build(mid+1,r,wd^1);
        pushup(k);
        return k;
    }
    inline bool check(reg int k,const Point& a){ //有交集
        return sqr(a.r+val[k].r)-(sqr(a.x[0]-val[k].x[0])+sqr(a.x[1]-val[k].x[1]))>=-eps;
    }
    inline bool far(reg int k,reg const Point& a){ //完全脱离
        return (a.x[0]+a.r<Min[k][0])||(a.x[0]-a.r>Max[k][0])||(a.x[1]+a.r<Min[k][1])||(a.x[1]-a.r>Max[k][1]);
    }
    inline void Query(reg int k,reg const Point& a){
        if(far(k,a))
            return;
        if(!ans[val[k].id]&&check(k,a))
            ans[val[k].id]=a.id;
        if(lson[k])
            Query(lson[k],a);
        if(rson[k])
            Query(rson[k],a);
        return;
    }
}

using namespace kD_Tree;

int main(void){
    n=read();
    srand(time(0));
    sina=sin((long double)rand()),cosa=sqrt(1-sina*sina);
    for(reg int i=1;i<=n;++i)
        p[i].Read(i);
    root=Build(1,n,0);
    sort(p+1,p+n+1);
    for(reg int i=1;i<=n;++i)
        if(!ans[p[i].id])
            Query(root,p[i]);
    for(reg int i=1;i<=n;++i)
        printf("%d%c",ans[i],i==n?'\n':' ');
    return 0;
}

方法二

思路

防止被卡的方法

考虑到坐标对应的投影区间有交是两圆相交的必要条件。不难发现，对于我们要查询的一个圆，它的半径必然小于之前的所有圆（我们以半径为降序排过序了）。

而且，前面的圆能够不被删去，它们必定两两不相交，且半径大于圆 $c_i$。

那么我们可以对前面的圆维护扫描线，每个圆和当前的直线 $x=x_0$ 相交两次，可以用括号表示。而且由于这些圆两两不相交，括号相对次序不会变。

由于前面的圆半径都比它大，因此和它有交的圆必然经过 $x=x_i-r_i$ 或 $x=x_i+r_i$ 或 $y=y_i−r_i$ 或 $y=y_i+r_i$。所以我们可以在做扫描线时，查询这四个位置的平衡树上，当前圆的前驱后继。

应对多重询问

有很多个询问，那就加上一个 CDQ 分治 咯。

时间复杂度为 $\Theta(n\log_2^2n)$。

代码

#include<bits/stdc++.h>
using namespace std;
#define reg register
#define INF 0X7F7F7F7F
typedef long long ll;
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++)
static char buf[100000],*p1=buf,*p2=buf;
inline int read(void){
    reg bool f=false;
    reg char ch=getchar();
    reg int res=0;
    while(ch<'0'||'9'<ch)f|=(ch=='-'),ch=getchar();
    while('0'<=ch&&ch<='9')res=10*res+ch-'0',ch=getchar();
    return f?-res:res;
}

const int MAXN=300010;

struct Circle{
    ll x,y,r;
    int id;
    inline void Read(reg int i){
        x=read(),y=read(),r=read();
        id=i;
        return;
    }
};

struct Event{
    ll t;
    int op,v;
    inline Event(reg ll t=0,reg int op=0,reg int v=0):t(t),op(op),v(v){
        return;
    }
};

inline bool cmp(const Circle& a,const Circle& b){
    return (a.r==b.r)?(a.id<b.id):(a.r>b.r);
}

inline bool cmp2(const Event& a,const Event& b){
    return a.t<b.t;
}

int n;
Circle a[MAXN];
int ans[MAXN];
int final[MAXN];
int b[MAXN];
Event E[2*MAXN];
set<pair<int,int>/**/> S;

inline ll sqr(reg ll x){
    return x*x;
}

inline bool inter(reg int x,reg int y){
    return sqr(a[x].x-a[y].x)+sqr(a[x].y-a[y].y)<=sqr(a[x].r+a[y].r);
}

inline void Solve(reg int l,reg int r){
    if(l==r){
        if(ans[l]==INF){
            ans[l]=l;
            b[l]=1;
        }
        return;
    }
    reg int mid=(l+r)>>1;
    Solve(l,mid);
    reg int cnt=0;
    for(reg int i=l;i<=mid;++i)
        if(b[i]){
            E[++cnt]=Event(3*(a[i].x-a[i].r)-2,1,i);
            E[++cnt]=Event(3*(a[i].x+a[i].r),2,i);
        }
    for(reg int i=mid+1;i<=r;++i){
        E[++cnt]=Event(3*(a[i].x-a[i].r)-1,3,i);
        E[++cnt]=Event(3*(a[i].x+a[i].r)-1,3,i);
    }
    sort(E+1,E+cnt+1,cmp2);
    for(reg int i=1;i<=cnt;++i)
        if(E[i].op==1)
            S.insert(make_pair(a[E[i].v].y,E[i].v));
    else if(E[i].op==2)
        S.erase(make_pair(a[E[i].v].y,E[i].v));
    else{
        auto it=S.lower_bound(make_pair(a[E[i].v].y,0));
        if(it!=S.end()){
            int x=it->second;
            if(inter(x,E[i].v))
                ans[E[i].v]=min(ans[E[i].v],x);
        }
        if(it!=S.begin()){
            --it;
            int x=it->second;
            if(inter(x,E[i].v))
                ans[E[i].v]=min(ans[E[i].v],x);
        }
    }
    cnt=0;
    for(reg int i=l;i<=mid;++i)
        if(b[i]){
            E[++cnt]=Event(3*(a[i].y-a[i].r)-2,1,i);
            E[++cnt]=Event(3*(a[i].y+a[i].r),2,i);
        }
    for(reg int i=mid+1;i<=r;++i){
        E[++cnt]=Event(3*(a[i].y-a[i].r)-1,3,i);
        E[++cnt]=Event(3*(a[i].y+a[i].r)-1,3,i);
    }
    sort(E+1,E+cnt+1,cmp2);
    for(reg int i=1;i<=cnt;++i)
        if(E[i].op==1)
            S.insert(make_pair(a[E[i].v].x,E[i].v));
    else if(E[i].op==2)
        S.erase(make_pair(a[E[i].v].x,E[i].v));
    else{
        auto it=S.lower_bound(make_pair(a[E[i].v].x,0));
        if(it!=S.end()){
            int x=it->second;
            if(inter(x,E[i].v))
                ans[E[i].v]=min(ans[E[i].v],x);
        }
        if(it!=S.begin()){
            --it;
            int x=it->second;
            if(inter(x,E[i].v))
                ans[E[i].v]=min(ans[E[i].v],x);
        }
    }
    Solve(mid+1,r);
    return;
}

int main(void){
    n=read();
    ll Minx=INF,Miny=INF;
    for(reg int i=1;i<=n;++i){
        a[i].Read(i);
        Minx=min(Minx,a[i].x),Miny=min(Miny,a[i].y);
    }
    for(reg int i=1;i<=n;++i){
        a[i].x-=Minx-1;
        a[i].y-=Miny-1;
    }
    sort(a+1,a+n+1,cmp);
    memset(ans,0X7F,sizeof(ans));
    Solve(1,n);
    for(reg int i=1;i<=n;++i)
        final[a[i].id]=a[ans[i]].id;
    for(reg int i=1;i<=n;++i)
        printf("%d%c",final[i],i==n?'\n':' ');
    return 0;
}

相关题目

Tags:

Lu-Anlai

江西赣州人

Share

You may also like...

发表评论

电子邮件地址不会被公开。 必填项已用*标注