Босоо задаргаа¶
Тойм¶
Босоо задаргаа бол геометрийн янз бүрийн бодлогод хэрэглэгддэг хүчирхэг арга юм. Ерөнхий санаа нь хавтгайг ямар нэг "сайн" шинж чанартай хэд хэдэн босоо судал болгон зүсэж, эдгээр судлын хувьд бодлогыг тусад нь бодох явдал юм. Бид санааг хэдэн жишээн дээр тайлбарлана.
Гурвалжнуудын нэгдлийн талбай¶
Хавтгай дээр $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$-ээс доош (доод хэрчмийн тохиолдолд дээш) хаалтын баланс тэг байх судлуудыг олох ёстойг анхаарна уу. Энэ нь судал бүрийн хувьд босоо шүүрдэлт ажиллуулахын оронд бид $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)$ болгон бууруулна.
Гүдгэр олон өнцөгтүүдийн огтлолцол¶
Босоо задаргааны өөр нэг хэрэглээ бол хоёр гүдгэр олон өнцөгтийн огтлолцлыг шугаман хугацаанд тооцоолох явдал юм. Хавтгайг олон өнцөгт бүрийн орой бүрээр дайрсан босоо шулуунуудаар босоо судал болгон хуваасан гэж үзье. Тэгвэл хэрэв бид оролтын олон өнцөгтүүдийн нэгийг ба ямар нэг судлыг авч үзвэл тэдгээрийн огтлолцол нь трапец, гурвалжин эсвэл цэг байна. Тиймээс бид эдгээр дүрсийг босоо судал бүрийн хувьд огтлолцуулаад эдгээр огтлолцлыг нэг олон өнцөгт болгон нийлүүлж болно.
Implementation¶
Below is the code that calculates area of the union of a set of triangles in $O(n^2\log n)$ time and $O(n)$ memory.
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;
}