Mergeselect

Quickselect všichni známe a milujeme – jde o algoritmus pro nalezení k-tého nejmenšího prvku neseřazeného pole v průměrném čase O(n). Jediný jeho nedostatek je, že v patologických případech může běžet v čase O(n2). Quicksort to má ale stejně a je to součást dohody, kterou jsme odsouhlasili, když jsme pustili do quick-věcí. Pro ty, kteří chtějí garance, existuje o něco rafinovanější algoritmus medián mediánů, který má garantovaný O(n) čas a nepříliš příznivé konstantní faktory.

To všechno víme už 55 let.

Když jsem před nějakou dobou připravoval přednášku o řazení, trávil jsem víc času, než by bylo záhodno, dumáním nad paralelami, podobnostmi a dualismy. Nedá moc práce a člověk si začne všímat analogií mezi algoritmy, často se lišící jen převrácenými šipkami. Bádání mě o pár týdnů později dovedlo k načrtnutí níže uvedeného diagramu2. Každá jeho osa představuje nějaký princip, způsob fungování nebo atribut. V místě, kde se osy protínají, existuje algoritmus, nemusím ho znát, ale vždycky tam je.

Quicksort a merge sort jsou dva nečekaní bratři, ležící na opačných koncích jedné osy.

A tady se dostávám k titulku článku. Pokud symetrie platí a existuje zrcadlový vztah mezi quicksortem a merge sortem, musí také existovat merge obdoba quickselectu – tedy hypotetický mergeselect, který stejně jako jeho quick-bratr běží v lineárním čase.

Nějakou dobu jsem nad touhle otázkou – zajímavou, i když z praktického hlediska nevýznamnou – přemýšlel, hnán touhou po uspořádaném vesmíru, který dává smysl a ale nemohl jsem přijít na řešení. Ale pak jednoho osamělého odpoledne jsem zažil heuréka moment. O(n) mergeselect funguje, pokud se dopustím všech hříchů asymptotické analýzy, budu předstírat, že konstanty nehrají žádnou roli1, a v jednom místě taktně přejdu nástrahy pravděpodobnosti.

Rekurzivní algoritmus mergeselectu pracuje ve čtyřech krocích znázorněných na následujícím schématu.

V prvním kroku rozdělím vstup na bloky konstantní délky k > 1 (např. 1024), kterých je tedy n/k.

Ve druhém začnu s merge sortem, nedotáhnu ho však do konce, ale sortuji jen uvnitř bloků. Výsledkem jsou seřazené bloky délky k, kdy každý z nich je zpracován v čase k * log(k). Složitost tohoto kroku je tedy: (n / k) * (k * log(k)) = n * log(k), ale protože k je konstanta, odpovídá to O(n).

Třetí krok je místo, kde se děje všechna magie. Chci najít m-tý nejmenší element, jinými slovy potřebuji najít prvek, kterému předchází m-1 menších bratří. Začnu tím, že si zvolím prostřední element x0 prvního seřazeného bloku a binárním hledání ve všech ostatních blocích najdu stejný nebo nejbližší větší prvek. V průběhu prohledávání sčítám pozice těchto prvků. Součet s0 udává kolik prvků je menších než mnou zvolený „pivot“ x0. Tady nastane rekurze třetího kroku. Pokud je s0 větší než m, zvolím x1 menší než x0, pokud je s0 větší než m, zvolím naopak. Provádím tedy další binární hledání.3

Složitost třetího kroku je (n / k) * log(k) * log(k) = (n / k) * log2(k). Stále platí, že k je konstanta, a výraz tedy odpovídá O(n).

Pokud platí, že sx = m, našel jsem m-tý nejmenší prvek. To však nebude příliš pravděpodobné, protože se hledaný element bude s největší pravděpodobností nacházet v ostatních blocích.

Ve čtvrtém kroku se nacházím v situaci, kdy jsem ve všech blocích binárním hledáním určil horní a dolní mez, kde by se m-tý element měl nacházet a tyto dílčí subsekvence jsou seřazené. Co udělám? Začnu je mergovat do bloků maximální délky k a celý algoritmus opakuji, nehledám však m-tý nejmenší prvek, ale ten na pozici (m-s), kde s je součet dolních mezí. Pokud je výsledkem jen jeden blok, můžu v něm přímo najít kýžený prvek.

Aby výsledná rekurze měla lineární čas, je třeba splnit jedinou podmínku – množství práce v další iteraci se musí zmenšit o konstantní násobek f, tedy ni+1 = f * ni. f nemusí být nijak veliké, stačí když je o něco málo větší než 1 a výsledný čas se poskládá do O(n).

Jaká je velké f v případě mergeselectu? Tady přijde na řadu zmiňovaný úkrok stranou kolem pravděpodobnosti. Sub-sekvence vymezená dolní a horní mezí v prvním bloku má délku 0 nebo 1 (buď obsahuje kýžený prvek nebo ne). Když budu uvažovat uniformní rozdělení dat, velikosti sub-sekvencí ve všech ostatních blocích by měly být podobně velké jako to v prvním bloku.4 Mělo by tedy dojít k redukci celkové velikosti na něco jako n/k. Pokud tohle platí, pak je očekávaná složitost celé rekurze mergeselectu na úrovni O(n), nejhorší případ má složitost O(n2) (stejně jako quickselect) a to je celkem fajn.

Dodatky:

  • Tenhle algoritmus jsem nikde neviděl. To ale moc neznamená, protože jsem jednak diletant, který nerozumí asymptotické analýze, a druhak jsem nehledal nijak důkladně. Kdybych si myslel, že jsem objevil něco zcela nového, tak by to nepříjemně zavánělo naivitou a arogancí. Už párkrát jsem byl přesvědčený, že jsem na stopě něčeho inovativního, jen abych o týden později našel v paperu z roku 1975 zmínku o té mojí novince. Navíc nešlo o stěžejní myšlenku publikace, ale jen drobnou noticku ve stylu: „Jo a tohle je tady věc, o které se ví.“
  • Mergeselect (stejně jako standardní merge sort) není in-place. Potřebuje jednu horní a jednu dolní mez pro každý blok, kterých je n/k. To odpovídá O(n) extra prostoru, což je víc než O(log(n)), který je považován za in-place. Parametr k může však být relativně velký a proto můžou být praktické nároky snesitelné.
  • Algoritmus nenaprogramoval, protože nechci ztratit pár dnů debugováním krajních případů.
  • Merge sort je možné nahradit Timsortem nebo dokonce zcela jiným O(n log(n)) řadícím algoritmem, takže vztah mergeselectu k merge sortu je spíš jen duchovní.
  • Po vzoru quickselectu se dá zkonstruovat radixselect, který místo podle pivotu dělí dle radixu.

Pseudokód:

def mergeselect[T](xs: Seq[T], m: Int): T = {
  require(m < xs.length)

  // #1
  val k = 1024
  val blocks = xs.grouped(k)

  // #2
  for (block <- blocks) {
    mergesort(block)
  }

  if (blocks.length == 1) {
    return blocks(m)
  }

  // #3
  val lowerBounds = Array.fill(blocks.length)(0)
  val upperBounds = Array.tabulate(blocks.length) { b => b.length }
  val positions = Array.fill(blocks.length)(0)

  // binary search for s == m
  var s = 0
  while (lowerBounds(0) <= upperBounds(0)) {
    val p = (lowerBounds(0) + upperBounds(0)) >>> 1 // mid pisition
    val x = blocks(0)(p) // mid value

    s = p

    for (b <- 1 until blocks.length) {
      val block = blocks(b)
      val pos = binarySearch(for = x, in = block, from = lowerBounds(b), to = upperBounds(b))
      positions(b) = pos
      s += pos
    }

    if (s < m) {
      lowerBounds(0) = p + 1
    } else if (s > m) {
      upperBounds(0) = p - 1
    } else {
      return x
    }

    if (s < m) {
      for (b <- 1 until blocks.length) {
        lowerBounds(b) = positions(b) + 1
      }
    } else {
      for (b <- 1 until blocks.length) {
        upperBounds(b) = positions(b) - 1
      }
    }
  }

  // #4
  val newXs: Seq[T] = 0 until blocks.length map { b => blocks(b).slice(from = lowerBounds(b), to = upperBounds(b)) }.flatten
  mergeselect(newXs, m-s)
}

Dále k tématu:

