【圆的反演变换+cdq分治】共点圆

最新推荐文章于 2022-01-20 22:22:20 发布
最新推荐文章于 2022-01-20 22:22:20 发布
阅读量1.2k 收藏
点赞数
分类专栏： 计算几何 数学 cdq分治 集训队互测
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/huyuncong/article/details/36895153
版权

http://www.tsinsen.com/ViewGProblem.page?gpid=A1381

因为所有的圆都过原点，因此可以以原点为反演中心，将所有的圆反演为一条直线，然后根据反演的性质，可以发现圆的内部会被反演到一个半平面上去，因此每次询问就相当于询问点是否在半平面的交之中，这个可以动态插半平面询问，也可以cdq分治来离线回答；询问的时候二分点在凸包的哪块三角形中，然后叉积一下就可以判断了。

我采用了第二种方法，然后开始了艰苦卓绝地卡常数的历程，历经各种常数优化，终于把最后两个点都卡在了3s整...

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
#include <cmath>
#include <ctime>
#define sqr(x) ((x)*(x))
const double pi=acos(-1.0);
const double eps=1e-8;
const double oo=1e15,LR=100000.0;
using namespace std;
struct point{
    double x,y,z,_,d;
    point () {}
    point (double X,double Y) {x=X,y=Y,z=1;}
    point operator -(const point &p) const {
        point e;
        e.x=x-p.x,e.y=y-p.y,e.z=1;
        return e;
    }
    point operator +(const point &p) const {
        point e;
        e.x=x+p.x,e.y=y+p.y,e.z=1;
        return e;
        }
        point operator *(const double a) const {
            point e;
            e.x=x*a,e.y=y*a,e.z=1;
            return e;
        }
        void print() {
            printf("X=%f Y=%f Z=%f\n",(double)x,(double)y,(double)z);
        }
        double dist() {return (sqr(x)+sqr(y));}
}p[500050],b[500050];
int ans[500050],a[500050],n,st[500050],u[500050],c[500050];
inline double cr(point e,point r) {return e.x*r.y-e.y*r.x;}
inline double dot(point &e,point &r) {return e.x*r.x+e.y*r.y+e.z*r.z;}
inline void cross(point &p,point &q,point &e)
{
    e.x=p.y*q.z-p.z*q.y;
    e.y=p.z*q.x-p.x*q.z;
    e.z=p.x*q.y-p.y*q.x;
}
inline void ori(point &a)
{
    a._=atan2(a.y,a.x);
    a.d=sqr(a.z)/(sqr(a.x)+sqr(a.y));
}
inline bool cmp(int i,int j)
{
    double tmp=p[i]._-p[j]._;
    if (fabs(tmp)>pi/2) return tmp<-eps;
    tmp=cr(p[i],p[j]);
    if (fabs(tmp)>eps) return tmp>eps;
    return p[i].d<p[j].d;
}
inline bool check(point L,point T,point I)
{
    point p;
    cross(L,T,p);
    if (dot(p,I)>-eps) return 1;
    return 0;
}
int search(point b[],int tot,point &P,point &Q,point &U)
{
        if (cr(P,U)<-eps || cr(Q,U)>eps) return 0;
        int l=2,r=tot-1;
				point pp; 
        for (;l<=r;) {
            int mid=((l+r)>>1);
            pp=b[mid]-b[1];
            if (cr(pp,U)>eps) l=mid+1;
            else r=mid-1;
        }
        return l-1;
}
point OB;
bool cmp2(point i,point j)
{
    double tmp=i._-j._;
        if (fabs(tmp)>pi/2) return tmp<-eps;
    tmp=cr(i-OB,j-OB);
    if (fabs(tmp)>eps) return tmp>eps;
    return i.d<j.d;
}
void modify(point b[],int tot,int a[],int sum) 
{
    for (int i=2;i<=tot;i++)
        if (b[i].y<b[1].y-eps || (fabs(b[i].y-b[1].y)<eps && b[i].x<b[1].x-eps)) swap(b[1],b[i]);
    for (int i=2;i<=tot;i++) 
        b[i]._=atan2(b[i].y-b[1].y,b[i].x-b[1].x);
    OB=b[1];
    sort(b+2,b+tot+1,cmp2);
    b[++tot]=b[1];
    point P=b[2]-b[1],Q=b[tot-1]-b[1];
    for (int i=1;i<=sum;i++) {
    		point U=p[a[i]]-b[1];
        int k=search(b,tot,P,Q,U);
        if (!k) ans[a[i]]=1;
        point e,r;
        e=b[k+1]-b[k],r=p[a[i]]-b[k];
        if (cr(e,r)<-eps) ans[a[i]]=1;
    }
}
void inter(int u[],int tot,int a[],int sum)
{
    for (int i=n+1;i<=n+4;i++) u[++tot]=i;
    sort(u+1,u+tot+1,cmp);
    int h,r;
    st[h=r=1]=u[1]; 
    for (int i=2;i<=tot;i++) {
        if (fabs(cr(p[u[i]],p[u[i-1]]))<eps) continue;
        for (;(h<r) && (!check(p[st[r-1]],p[st[r]],p[u[i]]));r--) ;
        for (;(h<r) && (!check(p[st[h]],p[st[h+1]],p[u[i]]));h++) ;
        st[++r]=u[i];
    }
    for (;(h<r) && (!check(p[st[r-1]],p[st[r]],p[st[h]]));r--) ;
    st[r+1]=st[h];
    int cnt=0;
    for (int i=h;i<=r;i++) {
        point e;
        cross(p[st[i]],p[st[i+1]],e);
        e.x/=e.z,e.y/=e.z,e.z=1;
        b[++cnt]=e;
    }
    modify(b,cnt,a,sum);
}
void calc(int l,int r)
{
    if (l==r) return ;
    int mid=(l+r)>>1;
    calc(l,mid);
    calc(mid+1,r);
    int tot=0,sum=0;
    for (int i=l;i<=mid;i++)
        if (!c[i]) u[++tot]=i;
    for (int i=mid+1;i<=r;i++)
        if (c[i] && !ans[i]) a[++sum]=i;
    if (sum && tot) inter(u,tot,a,sum);
}
void origin()
{
    int tot=n;
    ++tot,p[tot].x=1,p[tot].y=0,p[tot].z=oo;
    ++tot,p[tot].x=0,p[tot].y=1,p[tot].z=oo;
    ++tot,p[tot].x=-1,p[tot].y=0,p[tot].z=oo;
    ++tot,p[tot].x=0,p[tot].y=-1,p[tot].z=oo;
    for (int i=1;i<=tot;i++) ori(p[i]);    
}
void point_reverse(point &p)
{
        double len=LR/p.dist();
        p.x*=len,p.y*=len;
}
void circle_reverse(double x,double y,point &P)
{
        point p,q,O(x,y);
        p.x=x+x,p.y=y+y,p.z=1;
        q.x=x-y,q.y=y+x,q.z=1;
        point_reverse(p),point_reverse(q);
        cross(p,q,P);
        point_reverse(O);
        if (dot(P,O)<-eps) {
            P.x=-P.x,P.y=-P.y,P.z=-P.z;
        }
}
int main()
{
    scanf("%d",&n);
    for (int i=1;i<=n;i++) {
        double x,y;
        scanf("%d%lf%lf",&c[i],&x,&y);
        if (!c[i]) circle_reverse(x,y,p[i]);
        else {
            p[i].x=x,p[i].y=y;
            point_reverse(p[i]);
        }
    }
    origin();
    calc(1,n);
    for (int i=1;i<=n;i++)
        if (c[i]) ans[i]=1;
        else break;
    for (int i=1;i<=n;i++)
        if (c[i]) {
            if (ans[i]) printf("No\n");
            else printf("Yes\n");
        }
    return 0;
}


