计算几何-随机增量

随机增量法

随机增量法可以用来解决最小圆覆盖。

首先，我们先思考一下这个问题：

给定平面上\(n\)个点，求一个半径最小的圆去覆盖这\(n\)个点。

我们可以先设点集\(A\)的最小圆覆盖为\(c(A)\)，对于一个最小覆盖圆，它肯定满足以下性质：

• \(c(A)\) 是唯一的;

• 圆上有三个（或以上）\(A\)中的点;

• 圆上有两个点为一条直径的两端.

其中第二条和第三条必须满足其中之一。

我们先假设目前已经求出了\(i-1\)个点的最小覆盖圆，在加入第\(i\)个点后，最小覆盖圆一定是:

• 前\(i-1\)个点中的两个点与第\(i\)个点三点确定的圆;

• 前\(i-1\)个点中的一个点与第\(i\)个点为直径作的圆.

额。。易证！

主要的代码部分：

for(int i=2;i<=n;i++)
    {
        if(check(center,r,d[i]))continue;
        geometric now;double nr=0;
        for(int j=1;j<i;j++)
        {
            if(check(now,nr,d[j]))continue;
            now=(d[i]+d[j])/2.0;
            nr=opt.dis(d[i],d[j])/2.0;

            // 先以一条直径作圆
            for(int k=1;k<j;k++)
            {
                if(check(now,nr,d[k]))continue;
                now=Center(d[i],d[j],d[k]);
                nr=opt.dis(now,d[i]);

                // 找第三个点作圆
            }
        }
        center=now;r=nr;
    }

然后我们发现，我们似乎写了一个\(O(n^3)\)的暴力。。其实不然。

时间复杂度的证明

由于\(n\)个点最多由三个点确定的最小覆盖圆，因此每个点参与确定最下覆盖圆的概率不大于\(\frac{3}{n}\)。

所以每一层在第\(i\)个点处调用下一层的概率不大于\(\frac{3}{i}\)。

设三个循环的时间复杂度为\(T_1,T_2,T_3\):

可以解得\(T_1=T_2=T_3=O(n)\)。

证毕。

但是这只是在数据随机的情况下，但是出题人往往不这么做，所以我们需要用$ random $ _ \(shuffle\)函数进行扰动就很好了。

最后\(n^3\)过百万！！（迫真）

Code:

#include<bits/stdc++.h>
#define eps 1e-8
using namespace std;
const int maxn=1e5+10;
struct geometric{
    double x,y;
    geometric(double X=0,double Y=0):x(X),y(Y) {}
    friend geometric operator + (const geometric a,const geometric b){return geometric(a.x+b.x,a.y+b.y);} 
    friend geometric operator - (const geometric a,const geometric b){return geometric(a.x-b.x,a.y-b.y);} 
    friend geometric operator * (const geometric a,double p){return geometric(a.x*p,a.y*p);}
    friend geometric operator / (const geometric a,double p){return geometric(a.x/p,a.y/p);}
    double dis(geometric a,geometric b){return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));}
    double dot(geometric a1,geometric a2,geometric b1,geometric b2){return (a2.x-a1.x)*(b2.x-b1.x)+(a2.y-a1.y)*(b2.y-b1.y);}
    double cross(geometric a1,geometric a2,geometric b1,geometric b2){return (a2.x-a1.x)*(b2.y-b1.y)-(a2.y-a1.y)*(b2.x-b1.x);}
    double corner(geometric a1,geometric a2,geometric b1,geometric b2){return dot(a1,a1,b1,b2)/(dis(a1,a2)*dis(b1,b2));}
    double area(geometric a1,geometric a2,geometric b1,geometric b2){return fabs(cross(a1,a2,b1,b2));}
    double angle(geometric a){return atan2(a.y,a.x);}
    geometric rotate_clockwise(geometric a,double theta){return geometric(a.x*cos(theta)-a.y*sin(theta),a.x*sin(theta)+a.y*cos(theta));}
    geometric rotate_counterclockwise(geometric a,double theta){return geometric(a.x*cos(theta)+a.y*sin(theta),-a.x*sin(theta)+a.y*cos(theta));}
}opt,origin,d[maxn],st[maxn];
int n,top;double S=0,rx,ry;
int Sure(double x){return fabs(x)<eps?0:(x<0?-1:1);}
bool check(geometric a,double r,geometric b){return Sure(r-opt.dis(a,b))>=0;}
geometric Center(geometric a,geometric b,geometric c)
{
    geometric mid1,mid2,cen;
    double k1=0,k2=0,b1=0,b2=0;
    if(a.y!=b.y)k1=-(a.x-b.x)/(a.y-b.y);
    if(b.y!=c.y)k2=-(b.x-c.x)/(b.y-c.y);
    mid1=(a+b)/2.0;mid2=(b+c)/2.0;
    b1=mid1.y-mid1.x*k1;b2=mid2.y-mid2.x*k2;
    if(k1==k2)
        cen=(a+c)/2.0;
    else 
    {
        cen.x=(b2-b1)/(k1-k2);
        cen.y=k1*cen.x+b1;
    }
    return cen;
}
int main()
{
    srand(time(0));
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
        scanf("%lf%lf",&d[i].x,&d[i].y);
    random_shuffle(d+1,d+n+1);
    geometric center=d[1];double r=0;
    for(int i=2;i<=n;i++)
    {
        if(check(center,r,d[i]))continue;
        geometric now;double nr=0;
        for(int j=1;j<i;j++)
        {
            if(check(now,nr,d[j]))continue;
            now=(d[i]+d[j])/2.0;
            nr=opt.dis(d[i],d[j])/2.0;
            for(int k=1;k<j;k++)
            {
                if(check(now,nr,d[k]))continue;
                now=Center(d[i],d[j],d[k]);
                nr=opt.dis(now,d[i]);
            }
        }
        center=now;r=nr;
    }
    printf("%0.10lf\n%0.10lf %0.10lf",r,center.x,center.y);
    return 0;
}

一些例题

最小圆覆盖

[AHOI2012]信号塔