Poznámky:

  1. Překvapivě hrají. Quicksort byl publikován roku 1961 a největších pokroků bylo dosaženo zlepšováním právě konstantních faktorů. Původní Hoarova metoda zvolit jako pivot první (nebo náhodnou) položku má průměrný konstantní faktor 1.386, metoda nejlepší ze tří 1.188 a kvazi-nejlepší-z-devíti 1.094. Další pokrok byl three way quicksort, který je odolný proti patologickým vstupům.
  2. S největší pravděpodobnostní na něm chybí jedna až tři osy. Tři hlavní osy na diagramu jsou merge/diskrimi­nace, porovnání/radix, dělení na dva/široké. Další osy by měly být ještě efemérní proces/dynamická datová struktura, sekvenční/para­lelizovatelné a možná i líné/striktní. Ale když začnu zacházet to takových detailů, hranice se rozmazávají a víc než technický obor to začne připomínat teologické argumenty.
  3. Protože k je přijatelně velká konstanta je možné binární hledání provádět velice efektivně, jak je popsáno v Binary Search eliminates Branch Mispredictions.
  4. A to, podle několika simulací, které jsem udělal, vypadá, že platí.

Flattr this!

Posted in Uncategorized | 5 Comments

Persistentní datové struktury

Funkcionální programování je horké téma posledních let. Co to je, proč je to dobré, co to přináší a proč o tom všichni mluví právě teď? Odpovědi na tyto otázky jsou následující: Programováníreferenčně transparentními funkcemi6, zlepšuje to možnosti kompozice, porozumění programům, což v ideálním případě vede k lepším programům s menším počtem chyb.

Jedním z hlavních pilířů funkcionálního programování je fakt, že všechny datové struktury jsou neměnné. Když něco jednou vytvořím, už to nikdy nemůžu modifikovat. Jediný způsob, jak něco změnit, je vytvořit novou verzi struktury, která v sobě požadovanou změnu zahrnuje .

„Ale to bude příliš drahá operace…“ může vás těď napadat.

Není to tak docela pravda. Je fakt, že kdyby bylo nutné při každé změně udělat kompletní kopii celé datové struktury, funkcionální programování by, i přes všechny teoretické i praktické výhody, nestálo za tu námahu. Program by vždycky běžel příliš pomalu. Naštěstí odvážní lidé z FP gangů za poslední polovinu století přišli na mnoho způsobů, jak implementovat datové typy, které jsou jednak neměnné, okážou rychle vytvářet změněné verze, a všechny operace nad nimi mají asymptotickou složitost, která je velice blízká efemérním variantám.

Hlavní myšlenka těchto takzvaných persistentních datových struktur je jednoduchá: Sdílet co největí množství struktury s předchozí verzí a kopírovat jen nejmenší možnou část obsahující změnu (toto je také označováno jako structural sharing). Persistentní v tomto kontextu neznamená uložený na disku, ale trvalý, protože všechny verze dané struktury přežijí.

Jako příklad uvedu strom A, jehož poslední list změním z 4 na αω a výsledkem bude strom B.

Jak je vidět oba stromy sdílí většinu dat (modré uzly), zkopírována byla jen cesta od změněného listu ke kořeni stromu (takzvaná path copy8). Operace tedy nemá cenu O(n), ale jen O(log(n)) a to je něco s čím se už dá reálně pracovat. Vyvážený strom s 232 elementy (±4 miliardy) má hloubku 31 a je tedy potřeba zkopírovat jen 31 uzlů, což sice není konstantní množství práce, ale na druhou stranu ani není nijak gigantické.

Princip je to krásně jednoduchý, jak na něm však postavit efektivní persistentní datovou strukturu, která má ty správné asymptoty, je otázka na dlouhé zimní večery (například amortizace je v persitentním prostředí poněkud zrádná). Ale když se to povede, otevírá se před námi země nečekaných možností: S persistentními strukturami mám k dispozici všechny verze – update nikdy nezničí starý stav, vytvoří nový a já můžu bez starosti přistupovat k oběma, sdílet je, porovnávat je a číst a odvozovat z nich nové verze vrásek i ve vícevláknových programech. Protože data jsou po vytvoření zafixována a zmrazena v čase, nemusím se ničeho bát, nemusím zamykat, koordinovat nebo mezi sebou synchronizovat vlákna. Neměnnou datovou strukturu je vždycky bezpečné číst, nemusím se bát, že zatímco jsem přečetl jednu část, někdo jiný mi pod rukama změnil zbytek struktury a já teď nemám konzistentní stav. Neměnná struktura se chová jako hodnota, jako stabilní fakt, jako věc, o které se můžeme bavit a nehrozí nám, že na ní každý má jiný pohled.

Že to je všeobecně dobrý nápad ukazuje třeba kniha Effective Java1, kde Josh Bloch doporučuje i v Javě používat neměnné třídy, jak to je jen možné, právě proto, že to zjednoduší uvažování o kódu a odstraní jednu celou kategorii starostí. Když to není možné, tak radí přistoupit k defenzivnímu kopírování (které však není dokonalým řešením, protože jednak mi během kopírování někdo jiný může zdrojový objekt změnit pod rukama2 a druhak přináší určitou (a často zbytečkou) režii).

Intuitivně je jasné, že vytváření modifikovaných verzí persistentních DS je dražší než in-place modifikace běžných DS.5 Při každé změně musím něco nového alokovat, někdy jeden objekt (jako v případě Cons-Nil seznamů) jindy log<sub>b</sub>(n) objektů (jako v případě stromů, pro nějaký základ logaritmu b). I když je alokace většinou levná, bouře alokací se přesto může podepsat na výkonu a v nejhorším případě může způsobit, že výsledek je limitován propustností paměti. Existují určité techniky jak se alokací zbavit, nebo je aspoň částečně omezit, jako například escape analysis nebo deforestation, ale ani ty nejsou dokonalé.3

Ve funkcionálním prostředí se navíc prakticky nedá obejít bez garbage collectoru, který sleduje všechny reference a jakmile nějaká verze persistentní DS není používána, odstraní data, která jí patří a nejsou sdílena s ostatními živými verzemi. Díky GC je implementace persistentních DS velice jednoduchá.

Skoro4 všechny funkcionální struktury jsou implementované pomocí stromů právě proto, že stromy umožňují snadné a levné kopírování cest. Nejde jen o zjevné stromy jako Red-Black a ALV, ale také o rafinované vnitřní detaily komplikovaných struktur, které se chovají jako mapy, množiny, pole, fronty nebo oboustranné fronty.

Nechci tady však vysvětlovat, jak jsou různé persistentní datové struktury implementované, jak vypadají a jaké mají charakteristiky. Jde o příliš rozsáhlé téma, o němž se dá najít spousta dobrých zdrojů, které všechno vysvětlí lépe, než bych já kdy dovedl. Místo toho sem jen vysypu seznam odkazů s drobnými anotacemi a půjdu od toho.


Skvělá přednáška, ve které Daniel Spiewak vysvětlí implementaci několika základních persistetních DS (seznamy, fronty, deque, red-black stromy, vektory).

Chris Okasaki na prostoru dvou set stran této legendární knihy vysvětlí mnohem víc, než jste kdy chtěli vědět o funkcionálních datových strukturách – jak je napsat v čistě funkcionálním jazyku (autor používá ML, ale v dodatcích jsou ukázky přepsané do Haskellu) a jak analyzovat jejich asymptotickou složitost. Každá kapitola vždy přinese něco nového a nečekaného, co člověka naprosto ohromí (tedy aspoň mě to ohromilo): Persistence, líné vyhodnocování, amortizace, amortizace v persistetním prostředí, eliminace amortizace, lazy rebuilding, datové struktury postavené na různých numerických reprezentacích, bootstraping a další.

Na internetu se dá najít Okasakiho dizertační práce, z níž posléze malými úpravami vznikla samotná kniha. Autor deset let od prvního vydání napsal malou rekapitulaci.

Tato přednáška z kurzu Advanced data structures vyučovaného na MIT vysvětlí, že s persistentními datovými strukturami je to o něco komplikovanější, protože ty funkcionální nejsou jediné, které jsou. Kromě toho ještě existují partially (částečně), fully (plně) a confluent (souběžně?) persistentní datové struktury, které připouštějí možnost modifikovat existující data, ale stále takovým způsobem, že jsou všechny předchozí verze zachovány a můžu k nim přistupovat a číst je. Částečně persistentní DS dovolují modifikovat jen nejaktuálnější ver­zi.

Tyto dva články vysvětlí jak interně funguje Vector – stěžejní lineární datovou strukturu z Clojure a ze Scaly. Vektor je velice široký strom, založený na práci Phila Bagwella7 (Fast And Space Efficient Trie Searches a Ideal Hash Trees), do jeho současné podoby ho převedl Rich Hickey, ale poprvé ho implementoval Daniel Spiewak (pokud se nepletu).9

