包含點集所有點的最小圓的算法
?http://acm.zju.edu.cn/onlinejudge/showProblem.do?problemId=450
平面上有n個點,給定n個點的坐標,試找一個半徑最小的圓,將n
個點全部包圍,點可以在圓上。
1. 在點集中任取3點A,B,C。
2. 作一個包含A,B,C三點的最小圓,圓周可能通過這3點,也可能只通過
其中兩點,但包含第3點.后一種情況圓周上的兩點一定是位于圓的一條直
徑的兩端。
3. 在點集中找出距離第2步所建圓圓心最遠的D點,若D點已在圓內或圓周上,
則該圓即為所求的圓,算法結束.則,執行第4步。
4. 在A,B,C,D中選3個點,使由它們生成的一個包含這4個點的圓為最小,這3
點成為新的A,B,C,返回執行第2步。若在第4步生成的圓的圓周只通過A,B,C,D
中的兩點,則圓周上的兩點取成新的A和B,從另兩點中任取一點作為新的C。
?
程序設計題解上的解題報告:
對于一個給定的點集A,記MinCircle(A)為點集A的最小外接圓,顯然,對于所
有的點集情況A,MinCircle(A)都是存在且惟一的。需要特別說明的是,當A為空
集時,MinCircle(A)為空集,當A={a}時,MinCircle(A)圓心坐標為a,半徑為0;
?? 顯然,MinCircle(A)可以有A邊界上最多三個點確定(當點集A中點的個數大于
1時,有可能兩個點確定了MinCircle(A)),也就是說存在著一個點集B,|B|<=3
且B包含與A,有MinCircle(B)=MinCircle(A).所以,如果a不屬于B,則
MinCircle(A-{a})=MinCircle(A);如果MinCircle(A-{a})不等于MinCircle(A),則
a屬于B。
??? 所以我們可以從一個空集R開始,不斷的把題目中給定的點集中的點加入R,同
時維護R的外接圓最小,這樣就可以得到解決該題的算法。
代碼#include <stdio.h>
#include <math.h>
const int maxn = 1005;
const double eps = 1e-6;
struct TPoint {double x, y;TPoint operator-(TPoint & a) {TPoint p1;p1.x = x - a.x;p1.y = y - a.y;return p1;}
};
struct TCircle {double r;TPoint centre;
};
struct TTriangle {TPoint t[3];
};
TCircle c;
TPoint a[maxn];
double distance(TPoint p1, TPoint p2) {TPoint p3;p3.x = p2.x - p1.x;p3.y = p2.y - p1.y;return sqrt(p3.x * p3.x + p3.y * p3.y);
}
double triangleArea(TTriangle t) {TPoint p1, p2;p1 = t.t[1] - t.t[0];p2 = t.t[2] - t.t[0];return fabs(p1.x * p2.y - p1.y * p2.x) / 2;
}
TCircle circumcircleOfTriangle(TTriangle t) {//三角形的外接圓TCircle tmp;double a, b, c, c1, c2;double xA, yA, xB, yB, xC, yC;a = distance(t.t[0], t.t[1]);b = distance(t.t[1], t.t[2]);c = distance(t.t[2], t.t[0]);//根據S = a * b * c / R / 4;求半徑Rtmp.r = a * b * c / triangleArea(t) / 4;xA = t.t[0].x;yA = t.t[0].y;xB = t.t[1].x;yB = t.t[1].y;xC = t.t[2].x;yC = t.t[2].y;c1 = (xA * xA + yA * yA - xB * xB - yB * yB) / 2;c2 = (xA * xA + yA * yA - xC * xC - yC * yC) / 2;tmp.centre.x = (c1 * (yA - yC) - c2 * (yA - yB)) /((xA - xB) * (yA - yC) - (xA - xC) * (yA - yB));tmp.centre.y = (c1 * (xA - xC) - c2 * (xA - xB)) /((yA - yB) * (xA - xC) - (yA - yC) * (xA - xB));return tmp;
}TCircle MinCircle2(int tce, TTriangle ce) {TCircle tmp;if (tce == 0) tmp.r = -2;else if (tce == 1) {tmp.centre = ce.t[0];tmp.r = 0;} else if (tce == 2) {tmp.r = distance(ce.t[0], ce.t[1]) / 2;tmp.centre.x = (ce.t[0].x + ce.t[1].x) / 2;tmp.centre.y = (ce.t[0].y + ce.t[1].y) / 2;} else if (tce == 3) tmp = circumcircleOfTriangle(ce);return tmp;
}void MinCircle(int t, int tce, TTriangle ce) {int i, j;TPoint tmp;c = MinCircle2(tce, ce);if (tce == 3) return;for (i = 1; i <= t; i++) {if (distance(a[i], c.centre) > c.r) {ce.t[tce] = a[i];MinCircle(i - 1, tce + 1, ce);tmp = a[i];for (j = i; j >= 2; j--) {a[j] = a[j - 1];}a[1] = tmp;}}
}int main() {//freopen("circle.in", "r", stdin);//freopen("out.txt", "w", stdout);int n,i;while (scanf("%d", &n) != EOF && n) {for (i = 1; i <= n; i++)scanf("%lf%lf", &a[i].x, &a[i].y);TTriangle ce;MinCircle(n, 0, ce);printf("%.2lf %.2lf %.2lf\n", c.centre.x, c.centre.y, c.r);}return 0;
}
附上幾組測試數據
?
代碼
2 2
1 1
2
3
0 0
0
1
1
1 0
0
1
2
2
1
1
1
0 0
0
0
0
0
1
1
1
1
1
2
2
2
2
2
3
3
3
3
3
4
4
4
4
4
0 0
0
1
2
3
4
5
6
7
8
1
2
1
2
3
1
4
4
4
4
1
5
5
5
1
6
6
8
8
8
6
7
3
4
7
5
0 0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
0 0
1
2
3
4
5
6
7
8
9
1
2
3
4
5
6
7
8
2
3
4
5
6
7
3
4
5
6
4
5
0 0
1
2
3
0
1
2
3
0
1
2
3
4
0
0
0 3
2
1
0
1
2
3
2
1
0
1
2
3
?
總結
以上是生活随笔為你收集整理的ZOJ1450 Minimal Circle 最小圆覆盖的全部內容,希望文章能夠幫你解決所遇到的問題。
如果覺得生活随笔網站內容還不錯,歡迎將生活随笔推薦給好友。