2018wf Panda Preserve

计算几何好题

2018wf Panda Preserve

题意

给一个多边形在每个顶点画半径为r的圆，求覆盖整个多边形的最小半径。

分析

一开始看到很容易去想到二分+圆并，蛋感觉不好写，然后灵机一动把问题转化为：记多边形内的点，这个点到所有的顶点中距离的最小值为w§,求对于所有p的max(w§)。那么很容易想到，这个多边形可以被划分为若干个区域，每个区域中的任意一个点p到一个顶点q的距离为 w§。那么模拟一下题目发现，如果某个点到的最近点是顶点e，那么这个点肯定在e到与所有点的垂直平分线划分出来的靠e这侧的半平面的半平面交中。

那么对于每个顶点$E$去算其对于其他顶点的垂直平分线的半平面交，这个复杂度是 $O(n^2log (n)) $，算出半平面交之求解这个既半平面交也在多边形里的点到这个顶点$E$的最大距离。随便画个图可以发现其可能为半平面交的顶点和半平面交的边与多边形的边的交点到这个顶点$E$的距离。

学术上，这些垂直平分线构成的图叫Voronoi图，这个图的边数量是$O(n)$的所以枚举v图和多边形边的交点是$O(n^2)$的，所以总的时间复杂度是$O(n^2log(n))$的。

notice：判断点在多边形内的时候如果采用从p出发的一条直线经过几条边的话一定要注意这条“直线” 要足够长，本人因为这个wa了半年。

#include <vector>
#include <algorithm>
#include <iostream>
#include <array>
#include <numeric>
#include <cmath>
#include <functional>
#include <queue>
#include <complex>
#include <iomanip>
#include <cassert>

using ll = long long;
using T = double;
constexpr double eps = 1e-5;

using Point = std::complex<T>;
#define x real
#define y imag

T dot(const Point& a, const Point& b) {
	return (std::conj(a) * b).x();
}

T cross(const Point& a, const Point& b) {
	return (std::conj(a) * b).y();
}

struct Line {
	Point a;
	Point b;

	Line(const Point& a, const Point& b) : a(a), b(b) {}
};

Point rotate(const Point& a) {
	return Point(-a.y(), a.x());
}


Point Norm(const Point& a) {
	auto tmp = a / abs(a);
	tmp *= 10000;
	return tmp;
}

int sgn(const Point& a) {
	return a.y() > 0 || (a.y() == 0 && a.x() > 0) ? 1 : -1;
}

bool onLeft(const Point& a, const Line& l) {
	return cross(l.b - l.a, a - l.a) > 0;
}

Point intersection(const Line& l1, const Line& l2) {
	return l1.a + (l1.b - l1.a) * (cross(l2.b - l2.a, l1.a - l2.a) / cross(l2.b - l2.a, l1.a - l1.b));
}

using polygen = std::vector<Point>;

polygen hp(std::vector<Line>& lines) {
	std::sort(lines.begin(), lines.end(), [&](auto l1, auto l2) {
		auto d1 = l1.b - l1.a;
		auto d2 = l2.b - l2.a;

		if (sgn(d1) != sgn(d2)) {
			return sgn(d1) == 1;
		}

		return cross(d1, d2) > 0;
	});

	std::deque<Line> ls;
	std::deque<Point> ps;
	for (auto l: lines) {
		if (ls.empty()) {
			ls.push_back(l);
			continue;
		}

		while (!ps.empty() && !onLeft(ps.back(), l)) {
			ps.pop_back();
			ls.pop_back();
		}

		while (!ps.empty() && !onLeft(ps[0], l)) {
			ps.pop_front();
			ls.pop_front();
		}

		if (cross(l.b - l.a, ls.back().b - ls.back().a) == 0) {
			if (dot(l.b - l.a, ls.back().b - ls.back().a) > 0) {
				if (!onLeft(ls.back().a, l)) {
					assert(ls.size() == 1);
					ls[0] = l;
				}
				continue;
			}
			return {};
		}

		ps.push_back(intersection(ls.back(), l));
		ls.push_back(l);
	}

	while (!ps.empty() && !onLeft(ps.back(), ls[0])) {
		ps.pop_back();
		ls.pop_back();
	}
	if (ls.size() <= 2) {
		return {};
	}
	ps.push_back(intersection(ls[0], ls.back()));

	return {ps.begin(), ps.end()};
}


