BZOJ5316 : [Jsoi2018]绝地反击

若$R=0$，那么显然答案为离原点最远的点到原点的距离。

否则若所有点都在原点，那么显然答案为$R$。

否则考虑二分答案$mid$，检查$mid$是否可行。

那么每个点根据对应圆交，可以覆盖圆上的一部分，每个可行方案都可以通过平移使得刚好卡住某个交点。

枚举每个交点，算出圆上$n$个位置的坐标，然后匈牙利算法判断是否存在完美匹配，时间复杂度$O(n^4\log w)$，不能承受。

注意到这个图是个稠密图，所以可以用bitset对匈牙利进行加速，做到$O(\frac{n^3}{32})$每次匹配。

另一方面，可以先枚举一个点$x$，然后再二分答案$mid$，判断是否有可行方案使得$x$刚好匹配$x$和圆的交点。

在这里，显然只需要在之前答案$ans$的基础之上往下二分，如果$ans-eps$不可行那么就没有继续二分的必要。

即：设$f[x]$表示$x$得到的最优解，若$f[x]$不是$f[1,x]$的最小值，那么就没有继续二分的必要。

考虑将读入的$n$个点随机打乱，那么$f[x]$是$f[1,x]$的最小值的概率为$\frac{1}{x}$，一共只有期望$O(\log n)$个$x$有二分的必要。

检查次数骤降为$O(n+\log n\log w)$，时间复杂度$O(\frac{(n+\log n\log w)n^3}{32})$。

注意要特判$mid$过小或者过大导致$x$与圆没有交点的情况。

#include<cstdio>

#include<algorithm>

#include<cstdlib>

#include<cmath>

using namespace std;

typedef unsigned int U;

const int N=205,M=7;

const double eps=1e-9,pi=acos(-1.0);

int n,m,R,i,j;double lim,ans,rot[N][2],b[N][2];

int f[N];U v[M],g[N][M];

struct P{int x,y;}a[N];

inline int sgn(double x){

  if(x>eps)return 1;

  if(x<-eps)return -1;

  return 0;

}

bool find(int x){

  for(int i=0;i<=m;i++){

    U t=v[i]&g[x][i];

    while(t){

      int j=i<<5|__builtin_ctz(t);

      v[i]^=1U<<(j&31);

      if(f[j]<0||find(f[j]))return f[j]=x,1;

      t-=t&-t;

    }

  }

  return 0;

}

inline bool check(int A,int B,double C){

  double d=sqrt(a[A].x*a[A].x+a[A].y*a[A].y);

  double l=(a[A].x*a[A].x+a[A].y*a[A].y+C*C-R*R)/(2*d);

  double h=sqrt(max(C*C-l*l,0.0));

  double bx=-a[A].x/d,by=-a[A].y/d;

  double px=a[A].x+bx*l,py=a[A].y+by*l;

  by=-by;

  swap(bx,by);

  bx*=h,by*=h;

  if(B==0)px+=bx,py+=by;else px-=bx,py-=by;

  int i,j;

  b[0][0]=px,b[0][1]=py;

  for(i=1;i<n;i++){

    b[i][0]=px*rot[i][1]-py*rot[i][0];

    b[i][1]=px*rot[i][0]+py*rot[i][1];

  }

  C*=C;

  for(i=0;i<n;i++){

    for(j=0;j<=m;j++)g[i][j]=0;

    for(j=0;j<n;j++)if(sgn((a[i].x-b[j][0])*(a[i].x-b[j][0])+(a[i].y-b[j][1])*(a[i].y-b[j][1])-C)<=0)g[i][j>>5]|=1U<<(j&31);

  }

  for(i=0;i<n;i++)f[i]=-1;

  for(i=0;i<n;i++){

    for(j=0;j<=m;j++)v[j]=~0U;

    if(!find(i))return 0;

  }

  return 1;

}

int main(){

  scanf("%d%d",&n,&R);

  m=(n-1)>>5;

  for(i=0;i<n;i++)scanf("%d%d",&a[i].x,&a[i].y);

  if(!R){

    int ans=0;

    for(i=0;i<n;i++)ans=max(ans,a[i].x*a[i].x+a[i].y*a[i].y);

    double ret=sqrt(ans);

    return printf("%.15f",ret),0;

  }

  random_shuffle(a,a+n);

  for(i=1;i<n;i++){

    double o=pi*2*i/n;

    rot[i][0]=sin(o),rot[i][1]=cos(o);

  }

  for(i=0;i<n;i++){

    int t=a[i].x*a[i].x+a[i].y*a[i].y;

    double val;

    if(t<=R)val=R-sqrt(t);else val=sqrt(t)-R;

    lim=max(lim,val);

    ans=max(ans,sqrt(t)+R);

  }

  for(i=0;i<n;i++)if(a[i].x||a[i].y)for(j=0;j<2;j++){

    double l=lim,r=max(min(sqrt(a[i].x*a[i].x+a[i].y*a[i].y)+R,ans-eps),lim);

    if(!check(i,j,r))continue;

    while(l+eps<r){

      double mid=(l+r)/2;

      if(check(i,j,mid))r=ans=mid;else l=mid;

    }

  }

  return printf("%.15f",ans),0;

}

巴特西

BZOJ5316 : [Jsoi2018]绝地反击

最新文章

热门文章