Vstříc řazení v lineárním čase

Obecné řazení v lineárním čase je meta, které není možné dosáhnout1. Ve výsledné složitosti se vždycky vyskytne log(n) faktor. Může být maskovaný za něco jiného, ale vždycky tam bude. Třeba takový LSD radix sort má složitost O(kn), kde k je právě log(n) v převleku.

Ale to, že nemůžeme vyhrát na teoretické frontě, není důležité. Záleží jenom na tom, jestli můžeme vyhrát prakticky. Algoritmus nemusí být O(n), stačí jen když se bude tvářit jako O(n) pro všechny vstupy, které mu můžu realisticky předhodit. Dokonce ani O(n) čas není důležitý, záleží jen na rychlosti.

O tom je tenhle článek: Nesnažím se vyhrát, jen prohrávat co nejméně.


Jeden způsob velice rychlého řazení je least significan digit (LSD) radix sort. Pro vstupní pole libovolné délek je cena na prvek v postatě plochá (pro 4B hodnoty se pohybuje na úrovni 15ns/prvek). Jeho nevýhoda spočívá v tom, že nejlépe funguje jen pro krátké bitové stringy konstantní délky jako jsou 4B inty nebo 8B longy. Pracuje tedy na omezeném rozsahu hodnot (232 nebo 264) namísto nekonečné množiny obecných řetězců. Pro ty jsou jiné algoritmy: most significat digit (MSD) radix sort nebo burst sort. Ty můžou být rychlejší než klasické quick/merge/heap sorty, ale pořád jsou pomalejší než velice kompaktní LSD radix sort.

V ideálním světě by mělo být možné využít LSD radix sortu jako dílčí komponenty k řazení dlouhých stringů. Dlouhou dobu jsem přemýšlel jak to udělat co nejlépe, protože nejjednodušší řešení – řazení LSD radix sortem podle prefixů stringů – nefunguje zas tak dobře.

V pseudokódu, který se podobá Scale, to může vypadat nějak takhle:

val strings: Array[String] = ???

val packed = new Array[Long](strings.length)

for (i <- 0 until strings.length) {
  val s = strings(i)

  // zabalit 4B prefix a 4B index do jednoho 8B longu
  packed(i) = prefixOfString(s).toLong << 32 | i
}

// samotné řazení
val sorted = RadixSort(packed, sortBytes = 4 until 8)

sorted
  .groupBy(p => p >> 32) // rozlámat výstup na bloky se stejným prefixem
  .map { case (prefix, ps) =>
    ps.map(p => p & 0xffffffff) // extrahovat index
      .map(i => strings(i))     // převést na vstupní stringy
      .sorted                   // výsledné skupiny seřadit jiným algoritmem
  }.flatten

Problém tohoto přístupu je malá informační hodnota prefixů. Pokud bych takto řadil Javovské stringy obsahující především znaky [a-zA-Z0–9], bylo by to neslavné, protože Java string je v jádru pole dvoubajtových charů.

String "k47" je tvořen bajty 0 6B 0 34 0 37. Liché bajty nepřináší žádnou informaci a ty sudé v tomto případě obsahují jednu z 62 možných hodnot. Jeden char tedy má jen něco málo pod 6 bitů informace.

Ve skutečnosti je to ještě horší. Některé znaky jsou častější než jiné a tedy přináší ještě méně informace. 4 bajty prefixu můžou mít hodně pod 12 bitů informace, které můžou v nejlepším případě diskriminovat mezi 4096 různými prefixy. To není nic moc.5

Radixsortování podle prefixu nefunguje tak dobře jak by mělo navzdory tomu, že 4B int dokáže diskriminovat mezi více jak 4 miliardami hodnot. Tohle mi nějakou dobu leželo na mysli. Potřeboval bych funkci, která mapuje vstupní string na 4B int, který obsahuje víc informace než 4B prefix a tedy lépe rozdělí vstupní data.

