hdu 2215 Maple trees

Problem - 2215

　　题意是，给出一些直径为1的树的位置，要求求出最小的能够把所有的树围起来的圆的半径。

　　开始的时候打算用类似模拟退火的方法，假定一个位置是最小圆的圆心，然后向每一个方向尝试移动圆心。如果能得到更小的一个圆，那么就将圆心移动到这里。如果把所有的方向都找遍都不能找到更小的圆，这时就要将点移动的距离减少，继续移动。直到移动的距离足够小的时候停止搜索。不过这样做，精度不能保证，所以最后还是要换回最小包围圈算法来通过这题。

　　下面的是包围圈暴力算法的代码，复杂度O(n^3)：

  1 #include <cstdio>
  2 #include <cstring>
  3 #include <iostream>
  4 #include <algorithm>
  5 #include <cmath>
  6 
  7 using namespace std;
  8 
  9 const double PI = acos(-1.0);
 10 const double EPS = 1e-8;
 11 inline int sgn(double x) { return (EPS < x) - (x < -EPS);}
 12 struct Point {
 13     double x, y;
 14     Point() {}
 15     Point(double x, double y) : x(x), y(y) {}
 16     bool operator < (Point a) const { return sgn(x - a.x) < 0 || sgn(x - a.x) == 0 && sgn(y - a.y) < 0;}
 17     bool operator == (Point a) const { return sgn(x - a.x) == 0 && sgn(y - a.y) == 0;}
 18 } ;
 19 typedef Point Vec;
 20 Vec operator + (Vec a, Vec b) { return Vec(a.x + b.x, a.y + b.y);}
 21 Vec operator - (Vec a, Vec b) { return Vec(a.x - b.x, a.y - b.y);}
 22 Vec operator * (Vec a, double p) { return Vec(a.x * p, a.y * p);}
 23 Vec operator / (Vec a, double p) { return Vec(a.x / p, a.y / p);}
 24 inline double dotDet(Vec a, Vec b) { return a.x * b.x + a.y * b.y;}
 25 inline double crossDet(Vec a, Vec b) { return a.x * b.y - a.y * b.x;}
 26 inline double crossDet(Point o, Point a, Point b) { return crossDet(a - o, b - o);}
 27 inline double vecLen(Vec x) { return sqrt(dotDet(x, x));}
 28 inline double ptDis(Point a, Point b) { return vecLen(a - b);}
 29 inline Vec normal(Vec x) {
 30     double len = vecLen(x);
 31     return Vec(-x.y / len, x.x / len);
 32 }
 33 inline Vec vecUnit(Vec x) { return x / vecLen(x);}
 34 
 35 struct Line {
 36     Point s, t;
 37     Line() {}
 38     Line(Point s, Point t) : s(s), t(t) {}
 39     Point point(double x) { return s + (t - s) * x;}
 40 } ;
 41 typedef Line Seg;
 42 Point lineIntersect(Point P, Vec v, Point Q, Vec w) {
 43     Vec u = P - Q;
 44     double t = crossDet(w, u) / crossDet(v, w);
 45     return P + v * t;
 46 }
 47 inline Point lineIntersect(Line a, Line b) { return lineIntersect(a.s, a.t - a.s, b.s, b.t - b.s);}
 48 
 49 Line midLine(Point a, Point b) {
 50     Point mid = (a + b) / 2.0;
 51     return Line(mid, mid + normal(a - b));
 52 }
 53 
 54 inline Point getPoint(Point a, Point b, Point c) {    return lineIntersect(midLine(a, b), midLine(b, c));}
 55 
 56 double getDis(Point a, Point b, Point c) {
 57     if (sgn(crossDet(a, b, c)) == 0) return 0.0;
 58     if (sgn(dotDet(a - b, c - b)) <= 0) return ptDis(a, c) / 2.0;
 59     if (sgn(dotDet(a - c, b - c)) <= 0) return ptDis(a, b) / 2.0;
 60     if (sgn(dotDet(b - a, c - a)) <= 0) return ptDis(b, c) / 2.0;
 61     Point cen = getPoint(a, b, c);
 62     if (sgn(ptDis(cen, a) - ptDis(cen, b)) || sgn(ptDis(cen, b) - ptDis(cen, c))) {
 63         cout << "shit!" << endl;
 64         while (1) ;
 65     }
 66     return ptDis(cen, a);
 67 }
 68 
 69 int andrew(Point *pt, int n, Point *ch) {
 70     sort(pt, pt + n);
 71     int m = 0;
 72     for (int i = 0; i < n; i++) {
 73         while (m > 1 && crossDet(ch[m - 2], ch[m - 1], pt[i]) <= 0) m--;
 74         ch[m++] = pt[i];
 75     }
 76     int k = m;
 77     for (int i = n - 2; i >= 0; i--) {
 78         while (m > k && crossDet(ch[m - 2], ch[m - 1], pt[i]) <= 0) m--;
 79         ch[m++] = pt[i];
 80     }
 81     if (n > 1) m--;
 82     return m;
 83 }
 84 
 85 const int N = 111;
 86 
 87 int main() {
 88 //    freopen("in", "r", stdin);
 89     int n;
 90     Point pt[N], ch[N];
 91     while (cin >> n && n) {
 92         for (int i = 0; i < n; i++) cin >> pt[i].x >> pt[i].y;
 93         n = andrew(pt, n, ch);
 94 //        cout << n << endl;
 95 //        for (int i = 0; i < n; i++) cout << ch[i].x << ' ' << ch[i].y << endl;
 96         double dis = 0.0;
 97         if (n == 2) {
 98             dis = ptDis(ch[0], ch[1]) / 2.0;
 99         } else {
100             for (int i = 0; i < n; i++) {
101                 for (int j = i + 1; j < n; j++) {
102                     for (int k = j + 1; k < n; k++) {
103                         dis = max(dis, getDis(ch[i], ch[j], ch[k]));
104                     }
105                 }
106             }
107         }
108         printf("%.2f\n", dis + 0.5);
109     }
110     return 0;
111 }
View Code