template<typename T>
int sign(T rr) {
	return std::abs(rr) < eps ? 0 : rr < 0 ? -1 : 1;
}

bool crs(const Line& a, const Line& b) {
	if (sign(cross(a.a - b.a, b.b - b.a)) * sign(cross(a.b - b.a, b.b - b.a)) >= 0)return false;
	if (sign(cross(b.a - a.a, a.b - a.a)) * sign(cross(b.b - a.a, a.b - a.a)) >= 0)return false;
	return true;
}

int main() {
	std::ios::sync_with_stdio(false);
	std::cin.tie(nullptr);
	std::cout.tie(nullptr);

	int n;
	//	std::cin >> n;
	scanf("%d", &n);
	std::vector<Point> point(n);
	std::vector<Line> lines_max;
	lines_max.emplace_back(Point{1e4 + 10, -1e4 - 10}, Point{1e4 + 10, 1e4 + 10});
	lines_max.emplace_back(Point{1e4 + 10, 1e4 + 10}, Point{-1e4 - 10, 1e4 + 10});
	lines_max.emplace_back(Point{-1e4 - 10, 1e4 + 10}, Point{-1e4 - 10, -1e4 - 10});
	lines_max.emplace_back(Point{-1e4 - 10, -1e4 - 10}, Point{1e4 + 10, -1e4 - 10});
	for (int i = 0; i < n; ++i) {
		T xx, yy;
		//		std::cin >> xx >> yy;
		scanf("%lf%lf", &xx, &yy);
		point[i] = {xx, yy};
	}
	auto on_seg = [&](const Point& P, const Point& a, const Point& b) {
		return sign(dot((P - a), (P - b))) <= 0 && std::abs(cross(P - a, P - b)) < eps;
	};
	auto on_seg_line = [&](const Line& l, const Point& a, const Point& b) -> int {
		return sign(cross(l.b - l.a, l.a - a)) * sign(cross(l.b - l.a, l.a - b)) < 0 && sign(cross(b - a, a - l.a)) * sign(cross(b - a, a - l.b)) < 0;
	};
	auto in = [&](const Point& p) -> bool {
		for (int i = 0; i < n; ++i) {
			if (on_seg(p, point[i], point[(i + 1) % n]))return true;
		}
		Point dir = {43233.123422346L, 89939.12323486L};
		int ret = 0;
		for (int i = 0; i < n; ++i) {
			ret += on_seg_line(Line(p, p + dir), point[i], point[(i + 1) % n]);
		}
		return ret & 1;
	};

	T ans = -1e18;
	for (int i = 0; i < n; ++i) {
		std::vector<Line> lines(lines_max);
		for (int j = 0; j < n; ++j) {
			if (i == j)continue;
			Point mid = (point[i] + point[j]);
			mid /= 2;
			Point dir = point[j] - point[i];
			dir = rotate(dir);
//			dir = Norm(dir);
			lines.emplace_back(mid - dir, mid + dir);
		}
		auto nearest = hp(lines);
		int m = nearest.size();
		for (int j = 0; j < m; ++j) {
			auto p = nearest[j];
			if (in(p)) {
				ans = std::max(ans, abs(p - point[i]));
			}
			Line tmp = Line(nearest[(j + 1) % m], nearest[j]);
			for (int k = 0; k < n; k++) {
				Line tnp = Line(point[(k + 1) % n], point[k]);
				if (crs(tmp, tnp)) ans = std::max(ans, abs(intersection(tmp, tnp) - point[i]));
			}
		}
	}
	printf("%.4lf", ans);
	return 0;
}