LG1742 最小圆覆盖
给出(N)个点,让你画一个最小的包含所有点的圆。
(Nleq 100000)。
随机增量法
https://www.luogu.com.cn/blog/boshi/solution-p1742
定理:如果点(p)不在集合(S)的最小覆盖圆内,则(p)一定在(Scup{p})的最小覆盖圆上。
根据这个定理,我们可以分三次确定前(i)个点的最小覆盖圆。
- 令前(i-1)个点的最小覆盖圆为(C)
- 如果第(i)个点在(C)内,则前(i)个点的最小覆盖圆也是(C)
- 如果不在,那么第(i)个点一定在前(i)个点的最小覆盖圆上,接着确定前(i-1)个点中还有哪两个在最小覆盖圆上。因此,设当前圆心为(P_i),半径为0,做固定了第(i)个点的前(i)个点的最小圆覆盖。
- 固定了一个点:不停地在范围内找到第一个不在当前最小圆上的点(P_j),设当前圆心为((P_i+P_j)/2),半径为(frac{1}{2}|P_iP_j|),做固定了两个点的,前(j)个点外加第(i)个点的最小圆覆盖。
- 固定了两个点:不停地在范围内找到第一个不在当前最小圆上的点(P_k),设当前圆为(P_i,P_j,P_k)的外接圆。
复杂度证明
由于一堆点最多只有(3)个点确定了最小覆盖圆,因此(n)个点中每个点参与确定最小覆盖圆的概率不大于(frac{3}{n})。
所以,每一层循环在第(i)个点处调用下一层的概率不大于(frac{3}{i})。
那么设算法的三个循环的复杂度分别为(T_1(n),T_2(n),T_3(n)),则有:
不难解得,(T_1(n)=T_2(n)=T_3(n)=O(n))。
struct point {float128 x,y;};
IN point operator+(CO point&a,CO point&b){
return {a.x+b.x,a.y+b.y};
}
IN point operator-(CO point&a,CO point&b){
return {a.x-b.x,a.y-b.y};
}
IN point operator*(CO point&a,float128 k){
return {a.x*k,a.y*k};
}
IN point operator/(CO point&a,float128 k){
return {a.x/k,a.y/k};
}
IN float128 dot(CO point&a,CO point&b){
return a.x*b.x+a.y*b.y;
}
IN float128 cross(CO point&a,CO point&b){
return a.x*b.y-a.y*b.x;
}
IN float128 dist(CO point&a,CO point&b){
return sqrt(dot(b-a,b-a));
}
IN point rotate90(CO point&a){
return {a.y,-a.x};
}
IN point intersect(CO point&a,CO point&u,CO point&b,CO point&v){
return b+v*(cross(b-a,u)/cross(u,v));
}
IN point center(CO point&a,CO point&b,CO point&c){
return intersect((a+b)/2,rotate90(b-a),(a+c)/2,rotate90(c-a));
}
CO int N=1e5+10;
point p[N];
int main(){
srand(20030506);
int n=read<int>();
for(int i=1;i<=n;++i) scanf("%Lf%Lf",&p[i].x,&p[i].y);
random_shuffle(p+1,p+n+1);
point c={};
float128 r=0;
for(int i=1;i<=n;++i)if(dist(p[i],c)>r){
c=p[i],r=0;
for(int j=1;j<i;++j)if(dist(p[j],c)>r){
c=(p[i]+p[j])/2,r=dist(p[j],c);
for(int k=1;k<j;++k)if(dist(p[k],c)>r){
c=center(p[i],p[j],p[k]),r=dist(p[k],c);
}
}
}
printf("%.10Lf
%.10Lf %.10Lf
",r,c.x,c.y);
return 0;
}
LOJ6360 复燃「恋之埋火」
古明地恋(koishi)和小石子(koishi)是好朋友。旧地狱的空中散布着许多颗小石子。恋恋想找出一个位置,使得这个位置离最远的小石子的距离尽可能小。
需要注意的是,这里的空间可能是高维空间。
对于100%的数据,(nleq 20000,mleq 5,0leq)所有坐标(leq 10^4)
题解
http://jklover.hs-blog.cf/2020/04/13/Loj-6360-复燃「恋之埋火」/
高维空间的最小圆覆盖.
题意就是要求(m)维空间内(n)个点的最小圆覆盖。
把二维的随机增量算法直接拓展到高维,不难发现复杂度的分析是类似的。
当我们已经加入(k)个点时,需要在这(k)个点所在的(k−1)维平面上找一个圆心,使得它到这(k)个点距离相等.
从这(k)个点中选一个点作为原点,它到另外(k−1)个点的(k−1)个向量显然就是这个平面的一组基底。
只需要求出这些向量各自的系数即可确定圆心。
根据圆心到原点和到其他任意一个点的距离相同,可以列出(k−1)个方程,高斯消元求出系数,进而确定圆心位置.
时间复杂度(O(nm^3))。
CO int N=2e4+10,M=10;
int n,m;
typedef valarray<float128> point;
float128 dist(CO point&a,CO point&b){
float128 ans=0;
for(int i=0;i<m;++i) ans+=(a[i]-b[i])*(a[i]-b[i]);
return ans;
}
int id[M],k=0;
point p[N],center,vec[M];
float128 radius,a[M][M];
void find_center(){
for(int i=0;i<=k;++i) fill(a[i],a[i]+k+1,0);
for(int i=1;i<k;++i) vec[i]=p[id[i]]-p[id[0]];
for(int i=1;i<k;++i)for(int j=0;j<m;++j){
a[i][k]+=vec[i][j]*vec[i][j];
for(int l=1;l<k;++l) a[i][l]+=2*vec[i][j]*vec[l][j];
}
for(int i=1;i<k;++i){
int p=i;
for(int j=i+1;j<k;++j)if(abs(a[j][i])>abs(a[p][i])) p=j;
if(p!=i) swap(a[i],a[p]);
for(int j=1;j<k;++j)if(j!=i){
float128 t=a[j][i]/a[i][i];
for(int l=1;l<=k;++l) a[j][l]-=a[i][l]*t;
}
}
for(int i=0;i<m;++i) center[i]=0;
for(int i=1;i<k;++i) center=center+(vec[i]*a[i][k]/a[i][i]);
center=center+p[id[0]];
}
void dfs(int x,int bound){
if(x>m) return;
for(int i=0;i<bound;++i)if(dist(p[i],center)>radius){
id[k++]=i;
find_center();
radius=dist(center,p[i]);
dfs(x+1,i);
--k;
}
}
int main(){
srand(20030506);
read(n),read(m);
for(int i=0;i<n;++i){
p[i].resize(m);
for(int j=0;j<m;++j) scanf("%Lf",&p[i][j]);
}
random_shuffle(p,p+n);
center=p[0],radius=-1;
dfs(0,n);
for(int i=0;i<m;++i) printf("%.6Lf%c",center[i],"
"[i==m-1]);
return 0;
}
这个
valarray
到底是何方神圣?