Kdybych měl funkci, která přetransformuje string na jeho pořadové číslo, mám vyhráno. Kdyby ta funkce běžela v O(1), měl bych řazení v O(n). To ale není možné, protože taková funkce by záležela na všech n vstupních hodnotách.

Místo toho mi stačí funkce, která komprimuje prefix, tedy využije maximum možných informace z delšího prefixu a zkomprimuje je do 4B nebo 8B čísla, které pak použiji k řazení rychlým LSD radix sortem. Je potřeba jen aby funkce nezměnila vzájemné pořadí prvků. Může namapovat několik různých na jeden prefix, to je povolené a nevyhnutelné, ale nesmí změnit jejich vzájemné pořadí. Tedy například:

f(s) = Σ 1 / (s_i-96)^i

f("")       = 0.0
f("a")      = 0.0384
f("aa")     = 0.0399
f("aaaa")   = 0.03999
f("b")      = 0.0769
f("ggggga") = 0.2796
f("gggggz") = 0.2796 // kolize
f("z")      = 0.961

a = b  => f(a) = f(b)
a < b  => f(a) <= f(b)

Tuhle funkci můžu zkonstruovat ručně, pokud něco vím o řazených datech2 nebo ji můžu vytvořit statisticky ze vzorku vstupních stringů (předpokládám, že jednotlivé bajty jsou na sobě nezávislé).6 Ty nejčastější poskytují minimum informace a potřebuji tedy víc dalších bajtů, abych dostal dost informace pro přesné mapování. Ty nejméně časté poskytují víc informace a samy o sobě vedou k relativně přesné diskriminaci.7