Bagwellův paper popisující Relaxed Radix Balanced Trees – vylepšenou verzi vektoru, která dovoluje spojit dva RRB-vektory v logaritmickém čase (namísto lineárního, jak je běžné pro prostý Vector). RRB tak lookup, append, update pracující v efektivně konstantním čase, rozšíří o další efektivní operaci, bohužel za cenu drobného (konstantního) zpomalení všech ostatních operací. Lidé ze Scala komunity plánovali, že by zachovali současný vektor a RRB stromy by se používaly jen, když dojde ke spojení dvou vektorů nebo RRB stromů.

Následuje několik přednášek od Edwarda Kmetta. Jde o poněkud hutnější materiál, ale pokud jste se prokousali skrz Okasakiho, nemělo by být nemožné se zakousnout i do tohoto.


Neměnnost dat po jejich vytvoření značně zjednodušuje vícevláknové programování a synchronizaci (protože žádní není potřeba) a proto se těchto struktur využívá i v oblastech, kde by to člověk nečekal: V datech na disku.

CouchDB – používá sdílení struktury indexů na disku (data kopíruje celá). Když je nějaký dokument změněn, je na disk zapsána nová cesta ke kořeni dokumentu, která odkazuje na nezměněné části. Stará verze je pak nějakou dobu stále dostupná. Pokud se nemýlím, nejde však o hlavní myšlenku této databáze, ale jen o techniku, jak zrychlit zápisy, umožnit dočasné konzistentní okno a zlepšit paralelismus. Staré verze struktur jsou eventuálně smazány garbage collectorem.

Fulltext search engine Lucene, který tvoří základ pro Elastic a Solr během indexace také používá immutable přístup. Když zaindexuje nějaké dokumenty, Lucene vytvoří segment indexu, který je po vytvoření neměnný a ostatní procesy ho můžou použít k hledání. Na pozadí jsou tyto segmenty logaritmicky slučovány do větších segmentů (na principu Log-structured merge-tree, jako to třeba dělá LevelDB viz: Dominic Tarr: The database of the future: leveldb, Inside HyperLevelDB a Concurrency Improvements in HyperLevelDB), a i když nejsou persistentní, jsou stále neměnné. Jediné, co se mění, je seznam určující, které segmenty jsou dostupné a aktuální, a které už nikdo nepoužívá a budou smazány.

Analytická databáze Druid.io také vytváří neměnné (ale ne persistentní) segmenty, které poté distribuuje na dotazovací servery. O Druid.io dopodrobna píše Lukáš Havrlant: Druid.io: distribuovaná immutable databáze pro analytické výpočty, Architektura Druid.io a Jak Druid.io agreguje data.

Databáze Datomic je na myšlence neměnnosti a persistence postavená. To ji spolu s modelem subjekt-predikát-objekt-(transakce) ji umožňuje ukládat všechny verze všech faktů a efektivně tak pracovat s časem. Neměnnost vede k distribuovaným uzlům, které data jen čtou, a striktní požadavky na ACID pak mají za následek, že všechny transakce putují přes jeden uzel, tzv. trasactor. Více informací viz:

Bw-stromy, popsané v The Bw-Tree: A B-tree for New Hardware Platforms, také používají částečně neměnná data. V tomto případě proto, že je to ideální způsob, jak dosáhnout dobrého výkonu na mnohajádrových procesorech. Bw-strom je variantou B-stromu, která změny zaznamenává jako seznam delt. Když delta seznam dosáhne určité délky, dojde k vytvoření aktualizovaného uzlu Bw-stromu, který v sobě obsahuje všechny změny. Jediné, co se (atomicky) změní je pointer, který začne ukazovat na nový začátek delta seznamu. Protože delta nepřepisuje sdílená data uzlu, nedochází k invalidaci cache ostatních jader a eliminuje se tak množství zpráv cache coherency protokolu, což omezí nákladnou implicitní komunikaci mezi jádry procesoru a vede v lepšímu výkonu. K podobným závěrům došli i autoři Asynchronized Concurrency: The Secret to Scaling Concurrent Search Data Structures a Staring into the Abyss: An Evaluation of Concurrency Control with One Thousand Cores.

Persistentní datová struktura také leží v srdci rychlé a paralelní fronty SnapQueue, která dokáže rychle dělat snapshoty v čase a spojování dvou front. Má to však jeden háček – persistentní DS (v tomto případě Conc-Tree) není použita přímo, ale slouží jako kostra pro pole malých jednorázových front, které jsou pochopitelně měnitelné. SnapQueue přímo manipuluje pouze s první a poslední frontou (z první bere, do poslední ukládá, když je první prázdná, vytáhne si další, když je poslední plná, vloží ji do kostry a vytvoří další). Více detailů v SnapQueue: Lock-Free Queue with Constant Time Snapshots a Conc-Trees for Functional and Parallel Programming.

Persistentní struktury jsou také použity v Gitu. Btrfs a další systémy souborů, které interně používají Copy on Write B-stromy, udržují dočasnou formu persistence. Pokud jsou založené na append-only logu, musí řešit garbage collection a kompakci dat tím, že přesunou data z několika téměř prázdných stránek do jedné plné. Jde o nechvalně známý problém, protože vede k tzv. write amplification a může způsobit, že disk nebude stíhat uvolňovat stránky. Více v Log-structured file systems: There's one in every SSD a Stratified B-trees and versioning dictionaries.


Dále k tématu:


Pozn:

  1. Možná se za nějakou dobu dočkáme nové edice Effective Java aktualizované pro Javu 8 viz
  2. Technicky vzato na JVM není bez synchronizace bezpečné kopírovat ani long nebo double. Podle specifikace 64bitové primitivní typy nejsou atomické a zatímco jedno vlákno přečte prvních 32 bitů, nějaké jiné vlákno může změnit zbývajících 32 bitů. Prakticky vzato však všechna současná JVM implementují primitivní typy jako atomické. To se však může změnit v některé z následujících verzí, když přibudou value types.
  3. Jedna možná optimalizace, jak zrychlit program, je omezit alokace a kopírování cest prostřednictvím vkusného použití mutací. Pokud jsou všechny mutace stavu lokální nějaké funkci a nikdo jiný je nemůže pozorovat, vůbec nevadí, jestli byly provedeny jako několik čistě funkcionálních kroků nebo došlo k pochybným mutacím. Nějaká funkce zkonstruuje objekt in-place mutacemi, uvede ho do správného stavu,, finalizuje a zveřejní jako od té doby neměnný objekt. V Clojure tomuto mechanismu říkají transients.
  4. A ten zbytek je založený na spojových seznamech, které však nejsou nic víc než zdegenerované stromy.
  5. Nicméně čtení může být efektivnější. Jednak se nemusím strachovat se zamykáním a koordinací vláken při čtení, ale také sdílení je bezproblémové. Protože struktura je fixní, každý ji může sdílet a všechny procesory mají kopii ve své CPU cache a vesele si ji čtou. Kdyby docházelo ke změnám, procesory by mezi sebou museli posílat zprávy cache cohorenece protokolu a invalidovat změněné cache line v ostatních jádrech. – Navíc, když bych používal ke sdílení jednu atomickou referenci, která by ukazovala na neměnnou strukturu, mohlo by docházet ke CAS hell, kdy se mnoho vláken najednou snaží aktualizovat strukturu, ale vyhraje jen jedno, což vede k mnoha zprávám cache coherency protokolu nebo zbytečné práci. Řešením je stejně jako u sdílené paměti minimalizovat komunikaci mezi procesy a zbavit se míst, kde je běh programu serializován.
  6. Equational reasoning s referenčně transparentním programem není jen zbytečné akademické cvičení, ale skutečně užitečná věc, která se mi hodila, když jsem upravoval svojí knihovnu Matcher napsanou v PHP. Navzdory jazyku je Matcher čistě funkcionální a proto je možné s programem jednoduše manipulovat pomocí dosazovaní a eliminace. Právě takto – dosazováním za argumenty a vyhodnocováním – jsem upravil několik funkcí, které byly příliš složité a nepřehledné.
  7. Rich Hickey to sám zmiňuje v Clojure Concurrency
  8. Twigg et. al. v Stratified B-trees and versioning dictionaries zmiňují, že kopírování cest v roce 1986 Driscoll et al. (Making Data Structures Persistent) navrhl jako obecný způsob, jak z libovolné datové struktury založené na pointerech udělat persistentní verzi. Nešlo však o funkcionální persistenci, ale navrhovaný přístup povoluje modifikace již vytvořených objektů. Kopírování cest pro dosažení funkcionální persistence bylo použito už v roce 1982. Planar Point Location Using lan Munro Editor Persistent Search Trees obsahuje odkazy na několik zdrojů jejichž autoři více méně nezávisle došli ke stejným závěrům: Myers – AVL dags, Myers – Efficient applicative data types, Krijnen, Meertens – Making B-trees work for B, Reps, Teitelbaum, Demers – Incremental context dependent analysis for language-based editors a Swart – Efficient algorithms for computing geometric intersections.
  9. I když v Optimizing Hash-Array Mapped Tries for Fast and Lean Immutable JVM Collections se píše: „The first persistent HAMT implementation can be attributed to Rich Hickey, lead-developer of Clojure.“