updata：把交点存一下会快一点

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
#include <cmath>
#include <ctime>
#define sqr(x) ((x)*(x))
const double pi=acos(-1.0);
const double eps=1e-8;
const double oo=1e15,LR=100000.0;
using namespace std;
struct point{
    double x,y,z,_,d;
    point () {}
    point (double X,double Y) {x=X,y=Y,z=1;}
    point operator -(const point &p) const {
        point e;
        e.x=x-p.x,e.y=y-p.y,e.z=1;
        return e;
    }
    point operator +(const point &p) const {
        point e;
        e.x=x+p.x,e.y=y+p.y,e.z=1;
        return e;
        }
        point operator *(const double a) const {
            point e;
            e.x=x*a,e.y=y*a,e.z=1;
            return e;
        }
        void print() {
            printf("X=%f Y=%f Z=%f\n",(double)x,(double)y,(double)z);
        }
        double dist() {return (sqr(x)+sqr(y));}
}p[500050],b[500050],St[500050];
int ans[500050],a[500050],n,st[500050],u[500050],c[500050];
inline double cr(point e,point r) {return e.x*r.y-e.y*r.x;}
inline double dot(point &e,point &r) {return e.x*r.x+e.y*r.y+e.z*r.z;}
inline void cross(point &p,point &q,point &e)
{
    e.x=p.y*q.z-p.z*q.y;
    e.y=p.z*q.x-p.x*q.z;
    e.z=p.x*q.y-p.y*q.x;
}
inline void ori(point &a)
{
    a._=atan2(a.y,a.x);
    a.d=sqr(a.z)/(sqr(a.x)+sqr(a.y));
}
inline bool cmp(int i,int j)
{
    double tmp=p[i]._-p[j]._;
    if (fabs(tmp)>pi/2) return tmp<-eps;
    tmp=cr(p[i],p[j]);
    if (fabs(tmp)>eps) return tmp>eps;
    return p[i].d<p[j].d;
}
inline bool check(point &p,point I)
{
    if (dot(p,I)>eps) return 1;
    return 0;
}
inline bool equ(point &p) 
{
	return fabs(p.x)<eps && fabs(p.y)<eps && fabs(p.z)<eps;
}
int search(point b[],int tot,point &P,point &Q,point &U)
{
        if (cr(P,U)<-eps || cr(Q,U)>eps) return 0;
        int l=2,r=tot-1;
				point pp; 
        for (;l<=r;) {
            int mid=((l+r)>>1);
            pp=b[mid]-b[1];
            if (cr(pp,U)>eps) l=mid+1;
            else r=mid-1;
        }
        return l-1;
}
point OB;
bool cmp2(point i,point j)
{
    double tmp=i._-j._;
        if (fabs(tmp)>pi/2) return tmp<-eps;
    tmp=cr(i-OB,j-OB);
    if (fabs(tmp)>eps) return tmp>eps;
    return i.d<j.d;
}
void modify(point b[],int tot,int a[],int sum) 
{
    for (int i=2;i<=tot;i++)
        if (b[i].y<b[1].y-eps || (fabs(b[i].y-b[1].y)<eps && b[i].x<b[1].x-eps)) swap(b[1],b[i]);
    for (int i=2;i<=tot;i++) 
        b[i]._=atan2(b[i].y-b[1].y,b[i].x-b[1].x);
    OB=b[1];
    sort(b+2,b+tot+1,cmp2);
    b[++tot]=b[1];
    point P=b[2]-b[1],Q=b[tot-1]-b[1];
    for (int i=1;i<=sum;i++) {
    		point U=p[a[i]]-b[1];
        int k=search(b,tot,P,Q,U);
        if (!k) ans[a[i]]=1;
        point e,r;
        e=b[k+1]-b[k],r=p[a[i]]-b[k];
        if (cr(e,r)<-eps) ans[a[i]]=1;
    }
}
void inter(int u[],int tot,int a[],int sum)
{
    for (int i=n+1;i<=n+4;i++) u[++tot]=i;
    sort(u+1,u+tot+1,cmp);
    int h,r;
   st[h=r=1]=u[1]; 
   for (int i=2;i<=tot;i++) {
       if (fabs(cr(p[u[i]],p[u[i-1]]))<eps) continue;
       for (;(h<r) && (!check(St[r],p[u[i]]));r--) ;
       for (;(h<r) && (!check(St[h+1],p[u[i]]));h++) ;
       st[++r]=u[i];
       cross(p[st[r-1]],p[st[r]],St[r]);
       if (equ(St[r])) return ;
   }
   for (;(h<r) && (!check(St[r],p[st[h]]));r--) ;
   for (;(h<r) && (!check(St[h+1],p[st[r]]));h++) ;
   if (r-h+1<=2) return ;
    st[r+1]=st[h];
    int cnt=0;
    for (int i=h;i<=r;i++) {
        point e;
        cross(p[st[i]],p[st[i+1]],e);
        e.x/=e.z,e.y/=e.z,e.z=1;
        b[++cnt]=e;
    }
    modify(b,cnt,a,sum);
}
void calc(int l,int r)
{
    if (l==r) return ;
    int mid=(l+r)>>1;
    calc(l,mid);
    calc(mid+1,r);
    int tot=0,sum=0;
    for (int i=l;i<=mid;i++)
        if (!c[i]) u[++tot]=i;
    for (int i=mid+1;i<=r;i++)
        if (c[i] && !ans[i]) a[++sum]=i;
    if (sum && tot) inter(u,tot,a,sum);
}
void origin()
{
    int tot=n;
    ++tot,p[tot].x=1,p[tot].y=0,p[tot].z=oo;
    ++tot,p[tot].x=0,p[tot].y=1,p[tot].z=oo;
    ++tot,p[tot].x=-1,p[tot].y=0,p[tot].z=oo;
    ++tot,p[tot].x=0,p[tot].y=-1,p[tot].z=oo;
    for (int i=1;i<=tot;i++) ori(p[i]);    
}
void point_reverse(point &p)
{
        double len=LR/p.dist();
        p.x*=len,p.y*=len;
}
void circle_reverse(double x,double y,point &P)
{
        point p,q,O(x,y);
        p.x=x+x,p.y=y+y,p.z=1;
        q.x=x-y,q.y=y+x,q.z=1;
        point_reverse(p),point_reverse(q);
        cross(p,q,P);
        point_reverse(O);
        if (dot(P,O)<-eps) {
            P.x=-P.x,P.y=-P.y,P.z=-P.z;
        }
}
int main()
{
    scanf("%d",&n);
    for (int i=1;i<=n;i++) {
        double x,y;
        scanf("%d%lf%lf",&c[i],&x,&y);
        if (!c[i]) circle_reverse(x,y,p[i]);
        else {
            p[i].x=x,p[i].y=y;
            point_reverse(p[i]);
        }
    }
    origin();
    calc(1,n);
    for (int i=1;i<=n;i++)
        if (c[i]) ans[i]=1;
        else break;
    for (int i=1;i<=n;i++)
        if (c[i]) {
            if (ans[i]) printf("No\n");
            else printf("Yes\n");
        }
    return 0;
}



