6496. 【GDOI2020模拟03.08】圣痕

题目描述

题解

二分半径，每条直线变成圆上的一个区间，圆内交点数就是相交区间对数，只要区间不穿过x轴正半轴就不会算错

得到半径后暴力统计即可，10^7*log很稳

注意不要统计圆上的点（会被卡成n^2），圆上的点距离都为r

code

#include <bits/stdc++.h>
#define fo(a,b,c) for (a=b; a<=c; a++)
#define fd(a,b,c) for (a=b; a>=c; a--)
#define low(x) ((x)&-(x))
#define abs(x) ((x)>0?(x):-(x))
#define E 0.000000000001
#define ll long long
#define file
using namespace std;

struct type{int x,y,id;} b[50001];
struct Type{double x;int id;} c[100001];
double d[100001],a[50001][2],s,X,Y,x,y,ans,L,R,Mid;
int tr[100001],d2[100001],Id[100001],N,n,m,i,j,k,l,tot,Tot,sum;
vector<int> Tr[400001];

bool cmp(Type a,Type b) {return a.x<b.x;}
bool Cmp(type a,type b) {return a.x<b.x || a.x==b.x && a.y<b.y;}

void swap(double &x,double &y)
{
	double z=x;
	x=y;y=z;
}

void Change(int t,int s) {while (t<=N) {tr[t]+=s,t+=low(t);}}
void change(int x,int y) {Change(x,1);Change(y+1,-1);}
int find(int t) {int ans=0; while (t) {ans+=tr[t],t-=low(t);} return ans;}

int pd(double r)
{
	double A,B,C,D,S,x,y,t1,t2;
	int i,j,k,l,sum=0;
	
	tot=0;
	fo(i,1,n)
	{
		A=a[i][0]*a[i][0]+1;
		B=2*a[i][0]*a[i][1];
		C=a[i][1]*a[i][1]-r*r;
		
		D=B*B-4*A*C;
		
		if (D>=0)
		{
			A*=2;
			D=sqrt(D);
			x=(-B-D)/A;y=a[i][0]*x+a[i][1];t1=atan2(y,x);t1+=(t1<0)?M_PI*2:0;
			x=(-B+D)/A;y=a[i][0]*x+a[i][1];t2=atan2(y,x);t2+=(t2<0)?M_PI*2:0;
			
			if (t1>t2) swap(t1,t2);
			d[++tot]=t1,Id[tot]=i;
			d[++tot]=t2,Id[tot]=i;
		}
	}
	if (!tot) return 0;
	
	fo(i,1,tot) c[i]={d[i],i};
	stable_sort(c+1,c+tot+1,cmp);
	
	N=0;
	fo(i,1,tot)
	{
		N+=i==1 || abs(c[i].x-c[i-1].x)>E;
		d2[c[i].id]=N;
	}
	
	Tot=0;
	for (i=1; i<=tot; i+=2)
	b[++Tot]={d2[i+1],d2[i],Id[i]};
	
	stable_sort(b+1,b+Tot+1,Cmp);
	
	memset(tr,0,sizeof(tr));
	fo(i,1,Tot)
	{
		sum+=find(b[i].y);
		change(b[i].y,b[i].x);
	}
	
	return sum;
}

void change2(int t,int l,int r,int x,int y,int s)
{
	int mid=(l+r)/2;
	
	if (x<=l && r<=y)
	{
		Tr[t].push_back(s);
		return;
	}
	
	if (x<=mid)
	change2(t*2,l,mid,x,y,s);
	if (mid<y)
	change2(t*2+1,mid+1,r,x,y,s);
}

void find2(int t,int l,int r,int x,int id)
{
	int mid=(l+r)/2,id2,N=Tr[t].size(),i;
	double X,Y;
	
	fo(i,0,N-1)
	{
		id2=Tr[t][i];
		
		if (abs(a[id][0]-a[id2][0])>E)
		{
			X=(a[id2][1]-a[id][1])/(a[id][0]-a[id2][0]);
			Y=a[id][0]*X+a[id][1];
			
			ans+=sqrt(X*X+Y*Y);
		}
	}
	
	if (l==r) return;
	
	if (x<=mid)
	find2(t*2,l,mid,x,id);
	else
	find2(t*2+1,mid+1,r,x,id);
}

void work(double r)
{
	int i;
	
	fo(i,1,Tot)
	{
		find2(1,1,N,b[i].y,b[i].id);
		change2(1,1,N,b[i].y,b[i].x,b[i].id);
	}
	ans+=r*(m-sum);
}

int main()
{
	freopen("stigmata.in","r",stdin);
	#ifdef file
	freopen("stigmata.out","w",stdout);
	#endif
	
	scanf("%d%lf%lf%d",&n,&X,&Y,&m);X/=1000,Y/=1000;
	fo(i,1,n)
	scanf("%lf%lf",&a[i][0],&a[i][1]),a[i][0]/=1000,a[i][1]/=1000,a[i][1]-=Y-X*a[i][0];
	
	L=0;R=10000000;
	while (R-L>E)
	{
		Mid=(L+R)/2;
		
		sum=pd(Mid);
		if (sum<m)
		{
			L=Mid;
			if (R-L<E)
			break;
		}
		else
		R=Mid;
	}
	
	work(Mid);
	
	printf("%.9lf
",ans);
	
	fclose(stdin);
	fclose(stdout);
	
	return 0;
}