0xDEADBEEF (⇉english edition⇉)

[RSS] [odkazy]
««« »»»

Rychlý průnik množin a Jaccardův index přes SIMD instrukce

10. 8. 2020

Na funk­ci­o­nálně.cz jsem kdysi psal, jak rychle po­čí­tat Jac­car­dův index po­dob­nosti. Jac­card udává míru, do jaké jsou si dvě mno­žiny po­dobné jako ve­li­kost prů­niku vy­dě­le­nou ve­li­kostí sjed­no­cení. Ke své práci po­tře­buje co nej­rych­leji spo­čí­tat ve­li­kost prů­niku a to jsem v článku ro­ze­bí­ral. Od tri­vi­ál­ních pokusů jsem se pro­pra­co­val k mno­ži­nám ve formě se­řa­ze­ných polí, vek­to­ri­zo­vané verzi, minha­sho­vání a LSH.

Něco na ten způsob jsem zase po­tře­bo­val a proto je čas téma oprášit ještě jednou a na­po­sledy.


Intel do x86 v roz­ší­ření SSE 4.2 přidal pár in­strukcí za­mýš­le­ných pro zrych­lení práce se stringy (tak jsou po­u­žity např. v JVM). Umí toho ale mnohem víc. Jedna z nich – velice mocná ope­race pcmpestrm – se dá využít pro rychlý vý­po­čet průnik množin.

Jedna va­ri­anta pcmpestrm vezme dva 128bitové re­gis­try, obsah in­ter­pre­tuje jako pole osmi 16b intů a vy­pro­du­kuje bitmasku udá­va­jící, které hod­noty z prv­ního re­gis­tru se vy­sky­tují někde v tom dalším. Ne­po­rov­nává dva vek­tory pa­ra­lelně, první ele­ment s prvním, druhý s druhým a tak dál, ale každý ele­ment se každým na­jed­nou. Jedna in­strukce tak pro­vede 64 po­rov­nání a další verze, která vstupy bere jako šest­náct 8b intů, vykoná rovnou 256 po­rov­nání. A přesně tohle po­tře­bu­jeme.

Ska­lární verze ve­li­kosti prů­niku (v jazyku D) vypadá takhle:

int intersectionSize(T)(T[] a, T[] b) {
  int ai, bi, count;
  for (; ai < a.length && bi < b.length; ) {
    T aa = a[ai];
    T bb = b[bi];
    count += (aa == bb);
    ai += (aa <= bb);
    bi += (aa >= bb);
  }
  return count;
}

Vezme dva ele­menty, po­rovná je, když jsou stejné, in­kre­men­tuje počet a oba indexy, jinak posune index uka­zu­jící na menší z nich o jedno místo do­předu. Jde o velice jed­no­du­chý kód, bez pod­mí­ně­ných skoků ve smyčce a tak branch pre­dic­tor nemůže být zmaten z ne­před­ví­da­tel­ných vstupů.

Vek­to­rová verze z tohoto paperu (pozor, v ukázce kódu mají chybu) fun­guje na stej­ném prin­cipu, pra­cuje jen s celými vek­tory na­jed­nou, na­místo jed­not­li­vých ele­mentů.

while (ai < a.length-7 && bi < b.length-7) {
  ushort8 av = cast(ushort8) __builtin_ia32_loaddqu(cast(char*)&a.ptr[ai]);
  ushort8 bv = cast(ushort8) __builtin_ia32_loaddqu(cast(char*)&b.ptr[bi]);

  int4 res = cast(int4) __builtin_ia32_pcmpestrm128(
    cast(ubyte16)bv, 8, cast(ubyte16)av, 8,
    0x01 // _SIDD_UWORD_OPS | _SIDD_CMP_EQUAL_ANY | _SIDD_BIT_MASK
  );

  ushort a7 = av.array[7];
  ushort b7 = bv.array[7];
  ai += (a7 <= b7) * 8;
  bi += (a7 >= b7) * 8;
  count += popcnt(res.array[0]);
}

Kód je oš­k­livý a ma­toucí, jak bývá při po­u­žití in­trin­sics funkcí, fun­guje jen pro 2B ele­menty (pcmpestrm přímo ne­pod­po­ruje vek­tory 4B intů, ale exis­tují jiná řešení), ale jede rychle a jen na tom záleží. Vý­sledky v jedné úloze, která počítá Jac­carda pro velký počet množin, jsou tyto:

ska­lární průnik170.6 s
SSE + ska­lární113.6 s
pouze SSE78.1 s

Druhý řádek je v pod­statě úryvek kódu výše zpra­co­vá­va­jící bloky osmi 16 intů a po něm běží ska­lární verze, která dojede zbytek. Při pro­fi­lo­vání se uká­zalo, že právě tam pro­gram tráví zby­tečně dlouho. Třetí řádek v ta­bulce před­sta­vuje plně SSE verzi (zdro­ják tady), která všechno řeší přes pcmpestrm a rych­lostí značně vy­čnívá.

Ve vý­sledku je možné (v tomto pří­padě a na tomto stroji) do­sáh­nout víc než dvoj­ná­sob­ného zrych­lení jenom tím, že vy­u­ži­jeme všechny mož­nosti hard­waru. Navíc nejde o nic extra, SSE 4.2 Intel uvedl víc jak před de­ká­dou a proto se můžeme spo­leh­nout, že bude na každém pro­ce­soru, kde pro­gram poběží.


Do­da­tek: Jiný pro­gram za po­u­žití téhle metody zrych­lil asi 4.2×.

píše k47 (@kaja47, k47)