Flattr this!

Posted in DS, Funkcionální programování | 2 Comments

Úvod do podivností moderního hardwaru, které vás budou budit ze spaní

Nezasvěceným se může chování hardware a procesorů zvlášť zdát zcela bizarní a nepochopitelné, plné nesmyslných výjimek a nástrah čekajících na nepozorné, které můžou způsobit nevysvětlitelný propad výkonu. Optimalizace pak vypadá jako jakási prastará forma okultní magie, srozumitelná jen kruhu zasvěcených.

Pravdou je však naprostý opak. CPU jsou výjimečná soustrojí, kde má všechno do posledního detailu svůj důvod, a které neobsahují nic zbytečného. Všechny komponenty jsou vzájemně propojené a interagují spolu v delikátní souhře.

Některé z těchto důvodů se nám nemusí líbit. Například instrukční sada x86 i přes to, jak je barokní, redundantní a bizarní, má svůj důvod – zpětnou kompatibilitu. A Intel to ví nejlépe. V průběhu historie se snažil x86 zabít a nahradit něčím lepším čtyřikrát (iAPX 432, i960, i860 a naposledy v tandemu s HP velice ambiciózním Itaniem). Všechny tyto snahy ale selhaly, když AMD rozšířilo x86 o podporu 64 bitů. Instrukční sada už tak postrádala eleganci a rozšíření bylo ještě méně elegantní, ale ujalo se.


V tomto článku ukážu a vzápětí vysvětlím některé z těchto podivností na sérii testů a benchmarků napsaných v jazyce C. Proč zrovna v céčku? Proč ne ne třeba v Javě, když i v ní je možné dosáhnout low-level výkonu. To je pravda, ale když se chci bavit s hardwarem, vyberu si jazyk, který přede mě klade minimum překážek a je co nejpředvídatel­nější.

Všechny testy zkompiluji přes gcc -O3 a proženu perfem a budu měřit takové věci jako počet taktů procesoru, instrukcí, IPC (instructions-per-cycle/instrukce za takt), podíl branch-miss a podíl cache-miss. Každý z nich začíná tímto kusem kódu, který není měřený:

// argv[1] musí být mocnina dvou, ve všech testech používám hodnotu 16
const int len  = atoi(argv[1]) * 1024 * 1024;
const int iter = 100 * len;

long *seqarr  = sequential_order(len);
long *randarr = random_order(len);

Funkce sequential_order vygeneruje pole, kde platí, že na pozici n je hodnota (n+1) % len. Jde tedy o offset o jednu pozici vpravo. Dá se o tom taky uvažovat jako o pointeru, který ukazuje na vedlejší pozici.

Funkce random_order vygeneruje podobné pole offsetů/skoro-pointerů, které však mají náhodné pořadí. Když budu opakovaně následovat tyto pointery, v obou případech projdu všechny položky daných polí.

mod vs mask

Není instrukce jako instrukce a v případě superskalárních out-of-order procesorů to platí obzvlášť silně. Každá operace má různou latenci, která může být někdy zamaskována instrukční pipeline, a různou propustnost, kdy CPU dokáže vypravit několik instrukcí v jednom taktu.

Abych si tohle ověřil, budu testovat staré známé modulo a maskování. Jde o jednoduchý způsob jak zajistit, že hodnota i zůstane v určitém rozsahu len. Buď můžu spočítat i % len nebo, pokud je len mocnina dvou, můžu to samé napsat jako i & (len - 1). Část len - 1 vyrobí masku, která obsahuje samé jedničky na nejnižších pozicích a & vyfiltruje nejnižší bity čísla i.

Jaký je tedy v rozdíl v rychlosti těchto dvou fragmentů kódu:

// MOD
long sum = 0;
for (int i = 0; i < iter; i ++) {
  sum += seqarr[i % len];
}
// MASK
long sum = 0;
for (int i = 0; i < iter; i ++) {
  sum += seqarr[i & (len - 1)];
}

Výsledky jsou následující:

* MASK MOD
time 1.2 s 4.5 s
cycles 4103.9M 15655.8M
instructions 11745.5M 13426.9M
IPC 2.86 0.85
cache-references 104.0M 10.1M
cache-misses 43.1M 7.9M
cache-misses 41.4% 78.7%
branches 1677.9M 1678.6M
branch-misses 0.0M 0.0M
branch-misses 0.0% 0.0%

Jak se ukazuje modulo je v tomto případě skoro 4× pomalejší než odmaskování pár bitů. Proč? Celočíselné dělení2/modulo je pomalé, protože ho nikdo nepoužívá a protože ho nikdo nepoužívá, tak tyto operace není nutno optimalizovat (stejně jako např. instrukce INDEX na starých VAXech).

Tento fakt může mít vliv na návrh hashovacích tabulek a jiných datových struktur. Když použiji velikost, která je mocninou dvou, vyhnu se modulu/dělení3 a operace s tabulkou můžou být o něco rychlejší. Na druhou stranu to může mít neblahý dopad, pokud jsou (například) klíče tabulky násobky osmi nebo jiné malé mocniny dvou, hash funkce tyto korelace zachová, v důsledku toho bude obsazena jen 1/8 slotů tabulky, to povede k velkému počtu kolizí a to věci značně zpomalí. V tomto případě je třeba si na tohle dát pozor a používat dobrou hashovací funkci.

Pozn: Kompilátor dokáže triviální případy jako x % 64 optimalizovat na x & 63. Někdy si dokáže poradit i s mnohem záludnějšími případy, kdy to vypadá, že se nedělí/nemoduluje konstantou. V určitých případech může být modulo eliminováno nebo vytaženo před smyčku.1

Pozn 2: Větší počet cache-miss v případě maskující varianty může být způsobený tím, že kód je tak rychlý, že prefetch nestíhá data dodávat z paměti.

Závislé a nezávislé operace

Superskalární CPU mají pod kapotou každého vlákna značné množství paralelismu – v každém taktu dokážou načíst, dekódovat, naplánovat, vykonat několik instrukcí. Tato dech beroucí mašinérie může však fungovat jedině, když má k dispozici dostatečné množství nezávislých operací v jednom nebo v několika hyperthreadova­ných/SMT vláknech. Když nemá po ruce dostatečné množství práce, které by nasytilo OOO bestii, CPU musí čekat.4

Následující test testuje právě tohle. V prvním případě program iteruje polem sekvenčních offsetů. V druhém nepoužívá k navigaci offsety, ale čte hodnoty přímo. Oba programy prochází polem naprosto shodně a udělají stejné množství práce. Na konci má proměnná sum v obou případech stejnou hodnotu.

//závislý load

long sum = 0;
int x = arr[0];
for (int i = 0; i < iter; i ++) {
  sum += x;
  x = arr[x];
}
//nezávislý load

long sum = 0;
for (int i = 0; i < iter; i ++) {
  sum += arr[i & (len - 1)];
}

Výsledky:

* závislý load nezávislý load
time 2.61 s 1.23 s
cycles 8978.9M 4139.9M
instructions 8396.2M 11751.5M
IPC 0.93 2.83
cache-references 23.4M 104.9M
cache-misses 12.4M 43.2M
cache-misses 53.0% 41.2%
branches 1679.0M 1679.0M
branch-misses 0.0M 0.0M
branch-misses 0.0% 0.0%

Jak je vidět, i když nezávislá verze potřebuje více instrukcí, běží 2× rychleji.

