1. 程式人生 > 其它 >JZOJ 1082. 【GDOI2005】選址

JZOJ 1082. 【GDOI2005】選址

\(\text{Problem}\)

很久以前,在世界的某處有一個形狀為凸多邊形的小島,島上的居民們決定建一個祭壇,居民們任務祭壇的位置離島的頂點處越遠越好。
你的任務是求凸多邊形內一點,使其與各頂點的距離中最短的距離最遠,點在邊上也可以。 這樣的點可能有多個,你只需輸出這些點與各頂點的最短距離。

\(\text{Solution}\)

非常經典的題
以答案為半徑做圓,滿足圓心到所有點距離大於等於半徑
考慮二分半徑,判斷圓心存不存在
列舉任意兩點,考慮分別以此半徑做圓交出的點(選取內部的點,叉積判斷是否在內部)
如果這點沒有被所有圓覆蓋,即這點到凸邊形頂點的最短距離大於等於半徑,說明圓心存在
算交點用相似,如果沒有交點考慮其在邊上交出的點一樣判斷
精度很神奇,不要用 \(\text{long double}\)


本地過不了資料1 \(OJ\) 上卻過了?!

\(\text{Code}\)

#include <cstdio>
#include <algorithm>
#include <cmath>
#define RE register
#define IN inline
using namespace std;

const int N = 105;
const double eps = 1e-8;
double area;
int n;
struct Vector{
	double x, y;
	IN Vector(double xx = 0, double yy = 0){x = xx, y = yy;}
	IN Vector operator - (const Vector &B){return Vector(x - B.x, y - B.y);}
	IN double operator * (const Vector &B){return fabs(x * B.y - y * B.x);}
}p[N];

IN double sqr(double x){return x * x;}
IN double distance(Vector a, Vector b){return sqrt(sqr(a.x - b.x) + sqr(a.y - b.y));}
IN double Area(Vector a, Vector b, Vector c){return (a - b) * (c - b);}
IN double sum_Area(Vector a)
{
	double res = 0;
	for(RE int i = 1; i < n; i++) res += Area(a, p[i], p[i + 1]);
	return res + Area(a, p[n], p[1]);
}
IN int isIn(Vector a){if (fabs(sum_Area(a) - area) <= eps) return 1; return 0;}
IN int OK(Vector a, double r)
{
	double res = 1e18;
	for(RE int i = 1; i <= n; i++) res = min(res, distance(a, p[i]));
	return res >= r;
}

IN int check(double r)
{
	for(RE int i = 1; i <= n; i++)
		for(RE int j = 1; j < i; j++)
		{
			double d = distance(p[i], p[j]);
			if (d > r * 2)
			{
				double k = (p[i].y - p[j].y) / (p[i].x - p[j].x), x, y;
				y = min(p[i].y, p[j].y) + fabs(p[i].y - p[j].y) * r / d;
				if (k > 0)
				{
					x = min(p[i].x, p[j].x) + fabs(p[i].x - p[j].x) * r / d;
					if (OK(Vector{x, y}, r)) return 1;
					x = p[i].x + p[j].x - x, y = p[i].y + p[j].y - y;
					if (OK(Vector{x, y}, r)) return 1;
				}
				else{
					x = max(p[i].x, p[j].x) - fabs(p[i].x - p[j].x) * r / d;
					if (OK(Vector{x, y}, r)) return 1;
					x = p[i].x + p[j].x - x, y = p[i].y + p[j].y - y;
					if (OK(Vector{x, y}, r)) return 1;
				}
			}
			else{
				double w = sqrt(r * r - d * d / 4), k = w / d;
				double dx = (p[i].x + p[j].x) / 2, dy = (p[i].y + p[j].y) / 2;
				double x = dx - fabs(p[i].y - p[j].y) * k, y = dy + fabs(p[i].x - p[j].x) * k;
				if (isIn(Vector{x, y}) && OK(Vector{x, y}, r)) return 1;
				x = p[i].x + p[j].x - x, y = p[i].y + p[j].y - y;
				if (isIn(Vector{x, y}) && OK(Vector{x, y}, r)) return 1;
			}
		}
	return 0;
}

int main()
{
	scanf("%d", &n);
	for(RE int i = 1; i <= n; i++) scanf("%lf%lf", &p[i].x, &p[i].y);
	area = sum_Area(p[1]); double r = 0;
	for(RE int i = 1; i <= n; i++)
		for(RE int j = 1; j < i; j++) r = max(r, distance(p[i], p[j]));
	double l = 0, mid = (l + r) / 2, ans;
	for(RE int i = 0; i < 60; i++, mid = (l + r) / 2)
		if (check(mid)) ans = mid, l = mid; else r = mid;
	printf("%.3lf\n", ans);
}