pku3525 凸多边形的最大内切圆
2 评论132次阅读2008.09.06 17:36; 作者：Felicia
这个题可以用二分圆半径＋半平面交做，但是复杂度是`O(nlogn * logr)`的，应付这题足够，但是对于 sgu332`(n=10000)`来说，就需要更好的算法了。
我能想到的做法是，不断缩小多边形，每次把所有边向内移动同样的举例，保证恰好缩掉一条边，这样缩了`n-3`次后，剩下一个三角形，这个三角形的内心就是所求的内切圆圆心。
直接模拟这个做法是`O(n^2)`的，但是我们可以通过维护一个堆和一个双向循环链表，达到`O(nlogn)`的复杂度。
下面是我的核心代码，包括求角平分线和维护一个堆和一个双向循环链表。其中 dtype 是浮点数类型。
heap h;
poly p;
int n;
int L[10010];
int R[10010];
pair <point, point> get_bisector(const point &a, const point &b, const point &c, const point &d) {
if (sgn(cross(b - a, d - c)) == 0) {
if (b == c) {
return make_pair(c, c + rotate_left(d - c, 0));
}
else {
point dd = c + (a - b);
return make_pair((b + c) / 2, (a + dd) / 2);
}
}
point o;
inter_line_line(a, b, c, d, o);
point dd = o + ((d - o) / len(d - o) * len(a - o));
point e = (a + dd) / 2;
return make_pair(o, e);
}
void renew(const int &t) {
point La = p[L[t]], Lap1 = p[(L[t] + 1) % n];
point Ra = p[R[t]], Rap1 = p[(R[t] + 1) % n];
point a = p[t], ap1 = p[(t + 1) % n];
pair <point, point> p = get_bisector(La, Lap1, a, ap1);
pair <point, point> q = get_bisector(a, ap1, Ra, Rap1);
point pp;
inter_line_line(p.first, p.second, q.first, q.second, pp);
dtype new_dist = sqrt(dist2_point_line(pp, a, ap1));
h.decrease(t, new_dist);
}
void init() {
h.init(n);
for (int i = 0; i < n; i++)
L[i] = (i - 1 + n) % n, R[i] = (i + 1) % n;
for (int i = 0; i < n; i++)
renew(i);
}
dtype calc() {
while (h.size() > 3) {
int t = h.get_top();
h.pop();
R[L[t]] = R[t];
L[R[t]] = L[t];
renew(L[t]);
renew(R[t]);
}
return h.get_top_data();
}
poly p;
int n;
int L[10010];
int R[10010];
pair <point, point> get_bisector(const point &a, const point &b, const point &c, const point &d) {
if (sgn(cross(b - a, d - c)) == 0) {
if (b == c) {
return make_pair(c, c + rotate_left(d - c, 0));
}
else {
point dd = c + (a - b);
return make_pair((b + c) / 2, (a + dd) / 2);
}
}
point o;
inter_line_line(a, b, c, d, o);
point dd = o + ((d - o) / len(d - o) * len(a - o));
point e = (a + dd) / 2;
return make_pair(o, e);
}
void renew(const int &t) {
point La = p[L[t]], Lap1 = p[(L[t] + 1) % n];
point Ra = p[R[t]], Rap1 = p[(R[t] + 1) % n];
point a = p[t], ap1 = p[(t + 1) % n];
pair <point, point> p = get_bisector(La, Lap1, a, ap1);
pair <point, point> q = get_bisector(a, ap1, Ra, Rap1);
point pp;
inter_line_line(p.first, p.second, q.first, q.second, pp);
dtype new_dist = sqrt(dist2_point_line(pp, a, ap1));
h.decrease(t, new_dist);
}
void init() {
h.init(n);
for (int i = 0; i < n; i++)
L[i] = (i - 1 + n) % n, R[i] = (i + 1) % n;
for (int i = 0; i < n; i++)
renew(i);
}
dtype calc() {
while (h.size() > 3) {
int t = h.get_top();
h.pop();
R[L[t]] = R[t];
L[R[t]] = L[t];
renew(L[t]);
renew(R[t]);
}
return h.get_top_data();
}
发表回复
zuimeng | 1F
四月 21st, 2011 at 23:00
领悟中、、、、