　　之后将更新增量法的最小包围圈算法，平均复杂度O(n)。

增量法：

  1 #include <cstdio>
  2 #include <cmath>
  3 #include <iostream>
  4 #include <algorithm>
  5 #include <cstring>
  6 
  7 using namespace std;
  8 
  9 const double EPS = 1e-10;
 10 inline int sgn(double x) { return (x > EPS) - (x < -EPS);}
 11 struct Point {
 12     double x, y;
 13     Point() {}
 14     Point(double x, double y) : x(x),y(y) {}
 15     bool operator < (Point a) const { return sgn(x - a.x) < 0 || sgn(x - a.x) == 0 && sgn(y - a.y) < 0;}
 16     bool operator == (Point a) const { return sgn(x - a.x) == 0 && sgn(y - a.y) == 0;}
 17     Point operator + (Point a) const { return Point(x + a.x, y + a.y);}
 18     Point operator - (Point a) const { return Point(x - a.x, y - a.y);}
 19     Point operator * (double p) const { return Point(x * p, y * p);}
 20     Point operator / (double p) const { return Point(x / p, y / p);}
 21 } ;
 22 typedef Point Vec;
 23 inline double crossDet(Vec a, Vec b) { return a.x * b.y - a.y * b.x;}
 24 inline double crossDet(Point o, Point a, Point b) { return crossDet(a - o, b - o);}
 25 inline double dotDet(Vec a, Vec b) { return a.x * b.x + a.y * b.y;}
 26 inline double vecLen(Vec x) { return sqrt(dotDet(x, x));}
 27 inline Point normal(Vec x) { return Point(-x.y, x.x) / vecLen(x);}
 28 
 29 Point lineIntersect(Point P, Vec v, Point Q, Vec w) {
 30     Vec u = P - Q;
 31     double t = crossDet(w, u) / crossDet(v, w);
 32     return P + v * t;
 33 }
 34 inline Point getMid(Point a, Point b) { return (a + b) / 2.0;}
 35 
 36 struct Circle {
 37     Point c;
 38     double r;
 39     Circle() {}
 40     Circle(Point c, double r) : c(c), r(r) {}
 41 } ;
 42 
 43 Circle getCircle(Point a, Point b, Point c) {
 44     Vec v1 = b - a, v2 = c - a;
 45     if (sgn(dotDet(b - a, c - a)) <= 0) return Circle(getMid(b, c), vecLen(b - c) / 2.0);
 46     if (sgn(dotDet(a - b, c - b)) <= 0) return Circle(getMid(a, c), vecLen(a - c) / 2.0);
 47     if (sgn(dotDet(a - c, b - c)) <= 0) return Circle(getMid(a, b), vecLen(a - b) / 2.0);
 48     Point ip = lineIntersect(getMid(a, b), normal(v1), getMid(a, c), normal(v2));
 49     return Circle(ip, vecLen(ip - a));
 50 }
 51 
 52 int andrew(Point *pt, int n, Point *ch) {
 53     sort(pt, pt + n);
 54     int m = 0;
 55     for (int i = 0; i < n; i++) {
 56         while (m > 1 && sgn(crossDet(ch[m - 2], ch[m - 1], pt[i])) <= 0) m--;
 57         ch[m++] = pt[i];
 58     }
 59     int k = m;
 60     for (int i = n - 2; i >= 0; i--) {
 61         while (m > k && sgn(crossDet(ch[m - 2], ch[m - 1], pt[i])) <= 0) m--;
 62         ch[m++] = pt[i];
 63     }
 64     if (n > 1) m--;
 65     return m;
 66 }
 67 
 68 const int N = 555;
 69 Point pt[N], ch[N];
 70 int rnd[N];
 71 
 72 void randPoint(Point *pt, int n) {
 73     for (int i = 0; i < n; i++) rnd[i] = (rand() % n + n) % n;
 74     for (int i = 0; i < n; i++) swap(pt[i], pt[rnd[i]]);
 75 }
 76 
 77 inline bool inCircle(Point p, Circle C) { return sgn(vecLen(C.c - p) - C.r) <= 0;}
 78 
 79 int main() {
 80 //    freopen("in", "r", stdin);
 81     int n;
 82     while (cin >> n && n) {
 83         for (int i = 0; i < n; i++) scanf("%lf%lf", &pt[i].x, &pt[i].y);
 84         n = andrew(pt, n, ch);
 85         randPoint(ch, n);
 86         Circle ans = Circle(ch[0], 0.0), tmp;
 87         for (int i = 0; i < n; i++) {
 88             if (inCircle(ch[i], ans)) continue;
 89             ans = Circle(ch[i], 0.0);
 90             for (int j = 0; j < i; j++) {
 91                 if (inCircle(ch[j], ans)) continue;
 92                 ans = Circle(getMid(ch[i], ch[j]), vecLen(ch[i] - ch[j]) / 2.0);
 93                 for (int k = 0; k < j; k++) {
 94                     if (inCircle(ch[k], ans)) continue;
 95                     ans = getCircle(ch[i], ch[j], ch[k]);
 96                 }
 97             }
 98         }
 99 //        printf("%.2f %.2f %.2f\n", ans.c.x, ans.c.y, ans.r);
100         printf("%.2f\n", ans.r + 0.5);
101     }
102     return 0;
103 }
View Code

　　

——written by Lyon

【推广】 免费学中医，健康全家人
原文地址：https://www.cnblogs.com/LyonLys/p/hdu_2215_Lyon.html