V prvním případě jedna iterace plně závisí na té předchozí. Přiřazení x = arr[x] tvoří datovou závislost. Následující iterace nemůže pokračovat, pokud ještě neskončilo načítání arr[x], protože procesor neví, která hodnota bude načtena a tedy neví, na jakou adresu skočit. V tomto případě to bude x+1, ale to procesor nemůže vědět, ten vidí jen bajty naházené v paměti, a proto musí čekat, dokud načítání neskončí, než bude pokračovat.

Druhý kus kódu také obsahuje závislost, ale v tomto případě nejde o závislost datovou, ale o závislost řídící struktury (control flow). Procesor se musí jen rozhodnout, jestli i < iter bude vyhodnocena jako pravda nebo nepravda a podle toho skočit na další iteraci nebo smyčku zaříznout. Jde o jednoduchou binární volbu, kterou dokáže branch predictor odhadnout přesně i v komplikovaných případech. Tady jde o triviální smyčku, kterou procesor odhadne vždy, rozjede OOO mašinérii, začne spekulovat napříč iteracemi a najedou má k mání moře nezávislých operací (i za cenu, že se občas zmýlí ve spekulacích a bude muset zahodit některé výsledky).5

I když technicky vzato, procesor může odhadnout hodnotu před tím, než je načtena z paměti a na základě toho rozjet spekulativní exekuci. Pokud mají data jednoduše předvídatelný vzor jako v mém případě výše, procesor ho může velice snadno zachytit a využít. Obecně je to však problém, protože je třeba trefit správnou hodnotu z 264 možností (a navíc to má dopad na implementaci paměťového modelu). V případě control flow dependency se stačí rozhodnout mezi pouhýma dvěma možnostmi.6

Oba programy procházejí pamětí lineárně a proto prefetcher dodává všechna data s předstihem. Avšak stejně jako v minulém testu, se i tady ukazuje určitý počet cache miss. To může signalizovat, že RAM nestíhá dodávat data. Kdybych měl ve svém stroji rychlejší paměti, program by běžel o něco rychleji bez těchto cache-miss.

Rozdíl mezi závislými a nezávislými operacemi má překvapivě dopad úplně na všechno. Datové struktury založená na pointerech, jako jsou spojové seznamy nebo stromy, nedělají nic jiného než závislé operace – musím načítat jeden pointer za druhým než se prokoušu k listům stromu. Tohle platí pro všechny druhy stromů nebo hierarchických struktur, můžou být široké nebo mělké, ale v cestě za cílem vždycky stojí závislé operace. A shodou náhod všechny funkcionální persistentní datové struktury jsou implementovány právě jako stromy.


Pokud upravím předminulý mod/mask test tak, aby na sobě byly jednotlivé iterace závislé, situace je najednou hrozivě chmurná.

* mask mod
time 3.15 s 21.32 s
cycles 10805.6M 74504.3M
instructions 10076.0M 11768.4M
IPC 0.93 0.15
cache-references 15.1M 7.8M
cache-misses 7.4M 6.8M
cache-misses 49.1% 87.0%
branches 1679.4M 1681.8M
branch-misses 0.0M 0.0M
branch-misses 0.0% 0.0%

Najednou je závislé modulo 7× pomalejší než závislé maskování a 21× pomalejší než nezávislé plně pipelinované maskování. To (kromě jiného) není vůbec dobrá zpráva pro RRB stromy. IPC 0.15 mluví za své.

Závislé sekvenční čtení a závislé náhodné čtení

Je čas přitlačit na pilu. Co se stane, když porovnám závislé sekvenční čtení jako minule a závislé čtení z náhodných míst v paměti, tedy něčeho, co připomíná divoký hon na pointery napříč celou haldou.

// SEQ
long sum = 0, x = seqarr[0];
for (int i = 0; i < iter; i ++) {
  sum += x;
  x = seqarr[x];
}
// RAND
long sum = 0, x = randarr[0];
for (int i = 0; i < iter; i ++) {
  sum += x;
  x = randarr[x];
}

Jde o shodný kód, který se liší jen v tom, jakým polem prochází. V jenom případě polem lineárně rostoucích offsetů, v druhém případě offsetů rozházených bez ladu a skladu.

* SEQ RAND
time 2.45 s 130.98 s
cycles 8650.7M 445541.6M
instructions 8391.4M 8643.6M
IPC 0.97 0.02
cache-references 21.1M 1698.0M
cache-misses 10.1M 1662.8M
cache-misses 48.1% 97.9%
branches 1678.2M 1722.0M
branch-misses 3592 0.4M
branch-misses 0.0% 0.026%

Jestliže všechny ostatní testy představovaly jen drobné nepříjemnosti, tohle rozhodne o všem. Když můj program prochází data sekvenčně, prefetcher to pozná a postará se, že všechno potřebné bude včas k dispozici v rychlé cache, připravené ke čtení. Naproti tomu, když chci v každé iteraci číst z nepředvídatelné lokace v paměti, prefetcher nemůže pomoct, cache jsou k ničemu a výsledek je 53× pomalejší. Pokud jsou data větší než pár megabajtů, které má L3 cache, každý load jde přímo do RAM a RAM je zatraceně pomalá. Hlavní paměť má ucházející propustnost a proto nepředstavuje zas takový problém, když je program limitovaný právě průtokem bitů. Když však záleží na latenci každé operace, všechno jde do kytek, protože závislosti vytvoří kauzální řetězec.

Datové struktury, které jsou procházeny v předvídatelném lineárním pořadí nebo ty, jejichž horká část se vejde do některé úrovně cache procesoru budou rychlé. Ty, které tyto podmínky nesplňují, můžou trpět.

Například java.util.HashMap je implementovaná tak, že potřebuje 3 závislé load operace. Když je hashmapa velká a nevejde se do cache, každý přístup do ní, může stát několik cache-miss. Nedávno jsem napsal StringIntDicti­onary – mapu specializovanou pro stringové klíče a primitivní hodnoty, která je až 3× rychlejší než j.u.HashMap[String, Int] právě proto, že potřebuje načíst data jen z jednoho místa v paměti a proto vyvolá jen jeden cache-miss.7

Jako další příklad může sloužit scala.collection.immutable.Vector. Ten je implementován velice širokým stromem a v závislosti na velikosti potřebuje něco mezi 1 a 7 závislými load operacemi. Na druhou stranu, když k danému vektoru přistupuji často, prvních pár úrovní bude v cache a cena přístupu tak nemusí být přehnaná. V benchmarcích jsem zjistil, že Vector je 4–8× pomalejší než obyčejný plochý ArrayBuffer.

Za zmínku může stát také tento test specializovaných map, který ukazuje, že ty mapy které potřebují přistoupit jen k jednomu místu v paměti, jsou rychlejší než ty které potřebují přistoupit ke dvěma nebo třemi.

Nezávislé sekvenční čtení a nezávislé náhodné čtení

Další test: Jak je na tom nezávislé sekvenční čtení v porovnání s nezávislým náhodným čtením.

// SEQ
long sum = 0;
for (int i = 0; i < iter; i ++) {
  sum += seqarr[seqarr[i & (len - 1)]];
}
// RAND
long sum = 0;
for (int i = 0; i < iter; i ++) {
  sum += rndarr[rndarr[i & (len - 1)]];
}

Zase jde o stejný kód, který vede ke stejnému výsledku, liší se však jen pořadím v jakém přistupuje k elementům pole – jednou sekvenčním a podruhé náhodném.

* SEQ RAND
time 1.22 s 25.55 s
cycles 4073.1M 86769.1M
instructions 13431.2M 13718.2M
IPC 3.29 0.15
cache-references 91.8M 2895.3M
cache-misses 36.4M 1835.3M
cache-misses 39.7% 63.3%
branches 1679.3M 1727.9M
branch-misses 22135 386246
branch-misses 0.0% 0.022%

Rozdíl v rychlosti je tentokrát pouze 21× v neprospěch náhodného čtení. Jak je vidět z těchto výsledků, i když se čte z náhodných míst v DRAM, procesor dokáže spekulovat a začít pár čtení najednou a zamaskovat tak latence. Nezávislé sekvenční čtení je ještě o něco rychlejší než závislé sekvenční i když potřebuje o 50% více instrukcí.

Když porovnám všechny případy, začne se mi rýsovat následující obraz:

test čas 1 iterace 1 iterace IPC
SEQ INDEP L1 0.95 1.8 taktu 0.6 ns 4.45
SEQ INDEP 1.22 2.4 taktů 0.8 ns 3.29
SEQ DEP 2.45 4.9 taktů 1.5 ns 0.97
RAND INDEP 25.55 51.1 taktů 16.0 ns 0.15
RAND REP 130.98 262.0 taktů 81.6 ns 0.02

