Tenhle článek je pod licencí na způsob Public Domain uveřejněný na Addamově webu.

Hledání průsečíků mnohoúhelníka v rovině

Popis problému

Vstupními daty je posloupnost bodů v rovině, vrcholů mnohoúhelníka. Každé dva po sobě jdoucí body (bráno cyklicky, tedy včetně dvojice konec-začátek) určují úsečku. Mnohoúhelníkem je míněna uzavřená křivka, která vznikne navázáním těchto úseček na sebe. Mnohoúhelník protíná sám sebe, pokud určitým bodem roviny tato křivka prochází vícekrát. To odpovídá průsečíku dvou jeho úseček, vyjma každé dvojice sousedních úseček, které se smějí protínat v jednom koncovém bodě.

Výstupem je prostá odpověď, zda mnohoúhelník má, nebo nemá nějaký průsečík sám se sebou (TRUE/FALSE).

Datové struktury

Vstupem je pole objektů typu vektor. Ten obsahuje své souřadnice x, y v libovolném číselném formátu a dále metody:

Dále algoritmus pracuje s úsečkami. Ta má každá odkaz na své dva krajní body a na dvě jim příslušné sousední úsečky. Bod, který je více vlevo (tj. menší), označíme jako začátek, druhý bod jako konec. Navíc má metodu:

A konečně, algoritmus používá vyhledávací strom pro ukládání úseček protnutých jistou, řekněme aktuální svislou osou. Úsečky jsou v něm tedy všechny vzájemně porovnatelné. Tento strom má dvě metody:

Porovnání úseček

Každá úsečka určuje levou a pravou polorovinu, bráno při pohledu z jejího začátku do konce. Levá polorovina tedy obsahuje například body, které jsou přímo (po svislé ose) nad danou úsečkou. Přesněji vzato, levá polorovina úsečky AB obsahuje bod C právě tehdy, pokud úhel BAC je dutý - menší než přímý. Tento výrok je ekvivalentní s tím, že třetí složka vektorového součinu (AB) · (AC) je kladná. Pravá polorovina obdobně obsahuje právě body, pro které tato hodnota vychází záporná. Pro body na přímce AB vychází hodnota nulová.

Pokud oba krajní body porovnávané úsečky UV leží v téže, například pravé, polorovině úsečky AB, žádný bod úsečky UV nemůže být přímo nad úsečkou AB. Pokud jsou úsečky porovnatelné a nestejné, pak je v tomto případě AB nevyhnutelně výše než UV. Obdobnou otázku je potřeba vyřešit nejen pro případ, že U, V jsou v levé polorovině úsečky AB, ale také pro případy, že body A, B jsou v téže polorovině úsečky UV. Všechny tyto situace dávají jasný výsledek porovnání.

Pokud oba body U, V leží v různých polorovinách úsečky AB a zároveň A, B leží v různých polorovinách úsečky UV, průsečík přímek AB, UV leží i na obou úsečkách a tedy se jistě protínají.

V případě, kdy oba (a tedy všechny čtyři) součiny vyjdou nulové a úsečky jsou porovnatelné, leží na jedné přímce a mají tedy průsečík. Pokud jsou sousední, je potřeba ještě ověřit, zda se skutečně překrývají, k tomu ale stačí porovnání krajních bodů. Pokud sousední nejsou, průsečík mají jistě.

Pokud je nulový jen jeden z dvojice součinů, pro označení průsečíku je potřeba ověřit, že příslušný bod neleží jen na přímce, ale také na úsečce samotné. K tomu opět stačí porovnání s krajními body. Naopak, pokud příslušný bod leží mimo úsečku, poloha druhého bodu jednoznačně určuje výsledek porovnání.

