Алгоритмы C++
.pdf
Нахождение вписанной окружности в выпуклом многоугольнике с помощью тернарного поиска
Дан выпуклый многоугольник с N вершинами. Требуется найти координаты центра и радиус наибольшей вписанной окружности.
Здесь описывается простой метод решения этой задачи с помощью двух тернарных поисков, работающий за O (N log2 C), где C - коэффициент, определяемый величиной координат и требуемой точностью (см. ниже).
Алгоритм
Определим функцию Radius (X, Y), возвращающую радиус вписанной в данный многоугольник окружности с центром в точке (X;Y). Предполагается, что точки X и Y лежат внутри (или на границе) многоугольника. Очевидно, эту
функцию легко реализовать с асимптотикой O (N) - просто проходим по всем сторонам многоугольника, считаем для каждой расстояние до центра (причём расстояние можно брать как от прямой до точки, не обязательно рассматривать как отрезок), и возвращаем минимум из найденных расстояний - очевидно, он и будет наибольшим радиусом.
Итак, нам нужно максимизировать эту функцию. Заметим, что, поскольку многоугольник выпуклый, то эта функция будет пригодна для тернарного поиска по обоим аргументам: при фиксированном X0 (разумеется, таком, что
прямая X=X0 пересекает многоугольник) функция Radius(X0, Y) как функция одного аргумента Y будет сначала возрастать, затем убывать (опять же, мы рассматриваем только такие Y, что точка (X0, Y) принадлежит
многоугольнику). Более того, функция max (по Y) { Radius (X, Y) } как функция одного аргумента X будет сначала возрастать, затем убывать. Эти свойства ясны из геометрических соображений.
Таким образом, нам нужно сделать два тернарных поиска: по X и внутри него по Y, максимизируя значение функции Radius. Единственный особый момент - нужно правильно выбирать границы тернарных поисков, поскольку вычисление функции Radius за пределами многоугольника будет некорректным. Для поиска по X никаких сложностей нет, просто выбираем абсциссу самой левой и самой правой точки. Для поиска по Y находим те отрезки многоугольника, в которые попадает текущий X, и находим ординаты точек этих отрезков при абсциссе X (вертикальные отрезки не рассматриваем).
Осталось оценить асимптотику. Пусть максимальное значение, которое могут принимать координаты - это C1,
а требуемая точность - порядка 10-C2, и пусть C = C1 + C2. Тогда количество шагов, которые должен будет совершить каждый тернарный поиск, есть величина O (log C), и итоговая асимптотика получается: O (N log2 C).
Реализация
Константа steps определяет количество шагов обоих тернарных поисков.
В реализации стоит отметить, что для каждой стороны сразу предпосчитываются коэффициенты в уравнении прямой,
исразу же нормализуются (делятся на sqrt(A2+B2)), чтобы избежать лишних операций внутри тернарного поиска.
const double EPS = 1E-9; int steps = 60;
struct pt {
double x, y;
};
struct line {
double a, b, c;
};
double dist (double x, double y, line & l) { return abs (x * l.a + y * l.b + l.c);
}
double radius (double x, double y, vector<line> & l) { int n = (int) l.size();
double res = INF;
for (int i=0; i<n; ++i)
res = min (res, dist (x, y, l[i])); return res;
}
double y_radius (double x, vector<pt> & a, vector<line> & l) { int n = (int) a.size();
double ly = INF, ry = -INF; for (int i=0; i<n; ++i) {
int x1 = a[i].x, x2 = a[(i+1)%n].x, y1 = a[i].y, y2 = a
[(i+1)%n].y; |
== x2) |
continue; |
|
if (x1 |
|||
if (x1 |
> x2) |
swap |
(x1, x2), swap (y1, y2); |
if (x1 |
<= x+EPS && |
x-EPS <= x2) { |
|
double y |
= y1 |
+ |
(x - x1) * (y2 - y1) / (x2 - x1); |
ly = min |
(ly, |
y); |
|
ry = max |
(ry, |
y); |
|
} |
|
|
|
} |
++sy) { |
|
|
for (int sy=0; sy<steps; |
/ 3; |
||
double diff = (ry - ly) |
|||
double y1 = ly + |
diff, |
|
y2 = ry - diff; |
double f1 = radius (x, y1, l), f2 = radius (x, y2, l); if (f1 < f2)
ly = y1;
else
ry = y2;
}
return radius (x, ly, l);
}
int main() {
int n; vector<pt> a (n);
... чтение a ...
vector<line> l (n);
for (int i=0; i<n; ++i) {
l[i].a = a[i].y - a[(i+1)%n].y; l[i].b = a[(i+1)%n].x - a[i].x;
double sq = sqrt (l[i].a*l[i].a + l[i].b*l[i].b); l[i].a /= sq, l[i].b /= sq;
l[i].c = - (l[i].a * a[i].x + l[i].b * a[i].y);
}
double lx = INF, rx = -INF; for (int i=0; i<n; ++i) {
lx = min (lx, a[i].x); rx = max (rx, a[i].x);
}
for (int sx=0; sx<stepsx; ++sx) { double diff = (rx - lx) / 3;
double x1 = lx + diff, x2 = rx - diff;
double f1 = y_radius (x1, a, l), f2 = y_radius (x2, a, l); if (f1 < f2)
lx = x1;
else
rx = x2;
}
double ans = y_radius (lx, a, l); printf ("%.7lf", ans);
}
Нахождение вписанной окружности в выпуклом многоугольнике методом "сжатия сторон"
Дан выпуклый многоугольник с N вершинами. Требуется найти координаты центра и радиус наибольшей вписанной окружности.
В отличие от описанного здесь метода тернарного поиска, при данном методе решения время работы - O (N log N) - не зависит от ограничений на координаты и от точности, и поэтому этот метод проходит при значительно больших N.
Спасибо mf за описание этого красивого алгоритма.
Алгоритм
Итак, дан выпуклый многоугольник. Начнём одновременно и с одинаковой скоростью сдвигать все его стороны параллельно самим себе внутрь многоугольника:
Пусть, для удобства, это движение происходит со скоростью 1 координатная единица в секунду (т.е. время в данном случае - это расстояние от сторон до их новых положений).
Ясно, что в процессе этого движения стороны многоугольника будут постепенно исчезать (обращаться в точки).
const double EPS = 1E-9; const double INF = 1E+40;
struct pt { |
x, |
y; |
double |
||
pt() |
{ } |
|
pt (double x, double y) : x(x), y(y) { } |
||
}; |
|
|
double get_ang |
(pt & a, pt & b) { |
|
double |
ang1 = atan2 (a.y, a.x); |
|
double |
ang2 = atan2 (b.y, b.x); |
|
double |
ang = abs (ang1 - ang2); |
|
return |
min (ang, 2*M_PI-ang); |
|
} |
|
|
pt vec (pt & a, pt & b) { |
||
return |
pt |
(b.x-a.x, b.y-a.y); |
} |
|
|
double dist (pt & |
a, pt & b) { |
|
return |
sqrt ((a.x-b.x)*(a.x-b.x) + (a.y-b.y)*(a.y-b.y)); |
|
} |
|
|
double get_h (double a, double alpha, double beta) { |
||
return |
a * sin(alpha) * sin(beta) / sin(alpha+beta); |
|
} |
|
|
double det (double a, double b, double c, double d) { |
||
return |
a * d - b * c; |
|
} |
|
|
pt intersect_line |
(pt & p1, pt & p2, pt & q1, pt & q2) { |
|
double |
a1 |
= p1.y - p2.y; |
double |
b1 |
= p2.x - p1.x; |
double |
c1 |
= - a1 * p1.x - b1 * p1.y; |
double |
a2 |
= q1.y - q2.y; |
double |
b2 |
= q2.x - q1.x; |
double c2 = - a2 * q1.x - b2 * q1.y; return pt (
-det (c1, b1, c2, b2) / det (a1, b1, a2, b2),
-det (a1, c1, a2, c2) / det (a1, b1, a2, b2)
);
}
bool parallel (pt & p1, pt |
& p2, |
pt & q1, pt & q2) { |
|
double a1 = p1.y - |
p2.y; |
|
|
double b1 |
= p2.x - |
p1.x; |
|
double a2 |
= q1.y - |
q2.y; |
|
double b2 |
= q2.x - |
q1.x; |
a2, b2)) < EPS; |
return abs (det (a1, b1, |
|||
} |
|
|
|
double calc_val (pt & p1, pt & p2, pt & q1, pt & q2, pt & r1, pt & r2) { pt l1 = intersect_line (p1, p2, q1, q2);
pt l2 = intersect_line (p1, p2, r1, r2);
return get_h (dist (l1, l2), get_ang(vec(q1,q2),vec(p1,p2))/2, get_ang(vec(r1,r2),vec(p2,p1))/2);
}
int main() {
int n; vector<pt> a (n);
... чтение n и a ...
set < pair<double,int> > q; vector<double> val (n); for (int i=0; i<n; ++i) {
pt & p1 = a[i], & p2 = a[(i+1)%n], & q1 = a[(i-1+n)%n], &
q2 = a[(i+2)%n];
val[i] = calc_val (p1, p2, p1, q1, p2, q2); q.insert (make_pair (val[i], i));
}
vector<int> next (n), prev (n);
for (int i=0; i<n; ++i) { next[i] = (i + 1) % n;
prev[i] = (i - 1 + n) % n;
}
double last_time; while (q.size() > 2) {
last_time = q.begin()->first; int id = q.begin()->second; q.erase (q.begin());
val[id] = -1;
next[prev[id]] = next[id]; prev[next[id]] = prev[id];
int nxt = next[id], prv = prev[id];
if (parallel (a[nxt], a[(nxt+1)%n], a[prv], a[(prv+1)%n])) break;
q.erase (make_pair (val[nxt], nxt)); q.erase (make_pair (val[prv], prv));
val[nxt] = calc_val (a[nxt], a[(nxt+1)%n], a[(prv+1)%n], a[prv], a[next[nxt]], a[(next[nxt]+1)%n]);
val[prv] = calc_val (a[prv], a[(prv+1)%n], a[(prev[prv]+1)%n], a[prev[prv]], a[nxt], a[(nxt+1)%n]);
q.insert (make_pair (val[nxt], nxt)); q.insert (make_pair (val[prv], prv));
}
printf ("%.9lf", last_time);
}
Диаграмма Вороного в 2D
Определение
Даны точек |
на плоскости. Рассмотрим разбиение плоскости на областей |
(называемых многоугольниками Вороного или ячейками Вороного, иногда — многоугольниками близости, ячейками Дирихле, разбиением Тиссена), где — множество всех точек плоскости, которые находятся ближе к точке
, чем ко всем остальным точкам 
:
Само разбиение плоскости называется диаграммой Вороного данного набора точек 
. Здесь — заданная метрика, обычно это стандартная Евклидова
метрика: |
, однако ниже будет рассмотрен и случай так |
называемой манхэттенской метрики. Здесь и далее, если не оговорено иного, будет рассматриваться случай Евклидовой метрики
Ячейки Вороного представляют собой выпуклые многоугольники, некоторые являются бесконечными.
Точки, принадлежащие согласно определению сразу нескольким ячейкам Вороного, обычно так и относят сразу к нескольким ячейкам (в случае Евклидовой метрики множество таких точек имеет меру нуль; в случае манхэттенской метрики всё несколько сложнее).
Такие многоугольники впервые были глубоко изучены русским математиком Вороным (1868-1908 гг.).
Свойства
● |
Диаграмма Вороного является планарным графом, поэтому она имеет |
вершин и рёбер. |
||
● |
Зафиксируем любое |
. Тогда для каждого |
проведём прямую — |
|
|
серединный перпендикуляр отрезка |
; рассмотрим ту полуплоскость, образуемую этой прямой, в которой |
||
|
лежит точка . Тогда пересечение |
всех полуплоскостей для каждого |
даст ячейку Вороного . |
|
● Каждая вершина диаграммы Вороного является центром окружности, проведённой через какие-либо три точки