LOJ 2548 「JSOI2018」绝地反击 ——二分图匹配+网络流手动退流

题目：https://loj.ac/problem/2548

如果知道正多边形的顶点，就是二分答案、二分图匹配。于是写了个暴力枚举多边形顶点的，还很愚蠢地把第一个顶点枚举到 2*pi ，其实只要 ( frac{2*pi}{n} ) 就行了。

总之能得10分。

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define db double
using namespace std;
db Mx(db a,db b){return a>b?a:b;}
db Mn(db a,db b){return a<b?a:b;}
const int N=205,M=N*N;
const db eps=1e-7,Pls=1e-4,pi2=2*acos(-1);
int dcmp(db x)
{ if(x>eps)return 1;if(x<-eps)return -1;return 0;}

int n,R; db dis[N][N];
int hd[N],xnt,to[M],nxt[M],per[N],dfn[N],tim;
struct Node{ db x,y;}a[N],b[N];
db Sqr(db x){return x*x;}
db get_dis(Node u,Node v)
{return sqrt(Sqr(u.x-v.x)+Sqr(u.y-v.y));}
void add(int x,int y){to[++xnt]=y;nxt[xnt]=hd[x];hd[x]=xnt;}
bool xyl(int cr)
{
  for(int i=hd[cr],v;i;i=nxt[i])
    if(dfn[v=to[i]]!=tim)
      {
    dfn[v]=tim;
    if(!per[v]||xyl(per[v]))
      { per[v]=cr; return true;}
      }
  return false;
}
bool chk(db lm)
{
  memset(hd,0,sizeof hd); xnt=0;
  for(int i=1;i<=n;i++)
    for(int j=1;j<=n;j++)
      if(dcmp(dis[i][j]-lm)<=0)add(i,j);//<= not <
  memset(per,0,sizeof per);memset(dfn,0,sizeof dfn);//dfn!!
  for(int i=1;i<=n;i++)
    {tim=i; if(!xyl(i))return false;}
  return true;
}
int main()
{
  scanf("%d%d",&n,&R);
  for(int i=1;i<=n;i++)
    scanf("%lf%lf",&a[i].x,&a[i].y);
  db ans=300;
  for(db alp=0;alp<=pi2;alp+=Pls)
    {
      db bs=2*acos(-1)/n;
      for(int i=1;i<=n;i++)
    {
      alp+=bs; if(alp>pi2)alp-=pi2;
      b[i].x=R*cos(alp); b[i].y=R*sin(alp);
    }
      db l=0,r=0;
      for(int i=1;i<=n;i++)
    for(int j=1;j<=n;j++)
      dis[i][j]=get_dis(a[i],b[j]),r=Mx(r,dis[i][j]);
      r=Mn(r+eps,ans);
      while(r-l>=eps)
    {
      db mid=(l+r)/2;
      if(chk(mid))ans=mid,r=mid-eps;
      else l=mid+eps;
    }
    }
  printf("%.8f
",ans);
  return 0;
}

题解：https://www.cnblogs.com/cjyyb/p/10420259.html

首先，把二分放在外面，已知二分值，再考虑是否存在一个合法多边形。

已知二分值，一个点的可匹配范围是圆弧上一段区间。

合法方案可以转动多边形使得某个顶点卡在某个范围的边界上。一共有 O(n) 个边界，让多边形的第一个点分别卡上去即可。O( n⁴logn )。

多边形第一个点的转动角度可以对 ( frac{2*pi}{n} ) 取模。

因为多边形两点在圆弧上的间距是 ( frac{2*pi}{n} ) ，转动角度又小于 ( frac{2*pi}{n} ) ，所以不管怎么转，对一个点的匹配最多导致一条边被删除/加入。

每个点产生两个转动角度，可知一个使得可以多匹配一个顶点，另一个使得少匹配一个顶点；角度取模后排序，每次把一个角度的影响加入后，整个图多/少了一条边。

删掉一条边 ( u , v ) 的话，看看它如果有流量，就手动给源点到 u 的边、 v 到汇点的边改一下流量。注意给总流量减 1 。

每次删/加边之后跑一次网络流即可。注意删边也要跑。

找一个点在圆弧上的区间，用余弦公式可知 R , mid , dis 围成的三角形的一个角的 cos 。

要特判一个点可以走到圆上任一点，还要特判走不到圆上的任一点。不然 acos( ) 会有 nan 。

代码里是把圆的 x 轴上的那个点看作第 n 个点。转动看作顺时针转动。

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define db double
using namespace std;
int Mx(int a,int b){return a>b?a:b;}
int Mn(int a,int b){return a<b?a:b;}
db Sqr(db x){return x*x;}
const int N=205,N2=N<<1,M=N*N2;
const db jmp=1e-7,eps=1e-10,pi2=2*acos(-1);
int dcmp(db x)
{if(x>eps)return 1;if(x<-eps)return -1;return 0;}

int n,r2,nR,en; db x[N],y[N],alp[N],dis[N],d2[N],bs;
int hd[N2],cur[N2],xnt,to[M],nxt[M],cap[M],dfn[N2],q[N2];//N2
struct Node{
  db x;int u,v;bool fx;
  Node(db x=0,int u=0,int v=0,bool f=0):
    x(x),u(u),v(v),fx(f) {}
  bool operator< (const Node &b)const
  {return x<b.x;}
}t[N<<1];
void add(int x,int y)
{
  to[++xnt]=y;nxt[xnt]=hd[x];hd[x]=xnt;cap[xnt]=1;
  to[++xnt]=x;nxt[xnt]=hd[y];hd[y]=xnt;cap[xnt]=0;
}
void del(int x,int y,int &flow)
{
  int tp=0;
  for(int i=hd[x],pr;i;pr=i,i=nxt[i])
    if(to[i]==y)
      {if(i==hd[x])hd[x]=nxt[i]; else nxt[pr]=nxt[i]; break;}
  for(int i=hd[y],pr;i;pr=i,i=nxt[i])
    if(to[i]==x)
      { if(i==hd[y])hd[y]=nxt[i]; else nxt[pr]=nxt[i];
    tp=cap[i]; break;}
  if(!tp)return; flow--;
  for(int i=hd[x];i;i=nxt[i])
    if(to[i]==0){ cap[i]=0;cap[i^1]=1;break;}
  for(int i=hd[y];i;i=nxt[i])
    if(to[i]==en){ cap[i]=1;cap[i^1]=0;break;}
}
bool bfs()
{
  int he=0,tl=0; memset(dfn,0,sizeof dfn);
  dfn[0]=1; q[++tl]=0;
  while(he<tl)
    {
      int k=q[++he];
      for(int i=hd[k],v;i;i=nxt[i])
    if(cap[i]&&!dfn[v=to[i]])
      dfn[v]=dfn[k]+1, q[++tl]=v;
    }
  return dfn[en];
}
int dinic(int cr,int flow)
{
  if(cr==en)return flow;
  int use=0;
  for(int &i=cur[cr],v;i;i=nxt[i])
    if(cap[i]&&dfn[v=to[i]]==dfn[cr]+1)
      {
    int tmp=dinic(v,Mn(flow-use,cap[i]));
    if(!tmp)dfn[v]=0;
    use+=tmp; cap[i]-=tmp; cap[i^1]+=tmp;
    if(use==flow)return use;
      }
  return use;
}
bool chk(db mid)
{
  xnt=1;memset(hd,0,sizeof hd); int tot=0;
  db m2=Sqr(mid);
  for(int i=1;i<=n;i++)
    {
      if(mid>=dis[i]+nR)
    { for(int j=1;j<=n;j++)add(i,j+n); continue;}
      if(mid<dis[i]-nR)continue;////
      db fx=acos((r2+d2[i]-m2)/(2*nR*dis[i]));
      db l=alp[i]-fx, r=alp[i]+fx;
      if(l<0)l+=pi2; if(r<0)r+=pi2;//r<0 for alp<0
      int L=l/bs, R=r/bs;//bs not pi2
      //if(!L)L=n; if(!R)R=n;//not for l-L*bs
      t[++tot]=Node(l-L*bs,i,L?L:n,0);//bs not pi2
      t[++tot]=Node(r-R*bs,i,R?R:n,1);
      if(!L)L=n; if(!R)R=n;//
      if(L<=R)
    for(int j=L+1;j<=R;j++)add(i,j+n);
      else
    { for(int j=L+1;j<=n;j++)add(i,j+n);
      for(int j=1;j<=R;j++)add(i,j+n);}
    }
  int flow=0;
  for(int i=1;i<=n;i++)add(0,i),add(i+n,en);
  while(bfs())memcpy(cur,hd,sizeof hd),flow+=dinic(0,n);
  if(flow==n)return true;
  sort(t+1,t+tot+1);
  for(int i=1;i<=tot;i++)
    if(!t[i].fx)
      {
    add(t[i].u,t[i].v+n);
    if(bfs())memcpy(cur,hd,sizeof hd),flow+=dinic(0,n);
    if(flow==n)return true;
      }
    else
      {
    del(t[i].u,t[i].v+n,flow);
    if(bfs())memcpy(cur,hd,sizeof hd),flow+=dinic(0,n);//
      }
  return false;
}
int main()
{
  scanf("%d%d",&n,&nR); bs=pi2/n; en=(n<<1)+1; r2=Sqr(nR);
  for(int i=1;i<=n;i++)
    {
      scanf("%lf%lf",&x[i],&y[i]);
      alp[i]=atan2(y[i],x[i]);
      dis[i]=sqrt(Sqr(x[i])+Sqr(y[i]));
      d2[i]=Sqr(dis[i]);
    }
  db l=0,r=250,ans;
  while(r-l>jmp)
    {
      db mid=(l+r)/2;
      if(chk(mid))ans=mid,r=mid-jmp;
      else l=mid+jmp;
    }
  printf("%.8f
",ans);
  return 0;
}