První řádek reprezentuje případ, kdy se všechna data vejdou do L1 cache a program nepotřebuje číst a být limitován rychlostmi DRAM. Jak je vidět, v tomto případě procesor našel dost paralelismu, na to aby stihl vykonávat 4.45 instrukce každý takt.

Důležitý je rozdíl mezi nejrychlejším a nejpomalejší verzí – i když dělají to samé, jedna běží 140× pomaleji než ta druhá. To zhruba odpovídá rozdílu v rychlosti mezi dnešními počítači a těmi z roku 2003. Jde sice o vykonstruovaný případ, ale přesto ukazuje rozpětí, do kterého budou reálné programy spadat.

branch prediction

Na závěr ukážu neintuitivní chování podmínek.

Funkce generate_truly_random_array vygeneruje pole náhodných čísel z rozsahu INT_MIN až INT_MAX. Chci spočítat nějakou agregaci těchto hodnot, ale některé nevyhovující bych chtěl ignorovat. V prvním případě započítám jen hodnoty větší než 0 a v druhém všechny kromě nuly. Dva kusy kódu na rozdíl od předchozích povedou k různým součtům.

long* rndarr = generate_truly_random_array(len);

double d = 0.5;
long count = 0;
double sum = 0.0;
for (int i = 0; i < iter; i++) {
        if (rndarr[i & (len - 1)] > 0) { // tady se to liší
                sum += rndarr[i & (len - 1)] / d;
                count += 1;
        }
}
long* rndarr = generate_truly_random_array(len);

double d = 0.5;
long count = 0;
double sum = 0.0;
for (int i = 0; i < iter; i++) {
        if (rndarr[i & (len - 1)] != 0) { // tady se to liší
                sum += rndarr[i & (len - 1)] / d;
                count += 1;
        }
}

Výsledky:

* x > 0 x != 0
time 8.05s 1.98s
cycles 25693M 6308M
instructions 19304M 23491M
IPC 0.75 3.72
cache-references 19M 35M
cache-misses 15M 12M
cache-misses 79% 34%
branches 3357M 3355M
branch-misses 838M 0M
branch-misses 25% 0%

Jak je vidět, tak verze, která dělá dvakrát víc práce, je 4× rychlejší. Problém tentokrát spočívá v předvídatelnosti podmínek. Pokud je podmínka dobře předvídatelná, branch predictor ji uhodne, začne číst a dekódovat instrukce správné sekvence následujících instrukcí a všechno jede, jako kdyby se v kódu žádný skok nevyskytoval. Když však uhodne špatně a bude muset všechnu spekulativní práci zahodit a začít od správného konce. Nastane takzvaný pipeline stall nebo pipeline bubble. Množství ztracené práce odpovídá hloubce pipeline, která na mém Haswellu má délku 14–19 taktů. V benchmarku pomalá verze potřebuje 11.5 taktu na každou iteraci, což přibližně odpovídá tabulkovým hodnotám a faktu, že rychlá verze musí pořád udělat 2× víc užitečné práce.

Dále je vidět, že jen 25% instrukcí vedlo k branch-miss. Polovina podmíněných skoků tvoří smyčku, kterou hardware uhodl naprosto dokonale. Ve zbývající polovině se trefil na 50% – měl tedy stejnou úspěšnost, jako kdyby házel mincí.

K tomuto případu se sluší dodat, že kdyby bylo tělo smyčky o něco jednodušší (to dělení proměnnou d tam nebylo jen tak pro nic za nic), kompilátor by ji mohl přeložit na sérii cmov instrukcí. cmov neboli conditional move nepředstavuje větvení toku programu, ale operaci, která je vykonána tak jako tak, podmínka určí jestli má nějaký efekt. Kdyby se tak stalo, běžely by obě verze stejně rychle, protože by neobsahovaly žádný nepředvídatel­ný skok.

x86 má jedninou takovou instrukci – cmov, mnoho instrukcí na 32bitových ARMech je predikovatelných a celé neslavné Itanium bylo na tomto principu postaveno.

Takže co?

Doufám, že jsem ukázal, jaký může být rozdíl v rychlosti programů, když pracují s hardwarem a když pracují proti němu. Když jsou algoritmy a datové struktury malé, lokální, sekvenční a předvídatelné, na moderním hardwaru doslova poletí, pokud tyto vlastnosti chybí, s rychlostí to může drhnout.

To je všechno, jděte domů a použijte, co jste se dozvěděli.


Pozn:

  1. S modulo je samá zábava. Modulo negativního čísla vrátí záporný výsledek a ani funkce Math.abs nemusí pomoci, pokud není použita správně.
  2. Dělení v plovoucí čárce je o něco rychlejší, ale přesto pomalejší, než násobení převrácenou hodnotou.
  3. Dělení a modulo je (aspoň na x86) jedna a ta samá operace. DIV/IDIV vrátí jak výsledek dělení, tak i zbytek po dělení.
  4. V tomto kontextu stojí za zmínku Cray XMT/ThreadStorm, který dotáhl hyperthreading/SMT do extrému. ThreadStorm je procesor, určený pro úlohy, které vykazují jen mizivou lokalitu dat a kde cache nepomůžou (jako např. zpracování gigantických grafů). Zatímco běžné x86 procesory si poradí s 2 vlánky najednou, Xeon Phy, Power8 a nejnovější Sparcy od Oraclu zvládají 8 vláken, ThreadStorm dokáže zpracovat 64 vláken najednou. porcesor obrahuje všechny hardwarové prostředy všech 64 vláken (soubor registrů atd.) a v každém taktu najde jedno vlákno, které na nic nečeká a může běžet a spustí ho a to všechno bez zásahu operačního systému. Jde tedy o hardwarovou implementaci přepínání procesů, když je aktivně běžící vlákno blokuje. Viz: Overview of the Next Generation Cray XMT.
  5. A že může spekulovat zatraceně hluboko. Nový Skylake od Intelu obsahuje buffery pro 224 instrukcí za letu, 72 paralelních load instrukcí a 180 integer a 168 floating-point registrů, které spekulace potřebují k přejmenovávání.
  6. Navíc i když procesor nemá plnou podporu pro predikci hodnot, může obsahovat takzvaný branch target predictor, který se snaží odhadnout cíle nestatických skoků (jump register). Tohle může zrychlit volání virtuálních metod v C++ a podobných jazycích, které se nemůžou spolehnout na JIT a inline cache a oportunistickou optimalizaci. Ale i v tomto případě je počet možných cílů mnohem menší než počet možných hodnot.
  7. Tohoto jsem dosáhl tak, že krátké alfanumerické klíče neukládám jako String objekty, ale jako speciálně zakódovaných do 54 bitů nebo offset/délku od souvislého pole znaků.

Flattr this!

Posted in CPU, Paměť | 2 Comments

Escape analysis

Dneska jenom krátce: Nedávno jsem zakopl o paper popisující partial escape analysis a napadlo mě, že bych tu mohl napsat pár slov o tom, co je escape analysis vlastně zač (termín s dovolením nebudu překládat a vystačím si s kurzívovou metodou).

Escape analysis je optimalizační technika používaná v kompilátorech, JITech a virtuálních strojích k eliminaci určitých alokací, která povede ke zrychlení programu. Jde v podstatě o to, že pokud objekt neunikne z metody (odtud to escape), není ho třeba alokovat na haldě, ale postačí alokace na zásobníku. Když metoda skončí, její rámec je zahozen a s ním jsou automaticky dealokovány i všechny objekty na zásobníku. Díky tomu není třeba ztrácet čas alokací1, ulehčí se práce garbage collectoru a protože zásobník je skoro vždy horký v cache, je práce s ním velice rychlá. Kompilátor může jít ještě o krok dál a provést scalar replacement, kdy objekt není alokován ani na zásobníku, ale jeho členské proměnné skončí v registrech CPU.

Objekt/pole unikne pokud:

  • je návratovou hodnotou metody
  • je přiřazen do statické proměnné
  • je přiřazen do členské proměnné jiného objektu, který uniká
  • je předán jako argument metody (včetně implicitního argumentu this)

Jde o tranzitivní vlastnost. Kompilátor může alokovat na zásobníku jen pokud dokázal, že objekt za žádných okolností nemůže uniknout a z toho důvodu musí být konzervativní.

Jako v případě každé jiné optimalizace, je nutné k dobré funkci escape analysis inlinovat. Inlinování rozšíří kontext, na kterém JIT provádí analýzu a to vede k lepším výsledkům.