huyuncong
关注
专栏目录
【BZOJ4140】 共点圆加强版
CreationAugust is 14 years old forever
12-01 1911
Description在平面直角坐标系中，Wayne需要你完成n次操作，操作只有两种： 1.0 x y。表示在坐标系中加入一个以(x, y)为圆心且过原点的圆。 2.1 x y。表示询问点(x, y)是否在所有已加入的圆的内部（含圆周），且至少在一个圆内部（含圆周）。 为了减少你的工作量，题目保证圆心严格在x轴上方（纵坐标为正），且横坐标非零。 Input第1行一个整数n。 接下来n行，每
bzoj2961 共点圆
Elijahqi
02-10 341
http://www.elijahqi.win/2018/02/10/bzoj2961/ ‎ Description 　　在平面直角坐标系中，Wayne需要你完成n次操作，操作只有两种： 1.0 x y。表示在坐标系中加入一个以(x, y)为圆心且过原点的圆。 2.1 x y。表示询问点(x, y)是否在所有已加入的圆的内部（含圆周），且至少在一个圆内部（含圆周）。 为了减少你的工作量，...
【BZOJ2961】共点圆（CDQ分治）
weixin_30621959的博客
08-13 59
【BZOJ2961】共点圆（CDQ分治） 题面 BZOJ 题解 设询问点\((x,y)\)，圆心是\((X,Y)\) 那么如果点在园内的话就需要满足 \((X-x)^2+(Y-y)^2\le X^2+Y^2\) 拆开之后就变成了 \(x^2+y^2-2xX\le 2yY\) 除过去就是\(-\frac{x}{y}X+\frac{x^2+y^2}{2y}\le Y\) 显然左边是一个直线，那么，这个...
【BZOJ2961】共点圆（圆的反演）（半平面交）（CDQ分治）
zxyoi_dreamer的博客（不定期诈尸）
01-20 306
显而易见反演后上CDQ分治维护 半平面交 和 未定点 按照x排序即可
BZOJ 2961: 共点圆
月下沙茶树
05-11 619
突然想起来论文还没刷完QAQ简直智障 2013年xhr论文 卡斜率的时候用点积 #include #include #include #include #include using namespace std; #define rep(i,l,r) for(int i=l;i<=r;i++) #define per(i,r,l) for(int i=r;i>=l;i--) #define m
BZOJ 2961 共点圆 CDQ分治+凸包
世界
01-01 3036
题目大意：给定平面，多次插入点和圆，每次插入点时询问当前插入的点是否在之前插入的所有圆中并且至少在一个圆中 直接用数据结构维护这些点和圆不是很好写，我们考虑CDQ分治 对于每层分治，我们需要对于[mid+1,r]中的每个点求出[l,mid]中是否所有的圆都覆盖了这个点 设点的坐标为(x0,y0)，那么这个点在所有圆内必须满足对于所有的圆心(x,y)，(x-x0)^2+(y-y0)^2=x0^
精心整理2万字c++知识点
WangXinyu_090926的博客
11-30 4755
制作不易，求一键三连~ 文章目录1、基本数据结构数组链表队列、单调队列、双端队列栈2、中极数据结构堆并查集与带权并查集hash表自然溢出双hash3、高级数据结构树状数组线段树及其合并Zkw线段树Fhq线树超哥线段树平衡树Treap随机平衡二叉树Splay伸展树Scapegoat Tree替罪羊树后缀平衡树块状数组、块状链表树套树线段树套线段树线段树套平衡树平衡树套线段树可并堆左偏树配对堆KDTree、四分树舞蹈链（DLX）、二进制分组划分树4、可持久化数据结构可持久化线段树（主席树）可持久化平衡树可持久化
函数与表达式的反演变换
Recalling_Clouds的博客
01-20 944
下面提到的所有知识不一定同属一个范畴，但是它们的思想和作用相似，所以将其归为一篇文章。 莫比乌斯反演 请先学习狄利克雷卷积，否则您会一脸懵逼。 一般式： F(n)=∑d∣nG(d)→G(n)=∑d∣nμ(nd)F(d)F(n)=∑n∣dG(d)→G(n)=∑n∣dμ(dn)F(d) \begin{aligned} F(n) = \sum_{d | n} G(d) \rightarrow G(n) = \sum_{d | n} \mu(\frac{n}{d}) F(d) \\ F(n) = \sum_{n
【醒目】【业界偷懒】【Public】BZOJ题目一句话题解整理
热门推荐
CreationAugust is 14 years old forever
05-12 1万+
就当是复习一下自己做过的题,顺便提供一个简要题解给大家看. 做题时候实在想不出来看一下一句话题解,可以有一个提示的作用又不至于一下子知道了全部浪费了一道题吧.. 部分题目(如我A过得大部分奶牛题)是别人拿我的账号做的,不提供题解.已AC的题目(数学题均不提供分析过程,公式): 1000:A+B 1001:平面图最小割,转对偶图最短路 1002:矩阵树定理,也可以通过推矩阵的递推关系得到递推
[颓废史]蒟蒻的刷题记录
ws_fqk
01-05 3380
QAQ蒟蒻一枚，其实我就是来提供水题库的。 以下记录从2016年开始。1.1 1227: [SDOI2009]虔诚的墓主人 树状数组+离散化 3132: 上帝造题的七分钟 树状数组 二维区间加减+查询 3038: 上帝造题的七分钟2 线段树+剪枝1.2 1047: [HAOI2007]理想的正方形 二维单调队列维护最值1.4 2095: [Poi2010]Bridges 二分+混合图欧拉
信息学奥赛-省选及NOI课程表（2020.08.31）.pdf
11-25
4. **CDQ分治**：一种特殊的分治策略，用于处理复杂的问题。 5. **整体二分**：通过二分查找的方式优化搜索过程。 6. **块状链表**：一种特殊的链表数据结构，适用于某些特定的查询和更新操作。 7. **分治**：一种将...
【BZOJ4140】共点圆加强版（二进制分组）
小蒟蒻yyb的博客
08-13 277
题面 BZOJ 题解 我卡精度卡了一天。。。。 之前不强制在线的做法是CDQCDQCDQ分治，维护一个凸壳就好了。 现在改成二进制分组，每次重建凸壳就好了。。 #include&lt;iostream&gt; #include&lt;cstdio&gt; #include&lt;cstdlib&gt; #include&lt;cstring&gt; #include&lt;cm...
【BZOJ】2961: 共点圆 -圆的反演&cdq分治/二进制分组
远行客
02-21 328
圆的反演,调试巨坑
圆的反演
Nightmare_ak的博客
08-10 1137
首先给出参考的博客： http://wangzhpp.org/?p=106 http://blog.csdn.net/acdreamers/article/details/16966369 http://www.cnblogs.com/naturepengchen/articles/4298750.html 然后给出反演的动图： http://jingyan.baidu.com/artic
The Designer HDU - 6158 (圆的反演)
liexss的博客
08-28 180
题意：给你两个内切圆的半径，然后在两个圆之间放n个圆，这n个圆的面积是多少。 思路：学习完圆的反演后，这道题就变得很轻松了，设两个内切圆的切点为反演中心，这样两个圆就会反演成两条直线，n个小院反演成相同大小的圆切于两条直线中，然后依次计算每个小圆的面积即可，注意精度。 附上代码： #include <iostream> #include <algorithm&...
【科技】浅谈圆的反演
weixin_30849591的博客
02-05 374
一时兴起，就有了这篇博客。本人也学识浅薄，姑且讲一下我对于圆反演的一些皮毛之见。 首先我们要明白反演是什么： 反演是一种基本的几何变换。给定一个平面上的一个反演中心$O$和一个常数$k$，对于任意一个点$A(A \neq O)$，我们可以找到一个在直线$OA$上的点$A'$，使得线段$OA,OA'$的有向长度的乘积为$k$，那$A'$就是$A$关于$O$的反演点，可以证明这样的$A'$是唯一的...
bzoj2961 共点圆（cdq分治维护凸包，计算几何）
Icefox的博客
01-15 890
这题真是写的我有点绝望。菜不禁声。其实一开始是想拿KD tree树写的，听说能过，但是wa成狗？？至今改不对gg。一个点(x0,y0)在圆(x,y)内，要求 (x−x0)2+(y−y0)2<=x2+y2(x-x_0)^2+(y-y_0)^2<=x^2+y^2 化简一下就是 x20+y20<=2x0x+2y0yx_0^2+y_0^2<=2x_0x+2y_0y 假设y0>0y_0>0，我们有y>
圆的反演学习笔记
weixin_41597684的博客
10-16 1682
Part 1：相关性质： 一、反演的概念 设在平面内给定一点O和常数k（k不等于零），对于平面内任意一点A，确定A′，使A′为直线OA上一点，并且有向线段OA与OA′满足OA·OA′=k，我们称这种变换是以O为反演中心，以k为反演幂的反演变换，简称反演。称A′为A关于O（r）的互为反演点。 二、作已知点的反演点的方法 给出反演极O和反演幂k>0，作点A的反演点A′。 令k=r...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

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

抵扣说明：

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

余额充值