bzoj 2961 共点圆 cdq+凸包+三分

题目大意

两种操作
1）插入一个过原点的圆
2）询问一个点是否在所有的圆中

分析

在圆中则在半径范围内
设圆心 (x,y) 查询点(x_0,y_0)
则(sqrt{(x-x_0)^2+(y-y_0)^2} <= sqrt{x^2+y^2})
解得(2x_0 * x+2y_0 *y -(x_0^2+y_0^2)>=0)
x,y 为变量
是个半平面的式子
题意变成
1)插入一个点
2)询问是否所有点都在半平面内
插入互不干扰
点都在半平面内当且仅当凸包在半平面内
cdq，维护上下凸包，三分找出离直线最近的点
（凸包上点与直线的叉积是单峰的）

solution

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cctype>
#include <cmath>
#include <algorithm>
using namespace std;
typedef double db;
const int M=500007;

inline int rd(){
	int x=0;bool f=1;char c=getchar();
	for(;!isdigit(c);c=getchar()) if(c=='-') f=0;
	for(;isdigit(c);c=getchar()) x=x*10+c-48;
	return f?x:-x;
}

int m,n,nup,ndw;

struct pt{
	db x,y;
	pt(db X=0.0,db Y=0.0){x=X;y=Y;}
}p[M],up[M],dw[M];

struct node{
	int kd;
	pt d;
	node(int kk=0,pt pp=pt()){kd=kk;d=pp;}
}opr[M];

struct line{
	pt P,v;
	line(pt pp=pt(),pt vv=pt()){P=pp; v=vv;}
};

bool ans[M];

bool operator <(pt x,pt y){return (x.x!=y.x)?(x.x<y.x):(x.y<y.y);}
pt operator +(pt x,pt y){return pt(x.x+y.x,x.y+y.y);}
pt operator -(pt x,pt y){return pt(x.x-y.x,x.y-y.y);}
pt operator *(pt x,db d){return pt(x.x*d,x.y*d);}
pt operator /(pt x,db d){return pt(x.x/d,x.y/d);}

db dot(pt x,pt y){
	return x.x*y.x+x.y*y.y;
}

db cross(pt x,pt y){
	return x.x*y.y-x.y*y.x;
}

db length(pt x){
	return sqrt(dot(x,x));
}

db shadow(pt x,pt y,pt to){
	return dot(y-x,to-x)/length(to-x);
}

db area(pt x,pt y,pt z){
	return cross(y-x,z-x);
}

void convex(){
	sort(p+1,p+n+1);
	
	int i;
	nup=0;
	for(i=1;i<=n;i++){
		while(nup>1&&area(up[nup-1],up[nup],p[i])>=0) nup--;
		up[++nup]=p[i];
	}
	
	ndw=0;
	for(i=1;i<=n;i++){
		while(ndw>1&&area(dw[ndw-1],p[i],dw[ndw])>=0) ndw--;
		dw[++ndw]=p[i];
	}
}

bool check(pt x){
	if(x.y==0){
		return x.x <= 2*up[1].x;
	}
	else{
		line nw=line(pt(0,dot(x,x)/(2.0*x.y)),pt(1,-x.x/x.y));
		if(x.y>0){
			int l=1,r=ndw,m1,m2,len;
			db tp1,tp2;
			while(l+1<r){
				len=(r-l+1)/3;
				m1=l+len;
				m2=m1+len;
				tp1=cross(nw.v,dw[m1]-nw.P);
				tp2=cross(nw.v,dw[m2]-nw.P);
				if(tp1<0||tp2<0) return 0;
				if(tp1<tp2) r=m2-1;
				else l=m1+1;
			}
			tp1=cross(nw.v,dw[l]-nw.P);
			tp2=cross(nw.v,dw[r]-nw.P);
			if(tp1<0||tp2<0) return 0;
		}
		else{
			int l=1,r=nup,m1,m2,len;
			db tp1,tp2;
			while(l+1<r){
				len=(r-l+1)/3;
				m1=l+len;
				m2=m1+len;
				tp1=cross(nw.v,up[m1]-nw.P);
				tp2=cross(nw.v,up[m2]-nw.P);
				if(tp1<0||tp2<0) return 0;
				if(tp1<tp2) l=m1+1;
				else r=m2-1;
			}
			tp1=cross(nw.v,up[l]-nw.P);
			tp2=cross(nw.v,up[r]-nw.P);
			if(tp1<0||tp2<0) return 0;
		}
	}
	return 1;
}

void solve(int l,int r){
	if(l>=r) return;
	int mid=l+r>>1,i;
	solve(l,mid);
	solve(mid+1,r);
	for(n=0,i=l;i<=mid;i++)
		if(opr[i].kd==0) p[++n]=opr[i].d;
	if(n==0) return;
	convex();
	for(i=mid+1;i<=r;i++)
	if(opr[i].kd==1&&ans[i]){
		ans[i]=check(opr[i].d);		
	}
}

int main(){
	db x,y;
	int i,kd;
	
	m=rd();
	bool ok=0;
	for(i=1;i<=m;i++){
		kd=rd();
		if(kd==0) ok=1;
		scanf("%lf%lf",&x,&y);
		opr[i]=node(kd,pt(x,y));
		ans[i]=ok;
	}
	
	solve(1,m);
	
	for(i=1;i<=m;i++)
		if(opr[i].kd==1) puts(ans[i]?"Yes":"No");
	
	return 0;
}