class X {
  // metoda `doSomething` vrací Int, neuniká z ní `this`
  // pointer a je dost malá na to, aby byla inlinována
  def doSomething(): Int = ???
}

def f(): Int = {
  // objekt `x` neuniká z metody a může být alokován
  // na zásobníku
  val x = new X
  x.doSomething()
}

Efekt inlinování je vidět na následujícím kusu kódu.

// objekt x uniká z této metody
def f(): X = new X

def g(): Int = {
  // pokud je metoda `f` inlinována, JIT zjistí, že
  // objekt x neuniká z metody `g` a tedy může být
  // alokován na zásobníku pokud však `f` není
  // inlinována, alokuje objekt, vrátí ho a `g` se
  // s ním musí vypořádat
  val x = f()
  x.doSomething()
}

Přímočará escape analysis uspěje jen když si je jistá, že objekt nemůže uniknout nehledě na řídící struktury, podmínky a cykly. Pokud se však kód metody větví a v jedné větvi smyčky nic neuniká, ale v jiné, která se provede jen zřídka, uniká, JIT to vzdá a prohlásí, že objekt uniká.

Triviální příklad:

def f(): X = {
  val x = new X

  if (x.doSomething() != 0) {
    // v téhle větvi `x` neuniká
    null
  } else {
    // v téhle ano
    x
  }
}

Tohle je právě téma na začátku zmíněného paperu popisující partial escape analysis v GraalVM, která si dokáže poradit s podmínkami a control flow. Kompilátor generuje kód, který provádí alokaci na haldě až na poslední chvíli, kdy je jasné, že osud objektu je zpečetěn a nemůže se stát nic jiného, než že unikne z metody. Ale než tento okamžik nastane, objekt žije na haldě nebo v registrech.

Výsledek může odpovídat tomuto pseudokódu:

def f(): X = {
  val x1, x2, x3 = scalarReplacedFieldsOfX

  if (x.doSomething() != 0) {
    // nic se nealokuje
    null
  } else {
    // tady proběhne alokace, objekt je
    // rematerializován na haldě
    val x = new X
    setFields(x, x1, x2, x3)
    s
  }
}

V určitých benchmarcích použití partial escape analysis vede k redukci alokací o 58% a zrychlení programu o 33%.

Zmíněný GraalVM je projekt poskytující API, kterým je možné řídit chování JVM JITu. S touto kontrolou nad kompilátorem je možné na JVM jednoduše naroubovat podporu jazyků, které se chovají jinak než Java a které tedy standardní JIT nekompiluje ideálním způsobem. S GraalVM není třeba při implementaci nového jazyka začínat od nuly, protože má k dispozici všechny prostředky JVM jako je GC a efektivní JIT, který dokáže kód optimalizovat a (hlavně) deoptimalizovat.

O JRuby backendu založeném na Truffle a GraalVM píše Chris Seaton na svém blogu. Jeden z nejzajímavějších článků je ten, kde popisuje, jak překonali MRI Ruby s rozšířeními napsanými v céčku tím, že interpretovali C kód spolu Ruby kódem pomocí GraalVM.

Další přístup k eliminaci alokací je escape detection (letmo zmíněná Cliff Clickem), kterou používal Azul na jejich JVM běžící na jejich hardware.

Escape detection je optimistický přístup – všechny objekty jsou alokovány na zásobníku s tím, že hardware detekuje, jestli neunikl pointer na zásobník. Když unikl, procesor vyvolá výjimku, jejíž handler přesune objekt na haldu a přepíše pointer. Tento postup údajně vede k velice efektivní redukci alokací, protože je agresivní, nemusí mít 100% důkaz a sám o sobě dělá to, co partial escape analysis.


Dále k tématu:

Poznámky:

  1. Za alokaci se platí dvakrát v propustnosti paměti. Poprvé, když vytvářím objekt, musím jeho blok paměti načíst z RAM do cache, vynulovat (aspoň v případě JVM), přepsat užitečnými daty. Podruhé, když už není třeba, je vyhozen z cache a musí být zapsán zpátky do RAM na místo, kam patří na haldě. Tohle zabolí zvláště v případě objektů, které žijí jen deset nanosekund – po velice krátké době k nim už nikdy nepřistoupím, ale stále je třeba dodržet cache coherence protokol a zapsat je zpátky do RAM. I když propustnost pamětí je typicky velká, není nekonečná a obrovské tempo alokací na mnohajádrovém procesoru může trubky pořádně přiškrtit. Proto je alokace na zásobníku tak efektivní – opakovaně používá stejný blok paměti který je stále v cache a do RAM skoro vůbec nepřistupuje.

Flattr this!

Posted in Algo, JVM, VM | 5 Comments

Jak rychle řadit a šetřit čas

Když už tu mluvímřazení a řadící mašinerii, tak bych taky mohl stručně a srozumitelně shrnout co, jak a proč dělat a čemu se vyhnout.

Asi takhle: Obecné řadící algoritmy jako quicksort, mergesort nebo heapsort jsou super. Sami o sobě nejsou nejrychlejší, ale nejlepší hlavy strávily posledních padesát let jejich optimalizací a zlepšováním a proto dá se na ně spolehnout a výkonem nikoho nezahanbí.

Výsledkem oněch pěti dekád je například rafinovaný three-way quicksort, který je odolný proti zákeřným vstupům4 nebo celá paleta hybridů jako je stabilní Timsort (který se snaží využít seřazené běhy dat, které se přirozeně vyskytují ve vstupní kolekci, začíná práci insertion sortem a pak nastartuje merge sort) nebo nestabilní1 introsort (začne rychlým quicksortem, který v patologických případech přepadne do heapsortu a krátké sekvence řadí insertion sortem, který je sice asymptoticky pomalý, ale pro malé n je rychlejší než chytré algoritmy).

Samostatný k-ary heapsort také může být efektivní a rychlostí může konkurovat quicksortu, pokud je porovnání levné (třeba prosté porovnání intů) a přístup do paměti relativně drahý.


Pokud ale řadím krátké stringy pevné délky (jako třeba 32/64bit inty nebo floaty) a mám k dispozici velkou propustnost paměti least significant digit (LSD) radix sort je neporazitelný.

Klasická implementace řadící osm bitů v jedné iteraci potřebuje 5 iterací pro 32bit klíče (první spočítá frekvence jednotlivých bajtů, po které následují čtyři řadící iterace) a musí přenést všechna data třináctkrát (čtyři průchody potřebují data načíst a pak je zapsat na jiné místo, zapisování nejdřív natáhne z RAM stará data do cache, tam je přepíše a pak je později zapíše zpátky do RAM)2.

LSD radix sort data vždycky streamuje, vůbec nevyužívá cache a za všechno platí propustností pamětí. Quicksort nebo most significant digit (MSD) radix sort rekurzivně dělí vstup na menší a menší sekvence. Ty se eventuálně vejdou do cache a od toho okamžiku jsou všechny iterace limitované propustností cache, která bývá mnohem větší než propustnost RAM a je lokální každému jádru. Pokud je propustnosti málo, je možné nejdřív udělat quicksort, MSD radix sort nebo jinou diskriminaci do částí, které se vejdou do nějaké úrovně cache a na nich rozjet rychlý LSD radix sort. Vznikne tak několikaúrovňový hybrid, na který se dá narazit v literatuře.


Když mě zajímá řazení dlouhých stringů potenciálně proměnné délky, na současných cache CPU architekturách je nejefektivnějším algoritmem burstsort. Ten je založený na předpokladu, že největší dopad na rychlost řazení mají opakované cache-miss, kdy data nejsou připravena v cache, ale musí se pro ně do RAM. Quicksort nebo jiný diskriminační algoritmus musí na každé úrovni iterace projít všechny položky a pokud jich je hodně3, přístup ke každé z nich vyvolá jeden pomalý cache-miss5. A protože se dělá O(log n) iterací, bude třeba i tolik pomalých cache-miss na každý element vstupní sekvence. Jak rekurze postupuje, eventuálně se všechny položky jedné subsekvence vejdou do cache a od té doby jsou iterace rychlé a bez cache-miss. Pokud jsou elementy velké, náhodně rozhozené nebo kompozitní, tak každý z nich může zabrat několik cache-line a využití cache pak není nijak efektivní. Burstsort se snaží snížit počet cache-miss a zrychlit tak řazení. Nepracuje rekurzivně, jen jednou projede vstupní data a na základě prefixu stringového klíče je přidá do částečné trie. Nesnaží se umístit element na pozici určené plnou cestu, ale jen její částí, která vede k listu trie obsahujícímu kolekci neseřazených položek sdílejících daný prefix. Když algo dorazí k listu, vloží do něj nový element, v případě potřeby zvětší list, a když přesáhne určitou mez, kolekce v listu praskne (odtud burst) a vytvoří se další úroveň trie. V jedné iteraci je tedy vytvořena trie, kterou pak projdu v in-order pořadí, jiným algoritmem řadím kolekce, na které narazím, a získám tak seřazený výsledek. V ideálním případě jsou potřeba jen 2 cache-miss na položku: Jeden při vkládání do trie a jeden při závěrečném řazení9. Oproti diskriminačním algoritmům má burstsort tu výhodu, že se do cache nemusí vejít kolekce elementů, ale jen mnohem kompaktnější kostra trie stromu. Ta neobsahuje žádná data samotných elementů a navíc potřebuje jen jeden a něco pointer na každou neseřazenou kolekci v listu, která může obsahovat něco mezi 100 a 1000 elementy.

