zoukankan      html  css  js  c++  java
  • P1742 最小圆覆盖(计算几何)

    体验过(O(n^3))(10^5)吗?快来体验一波当(wys)的快感吧(QAQ)

    前置芝士1:二元一次方程组求解

    [egin{cases}a1 * x + b1*y=c1\a2 * x + b2*y=c2end{cases} ]

    (其中(a1,a2,b1,b2,c1,c2)为已知量)

    (②)式得:

    [x=frac{c2-b2*y}{a2} ]

    带入(①)式并化简得:

    [y=frac{c1-frac{a1*c2}{a2}}{b1-frac{a1*b2}{a2}} ]

    分子分母同时乘以(a2)得:

    [y=frac{a2*c1-a1*c2}{a2*b1-a1*b2} ]

    同理可得(把(a,b)互换即可):

    [x=frac{b2*c1-b1*c2}{b2*a1-b1*a2} ]

    前置芝士2:三点定圆

    给出三个点,求出圆心&半径

    [egin{cases}x1^2-2x1*x0+x0^2+y1^2-2y1*y0+y0^2=r^2\x2^2-2x2*x0+x0^2+y2^2-2y2*y0+y0^2=r^2\x3^2-2x3*x0+x0^2+y3^2-2y3*y0+y0^2=r^2end{cases} ]

    (②-①)(③-①),并化简得:

    [egin{cases}2*(x2-x1)x+2*(y2-y1)y=x2^2-x1^2+y2^2-y1^2\2*(x3-x1)x+2*(y3-y1)y=x3^2-x1^2+y3^2-y1^2end{cases} ]

    我们将三点定圆的柿子对应二元一次方程组中,可知:

    [a1=x2-x1,quad a2=x3-x1 ]

    [b1=y2-y1,quad b2=y3-y1 ]

    [c1=frac{x2^2-x1^2+y2^2-y1^2}{2},quad c2=frac{x3^2-x1^2+y3^2-y1^2}{2} ]

    然后就可以根据三个点求出圆心和半径了

    正文

    跟据前置芝士,我们知道对于任意三个不共线的点,我们可以求出三点定的圆,所以一个明显的想法就是枚举三个点

    我们先枚举第一个点,有两种情况

    ①:当前点在当前外面,即(dis()圆心,该点()>r)那么我们不管这个点

    ②:不是情况①的情况,那么我们就需要重新构造这个圆来包含所有的点了

    怎么构造呢?我们重新枚举两外两个已经遍历过的点,组成三个点。同理,若重新构造的圆包括了三个点,那么就不管,若有任意一个在圆外,那么我们根据前置芝士重新确定圆心和半径即可

    PS:本题出题人过于duliu,故意构造数据卡掉了上述解法,所以我们需要一个神奇的东西:随((da))((luan))((shu))((ju))法,来防止掉精度

    代码如下:

    #include<bits/stdc++.h>
    using namespace std;
    #define il inline
    #define re register
    #define D double
    il int read()
    {
        re int x = 0, f = 1; re char c = getchar();
        while(c < '0' || c > '9') { if(c == '-') f = -1; c = getchar();}
        while(c >= '0' && c <= '9') x = x * 10 + c - 48, c = getchar();
        return x * f;
    }
    #define rep(i, s, t) for(re int i = s; i <= t; ++ i)
    #define eps 1e-12
    #define maxn 100005
    #define ff(x) (x) * (x)
    int n, m;
    D r;
    struct node
    {
    	D x, y;
    }o, e[maxn];
    il D dis(node a, node b){return sqrt(ff(a.x - b.x) + ff(a.y - b.y));}
    il void get(node a, node b, node c)
    {
    	D a1 = b.x - a.x, a2 = c.x - a.x, b1 = b.y - a.y, b2 = c.y - a.y;
    	D c1 = (ff(b.x) - ff(a.x) + ff(b.y) - ff(a.y));
    	D c2 = (ff(c.x) - ff(a.x) + ff(c.y) - ff(a.y));
    	o = (node){(b2 * c1 - b1 * c2) / (b2 * a1 * 2 - b1 * a2 * 2), 
    			   (a2 * c1 - a1 * c2) / (a2 * b1 * 2 - a1 * b2 * 2)};
    	r = dis(a, o);
    }
    il void work()
    {
    	o = e[1], r = 0;
    	rep(i, 2, n)
    	{
    		if(dis(o, e[i]) > r + eps)
    		{
    			o = e[i], r = 0;
    			rep(j, 1, i - 1)
    			{
    				if(dis(o, e[j]) > r + eps)
    				{
    					o.x = (e[i].x + e[j].x) / 2, o.y = (e[i].y + e[j].y) / 2;
    					r = dis(o, e[j]);
    					rep(k, 1, j - 1) if(dis(o, e[k]) > r + eps) get(e[i], e[j], e[k]);
     				}
    			}
    		}
    	//	printf("%.10lf
    %.10lf %.10lf
    ", r, o.x, o.y);
    	}
    }
    int main()
    	n = read();
    	rep(i, 1, n) scanf("%lf%lf", &e[i].x, &e[i].y);
        random_shuffle(e + 1, e + n + 1);
    	work();
    	printf("%.10lf
    %.10lf %.10lf", r, o.x, o.y);
    	return 0;
    }
    
  • 相关阅读:
    《作业二》总结
    《作业一》总结
    团队项目-需求分析报告
    团队项目-选题报告
    第一次结对编程作业
    第一次个人编程作业
    第一次博客作业
    第12组 团队项目-需求分析报告
    团队项目-选题报告
    第二次结对编程作业
  • 原文地址:https://www.cnblogs.com/bcoier/p/10515731.html
Copyright © 2011-2022 走看看