수직 분할¶
개요¶
수직 분할은 다양한 기하 문제에서 쓰이는 강력한 기법입니다. 기본 아이디어는 몇 가지 "좋은" 성질을 가진 여러 개의 수직 띠로 평면을 자르고, 각 띠에 대해 독립적으로 문제를 푸는 것입니다. 몇 가지 예시를 통해 이 아이디어를 설명하겠습니다.
삼각형 합집합의 넓이¶
평면 위에 $n$개의 삼각형이 있고, 이들의 합집합의 넓이를 구해야 한다고 가정합시다. 삼각형들이 서로 교차하지 않는다면 문제는 쉬울 것이므로, 모든 정점과 서로 다른 삼각형의 변들이 만나는 모든 교점을 지나는 수직선을 그어 평면을 수직 띠로 나눔으로써 이러한 교차를 없애 봅시다. 이러한 직선은 최대 $O(n^2)$개 있을 수 있으므로 $O(n^2)$개의 띠를 얻게 됩니다. 이제 어떤 수직 띠 하나를 생각해 봅시다. 수직이 아닌 각 선분은 이 띠를 왼쪽에서 오른쪽으로 가로지르거나, 전혀 가로지르지 않습니다. 또한 띠의 내부에서 엄밀하게 교차하는 두 선분은 존재하지 않습니다. 이는 삼각형 합집합 중 이 띠 내부에 놓인 부분이, 밑변이 띠의 양쪽 경계 위에 놓인 서로소인 사다리꼴들로 이루어진다는 뜻입니다. 이 성질 덕분에 다음과 같은 스캔라인 알고리즘으로 각 띠 내부의 넓이를 계산할 수 있습니다. 띠를 가로지르는 각 선분은 해당 삼각형의 내부가 선분의 위쪽에 있는지 아래쪽에 있는지에 따라 위쪽 선분 또는 아래쪽 선분이 됩니다. 각 위쪽 선분을 여는 괄호로, 각 아래쪽 선분을 닫는 괄호로 생각하면, 괄호 수열을 더 작은 올바른 괄호 문자열로 분해하는 방식으로 띠를 사다리꼴들로 분해할 수 있습니다. 이 알고리즘은 $O(n^3\log n)$의 시간과 $O(n^2)$의 메모리를 필요로 합니다.
최적화 1¶
먼저 실행 시간을 $O(n^2\log n)$으로 줄이겠습니다. 각 띠마다 사다리꼴을 생성하는 대신, 삼각형의 한 변(선분 $s = (s_0, s_1)$)을 고정하고 이 선분이 어떤 사다리꼴의 변이 되는 띠들의 집합을 구해 봅시다. 이 경우 $s$ 아래쪽(아래쪽 선분인 경우에는 위쪽)의 괄호 균형이 0인 띠들만 찾으면 된다는 점에 유의하십시오. 즉, 각 띠마다 수직 스캔라인을 실행하는 대신, $s$에 대한 괄호 균형에 영향을 주는 다른 선분들의 모든 부분에 대해 수평 스캔라인을 실행할 수 있습니다. 간단히 하기 위해 위쪽 선분에 대해 이를 수행하는 방법을 보이겠으며, 아래쪽 선분에 대한 알고리즘도 비슷합니다. 수직이 아닌 다른 선분 $t = (t_0, t_1)$을 생각하고, $s$와 $t$를 $Ox$에 사영한 구간들의 교집합 $[x_1, x_2]$를 구합니다. 이 교집합이 비어 있거나 한 점으로만 이루어져 있다면, $s$와 $t$는 같은 띠의 내부를 지나지 않으므로 $t$를 버릴 수 있습니다. 그렇지 않다면 $s$와 $t$의 교집합 $I$를 생각해 봅시다. 세 가지 경우가 있습니다.
-
$I = \varnothing$
이 경우 $t$는 $[x_1, x_2]$에서 $s$의 위쪽 또는 아래쪽에 있습니다. $t$가 위쪽에 있다면, $s$가 어떤 사다리꼴의 변인지 여부에 영향을 주지 않습니다. $t$가 $s$의 아래쪽에 있다면, $t$가 위쪽 선분인지 아래쪽 선분인지에 따라 $[x_1, x_2]$에 있는 모든 띠의 괄호 수열 균형에 $1$ 또는 $-1$을 더해야 합니다.
-
$I$가 한 점 $p$로 이루어진 경우
이 경우는 $[x_1, x_2]$를 $[x_1, p_x]$와 $[p_x, x_2]$로 나누어 앞의 경우로 환원할 수 있습니다.
-
$I$는 어떤 선분 $l$입니다
이 경우는 $x\in[x_1, x_2]$ 범위에서 $s$와 $t$의 해당 부분이 일치함을 의미합니다. $t$가 더 아래에 있다면 $s$는 분명히 사다리꼴의 변이 아닙니다. 그렇지 않다면 $s$와 $t$ 모두 어떤 사다리꼴의 변으로 간주될 수 있는 상황이 생길 수 있습니다. 이 모호함을 해소하기 위해, 인덱스가 가장 작은 선분만 변으로 간주하기로 정할 수 있습니다(여기서는 삼각형의 변들에 어떤 방식으로든 번호가 매겨져 있다고 가정합니다). 따라서 $index(s) < index(t)$이면 이 경우를 무시하고, 그렇지 않으면 $s$가 $[x_1, x_2]$에서 결코 변이 될 수 없음을 표시해야 합니다(예를 들어, 균형값이 $-2$인 대응 이벤트를 추가하는 방식으로).
다음은 세 가지 경우를 그림으로 나타낸 것입니다.
마지막으로 $[x_1, x_2]$에 속하는 모든 띠에 $1$ 또는 $-1$을 더하는 연산들을 처리하는 방법을 짚고 넘어가야 합니다. $[x_1, x_2]$에 $w$를 더하는 각 연산마다 이벤트 $(x_1, w),\ (x_2, -w)$를 만들고, 이 모든 이벤트를 스위프 라인으로 처리할 수 있습니다.
최적화 2¶
앞의 최적화를 적용하면 더 이상 모든 띠를 명시적으로 찾을 필요가 없다는 점에 유의하십시오. 이로써 메모리 사용량이 $O(n)$으로 줄어듭니다.
볼록 다각형의 교집합¶
수직 분해의 또 다른 활용은 두 볼록 다각형의 교집합을 선형 시간에 계산하는 것입니다. 각 다각형의 각 정점을 지나는 수직선들에 의해 평면이 수직 띠들로 나뉜다고 가정합시다. 그러면 입력 다각형 중 하나와 어떤 띠를 생각할 때, 이들의 교집합은 사다리꼴, 삼각형 또는 점 중 하나입니다. 따라서 각 수직 띠마다 이 도형들의 교집합을 구하고, 이 교집합들을 하나의 다각형으로 합치기만 하면 됩니다.
구현¶
아래는 삼각형 집합의 합집합 넓이를 $O(n^2\log n)$ 시간과 $O(n)$ 메모리로 계산하는 코드입니다.
typedef double dbl;
const dbl eps = 1e-9;
inline bool eq(dbl x, dbl y){
return fabs(x - y) < eps;
}
inline bool lt(dbl x, dbl y){
return x < y - eps;
}
inline bool gt(dbl x, dbl y){
return x > y + eps;
}
inline bool le(dbl x, dbl y){
return x < y + eps;
}
inline bool ge(dbl x, dbl y){
return x > y - eps;
}
struct pt{
dbl x, y;
inline pt operator - (const pt & p)const{
return pt{x - p.x, y - p.y};
}
inline pt operator + (const pt & p)const{
return pt{x + p.x, y + p.y};
}
inline pt operator * (dbl a)const{
return pt{x * a, y * a};
}
inline dbl cross(const pt & p)const{
return x * p.y - y * p.x;
}
inline dbl dot(const pt & p)const{
return x * p.x + y * p.y;
}
inline bool operator == (const pt & p)const{
return eq(x, p.x) && eq(y, p.y);
}
};
struct Line{
pt p[2];
Line(){}
Line(pt a, pt b):p{a, b}{}
pt vec()const{
return p[1] - p[0];
}
pt& operator [](size_t i){
return p[i];
}
};
inline bool lexComp(const pt & l, const pt & r){
if(fabs(l.x - r.x) > eps){
return l.x < r.x;
}
else return l.y < r.y;
}
vector<pt> interSegSeg(Line l1, Line l2){
if(eq(l1.vec().cross(l2.vec()), 0)){
if(!eq(l1.vec().cross(l2[0] - l1[0]), 0))
return {};
if(!lexComp(l1[0], l1[1]))
swap(l1[0], l1[1]);
if(!lexComp(l2[0], l2[1]))
swap(l2[0], l2[1]);
pt l = lexComp(l1[0], l2[0]) ? l2[0] : l1[0];
pt r = lexComp(l1[1], l2[1]) ? l1[1] : l2[1];
if(l == r)
return {l};
else return lexComp(l, r) ? vector<pt>{l, r} : vector<pt>();
}
else{
dbl s = (l2[0] - l1[0]).cross(l2.vec()) / l1.vec().cross(l2.vec());
pt inter = l1[0] + l1.vec() * s;
if(ge(s, 0) && le(s, 1) && le((l2[0] - inter).dot(l2[1] - inter), 0))
return {inter};
else
return {};
}
}
inline char get_segtype(Line segment, pt other_point){
if(eq(segment[0].x, segment[1].x))
return 0;
if(!lexComp(segment[0], segment[1]))
swap(segment[0], segment[1]);
return (segment[1] - segment[0]).cross(other_point - segment[0]) > 0 ? 1 : -1;
}
dbl union_area(vector<tuple<pt, pt, pt> > triangles){
vector<Line> segments(3 * triangles.size());
vector<char> segtype(segments.size());
for(size_t i = 0; i < triangles.size(); i++){
pt a, b, c;
tie(a, b, c) = triangles[i];
segments[3 * i] = lexComp(a, b) ? Line(a, b) : Line(b, a);
segtype[3 * i] = get_segtype(segments[3 * i], c);
segments[3 * i + 1] = lexComp(b, c) ? Line(b, c) : Line(c, b);
segtype[3 * i + 1] = get_segtype(segments[3 * i + 1], a);
segments[3 * i + 2] = lexComp(c, a) ? Line(c, a) : Line(a, c);
segtype[3 * i + 2] = get_segtype(segments[3 * i + 2], b);
}
vector<dbl> k(segments.size()), b(segments.size());
for(size_t i = 0; i < segments.size(); i++){
if(segtype[i]){
k[i] = (segments[i][1].y - segments[i][0].y) / (segments[i][1].x - segments[i][0].x);
b[i] = segments[i][0].y - k[i] * segments[i][0].x;
}
}
dbl ans = 0;
for(size_t i = 0; i < segments.size(); i++){
if(!segtype[i])
continue;
dbl l = segments[i][0].x, r = segments[i][1].x;
vector<pair<dbl, int> > evts;
for(size_t j = 0; j < segments.size(); j++){
if(!segtype[j] || i == j)
continue;
dbl l1 = segments[j][0].x, r1 = segments[j][1].x;
if(ge(l1, r) || ge(l, r1))
continue;
dbl common_l = max(l, l1), common_r = min(r, r1);
auto pts = interSegSeg(segments[i], segments[j]);
if(pts.empty()){
dbl yl1 = k[j] * common_l + b[j];
dbl yl = k[i] * common_l + b[i];
if(lt(yl1, yl) == (segtype[i] == 1)){
int evt_type = -segtype[i] * segtype[j];
evts.emplace_back(common_l, evt_type);
evts.emplace_back(common_r, -evt_type);
}
}
else if(pts.size() == 1u){
dbl yl = k[i] * common_l + b[i], yl1 = k[j] * common_l + b[j];
int evt_type = -segtype[i] * segtype[j];
if(lt(yl1, yl) == (segtype[i] == 1)){
evts.emplace_back(common_l, evt_type);
evts.emplace_back(pts[0].x, -evt_type);
}
yl = k[i] * common_r + b[i], yl1 = k[j] * common_r + b[j];
if(lt(yl1, yl) == (segtype[i] == 1)){
evts.emplace_back(pts[0].x, evt_type);
evts.emplace_back(common_r, -evt_type);
}
}
else{
if(segtype[j] != segtype[i] || j > i){
evts.emplace_back(common_l, -2);
evts.emplace_back(common_r, 2);
}
}
}
evts.emplace_back(l, 0);
sort(evts.begin(), evts.end());
size_t j = 0;
int balance = 0;
while(j < evts.size()){
size_t ptr = j;
while(ptr < evts.size() && eq(evts[j].first, evts[ptr].first)){
balance += evts[ptr].second;
++ptr;
}
if(!balance && !eq(evts[j].first, r)){
dbl next_x = ptr == evts.size() ? r : evts[ptr].first;
ans -= segtype[i] * (k[i] * (next_x + evts[j].first) + 2 * b[i]) * (next_x - evts[j].first);
}
j = ptr;
}
}
return ans/2;
}