Celá tahle záležitost zapsaná v Pythonu:

 def is_a_below_b(A: Segment, B: Segment) -> bool:
   assert (A.end >= B.start) and (B.end >= A.start)
   
   cross_b1 = (A.end - A.start) · (B.start - A.start)
   cross_b2 = (A.end - A.start) · (B.end - A.start)
   
   if (cross_b1 == 0 and B.start >= A.start) or (cross_b2 == 0 and B.end <= A.end):
     if (A.next is B) or (B.next is A):
       if (cross_a == cross_b == 0) and (A.end > B.start and B.end > A.start):
         raise Intersection()
       return False
     raise Intersection()
   if (cross_b1 >= 0) and (cross_b2 >= 0):
     return True
   if (cross_b1 <= 0) and (cross_b2 <= 0):
     return False
   
   cross_a1 = (B.end - B.start) · (A.start - B.start)
   cross_a2 = (B.end - B.start) · (A.end - B.start)
   if (cross_a1 >= 0) and (cross_a2 >= 0):
     return False
   if (cross_a1 <= 0) and (cross_a2 <= 0):
     return True
   
   raise Intersection()

Třída Intersection je výjimka, která je případně zpracována volajícím kódem.

Algoritmus

Samotný algoritmus je pak stručný a zakládá se na myšlené svislé ose (sweepline), kterou posouváme zleva doprava přes mnohoúhelník.

Nejprve mnohoúhelník převedeme na pole úseček. Dále vytvoříme dvě pole, obě obsahující dvojice vektor a odkaz na úsečku, kde vektor slouží jako klíčová hodnota pro třídění. První pole představuje události, kdy svislá osa protne začátek některé úsečky: obsahuje tedy všechny úsečky utříděné vzestupně podle začátků. Obdobně druhé představuje události, kdy osa protne konec úsečky a tedy obsahuje úsečky utříděné podle konců.

Události projdeme slévané zleva doprava. Tedy, jedním cyklem vytáhneme vždy další úsečku pro odebrání a dalším vnořeným cyklem přidáme všechny ještě nezpracované úsečky pro přidání, které ji svým začátkem (neostře) předcházejí. Potom odebereme aktuální úsečku. Pokud nastane výjimka, některé úsečky se protínají, program tedy může navrátit TRUE a ukončit se. Pokud zpracujeme všechny události a žádná výjimka nenastane, program nakonec navrátí FALSE.

Funkčnost

Samotné funkci is_a_below_b jsem se věnoval snad už dost podrobně: za předpokladu, že jsou obě vstupní úsečky porovnatelné, výjimku vyvolá právě pokud se protínají a nesousedí spolu, anebo sousedí a protínají se ve více bodech. Když program navrací TRUE, jistě tedy našel průsečík dvou úseček odpovídající zadání problému. A navíc, jakmile se dvojice protínajících se úseček objeví ve stromu těsně vedle sebe, program najde jejich průsečík a skončí. Zbývá tedy ověřit, že v každém průsečíku několika přímek se alespoň jedna dvojice dostane ve stromu vedle sebe.

Sweepline si můžeme představit jako svislou čáru nepatrně nakloněnou proti směru hodinových ručiček, protože události jsou řazeny lexikograficky. Všechny úsečky ve stromu jsou vzhledem k aktuální sweepline (procházející bodem aktuální události) porovnatelné - to by šlo formálně podepřít indukcí, snad není potřeba. Když je ve stromu dvojice protínajících se úseček a ještě se na to nepřišlo, musí to být tím, že jsou nějaké úsečky mezi nimi. Tyto další úsečky ale sweepline protínají a porovnatelné jsou, tedy pořadí úseček ve stromu vždy odpovídá tomu, v jakém svislém pořadí jsou jejich průsečíky se sweepline.

Průsečík dvojice úseček může mít speciální podobu, kdy se koncem a začátkem dotýkají dvě úsečky, které v mnohoúhelníku nebyly za sebou. Díky neostrému porovnávání ale je před odebráním levé z těchto úseček přidána pravá. Do stromu bude přidána přímo vedle levé a jejich průsečík tedy bude nalezen.