Javovská verze7 burstsortu používající Timsort (kterou jsem opilý napsal ve tři ráno a zkompiloval jsem ji až druhý den) je 2.5–3× rychlejší než Timsort ze standardní knihovny Javy. Je stejně rychlá jako three-way radix quicksort* z má zahrádky (který se výborně hodí pro řazení stringů se společným prefixem). Hybrid burstsort + *three-way radix quicksort je v mých testech 3–4× rychlejší než samotný Timsort. To podle mě není vůbec špatný výsledek pár hodin programování.

Dva dny po přednášce mě napadl způsob, jak vylepšit burstsort pomocí schwartzian transform a v ideálním případě snížit počet cache-miss na polovinu. O tom se však rozepíšu až někdy příště.


Pro problémy, kdy je počet možných hodnot srovnatelný nebo menší než velikost vstupní kolekce, dá se řadit v lineárním čase pomocí counting sortu nebo pigeonhole sortu. A když dopředu znám rozložení vstupních dat, můžu řadit v očekávaném lineárním čase pomocí bucket sortu nebo proxmap sortu.

Ale to není všechno. Co když chci řadit líně? Heap sort je klasická volba, ale líně se dá řadit i pomocí quicksortu, radix sortu nebo quicksortu reifikovaného do podoby stromu (na tohle jsem nikde v literatuře nenarazil a tak jsem to skromně pojmenoval YOLO tree, relevantní informace v Implementing sets efficiently in a functional language).

Když mě zajímají percentily, klasickou cestou je quickselect. Ten se dá také naroubovat na radixsort nebo jiný diskriminační algoritmus. *Median of medians* zaručí lineární čas, ale konstantní faktory nejsou nijak příznivé.6 Po­dobného výsledku se dá dosáhnout adaptací merge sortu a trochu pravděpodobnostní magie. Tento algoritmus jsem také nikde nenašel a proto jsem mu promptně dal dočasné jméno mergeselect. Pokusím se o něm v budoucnu něco napsat.

A co paralelní řazení? Quicksort se dá triviálně paralelizovat pomocí součtů prefixu (prefix sum) jak je popsáno v Programming on Parallel Machines, stejně tak radix sorty nebo merge sort. Určité verze Radixsortu můžou být nejen paralelní, ale také in-place.

If a thing like this is worth doing at all, it’s worth doing right

Klasický přístup k řazení quicksort and chill funguje, rychlost nikoho neurazí a většinou připraven k okamžitému použití ve standardních knihovnách většiny jazyků. Tohle funguje, ale vždycky se spokjit jen s tím, co pouze funguje, může být rychlé a pohodlné řešení, ale jen stěží je intelektuálně uspokojivé. Pokud na něčem záleží, mohli bychom se poohlédnou po způsobech, jak to udělat co nejlépe.

Já v současné době přemýšlím, jak vzít Hengleinovy nápady, z kombinátorů popisujících řazení makry vygenerovat kompaktní a efektivní kód, který nealokuje dočasné datové struktury, a zapojit to do vylepšeného burstsortu. Výsledkem by mohla být knihovna, která je nejen generická, ale také rychlejší než všechno ostatní. To je podle mě plán dostatečně odvážný na to, aby se do něj někdo pustil.


Dále k tématu:


Pozn:

  1. Z každého nestabilního sortu se dá snadno udělat stabilní, když se neřadí podle klíče, ale podle dvojice [klíč, pořadové číslo elementu].
  2. Situace se dá v tomto případě vylepšit použitím tzv. non-temporal store instrukcí, které nezpůsobí natažení cache-line do cache, ale rovnou data zapíšou do RAM. Non-temporal instrukce jsou výhodné, když vím, že data po zápisu nebudu dlouho potřebovat a nechci jimi špinit cache a tahat je do věcí cache coherence protokolu.
  3. A jsou to pointery na dynamicky alokované objekty a ne pole structů, se kterými pomůže prefetcher, viz Intel optimization manual, který poskytuje představu o tom, co prefetch dokáže.
  4. Problematický vstup je například sestupně nebo vzestupně seřazená sekvence. Klasický Hoarův quicksort, který jako pivot vybírá první element si na ní vyláme zuby. Změnou volby pivotu na best-of-3, quasi-best-of-9 nebo náhodný výběr se tento patologický případ dá eliminovat. Ale změna pivotu nepomůže se vstupem obsahujícím všechny stejné hodnoty. Na ty je třeba pozměnit algoritmus, aby dělil položky do tří skupin: menší než pivot, rovné pivotu a menší než pivot. S touto změnou quicksort pole identických položek trvá lineární čas místo katastrofického kvadratického času.
  5. Moderní OOO hardware dokáže spekulovat dopředu a zatímco čeká na data z RAM, spekulativně vykoná kód následujících iterací, které můžou obsahovat další load instrukce, které skončí jako cache-miss. Takto může například načítat čtyři věci najednou a efektivně taz zredukovat latenci na čtvrtinu. I když je to značné zrychlení, může to být až 100× pomalejší než sekvenční čtení z paměti nebo přístup do cache. (viz MLP)
  6. Když jsem připravoval výše zmíněnou přednášku, na mysl se mi stále vracela slova z knihy/filmu No Country For Old Men, která řekl Anton Chigurh lovci odměn Carsonovi chvíli před tím, než ho zabil: „If the rule you followed brought you to this, of what use was the rule?“ K čemu jsou nám záruky dobré asymptotické komplexity, když výsledek bude pomalý. Stejně tak není žádná hanba použít teoreticky neefektivní algo, který je ale v praxi na reálných datech rychlý.
  7. Někteří jistě budou namítat, že JVM je špatná volba pro kód, pro který je kritický výkon. Těm lidem vzkazuji, aby si šli hrát na pískoviště. JVM má jen problém v tom, že není možné ovlivnit rozložení paměti – buď mám objekty nebo pole primitivních typů/referencí a to je všechno. To člověku svazuje ruce v návrhu datových struktur, protože nemůže požít například objekt, který má na konci inlinované pole nebo struct tam, kde by se přirozeně hodily. Například burstsort řazení Java string objektů je o 50% pomalejší než řazení nahých polí charů/bajtů, protože string objekt obsahuje pointer na pole znaků, který je třeba dereferencovat, což vede k datovým závislostem a (potenciálně8) dalším zbytečným cache-miss. Kdyby string objekt obsahoval pole znaků na konci, bylo by to mnohem lepší.
  8. Tady situaci značně komplikuje fakt, že JVM je dynamické prostředí s garbage collectorem, který přesouvá data z místa na místo a může (ale tady nemusí) data, která jsou použita spolu (např. obsahují datově závislé pointery), skončí vedle sebe. Někdy se může stát, že výsledek benchmarku bude záviset jestli proběhl GC, v jakém momentu proběhl a v jakém pořadí zpracovává kořenové reference. Když je kód velice těsný a jde o nanosekundy, takováto šaráda může mít dopad v řádech desítek procent.
  9. Toto platí pokud se celý obsah listu vejde do cache, což není problém, když listy mají maxilální velikost v rozmezí od 1k do 8k. Dále bude pár cache-miss potřeba během praskání listů, ale to se děje jen zřídka. Když se kostra stromu nevejde do cache, budou třeba další cache miss. Naštěstí se burst trie široce větví (v literatuře se uvádí praktická hodnota 256) a proto, když se do cache vejde burst trie např. pro 1M položek, jeden extra cache-miss je třeba pro vstupní kolekci do 256M položek.

Flattr this!

Posted in Algo, CPU, Paměť, řazení, Uncategorized | Leave a comment