Obecn´e v´ypoˇcty na grafick´ych kart´ach – pouˇzitelnost, vlastnosti, praktick´e zkuˇsenosti Martin Kruliˇs, Jakub Yaghob KSI MFF UK Malostransk´e n´am. 25, Praha {krulis,yaghob}@ksi.mff.cuni.cz
Abstrakt. Ned´ avn´ y pokrok ve v´ yvoji grafick´ ych ˇcip˚ u umoˇznil prov´adˇet na bˇeˇzn´ ych grafick´ ych kart´ ach obecn´e v´ ypoˇcty. V kombinaci s univerz´aln´ı platformou OpenCL, jej´ıˇz prvn´ı implementace pro GPU byly uvedeny v druh´e polovinˇe roku 2009, slibuj´ı tyto technologie znaˇcn´e zrychlen´ı paralelizovateln´ ych u ´loh. Tento ˇcl´anek se zab´ yv´a pouˇzitelnost´ı tˇechto technologi´ı na r˚ uzn´e praktick´e probl´emy. Z´aroveˇ n nab´ıdne v´ ykonnostn´ı srovn´ an´ı prototypov´ ych implementac´ı na soudob´em hardware. Kl´ıˇ cov´ a slova: GPGPU, grafick´e karty, paralelizace
1
´ Uvod
Jedn´ım z hlavn´ıch c´ıl˚ u informatiky jak na poli teoretick´em tak praktick´em je snaha optimalizovat ˇreˇsen´ı nejr˚ uznˇejˇs´ıch probl´em˚ u za u ´ˇcelem jejich rychlejˇs´ıho zpracov´an´ı. Efektivnˇejˇs´ı algoritmy zvl´ adnou ˇreˇsit probl´emy rychleji, pˇr´ıpadnˇe vyˇreˇsit vˇetˇs´ı probl´emy ve stejn´em ˇcase. Paralelizace u ´loh byla vˇzdy jedn´ım ze slibn´ ych smˇer˚ u optimalizace ˇreˇsen´ı. V minulosti byl tento pˇr´ıstup v´ ysadou specializovan´ ych v´ ypoˇcetn´ıch stanic a s´ alov´ ych poˇc´ıtaˇc˚ u. D´ıky ned´avn´e revoluci na poli ˇcip˚ u grafick´ ych karet je dnes moˇzn´e z bˇeˇznˇe dostupn´ ych komponent postavit dom´ac´ı poˇc´ıtaˇc, kter´ y bude m´ıt v´ ypoˇcetn´ı v´ ykon srovnateln´ y s o nˇekolik let starˇs´ımi superpoˇc´ıtaˇci. Tento ˇcl´ anek shrnuje v´ ysledky zkoum´an´ı pouˇzitelnosti grafick´ ych karet pro obecn´e v´ ypoˇcty nejr˚ uznˇejˇs´ıch typ˚ u. U vˇsech pˇr´ıklad˚ u uv´ ad´ıme tak´e pouˇzit´ y algoritmus a popisujeme nutn´e u ´pravy pro jeho nasazen´ı na GPU. Kaˇzd´ y pˇr´ıklad tak´e demonstruje urˇcit´ y probl´em, kter´ y se m˚ uˇze pˇri tˇechto u ´prav´ach projevit, zejm´ena pak ot´azka latence a propustnosti pamˇet´ı.
1.1
Hardware a metodika mˇ eˇ ren´ı
Vˇsechny pˇr´ıklady byly testov´ any na grafick´e kartˇe ATI Radeon HD 5870 s ˇcipem RV870 a 1GB intern´ı pamˇeti DDR5 taktovan´e na 4800MHz s 256 bitovou sbˇernic´ı. J´adro pracuje na frekvenci 850MHz a obsahuje 1600 5-cestn´ ych shader jednotek (coˇz 1
odpov´ıd´ a 320 v´ ypoˇcetn´ım j´ adr˚ um), kter´e bˇeˇz´ı na stejn´e frekvenci [1]. Testovac´ı sestavu poh´an´ı procesor Intel Core i7 860 obsahuj´ıc´ı 4 fyzick´a j´adra s technologi´ı Hyperthreading (tzn. 8 logick´ ych jader) taktovan´a na 2.8GHz. Sestava obsahovala 4GB operaˇcn´ı pamˇeti DDR3 na frekvenci 1333MHz a 64-bitov´ y operaˇcn´ı syst´em Windows 7. Kaˇzd´ y test byl opakov´ an desetkr´at a v´ ysledn´e namˇeˇren´e ˇcasy jsou aritmetick´ ym pr˚ umˇerem tˇechto mˇeˇren´ı. Jednotliv´a mˇeˇren´ı kaˇzd´eho testu se od sebe neliˇsila o v´ıce neˇz 5% od v´ ysledn´eho pr˚ umˇeru.
1.2
Osnova
Sekce 2 struˇcnˇe shrnuje nejd˚ uleˇzitˇejˇs´ı aspekty architektury GPU a jej´ı abstrakce, kterou poskytuje framework OpenCL. Sekce 3 pˇredstavuje nejednoduˇsˇs´ı moˇzn´ y probl´em – jednoduch´e SIMD v´ ypoˇcty nad jednorozmˇern´ ym vektorem ˇc´ısel. Sekce 4 rozeb´ır´a velmi zn´ am´ y probl´em n´asoben´ı velk´ ych matic, na kter´em demonstruje probl´emy s pˇr´ıstupem do glob´ aln´ı pamˇeti GPU. Sekce 5 navazuje Floyd-Warshallov´ ym algoritmem, kter´ y je do jist´e m´ıry podobn´ y probl´emu n´asoben´ı matic. Sekce 6 se vˇenuje ot´azce NP probl´em˚ u a jejich ˇreˇsen´ı pomoc´ı backtrackingu. Sekce 7 shrnuje poznatky z cel´eho ˇcl´anku.
2
Architektura GPU a OpenCL framework
V t´eto sekci se budeme vˇenovat nejd˚ uleˇzitˇejˇs´ım vlastnostem architektury GPU [1, 2] a frameworku OpenCL [3] pro paraleln´ı v´ ypoˇcty. Neklademe pˇritom d˚ uraz na u ´plnost, sp´ıˇse se snaˇz´ıme zd˚ uraznit aspekty, na kter´e je tˇreba br´at zvl´aˇstn´ı ohled pˇri paralelizaci algoritm˚ u.
2.1
Architektura GPU
Architektury GPU a CPU se znaˇcnˇe liˇs´ı v cel´e ˇradˇe ohled˚ u, zejm´ena • v poˇctu v´ ypoˇcetn´ıch jader, • ve v´ ypoˇcetn´ım v´ ykonu (zejm´ena v plovouc´ı desetin´e ˇc´arce) a • v pˇr´ıstupu k ˇreˇsen´ı probl´emu latence pamˇet´ı. GPU m´a mnohem vyˇsˇs´ı poˇcet jader neˇz soudob´a CPU, avˇsak tato j´adra maj´ı celou ˇradu omezen´ı. J´adra jsou seskupena do SMP jednotek1 , kter´e aplikuj´ı SIMD2 v´ ypoˇcetn´ı model. Vˇsechna j´ adra jedn´e SMP jednotky maj´ı tedy stejn´ y program a vykon´ avaj´ı vˇzdy tut´eˇz instrukc´ı (pouze nad jin´ ymi daty). Napˇr´ıklad ATI Radeon pouˇzit´ y pro testy obsahuje 320 jader, kter´a jsou seskupena do 20 SMP jednotek po 16 j´adrech. Na jednu SMP jednotku b´ yv´ a typicky napl´anov´ano v´ıce vl´aken, neˇz kolik m´a jednotka procesor˚ u, pˇriˇcemˇz poˇcet tˇechto vl´aken je n´asobkem poˇctu procesor˚ u. Vˇsechna tato vl´ akna (jak bˇeˇz´ıc´ı tak napl´ anovan´a) bˇeˇz´ı ve stˇr´ıdav´em SIMD reˇzimu – tedy mus´ı 1 Streaming 2 Single
MultiProcessor Units Instruction Multiple Data
2
m´ıt stejn´ y program, kter´ y ale nevykon´avaj´ı zcela soubˇeˇznˇe, n´ ybrˇz se po skupin´ach stˇr´ıdaj´ı na dostupn´ ych j´ adrech. V pˇr´ıpadˇe, ˇze pr´avˇe bˇeˇz´ıc´ı sada vl´aken je nucena z nˇejak´eho d˚ uvodu ˇcekat (napˇr. na naˇcten´ı dat z glob´aln´ı pamˇeti), pl´anovaˇc tuto sadu odstav´ı a mezi t´ım nech´ a poˇc´ıtat jinou sadu napl´anovan´ ych vl´aken, aby byly procesory co nejv´ıce vyt´ıˇzeny. Touto technikou se GPU snaˇz´ı minimalizovat efekt latence pamˇet´ı. Grafick´ a karta obsahuje nˇekolik druh˚ u pamˇet´ı: • Registry – velmi rychl´ a pamˇet’, kterou m´a k dispozici kaˇzd´e j´adro. Nejnovˇejˇs´ı GPU nab´ız´ı aˇz 1024 registr˚ u po 32 bitech. • Sd´ılen´ a pamˇet’ – vyrovn´ avac´ı pamˇet’ SMP jednotky, do kter´e mohou pˇristupovat vˇsechna vl´ akna. Pro grafick´e v´ ypoˇcty funguje zpravidla podobnˇe jako L1 cache na CPU. Tato pamˇet’ je stejnˇe rychl´a jako registry, ale pˇr´ıstup do n´ı m´a jist´a omezen´ı (viz d´ ale). Souˇcasn´ a GPU maj´ı ˇr´adovˇe des´ıtky kB pamˇeti na SMP jednotku (typicky 32 nebo 64 kB). • Glob´ aln´ı pamˇet’ – je spoleˇcn´ a pro cel´ y ˇcip (vˇsechny SMP jednotky). Pˇr´ıstup do n´ı je v´ yraznˇe pomalejˇs´ı nebot’ vˇsechny SMP sd´ıl´ı pamˇet’ovou sbˇernici a data z n´ı nejsou vˇetˇsinou ukl´ ad´ ana do ˇz´adn´e (nebo jen velmi mal´e) vyrovn´avac´ı pamˇeti. Souˇcasn´e karty maj´ı stovky MB aˇz jednotky GB t´eto pamˇeti. Kromˇe tˇechto pamˇet´ı je potˇreba jeˇstˇe zd˚ uraznit, ˇze pˇri kombinov´an´ı v´ ypoˇct˚ u na CPU a GPU jsou data, kter´ a chceme pouˇz´ıt, typicky um´ıstˇena v pamˇeti hostuj´ıc´ıho PC (v operaˇcn´ı pamˇeti CPU). Proto je potˇreba nejprve data pˇren´est do glob´aln´ı pamˇeti grafick´e karty, neˇz je v˚ ubec moˇzn´e zaˇc´ıt samotn´ y v´ ypoˇcet.
2.2
Spouˇ stˇ en´ı vl´ aken
Framework OpenCL nab´ız´ı ˇradu mechanism˚ u, pomoc´ı kter´ ych je schopen detekovat paraleln´ı hardware a spustit na nˇem poˇzadovan´ y k´od. K´od je v programu reprezentov´an zpravidla ve zdrojov´em tvaru a pˇred spuˇstˇen´ım je zkompilov´an pˇr´ımo pro c´ılovou architekturu. Funkce, kter´e slouˇz´ı jako vstupn´ı body vl´aken, se naz´ yvaj´ı kernely. Kernel m˚ uˇze b´ yt spuˇstˇen bud’ klasick´ ym zp˚ usobem (jako jedna instance vl´akna), nebo paralelnˇe na v´ıce dat. Pˇri paraleln´ım spuˇstˇen´ı jednoho kernelu jsou instance vl´aken organizov´any do jedno aˇz tˇr´ı rozmˇern´eho pole. Svou pozici v poli m˚ uˇze kernel zjistit pomoc´ı vestavˇen´ ych funkc´ı a tato pozice tak´e jednoznaˇcnˇe urˇcuje, jakou ˇc´ast pr´ace m´a instance vl´akna vykonat. Kromˇe toho se instance vl´aken sdruˇzuj´ı do skupin. Velikost skupiny mus´ı b´ yt v kaˇzd´em rozmˇeru soudˇeln´a s velikost´ı pole vl´aken. Vl´akna v jedn´e skupinˇe jsou fyzicky mapov´ ana na jednu SMP jednotku se vˇsemi v´ yhodami (napˇr. sd´ılen´a pamˇet’) a omezen´ımi (napˇr. SIMD reˇzim), kter´a z toho vypl´ yvaj´ı. Vˇsechny operace prov´ adˇen´e na zaˇr´ızen´ı jsou spravov´any tzv. frontou poˇzadavk˚ u. Do t´eto fronty se vkl´ adaj´ı poˇzadavky na spuˇstˇen´ı kernel˚ u, pˇresuny dat z/do grafick´e karty a tak´e synchronizaˇcn´ı znaˇcky a bari´ery. Technick´e detaily si dovol´ıme z d˚ uvodu u ´spory m´ısta vynechat.
3
2.3
Spr´ ava pamˇ eti
Nyn´ı se pod´ıv´ ame na pamˇet’ jeˇstˇe jednou, tentokr´at z pohledu architektury frameworku, a zm´ın´ıme tak´e nˇekter´ a podstatn´a omezen´ı. Z hlediska OpenCL dˇel´ıme pamˇet’ na: • Priv´ atn´ı – pamˇet’ vlastn´ı jednotliv´ ym instanc´ım kernel˚ u. Tato pamˇet’ v podstatˇe odpov´ıd´ a registr˚ um jednotliv´ ych jader. • Lok´ aln´ı – pamˇet’ sd´ılen´ a mezi vl´akny v jedn´e skupinˇe. Tato pamˇet’ odpov´ıd´a sd´ılen´e pamˇeti SMP jednotek. • Glob´ aln´ı – pamˇet’ sd´ılen´ a vˇsemi instancemi kernelu, coˇz odpov´ıd´a glob´aln´ı pamˇeti GPU. • Konstantn´ı – speci´ aln´ı pˇr´ıpad glob´aln´ı pamˇeti urˇcen´e pouze pro ˇcten´ı. D´ıky tomu, ˇze je dopˇredu deklarov´ano, ˇze se jedn´a o konstantn´ı pamˇet’, m˚ uˇze b´ yt pˇr´ıstup k dat˚ um v t´eto pamˇeti optimalizov´an. Bˇeˇz´ıc´ı vl´ akna si nemohou ˇz´ adnou pamˇet’ alokovat. Veˇsker´e alokace mus´ı b´ yt deklarov´ any pˇred spuˇstˇen´ım kernelu a jednotliv´e instance mohou manipulovat pouze s pamˇet´ı, na kterou dostanou ukazatele v argumentech kernelu. Z toho vypl´ yv´a pomˇernˇe z´asadn´ı omezen´ı, ˇze jednotliv´e instance mus´ı vracet v´ ysledky konstantn´ı velikosti nebo je potˇreba velikost v´ ysledku nejprve dopˇredu spoˇc´ıtat. Pˇri pr´ aci s pamˇet´ı GPU mus´ıme m´ıt na pamˇeti jeˇstˇe dvˇe d˚ uleˇzit´a omezen´ı [1, 2]. Prvn´ı omezen´ı se t´ yk´ a naˇc´ıt´ an´ı dat z glob´aln´ı pamˇeti. Glob´aln´ı pamˇet’ je pˇr´ıstupn´a skrze relativnˇe ˇsirokou sbˇernici, avˇsak latence poˇzadavk˚ u je pomˇernˇe v´ yrazn´a. Z tohoto d˚ uvodu se optimalizuj´ı pˇrenosy, pˇri kter´ ych sousedn´ı vl´akna pˇren´aˇs´ı sousedn´ı bloky pamˇeti (tzv. coalesced load). Pokud alespoˇ n 16 sousedn´ıch vl´aken z jedn´e skupiny zaˇcne ve stejn´em okamˇziku ˇc´ıst nebo zapisovat 16 sousedn´ıch blok˚ u pamˇeti (o velikosti 32 bit˚ u) a cel´ y tento blok (64B) je spr´avnˇe zarovnan´ y, hardware zpracuje tento pˇresun v jedin´em poˇzadavku. Druh´e omezen´ı se t´ yk´ a pˇr´ıstupu do lok´aln´ı (sd´ılen´e) pamˇeti. Lok´aln´ı pamˇet’ je rozdˇelena do tzv. bank (na nejnovˇejˇs´ı architektur´ach od je tˇechto bank zpravidla 16). Dvˇe sousedn´ı 32-bit. buˇ nky pamˇeti se nach´az´ı v n´asleduj´ıc´ıch dvou bank´ach (modulo 16). Pokud dvˇe j´ adra v jedn´e instrukci pˇristoup´ı do stejn´e banky, je tento pˇr´ıstup serializov´ an, ˇc´ımˇz dojde ke zpomalen´ı v´ ypoˇctu. V´ yjimkou je, pokud vˇsechna vl´akna ˇctou stejnou hodnotu z lok´ aln´ı pamˇeti. Takov´ y pˇr´ıstup je detekov´an v hardware a data jsou k j´ adr˚ um pˇrenesena pomoc´ı broadcastu.
3
Jednoduch´ e SIMD v´ ypoˇ cty
Jak jsme zm´ınili v pˇredchoz´ı sekci, architektura grafick´ ych procesor˚ u je od z´akladu postavena na SIMD3 paralelizaci. Prvn´ı pˇr´ıklad proto otestuje v´ ykonnost grafick´e karty pˇri jednoduch´ ych vektorov´ ych v´ ypoˇctech, kter´e se nejl´epe zpracov´avaj´ı pr´avˇe na SIMD v´ ypoˇcetn´ım modelu. 3 Single
Instruction Multiple Data
4
V pˇr´ıkladu vyzkouˇs´ıme dvˇe vektorov´e operace. Obˇe maj´ı dvˇe vstupn´ı pole x a y a v´ ystupn´ı pole z. V´ ypoˇcet prob´ıh´ a nez´avisle nad vˇsemi prvky pole, tedy pro vˇsechna i z rozsahu pole provedeme z[i] = op(x[i], y[i]). Prvn´ı operace pouze vyn´asob´ı oba prvky. Druh´ a operace je o trochu sloˇzitˇejˇs´ı a pˇri v´ ypoˇctu uplatn´ı tak´e druhou odmocninu a goniometrick´e funkce: √ y[i] x[i] op2 (x[i], y[i]) = + x[i] cos(y[i]) x[i] Tato operace nem´ a ˇz´ adn´ y praktick´ y z´aklad. Byla pouze umˇele vytvoˇrena k otestov´ an´ı v´ıce druh˚ u operac´ı a sloˇzitˇejˇs´ıch matematick´ ych funkc´ı.
3.1
Namˇ eˇ ren´ e v´ ysledky
Pro u ´ˇcely n´asleduj´ıc´ıch test˚ u byly pouˇzity vektory ˇc´ısel x a y o velikosti 16M (t.j. 224 ) 32-bitov´ ych re´ aln´ ych ˇc´ısel (float˚ u). Kernely obou operac´ı byly spuˇstˇeny paralelnˇe pro kaˇzd´ y prvek pole, pˇriˇcemˇz velikost skupiny byla nastavena na 256 (coˇz je maximum na pouˇzit´em GPU). Tato hodnota byla vybr´ana na z´akladˇe empirick´ ych pozorov´an´ı. V tabulce 1 jsou shrnuty v´ ysledky mˇeˇren´ı. Jednotliv´e sloupce obsahuj´ı srovn´an´ı s´eriov´eho algoritmu, paraleln´ı implementace pouˇz´ıvaj´ıc´ı knihovnu Threading Building Blocks [6], implementace OpenCL bˇeˇz´ıc´ı na CPU a na GPU. Sloupec GPU obsahuje ˇcasy cel´eho v´ ypoˇctu vˇcetnˇe pˇresunu dat na grafickou kartu a zpˇet, zat´ımco GPU∗ obsahuje pouze ˇcasy samotn´eho v´ ypoˇctu.
op1 op2
s´eriov´ y 38 505
TBB 22 130
CPU 65 165
GPU 189 196
GPU∗ 15 21
Tabulka 1: Namˇeˇren´e ˇcasy SIMD v´ ypoˇct˚ u v milisekund´ach Z namˇeˇren´ ych v´ ysledk˚ u plyne, ˇze GPU zvl´adne v´ ypoˇcty prov´est pomˇernˇe rychle, avˇsak u ´zk´ ym hrdlem cel´e operace je pˇrenos dat z hlavn´ı pamˇeti do pamˇeti grafick´e karty a zpˇet. U vˇsech v´ ypoˇct˚ u prov´adˇen´ ych na GPU mus´ıme tedy nutnˇe zapoˇc´ıtat dobu potˇrebnou na pˇrenos dat, pˇr´ıpadnˇe pl´anovat v´ıce operac´ı za sebou, jinak se trivi´aln´ı vektorov´e operace nevyplat´ı paralelizovat na GPU, nebot’ procesor je zvl´adne prov´est ve srovnateln´em ˇcase.
4
N´ asoben´ı matic
Druh´ ym pˇr´ıkladem je ˇcasto pouˇz´ıvan´ y matematick´ y probl´em – n´asoben´ı matic. Pˇrestoˇze tato operace nen´ı pˇr´ıliˇs ˇcast´a v oblasti zpracov´an´ı dat, pom˚ uˇze n´am demonstrovat dalˇs´ı aspekty pˇrenos˚ u dat mezi pamˇetmi. Tentokr´at se zamˇeˇr´ıme na pˇrenosy dat mezi glob´ aln´ı pamˇet´ı grafick´e karty a lok´aln´ı (a priv´atn´ı) pamˇet´ı jednotliv´ ych jader. Pro jednoduchost se omez´ıme na n´asoben´ı dvou ˇctvercov´ ych matic, jejichˇz strany maj´ı d´elku mocniny 2. Rovnˇeˇz pro jednoduchost pˇredpokl´ad´ame, ˇze druh´a matice je jiˇz transponovan´ a, abychom l´epe vyuˇzili cache pˇri v´ ypoˇctu na CPU. 5
Pˇrestoˇze existuji lepˇs´ı algoritmy [4], v naˇsem pˇr´ıkladu pouˇzijeme klasick´ y algoritmus pracuj´ıc´ı v ˇcase O(N 3 ), kde N je d´elka strany matice.
4.1
Naivn´ı implementace
Kernel naivn´ı implementace obsahuje pouze nejv´ıce vnoˇren´ y cyklus. Zb´ yvaj´ıc´ı dvˇe u ´rovnˇe cykl˚ u jsou nahrazeny vytvoˇren´ım pˇr´ısluˇsn´eho poˇctu vl´aken. Kaˇzd´ y prvek v´ ysledn´e matice je tedy poˇc´ıt´ an paralelnˇe a v pˇr´ıpadˇe, ˇze bychom mˇeli k dispozici dostatek jader (ˇr´ adovˇe O(N 2 )), mohli bychom dos´ahnout asymptotick´e sloˇzitosti O(N ). __kernel void mul_matrix (__global const float *m1, __global const float *m2, __global float *mRes) { int n = get_global_size(0); int r = get_global_id(0); int c = get_global_id(1); float sum = 0; for (int i = 0; i < n; ++i) sum += m1[r*n + i] * m2[c*n + i]; mRes[r*n + c] = sum; } Pˇrestoˇze zvolen´ y algoritmus je velmi dobˇre paralelizovateln´ y a na prvn´ı pohled ˇ slibuje v´ yrazn´e zrychlen´ı, experiment´aln´ı v´ ysledky naznaˇcuj´ı opak. Casy uveden´e v tabulce 2 byly namˇeˇreny pˇri n´ asoben´ı matic 32-bitov´ ych re´aln´ ych ˇc´ısel pro N rovno 1024 a 2048. Vl´ akna uspoˇr´ ad´ ana do dvourozmˇern´eho pole, kter´e pˇresnˇe kop´ıruje tvar matice, a spojena do skupin velikosti 16 × 16. Matice 1024 × 1024 2048 × 2048
s´eriov´ y 3630 29370
TBB 542 4060
CPU 319 2311
GPU 392 3400
ˇ Tabulka 2: Casy n´asoben´ı matic v milisekund´ach V´ ypoˇcet na GPU, na kter´em spolupracovalo 320 jader, byl znatelnˇe pomalejˇs´ı neˇz v´ ypoˇcet na ˇctyˇrj´ adrov´em CPU s hyperthreadingem. V tomto pˇr´ıpadˇe ale neleˇz´ı probl´em v pˇrenosu dat z operaˇcn´ı pamˇeti do pamˇeti grafick´e karty, nebot’ tento pˇrenos zabere pouze 20ms v pˇr´ıpadˇe N = 1024, resp. 70ms v pˇr´ıpadˇe N = 2048. Ot´azku neuspokojiv´eho v´ ykonu naˇseho ˇreˇsen´ı zodpov´ı profilovac´ı n´astroj. V n´asleduj´ıc´ı tabulce jsou uvedeny nejd˚ uleˇzitˇejˇs´ı hodnoty namˇeˇren´e pˇri n´asoben´ı dvou matic velikosti 1024 × 1024.
6
Fetch (kolikr´ at ˇcetlo kaˇzd´e vl´akno z glob´aln´ı pamˇeti) ALUFetchRatio (pomˇer operac´ı ALU v˚ uˇci ˇcten´ı) ALUFetchBusy (pod´ıl operace ˇcten´ı na celkov´em ˇcasu) FetchUnitStalled (kolik ˇcasu ˇcekaly fetch jednotky na data)
2048× 2.51% 94.37% 83.15%
Z namˇeˇren´ ych u ´daj˚ u je jasn´e, ˇze vˇetˇsinu ˇcasu cel´eho v´ ypoˇctu se pouze pˇren´aˇsela data z hlavn´ı pamˇeti k j´adr˚ um, pˇriˇcemˇz se t´emˇeˇr v´ yhradnˇe na data ˇcekalo. Pˇri optimalizaci se tedy zamˇeˇr´ıme hlavnˇe na pˇrenosy dat v r´amci grafick´e karty.
4.2
Cache-aware implementace
Vylepˇsen´ a implementace vyuˇz´ıv´ a faktu, ˇze vl´akna jsou spouˇstˇena v SIMD reˇzimu po skupin´ ach. Kaˇzd´ a tato skupina m´a k dispozici lok´aln´ı pamˇet’4 , kterou m˚ uˇze vyuˇz´ıt jako cache pro data z glob´ aln´ı pamˇeti. Pˇri vhodn´e u ´pravˇe algoritmu m˚ uˇzeme doc´ılit v´ yrazn´eho zrychlen´ı, nebot’ • pˇr´ıstup do lok´ aln´ı pamˇeti je stejnˇe rychl´ y jako pˇr´ıstup do priv´atn´ı pamˇeti vl´akna, • vl´ akna mohou spolupracovat pˇri naˇc´ıt´an´ı dat (tzv. coalesced load) a • vˇetˇsina dat je sd´ılena mezi vl´akny, takˇze postaˇc´ı jedno naˇcten´ı tˇechto dat do lok´ aln´ı pamˇeti m´ısto aby si kaˇzd´e vl´akno tato data stahovalo zvl´aˇst’. Vl´ akna jsou uspoˇr´ ad´ ana do skupin velikosti 16 × 16. Pˇri zpracov´an´ı matice provedou vl´ akna v kaˇzd´em pr˚ uchodu hlavn´ım cyklem dva kroky. V prvn´ım kroku kooperativnˇe naˇctou ˇctvercov´e v´ yˇrezy velikosti 16 × 16 prvk˚ u z n´asoben´ ych matic. Tyto v´ yˇrezy obsahuj´ı prvky ˇr´ adk˚ u resp. sloupc˚ u n´asoben´ ych matic, takˇze vˇzdy jeden ˇr´adek resp. sloupec je sd´ılen 16 vl´ akny ve skupinˇe. V druh´em kroku si kaˇzd´e vl´akno pˇriˇcte n´asleduj´ıc´ıch 16 souˇcin˚ u prvk˚ u k ˇc´asteˇcn´emu souˇctu. Za prvn´ım i druh´ ym krokem n´asleduje synchronizaˇcn´ı bari´era, kter´a zajist´ı, ˇze vˇsechna vl´akna vykon´avaj´ı stejn´ y krok. Dalˇs´ı detaily jsou patrn´e z n´asleduj´ıc´ıho k´odu kernelu. __kernel void mul_matrix_opt(__global const float *m1, __global const float *m2, __global float *mRes, __local float *tmp1, __local float *tmp2) { int size = get_global_size(0); int lsize_x = get_local_size(0); int lsize_y = get_local_size(1); int block_size = lsize_x * lsize_y; int gid_x = get_global_id(0); int gid_y = get_global_id(1); int lid_x = get_local_id(0); int lid_y = get_local_id(1); 4V
pˇr´ıpadˇ e pouˇ zit´ e grafick´ e karty AMD Radeon HD5870 je to 32kB.
7
int offset = lid_y*lsize_x + lid_x; float sum = 0; for (int i = 0; i < size; i += lsize_x) { // Load data to local memory tmp1[offset] = m1[gid_y*size + i + lid_x]; for (int j = 0; j < lsize_x / lsize_y; ++j) tmp2[offset + j*block_size] = m2[(gid_x + lsize_y*j)*size + i + lid_x]; barrier(CLK_LOCAL_MEM_FENCE); // Add data from block to the sum for (int k = 0; k < lsize_x; ++k) sum += tmp1[lid_y*lsize_x + k] * tmp2[lid_x*lsize_x + k]; barrier(CLK_LOCAL_MEM_FENCE); } mRes[gid_y*size + gid_x] = sum; } V tabulce 3 jsou opˇet namˇeˇren´e v´ ysledky, vˇcetnˇe optimalizovan´e verze pro GPU, kter´a je oznaˇcena GPU+ . Toto ˇreˇsen´ı jiˇz vykazuje v´ yrazn´e zrychlen´ı – 45× resp. 52× proti s´eriov´e verzi a pˇribliˇznˇe 4× proti OpenCL verzi spuˇstˇen´e na vˇsech j´adrech CPU. Matice 1024 × 1024 2048 × 2048
s´eriov´ y 3630 29370
TBB 542 4060
CPU 319 2311
GPU 392 3400
GPU+ 81 564
ˇ Tabulka 3: Casy n´asoben´ı matic v milisekund´ach
5
Floyd-Warshall˚ uv algoritmus
Podobnou charakteristiku z hlediska ˇcasov´e sloˇzitosti a pˇr´ıstupu do pamˇeti jako mˇelo n´asoben´ı matic vykazuje tak´e zn´ am´ y algoritmus na hled´an´ı nejkratˇs´ıch cest v grafu. N´asoben´ı matic vˇsak mohlo paralelizovat vnˇejˇs´ı dva cykly, d´ıky ˇcemuˇz se vnitˇrn´ı cyklus vykon´ aval pˇr´ımo v kernelu. Floyd-Warshall˚ uv algoritmus proti tomu dok´aˇze paralelismu vyuˇz´ıt pouze u vnitˇrn´ıch dvou cykl˚ u, protoˇze po kaˇzd´em pr˚ uchodu vnˇejˇs´ım cyklem mus´ı doj´ıt k synchronizaci, abychom udrˇzeli integritu dat. OpenCL bohuˇzel neobsahuje glob´ aln´ı bari´eru pro vˇsechny instance jednoho kernelu a ani jej´ı implementace by nebyla pˇr´ıliˇs efektivn´ı. Proto pˇresuneme vnˇejˇs´ı cyklus mimo grafickou kartu a v kaˇzd´em jeho kroku spust´ıme paralelnˇe kernel ˇreˇs´ıc´ı vnitˇrn´ı dva cykly. Tabulka 4 shrnuje namˇeˇren´e v´ ysledky, pˇriˇcemˇz verze pro GPU jiˇz pouˇz´ıv´a obdobnou optimalizaci, jako algoritmus n´asoben´ı matic. Z v´ ysledk˚ u je patrn´e, ˇze cyklick´a invokace kernel˚ u je v´ yraznˇe m´enˇe efektivn´ı neˇz proveden´ı cyklu uvnitˇr kernelu. I pˇresto 8
Vrchol˚ u 1024 2048
s´eriov´ y 1262 10600
TBB 1006 7775
CPU 5200 35800
GPU 289 3285
ˇ Tabulka 4: Casy Floyd-Warshallova algoritmu v ms vˇsak pod´av´ a naˇse ˇreˇsen´ı v´ıce neˇz uspokojiv´e v´ ysledky a pouˇzit´ı GPU se i v tomto pˇr´ıpadˇe v´ıce neˇz vyplat´ı.
6
NP probl´ emy a backtracking
Dalˇs´ım uk´ azkov´ ym probl´emem je backtracking, kter´ y zde zastupuje exaktn´ı zp˚ usob ˇreˇsen´ı NP probl´em˚ u. NP probl´emy nen´ı5 v souˇcasn´e dobˇe moˇzn´e ˇreˇsit v polynomi´aln´ım ˇcase. Masivn´ı paralelismus je proto jednou z cest, jak alespoˇ n o trochu posunout hranice velikost´ı probl´em˚ u, kter´e je moˇzn´e v rozumn´e dobˇe spoˇc´ıtat. Jako reprezentanta jsme vybrali NP-´ uplnou u ´lohu zn´amou pod n´azvem Souˇcet podmnoˇziny, jej´ıˇz zad´an´ı je n´ asleduj´ıc´ı. Je d´ana mnoˇzina M cel´ ych ˇc´ısel a cel´e ˇc´ıslo s. Pt´ame se, zdali existuje vybran´ a podmnoˇzina M′ ⊂ M takov´a, ˇze souˇcet vˇsech jej´ıch prvk˚ u je pr´ avˇe s. Princip ˇreˇsen´ı je velice snadn´ y. Vybereme vˇsechny existuj´ıc´ı podmnoˇziny M (kter´ ych je 2|M| ) a u kaˇzd´e z nich ovˇeˇr´ıme, zda nem´a souˇcet s. Testy budeme prov´adˇet na mnoˇzinˇe velikosti |M| = 30. Kaˇzdou podmnoˇzinu identifikujeme 30bitov´ ym ˇc´ıslem, kde bity odpov´ıdaj´ı jednotliv´ ym prvk˚ um, pˇriˇcemˇz 1 znamen´a, ˇze je dan´ y prvek v podmnoˇzinˇe pˇr´ıtomen a 0, ˇze v n´ı pˇr´ıtomen nen´ı. Kaˇzd´ y kernel pak dostane prefix d´elky 24 bit˚ u, kter´e pouˇzije jako pevn´ y z´aklad a pro zb´ yvaj´ıc´ıch 8 bit˚ u vyzkouˇs´ı vˇsechny moˇznosti. Tabulka 5 shrnuje namˇeˇren´e v´ ysledky. |M| 30
s´eriov´ y 6625
TBB 1865
CPU 3820
GPU 595
ˇ Tabulka 5: Casy hled´ an´ı podmnoˇziny s dan´ ym souˇctem v ms Pˇri testov´ an´ı jsme si dovolili mal´e zjednoduˇsen´ı – vˇsechna ˇc´ısla v mnoˇzinˇe byla sud´a, zat´ım co hledan´ y souˇcet s byl lich´ y, abychom donutili algoritmus proj´ıt vˇsechny podmnoˇziny. Tento pˇr´ıstup ponech´av´a ve vzduchu ot´azku jak zastavit bˇeˇz´ıc´ı v´ ypoˇcet, kdyˇz jedno z vl´ aken nalezne ˇreˇsen´ı. Vzhledem k mnoˇzstv´ı vl´aken a nemoˇznosti efektivn´ı komunikace mezi nimi se jako nejlepˇs´ı zp˚ usob jev´ı spouˇstˇen´ı kernel˚ u po vhodnˇe velk´ ych ˇc´ astech tak, aby se pˇr´ıliˇs nezv´ yˇsila reˇzie, ale z´aroveˇ n aby se jednotliv´e ˇc´asti nepoˇc´ıtali pˇr´ıliˇs dlouho. Po skonˇcen´ı v´ ypoˇctu kaˇzd´e ˇc´asti se provede kontrola, zda byl jiˇz v´ ysledek nalezen, a pokud ano, v´ ypoˇcet skonˇc´ı. 5 Autor
se zde pˇrikl´ an´ı k vˇseobecnˇ e rozˇs´ıˇren´ e domnˇ ence, ˇ ze P̸=NP.
9
7
Z´ avˇ er
V tomto ˇcl´ anku jsme pˇredstavili relativnˇe novou architekturu na poli paraleln´ıch v´ ypoˇct˚ u a ovˇeˇrili jej´ı pouˇzitelnost na ˇradˇe r˚ uznorod´ ych probl´em˚ u. V´ ysledky naznaˇcuj´ı, ˇze se jedn´ a o slibnou technologii, avˇsak jej´ı pouˇzitelnost je do znaˇcn´e m´ıry omezov´ana nutnost´ı v nˇekter´ ych pˇr´ıpadech pomˇernˇe rozs´ahl´ ych u ´prav algoritm˚ u a programovac´ıch technik. Zat´ım co na CPU je program´ator ˇcasto zachr´anˇen pˇr´ıtomnost´ı velk´e cache, na GPU je potˇreba pˇr´ıstupy do pamˇeti peˇclivˇe pl´anovat a optimalizovat. V budouc´ı pr´ aci bychom se r´ adi zamˇeˇrili na vyuˇzit´ı GPU pˇri zpracov´an´ı datab´azov´ ych dotaz˚ u. Zejm´ena n´ as zaj´ımaj´ı moˇznosti nasazen´ı na relaˇcn´ı a s´emantick´a data (RDF [5], linked-data). V t´eto souvislosti bychom tak´e r´adi prozkoumali moˇznosti zpracov´ an´ı grafov´ ych algoritm˚ u, kter´e by mohli b´ yt vyuˇzity pˇri zpracov´an´ı index˚ u s´ıt’ov´ ych dat.
Literatura [1] ATI Stream Computing – OpenCL Programming Guide. http://developer. amd.com/gpu/ATIStreamSDK/assets/ATI_Stream_SDK_OpenCL_Programming_ Guide.pdf. [2] NVIDIA OpenCL Programming Guide. http://developer.download.nvidia. com/compute/cuda/3_0/toolkit/docs/NVIDIA_OpenCL_ProgrammingGuide. pdf. [3] OpenCL Framework Manual. http://www.khronos.org/opencl/. [4] R.P. Brent and STANFORD UNIV CALIF DEPT OF COMPUTER SCIENCE. Algorithms for matrix multiplication. Citeseer, 1970. [5] Frank Manola and Eric Miller. RDF Primer, W3C Recommendation, February 2004. http://www.w3.org/TR/2004/REC-rdf-primer-20040210/. [6] J. Reinders. Intel threading building blocks. O’Reilly, 2007.
10