V ostatních případech, pokud se několik úseček protíná v jednom bodě, pak vždy je možné vést sweepline tak, aby mezi jejich průsečíky s ní nebyl průsečík s žádnou další úsečkou. Nejpozději po zpracování poslední události před určitým průsečíkem několika úseček musí být všechny tyto úsečky ve stromu přímo za sebou. Jedna dvojice z nich tedy bude porovnána a program se ukončí.

Se svislými úsečkami algoritmus pracuje bez potíží díky lexikografickému porovnávání bodů.

Složitost

Převést mnohoúhelník na pole úseček lze přímo, v lineárním čase. Utřídit n-prvkové pole podle vektorů lze vždy například quicksortem v čase O(n·log(n)). Přidání nebo odebrání prvku vhodnému vyhledávacímu stromu trvá O(log(n)) podle počtu prvků, které strom právě obsahuje, a v nejhorším případě může obsahovat skutečně všechny úsečky. Nalezení obou sousedů ve stromu může trvat také až log(n). Událostí je celkem 2·n, jejich zpracování tedy trvá O(n·log(n)). Celý algoritmus skončí v čase O(n·log(n)).

Paměťová složitost všech použitých struktur je přinejhorším lineární podle počtu úseček.

Vylepšení

Nenašel jsem žádný algoritmus řešící tenhle problém, který by byl v nejhorším případě asymptoticky rychlejší. Určitého zrychlení pro obvyklé mnohoúhelníky ale lze dosáhnout tím, že označíme úsečky, na které navazuje sousední úsečka směřující doprava. Jinými slovy, rozdělíme mnohoúhelník na řetězy z úseček začínající a končící v lokálních extrémech mnohoúhelníka. V případě dvou takhle navazujících úseček lze událost odebrání a přidání nahradit událostí nahrazení úsečky, kterou lze zpracovat v konstantním čase. Příslušný vrchol stromu, ve kterém je potřeba jednu úsečku nahradit druhou, i její oba sousedy ve stromu si můžeme ukládat jako další vlastnosti úsečky a nastavit je stačí při událostech přidávání a odebírání jí nebo jejích sousedů ve stromu. Potom metoda nahrazení prostě z levé úsečky zjistí oba sousedy, novou úsečku s nimi porovná kvůli možným průsečíkům a ve stromu nakonec nahradí jednu úsečku druhou.

Převést mnohoúhelník na takovéhle řetězy lze prostě v lineárním čase. V rámci každého řetězu jsou navíc už úsečky utříděné podle konců i začátků zároveň, sléváním tedy jdou všechny události utřídit v O(n·log(m)), kde m je počet řetězů. Každý řetěz je potřeba jen jednou přidat a jednou odebrat, ostatní události se stanou událostmi nahrazení, konkrétně jich bude n-m. Ve stromu bude od každého řetězu vždy nejvýše jedna úsečka, takže všechny události zpracujeme v čase O(n + m·log(m)). Celý algoritmus tedy proběhne v čase O(n·log(m)). Počet řetězů je číslo v rozsahu 2 až n, navíc je závislý na otočení mnohoúhelníka. Právě vhodným otočením jej lze ještě před samotným začátkem minimalizovat - také v čase O(n·log(m)) - ale to už nepovažuju za významný zisk, abych se mu podrobně věnoval.

Výsledky porovnávání úseček je asi rozumné cachovat, protože je to poměrně výpočetně náročná záležitost a pokud je nějaká dvojice úseček porovnávaná jednou při přidání, často je porovnávaná i podruhé při odebírání.

Související problémy

Časté rozšíření této úlohy je otázka, které všechny průsečíky mnohoúhelník obsahuje. To vyžaduje trochu pečlivého ošetření speciálních případů, ale podstata algoritmu zůstává stejná.

Původně jsem s mnohoúhelníkem zacházel jako s plochou ostře uvnitř křivky, takže průsečíky okraj na okraj by se nepočítaly. Řešení takového problému by ale vyžadovalo poměrně složité dodatečné zpracování výsledků (hlavně kvůli průsečíku více úseček v jednom bodě a překryvu souběžných úseček) a v takovém světle se mi to nezdá jako moc užitečná úloha.