Spočítám relativní frekvence jednotlivých bajtů, ty využiji, abych převedl každý string na pozici na intervalu 0.0 – 1.0 a tu pak převedu na plný interval longu` (jen se bojím, že se v tomto kroku můžou vyskytnout chyby zaokrouhlování, které poruší požadavek na zachování pořadí).

Ukázka mapování se čtyřmi možnými hodnotami – a je velice časté, b a c jsou méně časté a d je raritní:

Když mám tuhle funkci, aplikuji ji na každý řazený prvek (v ukázce pseudokódu to nahradí volání prefixOfString), sestavím pole dvojic (prefix, index) zabalené do longů, seřadím radix sortem a pak rozsekám remízy jiným řadícím algoritmem3 .

A to je všechno. Nic víc v tom není. Používám starý známý radix sort, změním jen data na kterých ho pouštím.

Dál se dá aplikovat jen pár triků. Například můžu využít co nejvíce bitů v longu do kterého je zakódován pár (prefix, index). Index má omezené rozmezí a na jeho reprezentaci je třeba méně místa než plných 32 bitů. Například pro 1 milion indexů stačí 20 bitů.

Místo

|*******|*******|*******|*******|*******|*******|*******|*******|
[ prefix 32 bitů                |                 index 32 bitů ]

pak můžu použít

|*******|*******|*******|*******|*******|*******|*******|*******|
[ prefix 44 bitů                            |     index 20 bitů ]

Ve stejných 8 bajtech dostanu delší prefix, který bude lépe diskriminovat. Radix sort sice udělá víc průchodů (v mém případě řadí v jenom průchodu jeden bajt, ale je možné napsat verzi, která řadí po 10 bitech), ale to není problém, protože RS je velice rychlý a v porovnání s ostatními fázemi běží v zanedbatelném ča­se.

Ale pořád je třeba alokovat 24 bajtů na každý řazený prvek. To se nedá nijak obejít a je s tím třeba počítat. Samotná rychlost alokace ale nepředstavuje problém.

Výsledky jsou uspokojivé. f+sort je první jednodušší verze, která tohle používá jen jako preprocessing. Vezme celý výsledek LSD radix sortu tak jak je a seřadí ho. Javovský TimSort doslova letí na už téměř seřazených vstupech. f+s*rt je to, co popisuji výše.

Řazení slov z českých textů:

Rychlost na prvek byla téměř shodná pro 900000 slov i 4 miliony.

Řazení náhodně generovaných stringů [a-zA-Z0–9]+:

V tomto případě mapování funguje fantasticky a téměř vůbec není třeba dělat druhotné řazení. To je způsobeno tím, že stringy jsou dost dlouhé, náhodně rozložené a téměř bez duplicit.

Úzká hrdla výkonu překvapivě představuje procházení prvků při komprimaci prefixu, přerovnání pointerů podle výsledného pořadí radix sortu a inspekce samotných prvků při řazení (V ukázce kódu výše to odpovídá funkci prefixOfString, řádku „převést na vstupní stringy“ a pak závěrečné metodě sorted). Každý tento přesun způsobí jeden nevyhnutelný cache-miss, který na mém stroji v nejlepším případě trvá 16ns. Tedy něco kolem 32–48 nanosekund je promrháno na dvou nebo třech mov instrukcích4, zbytek trvá komprimace prefixů, příprava dat pro radix sort, samotný radix sort a závěrečné řazení.

Prostoru pro zlepšení už moc není a tohle může být velice blízko maximální rychlosti, které je možné vůbec dosáhnout. A to podle mě docela dobře odpovídá původnímu požadavku prohrávat co nejméně.


Dále k tématu:


Pozn:

  1. Řadit v O(n) se dá jen ve speciálních případech, kdy vstupní data mají omezený rozsah hodnot (Pigeonhole sort, Counting sort) nebo známé rozložení (Bucket sort).
  2. Můžu třeba stringy převést na číselníky, kdy číselná reprezentace pořadím odpovídá pořadí stringů jako v některých sloupcových databázích.
  3. Tento případ lehce komplikuje fakt, že nevím kolik bajtů prefixu jsem použil a kolik tedy můžu bezpečně přeskočit, kdybych chtěl druhé kolo řazení provést algoritmem založeném na radixech. Můžu určit jenom dolní mez.
  4. Tady by se hodilo mít nějaký Céčkovitější jazyk nebo jazyk, který má structy. Pak bych mohl index nahradit přímo pointerem a ušetřit si jeden cache-miss.
  5. Ale i když to nebude fungovat ideálně, LSD RS nezabere tolik času a udělá aspoň nějakou diskriminaci, která aspoň trochu urychlí následující řazení.
  6. Logický další krok je použít strojové učení pro konstrukci kompresní funkce.
  7. Pro stringy se vyplatí zjišťovat frekvence sudých a lichých bajtů zvlášť, protože mají velice rozdílné frekvence.

Flattr this!

This entry was posted in Algo, řazení. Bookmark the permalink.

2 Responses to Vstříc řazení v lineárním čase

  1. Filip Jirsák says:

    Nebylo by lepší pro tu kompresní funkci použít znalost toho, že se (například) jedná jen o číslice a písmena anglické abecedy a použít jen daný počet bitů? Protože při tom mapování podle četností se v té kompresní funkci zase musí vyhledávat. Jasně, pro každý řetězec se to vyhledávání provede jen jednou, ale stejně…

    • To je jedna z možností, pokud si můžu být zcela jistý, co je ve vstupu za data. Ale když jde o nějaký neznámý vstup jistoty jsou pryč a kompresní funkci musím vytvořit nějak mechanicky + to nebere v potaz relativní frekvence jednotlivých znaků, které se můžou dramaticky lišit.

      Nic se vyhledávat nemusí. Když si předem spočítám tabulky relativní frekvence a kumulativní součet, funkce může vypadat nějak takhle.

      var res = 0.0
      var m = 1.0
      var i = 0
      while (i < str.length) {
        val char = str.charAt(i)
        val b0 = (char >> 8) & 0xff
        val b1 = char & 0xff
      
        res += cumsum0(b0) * m
        m *= relfreq0(b0)
        res += cumsum1(b1) * m
        m *= relfreq1(b1)
      
        i += 1
      }

      Nejhorší, co dělám je pár floating-point násobení.

Leave a Reply

Your email address will not be published. Required fields are marked *