1. 程式人生 >實用技巧 >其它 >ZOJ 1450 Minimal Circle

ZOJ 1450 Minimal Circle


問題連結:https://zoj.pintia.cn/problem-sets/91827364500/problems/91827364949

隨機增量演算法主要思想:

  1. 隨機打亂點集中所有點的順序;

  2. 以序列中前2個點構造最小包圍圓;

  3. 一次搜尋其餘的點,若某個點v在D外,則重新構造最小包圍圓D',使D'包含D中的所有點且以v為邊界點。重複執行第三步,直到求出最小包圍圓。

這個演算法的關鍵在於第 3 步,以 v 作為邊界點構造最小包圍圓,從而增加了構造約束,減少了自由度。該演算法的記憶體耗用少,速度快。在隨機點集的條件下,期望執行時間為 O(n)。但是,若原始輸入序列是有序的,則打亂點集的順序也是需要耗費時間的。

核心程式碼:

void min_cover_circle(point &c, double &r, int n) {
	random_shuffle(p + 1, p + 1 + n);
	c = p[1];
	r = 0;
	for (int i = 2; i <= n; ++i) {
		if (dis(c, p[i] ) > r + eps) {
			c = p[i];
			r = 0;
			for (int j = 1; j < i; ++j) {
				if (dis(c, p[j]) > r + eps) {
					c.x = (p[i].x + p[j].x) / 2;
					c.y = (p[i].y + p[j].y) / 2;
					r = dis(c, p[i]);
					for (int k = 1; k < j; ++k) {
						if (dis(c, p[k]) > r + eps) {
							// 根據點p[i],p[j],p[k]計算包圍圓圓心
							c = circumcenter(p[i], p[j], p[k]);
							r = dis(c, p[i]);
						}
					}
				}
			}
		}
	}
}

計算圓心演算法:

設圓座標為

ZOJ 1450 Minimal Circle

,且過點

ZOJ 1450 Minimal Circle

,則圓心c座標為

ZOJ 1450 Minimal Circle

,

ZOJ 1450 Minimal Circle

其中

ZOJ 1450 Minimal Circle

,

ZOJ 1450 Minimal Circle

,

ZOJ 1450 Minimal Circle

ZOJ 1450 Minimal Circle

,

ZOJ 1450 Minimal Circle

,

ZOJ 1450 Minimal Circle

point circumcenter(point A, point B, point C) {
	point ans;

	//double a1 = 2 * (B.x - A.x), b1 = 2 * (B.y - A.y), c1 = pow(B.x, 2) + pow(B.y, 2) - pow(A.x, 2) - pow(A.y, 2);
	//double a2 = 2 * (C.x - B.x), b2 = 2 * (C.y - B.y), c2 = pow(C.x, 2) + pow(C.y, 2) - pow(B.x, 2) - pow(B.y, 2);
	//double d = a1 * b2 - a2 * b1;
	//ans.x = (c1 * b2 - c2 * b1) / d;
	//ans.y = (a1 * c2 - a2 * c1) / d;

	double a1 = B.x - A.x, b1 = B.y - A.y, c1 = (pow(a1, 2) + pow(b1, 2)) / 2;
	double a2 = C.x - A.x, b2 = C.y - A.y, c2 = (pow(a2, 2) + pow(b2, 2)) / 2;
	double d = a1 * b2 - a2 * b1;
	ans.x = A.x + (c1 * b2 - c2 * b1) / d;
	ans.y = A.y + (a1 * c2 - a2 * c1) / d;
	
	return ans;
}

完整的程式碼:

#include <iostream>
#include <algorithm>
#include <math.h>

using namespace std;

const int maxn = 500;
const double eps = 1e-8;

struct point {
	double x, y;
	point(double _x = 0, double _y = 0) {
		x = _x; y = _y;
	}
	point operator + (const point &p) {
		return point(x + p.x, y + p.y);
	}
	point operator - (const point &p) {
		return point(x - p.x, y - p.y);
	}
};

double dis(point p, point q) {
	return pow(pow((p.x - q.x), 2) + pow((p.y - q.y), 2), 0.5);
}

point p[maxn];

point circumcenter(point A, point B, point C) {
	point ans;

	double a1 = B.x - A.x, b1 = B.y - A.y, c1 = (pow(a1, 2) + pow(b1, 2)) / 2;
	double a2 = C.x - A.x, b2 = C.y - A.y, c2 = (pow(a2, 2) + pow(b2, 2)) / 2;
	double d = a1 * b2 - a2 * b1;
	ans.x = A.x + (c1 * b2 - c2 * b1) / d;
	ans.y = A.y + (a1 * c2 - a2 * c1) / d;

	return ans;
}

void min_cover_circle(point &c, double &r, int n) {
	random_shuffle(p + 1, p + 1 + n);
	c = p[1];
	r = 0;
	for (int i = 2; i <= n; ++i) {
		if (dis(c, p[i] ) > r + eps) {
			c = p[i];
			r = 0;
			for (int j = 1; j < i; ++j) {
				if (dis(c, p[j]) > r + eps) {
					c.x = (p[i].x + p[j].x) / 2;
					c.y = (p[i].y + p[j].y) / 2;
					r = dis(c, p[i]);
					for (int k = 1; k < j; ++k) {
						if (dis(c, p[k]) > r + eps) {
							// 由p[i],p[j],p[k]三點計算包圍圓圓心
							c = circumcenter(p[i], p[j], p[k]);
							r = dis(c, p[i]);
						}
					}
				}
			}
		}
	}
}

int main() {
	int n;
	point c;
	double r;

	while (scanf("%d", &n) != EOF && n) {
		for (int i = 1; i <= n; i++) {
			scanf("%lf %lf", &p[i].x, &p[i].y);
		}
		min_cover_circle(c, r, n);
		printf("%.2f %.2f %.2f\n", c.x, c.y, r);
	}
	return 0;
}

參考:《離散點集最小包圍圓演算法分析與改進》