【HDOJ】3007 Buried memory

1. 题目描述
有n个点，求能覆盖这n个点的半径最小的圆的圆心及半径。

2. 基本思路
算法模板http://soft.cs.tsinghua.edu.cn/blog/?q=node/1066
定义Di表示相对于P[1]和P[i]组成的最小覆盖圆，如果P[2..i-1]都在这个圆内，那么当前的圆心和半径即为最优解。
如果P[j]不在这个圆内，那么P[j]一定在新的最小覆盖圆的边界上即P[1]、P[j]、P[i]组成的圆。
因为三点可以确定一个圆，因此只需要不断的找到不满足的P[j]，进而更新最优解即可。
其实就是三层循环，不断更新最优解。然而，这个算法的期望复杂度是O（n）。这个比较难以理解。

3. 代码

  1 /* 3007 */
  2 #include <iostream>
  3 #include <sstream>
  4 #include <string>
  5 #include <map>
  6 #include <queue>
  7 #include <set>
  8 #include <stack>
  9 #include <vector>
 10 #include <deque>
 11 #include <bitset>
 12 #include <algorithm>
 13 #include <cstdio>
 14 #include <cmath>
 15 #include <ctime>
 16 #include <cstring>
 17 #include <climits>
 18 #include <cctype>
 19 #include <cassert>
 20 #include <functional>
 21 #include <iterator>
 22 #include <iomanip>
 23 using namespace std;
 24 //#pragma comment(linker,"/STACK:102400000,1024000")
 25 
 26 #define sti                set<int>
 27 #define stpii            set<pair<int, int> >
 28 #define mpii            map<int,int>
 29 #define vi                vector<int>
 30 #define pii                pair<int,int>
 31 #define vpii            vector<pair<int,int> >
 32 #define rep(i, a, n)     for (int i=a;i<n;++i)
 33 #define per(i, a, n)     for (int i=n-1;i>=a;--i)
 34 #define clr                clear
 35 #define pb                 push_back
 36 #define mp                 make_pair
 37 #define fir                first
 38 #define sec                second
 39 #define all(x)             (x).begin(),(x).end()
 40 #define SZ(x)             ((int)(x).size())
 41 #define lson            l, mid, rt<<1
 42 #define rson            mid+1, r, rt<<1|1
 43 
 44 typedef struct {
 45     double x, y;
 46 } Point;
 47 
 48 const double eps = 1e-6;
 49 const int maxn = 505;
 50 Point P[maxn];
 51 Point o;
 52 double r;
 53 int n;
 54 
 55 double Length(Point a, Point b) {
 56     return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
 57 }
 58 
 59 double Cross(Point a, Point b, Point c) {
 60     return (c.x-a.x)*(b.y-a.y) - (c.y-a.y)*(b.x-a.x);
 61 }
 62 
 63 Point Intersect(Point a, Point b, Point c, Point d) {
 64     Point ret = a;
 65     double t = ((a.x - c.x) * (c.y - d.y) - (a.y - c.y) * (c.x - d.x)) / 
 66                ((a.x - b.x) * (c.y - d.y) - (a.y - b.y) * (c.x - d.x)); 
 67     ret.x += (b.x - a.x) * t;
 68     ret.y += (b.y - a.y) * t;
 69     return ret;
 70 }
 71 
 72 Point circumcenter(Point a, Point b, Point c) {
 73     Point ua, ub, va, vb;
 74     
 75     ua.x = (a.x + b.x) / 2.0;
 76     ua.y = (a.y + b.y) / 2.0;
 77     ub.x = ua.x - a.y + b.y;
 78     ub.y = ua.y + a.x - b.x;
 79     
 80     va.x = (a.x + c.x) / 2.0;
 81     va.y = (a.y + c.y) / 2.0;
 82     vb.x = va.x - a.y + c.y;
 83     vb.y = va.y + a.x - c.x;
 84     return Intersect(ua, ub, va, vb);
 85 }
 86 
 87 void min_center() {
 88     o = P[0];
 89     r = 0;
 90     
 91     rep(i, 1, n) {
 92         if (Length(P[i], o)-r > eps) {
 93             o = P[i];
 94             r = 0;
 95             rep(j, 0, i) {
 96                 if (Length(P[j], o)-r > eps) {
 97                     o.x = (P[i].x + P[j].x) / 2;
 98                     o.y = (P[i].y + P[j].y) / 2;
 99                     r = Length(o, P[j]);
100                     
101                     rep(k, 0, j) {
102                         if (Length(P[k], o)-r > eps) {
103                             o = circumcenter(P[i], P[j], P[k]);
104                             r = Length(o, P[k]);
105                         }
106                     }
107                 }
108             }
109         }
110     }
111 }
112 
113 void solve() {
114     min_center();
115     printf("%.2lf %.2lf %.2lf
", o.x, o.y, r);
116 }
117 
118 int main() {
119     ios::sync_with_stdio(false);
120     #ifndef ONLINE_JUDGE
121         freopen("data.in", "r", stdin);
122         freopen("data.out", "w", stdout);
123     #endif
124     
125     while (scanf("%d", &n)!=EOF && n) {
126         rep(i, 0, n)
127             scanf("%lf%lf", &P[i].x, &P[i].y);
128         solve();
129     }
130     
131     #ifndef ONLINE_JUDGE
132         printf("time = %d.
", (int)clock());
133     #endif
134     
135     return 0;
136 }