Skip to Content
CTRL K

Marshland Rescues

Explicación

Cada región anegada se da como un polígono convexo. MAPS quiere saber hasta dónde puede tener que vadear un rescatista dentro de una región así.

Formalmente, para un polígono convexo dado, debemos hallar un punto en su interior que esté lo más lejos posible de la frontera, y reportar esa distancia máxima. Esto equivale a hallar el radio del círculo más grande que se puede inscribir por completo dentro del polígono.

Así, la tarea se reduce a calcular la máxima distancia mínima posible de un punto interior del polígono a todos sus lados.

Idea de semiplanos

Para cualquier distancia fija d, queremos comprobar si existe un punto dentro del polígono que esté al menos a distancia d de cada lado.

Para cada lado del polígono convexo, este conjunto de puntos a distancia al menos d del lado se puede representar como un semiplano cuya recta frontera está a distancia d del lado (hacia el interior).

Si existe una intersección común de todos estos planos, entonces podemos estar seguros de que hay un punto que está al menos a distancia d de TODOS los lados. (Es más, todos los puntos de esta intersección lo cumplen.)

Solución

Para cualquier d sabemos comprobar si existe una solución. Ahora podemos hacer búsqueda binaria sobre todos los valores posibles de d.

Implementación

Complejidad temporal: O(N∗log⁡N∗log⁡(R/ε))\mathcal{O}(N * \log N * \log(R / \varepsilon))

donde

  • NN es el número de vértices,
  • RR es el rango de búsqueda
  • ε\varepsilon es la precisión requerida (aquí, ε=10−12\varepsilon = 10^{-12}).
 
 
#include <bits/stdc++.h>
using namespace std;
 
using ld = long double;
const ld EPS = 1e-9, INF = 1e18;
 
vector<Point> poly;
int n;
 
// BeginCodeSnip{Point Template}
struct Point {
	ld x, y;
	Point(ld x = 0, ld y = 0) : x(x), y(y) {}
	Point operator+(const Point &o) const { return {x + o.x, y + o.y}; }
	Point operator-(const Point &o) const { return {x - o.x, y - o.y}; }
	Point operator*(ld k) const { return {x * k, y * k}; }
	ld dot(const Point &o) const { return x * o.x + y * o.y; }
	ld cross(const Point &o) const { return x * o.y - y * o.x; }
	ld norm() const { return sqrtl(x * x + y * y); }
	Point rot90() const { return {-y, x}; }
};
// EndCodeSnip
 
// BeginCodeSnip{Half Plane Template}
struct HalfPlane {
	Point p, dir;
	ld ang;
	HalfPlane() {}
	HalfPlane(Point a, Point b) : p(a), dir(b - a) { ang = atan2l(dir.y, dir.x); }
	bool outside(const Point &r) const { return dir.cross(r - p) < -EPS; }
};
 
Point intersect(const HalfPlane &a, const HalfPlane &b) {
	ld t = (b.p - a.p).cross(b.dir) / a.dir.cross(b.dir);
	return a.p + a.dir * t;
}
 
vector<Point> halfPlaneIntersection(vector<HalfPlane> &h) {
	vector<Point> box = {{INF, INF}, {-INF, INF}, {-INF, -INF}, {INF, -INF}};
	for (int i = 0; i < 4; i++) h.emplace_back(box[i], box[(i + 1) % 4]);
 
	sort(h.begin(), h.end(), [](auto &a, auto &b) { return a.ang < b.ang; });
	deque<HalfPlane> dq;
 
	for (auto &hp : h) {
		while (dq.size() > 1 && hp.outside(intersect(dq.back(), dq[dq.size() - 2])))
			dq.pop_back();
		while (dq.size() > 1 && hp.outside(intersect(dq[0], dq[1]))) dq.pop_front();
 
		if (!dq.empty() && fabsl(hp.dir.cross(dq.back().dir)) < EPS) {
			if (hp.dir.dot(dq.back().dir) < 0) return {};
 
			if (hp.outside(dq.back().p)) dq.pop_back();
			else continue;
		}
		dq.push_back(hp);
	}
 
	while (dq.size() > 2 && dq[0].outside(intersect(dq.back(), dq[dq.size() - 2])))
		dq.pop_back();
	while (dq.size() > 2 && dq.back().outside(intersect(dq[0], dq[1]))) dq.pop_front();
	if (dq.size() < 3) return {};
 
	vector<Point> poly(dq.size());
	for (int i = 0; i + 1 < (int)dq.size(); i++) poly[i] = intersect(dq[i], dq[i + 1]);
	poly.back() = intersect(dq.back(), dq[0]);
	return poly;
}
// EndCodeSnip
 
// BeginCodeSnip{Check function for Binary Seach}
bool ok(ld d) {
	vector<HalfPlane> h;
	for (int i = 0; i < n; i++) {
		Point a = poly[i], b = poly[(i + 1) % n];
		Point nrm = (b - a).rot90();
		nrm = nrm * (d / nrm.norm());
		h.emplace_back(a + nrm, b + nrm);
	}
	return !halfPlaneIntersection(h).empty();
}
// EndCodeSnip
 
int main() {
	ios::sync_with_stdio(false);
	cin.tie(nullptr);
 
	cin >> n;
	poly.resize(n);
	for (auto &p : poly) cin >> p.x >> p.y;
 
	ld lo = 0, hi = 1e5;
	for (int i = 0; i < 300; i++) {
		ld mid = (lo + hi) / 2;
		if (ok(mid)) lo = mid;
		else hi = mid;
	}
 
	cout << fixed << setprecision(12) << lo << '\n';
}
Fundamentos