Obsah
1. Instrukční sady SIMD a automatické vektorizace prováděné překladačem GCC
2. První příklad: vynulování obsahu pole o předem známém počtu prvků
3. Výsledek překladu se zákazem i povolením SSE instrukcí: neoptimalizované varianty
4. Výsledek překladu se zákazem i povolením SSE instrukcí: optimalizované varianty a vliv zarovnání
5. Vliv změny velikosti nulovaného pole na kód vytvořený překladačem
6. Druhý příklad: přičtení konstantní hodnoty ke všem prvkům pole
7. Výsledek překladu se zákazem a povolením vektorizace
8. Vliv změny velikosti pole na kód vytvořený překladačem
9. Třetí příklad: přičtení obsahu jednoho pole k poli druhému
10. Výsledek překladu se zákazem a povolením vektorizace
11. Pomáháme překladači: klíčové slovo restrict
12. Kombinace obou nápověd: zarovnání i aliasu ukazatelů
13. Vliv změny velikosti pole na kód vytvořený překladačem
14. Výpočet druhé odmocniny všech prvků polí
15. Výsledek překladu se zákazem a povolením vektorizace
16. Vliv přepínače -ffast-math
17. Automatická vektorizace u operací typu reduce
18. Repositář s demonstračními příklady
19. Seznam všech předchozích částí tohoto seriálu a článků o SIMD instrukcích
1. Instrukční sady SIMD a automatické vektorizace prováděné překladačem GCC
Na stránkách Roota jsme se již několikrát setkali s instrukcemi typu SIMD, které dokážou provádět zvolené operace paralelně, například se čtyřmi prvky vektorů typu float nebo se šestnácti jednobajtovými prvky. Prozatím jsme si ukázali tři možnosti využití těchto instrukcí:
- Tzv. vektorové rozšíření GCC o další datové typyv16us atd. Základní operace s těmito typy (což jsou vektory) jsou totiž prováděny právě SIMD instrukcemi.
- Voláním takzvaných intrinsic zabudovaných přímo do překladače, které opět generují SIMD instrukce.
- Naprogramováním příslušných subrutin v assembleru
Ovšem existuje ještě jeden způsob nabízený moderními překladači (a to nejenom překladači jazyka C). Jedná se o takzvanou autovektorizaci, což je sada optimalizací, které dokážou například převést některé formy programových smyček do podoby se SIMD instrukcemi. Ovšem v mnoha případech je nutné překladači vhodným způsobem pomoci, například specifikací modifikátoru restrict, uvedením, že pole jsou zarovnána na 8, 16 atd. bajtů atd.
V dnešním článku si některé tyto techniky ukážeme. Přitom se (prozatím) zaměříme na překladač GCC a jeho podporu autovektorizace. Ta se povoluje přepínačem -ftree-vectorize. Navíc se dnes omezíme pouze na SSE instrukce (ty jsou totiž poměrně dobře pochopitelné) a pole/vektory typu float/single, takže použijeme i přepínač -msse.
2. První příklad: vynulování obsahu pole o předem známém počtu prvků
Pokusy s autovektorizací začneme na triviálním příkladu. Jedná se o funkci, která vynuluje obsah pole hodnot typu float (v SIMD označovaných jako single podle IEEE 754). Počet prvků nulovaného pole je pro jednoduchost dopředu známý a je uložen v konstantě nazvané SIZE. Zdrojový kód, který vůbec nebere v úvahu případné optimalizace, vypadá následovně:
void clear(float *a) {
#define SIZE 8
int i;
for (i=0; i<SIZE; i++) {
a[i] = 0.0;
}
}
3. Výsledek překladu se zákazem i povolením SSE instrukcí: neoptimalizované varianty
Jak bude vypadat překlad této jednoduché funkce do strojového kódu mikroprocesorů s architekturou x86 (resp. přesněji řečeno x86–64)? Výsledky se v tomto případě budou zásadně lišit podle toho, zda je povoleno používat instrukce SSE nebo je to naopak zakázáno.
V případě, že jsou optimalizace zakázány a současně je zakázáno i použití SSE, bude výsledný kód používat instrukce FPU a vše bude realizováno v neoptimalizované smyčce:
clear(float*):
push rbp
mov rbp, rsp
mov QWORD PTR [rbp-24], rdi
mov DWORD PTR [rbp-4], 0
jmp .L2
.L3:
mov eax, DWORD PTR [rbp-4]
cdqe
lea rdx, [0+rax*4]
mov rax, QWORD PTR [rbp-24]
add rax, rdx
fldz
fstp DWORD PTR [rax]
add DWORD PTR [rbp-4], 1
.L2:
cmp DWORD PTR [rbp-4], 7
jle .L3
nop
nop
pop rbp
ret
Při zákazu optimalizací, ale povolení SSE je opět výsledkem neoptimalizovaná smyčka, v níž ale nalezneme registr XMM0 a skalární instrukci MOVSS:
clear(float*):
push rbp
mov rbp, rsp
mov QWORD PTR [rbp-24], rdi
mov DWORD PTR [rbp-4], 0
jmp .L2
.L3:
mov eax, DWORD PTR [rbp-4]
cdqe
lea rdx, [0+rax*4]
mov rax, QWORD PTR [rbp-24]
add rax, rdx
pxor xmm0, xmm0
movss DWORD PTR [rax], xmm0
add DWORD PTR [rbp-4], 1
.L2:
cmp DWORD PTR [rbp-4], 7
jle .L3
nop
nop
pop rbp
ret
4. Výsledek překladu se zákazem i povolením SSE instrukcí: optimalizované varianty a vliv zarovnání
V dalším kroku se pokusme o povolení optimalizací, prozatím s využitím přepínače -O2. V případě, že jsou instrukce SSE zakázány, bude zápis (na platformě x86–64) prováděn po 64bitových slovech, tj. po dvojici hodnot typu float. To je již mnohem rychlejší a není zapotřebí realizovat pomalou počítanou programovou smyčku:
clear(float*):
mov QWORD PTR [rdi], 0
mov QWORD PTR [rdi+8], 0
mov QWORD PTR [rdi+16], 0
mov QWORD PTR [rdi+24], 0
ret
Při povolení instrukcí SSE je kód odlišný. Nejdříve se smaže obsah 128bitového registru XMM0 a následně se provedou pouze dva zápisy, z nichž každý uloží do paměti čtyři hodnoty typu float:
clear(float*):
pxor xmm0, xmm0
movups XMMWORD PTR [rdi], xmm0
movups XMMWORD PTR [rdi+16], xmm0
ret
Překladač jazyka C v tomto případě neví, jak budou prvky pole (ne)zarovnány a proto musí použít instrukci MOVUPS a nikoli MOVAPS. Ovšem pokud existuje jistota, že pole začíná na zarovnané adrese, je možné překladači „napovědět“:
void clear(float *a) {
#define SIZE 8
a = __builtin_assume_aligned (a, 16);
int i;
for (i=0; i<SIZE; i++) {
a[i] = 0.0;
}
}
Výsledkem nyní bude obecně rychlejší kód s instrukcemi MOVAPS namísto MOVUPS:
clear(float*):
pxor xmm0, xmm0
movaps XMMWORD PTR [rdi], xmm0
movaps XMMWORD PTR [rdi+16], xmm0
ret
5. Vliv změny velikosti nulovaného pole na kód vytvořený překladačem
Předchozí zdrojový kód byl pro překladač snadným úkolem, protože délka nulovaného pole byla rovna osmi prvkům. Co se však stane ve chvíli, kdy tuto délku pozměníme? Nejprve nastavíme SIZE na hodnotu 16 a provedeme překlad s využitím SSE instrukcí. Výsledkem je stále „rozbalená“ smyčka
clear(float*):
pxor xmm0, xmm0
movups XMMWORD PTR [rdi], xmm0
movups XMMWORD PTR [rdi+16], xmm0
movups XMMWORD PTR [rdi+32], xmm0
movups XMMWORD PTR [rdi+48], xmm0
ret
Pro SIZE nastavené na hodnotu 18 již není možné (alespoň ne jednoduše) všechny prvky zapsat instrukcemi MOVUPS nebo MOVAPS, protože 18 není dělitelná čtyřmi, Překladač si poradí tak, že 16 prvků zapíše po čtveřicích a další dva prvky jako jediné čtyřslovo:
clear(float*):
pxor xmm0, xmm0
mov QWORD PTR [rdi+64], 0
movups XMMWORD PTR [rdi], xmm0
movups XMMWORD PTR [rdi+16], xmm0
movups XMMWORD PTR [rdi+32], xmm0
movups XMMWORD PTR [rdi+48], xmm0
ret
A pro SIZE nastavené na 17 prvků dopadne výsledek takto (poslední prvek je zapsán jako dvojslovo):
clear(float*):
pxor xmm0, xmm0
mov DWORD PTR [rdi+64], 0
movaps XMMWORD PTR [rdi], xmm0
movaps XMMWORD PTR [rdi+16], xmm0
movaps XMMWORD PTR [rdi+32], xmm0
movaps XMMWORD PTR [rdi+48], xmm0
ret
6. Druhý příklad: přičtení konstantní hodnoty ke všem prvkům pole
V dnešním druhém ukázkovém příkladu je realizována funkce, která ke všem prvkům pole o známé délce přičítá konstantu předanou formou argumentu. Celý příklad je naprogramován naivním způsobem, bez jakékoli snahy o optimalizace na úrovni céčkovského zdrojového kódu:
void add_delta(float *a, float delta) {
#define SIZE 16
int i;
for (i=0; i<SIZE; i++) {
a[i] += delta;
}
}
7. Výsledek překladu se zákazem a povolením vektorizace
Nejprve se podívejme na výsledek překladu v situaci, kdy není automatická vektorizace povolena. V tomto případě se sice používají SSE instrukce a registry SSE (XMM0, XMM1), ovšem realizace přičtení konstanty používá programovou smyčku se skalárními instrukcemi ADDSS a MOVSS (ty dobře známe):
add_delta(float*, float):
lea rax, [rdi+64]
.L2:
movss xmm1, DWORD PTR [rdi]
add rdi, 4
addss xmm1, xmm0
movss DWORD PTR [rdi-4], xmm1
cmp rdi, rax
jne .L2
ret
Vektorizovaná varianta je opět realizována programovou smyčkou, ovšem nyní se v každé iteraci provede přičtení konstanty ke čtyřem prvkům pole a počet opakování smyčky je čtvrtinový. Díky tomu, že je délka pole dělitelná čtyřmi, je tato optimalizace zcela legální. Za povšimnutí stojí především použití instrukce SHUFPS, v níž je zdrojový i cílový registr totožný. Díky tomu se nejnižší prvek registru XMM0 rozkopíruje do ostatních tří prvků, což je trik, který je užitečné znát:
add_delta(float*, float):
lea rax, [rdi+64]
shufps xmm0, xmm0, 0
.L2:
movups xmm1, XMMWORD PTR [rdi]
add rdi, 16
addps xmm1, xmm0
movups XMMWORD PTR [rdi-16], xmm1
cmp rax, rdi
jne .L2
ret
8. Vliv změny velikosti pole na kód vytvořený překladačem
A jak bude vypadat výsledný kód při změně velikosti pole? Nejprve si ukažme výsledek překladu pro SIZE nastavenou na hodnotu 24 (což je stále dělitelné čtyřmi). Opět se použije plně optimalizovaná varianta (nerozbalené) programové smyčky se součtem vektorů:
add_delta(float*, float):
lea rax, [rdi+96]
shufps xmm0, xmm0, 0
.L2:
movups xmm1, XMMWORD PTR [rdi]
add rdi, 16
addps xmm1, xmm0
movups XMMWORD PTR [rdi-16], xmm1
cmp rax, rdi
jne .L2
ret
Zajímavá situace nastane pro SIZE nastavenou na hodnotu 17, což v žádném případě není dělitelné čtyřmi. V takovém případě se překladač bude snažit o rozdělení zápisu s využitím součtu vektorů ve smyčce. A zbývající jeden až tři součty budou provedeny až za touto smyčkou:
add_delta(float*, float):
movaps xmm2, xmm0
mov rax, rdi
lea rdx, [rdi+64]
shufps xmm2, xmm2, 0
.L2:
movups xmm1, XMMWORD PTR [rax]
add rax, 16
addps xmm1, xmm2
movups XMMWORD PTR [rax-16], xmm1
cmp rdx, rax
jne .L2
addss xmm0, DWORD PTR [rdi+64]
movss DWORD PTR [rdi+64], xmm0
ret
Pro SIZE nastavené na hodnotu 18 se po dokončení smyčky provede součet dvou prvků s využitím součtu vektorů:
add_delta(float*, float):
movaps xmm2, xmm0
mov rax, rdi
lea rdx, [rdi+64]
shufps xmm2, xmm2, 0
.L2:
movups xmm1, XMMWORD PTR [rax]
add rax, 16
addps xmm1, xmm2
movups XMMWORD PTR [rax-16], xmm1
cmp rdx, rax
jne .L2
movq xmm1, QWORD PTR [rdi+64]
shufps xmm0, xmm0, 0xe0
movq xmm0, xmm0
addps xmm0, xmm1
movlps QWORD PTR [rdi+64], xmm0
ret
Dtto pro velikost 19:
add_delta(float*, float):
movaps xmm2, xmm0
mov rax, rdi
lea rdx, [rdi+64]
shufps xmm2, xmm2, 0
.L2:
movups xmm1, XMMWORD PTR [rax]
add rax, 16
addps xmm1, xmm2
movups XMMWORD PTR [rax-16], xmm1
cmp rdx, rax
jne .L2
movq xmm1, QWORD PTR [rdi+64]
movaps xmm2, xmm0
addss xmm0, DWORD PTR [rdi+72]
shufps xmm2, xmm2, 0xe0
movq xmm2, xmm2
addps xmm1, xmm2
movss DWORD PTR [rdi+72], xmm0
movlps QWORD PTR [rdi+64], xmm1
ret
9. Třetí příklad: přičtení obsahu jednoho pole k poli druhému
Třetí příklad, který se pokusíme přeložit do strojového kódu ve vektorizované podobě, provádí přičtení obsahu jednoho pole typu float ke druhému poli typu float. Sčítají se odpovídající si prvky pole, tedy nultý s nultým, první s prvním atd. Opět si nejdříve ukažme naivní variantu zápisu tohoto algoritmu, ve které se nijak nesnažíme překladači pomoci přidáním dalších informací:
void add_arrays(float *a, float *b) {
#define SIZE 4
int i;
for (i=0; i<SIZE; i++) {
a[i] += b[i];
}
}
10. Výsledek překladu se zákazem a povolením vektorizace
Překlad se zákazem automatické vektorizace vede k vytvoření programové smyčky, ve které se postupně provádí skalární součet s využitím instrukce addss. Jedná se sice o poměrně primitivní způsob překladu, ovšem jeho výhodou je, že není zapotřebí zjišťovat překryvy polí atd.:
add_arrays(float*, float*):
xor eax, eax
.L2:
movss xmm0, DWORD PTR [rdi+rax]
addss xmm0, DWORD PTR [rsi+rax]
movss DWORD PTR [rdi+rax], xmm0
add rax, 4
cmp rax, 96
jne .L2
ret
Při povolení vektorizace se vygeneruje poměrně dlouhý programový kód:
add_arrays(float*, float*):
lea rdx, [rsi+4]
mov rax, rdi
sub rax, rdx
cmp rax, 8
mov eax, 0
jbe .L2
.L3:
movups xmm0, XMMWORD PTR [rdi+rax]
movups xmm1, XMMWORD PTR [rsi+rax]
addps xmm0, xmm1
movups XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 96
jne .L3
ret
.L2:
movss xmm0, DWORD PTR [rdi+rax]
addss xmm0, DWORD PTR [rsi+rax]
movss DWORD PTR [rdi+rax], xmm0
add rax, 4
cmp rax, 96
jne .L2
ret
Ovšem proč tomu tak je? Překladač si nemůže být jistý, že se předaná pole nepřekrývají. Musí tedy provést kontrolu, jestli nedochází k překryvu (o čtyři prvky) a pokud ano, provede se neoptimalizovaná „skalární“ část algoritmu. Tyto části jsou dobře viditelné. „Vektorová část“ provádějící součty po čtveřicích:
.L3:
movups xmm0, XMMWORD PTR [rdi+rax]
movups xmm1, XMMWORD PTR [rsi+rax]
addps xmm0, xmm1
movups XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 96
jne .L3
ret
Skalární část:
.L2:
movss xmm0, DWORD PTR [rdi+rax]
addss xmm0, DWORD PTR [rsi+rax]
movss DWORD PTR [rdi+rax], xmm0
add rax, 4
cmp rax, 96
jne .L2
ret
11. Pomáháme překladači: klíčové slovo restrict
V případě, že sčítáme dvě rozdílná pole, by bylo vhodné překladači nějakým způsobem napovědět, že nemusí do výsledného kódu vkládat i variantu s překrývajícími se poli. Pokud by byla pole odlišného typu (resp. pokud by ukazatele byly odlišného typu), je řešení snadné – podle normy je chování překladače nedefinované :-). Ovšem naše pole mají shodný typ (přesněji řečeno oba ukazatele mají totožný typ float *), takže překladač musí předpokládat, že se může jednat o stejné ukazatele nebo že se budou oblasti překrývat. Jak překladači napovědět, že k této situaci nedojde a nemusí se jí zabývat? V jazyku C, konkrétně ve verzi C99 a vyšší, je pro tento účel rezervováno slovo restrict, které se používá následujícím způsobem:
void add_arrays(float *restrict a, float * restrict b) {
#define SIZE 4
int i;
for (i=0; i&tl;SIZE; i++) {
a[i] += b[i];
}
}
Překladač v tomto případě může provést optimalizace předpokládající, že se sčítají odlišná pole a výsledek bude mnohem kratší:
add_arrays(float* restrict, float* restrict):
movups xmm1, XMMWORD PTR [rdi]
movups xmm0, XMMWORD PTR [rsi]
addps xmm0, xmm1
movups XMMWORD PTR [rdi], xmm0
ret
12. Kombinace obou nápověd: zarovnání i aliasu ukazatelů
Samozřejmě můžeme zkombinovat obě nápovědy, které předáváme překladači, tj. informaci o zarovnání polí i informaci o tom, že ukazatele obsahují adresu rozdílných polí:
void add_arrays(float *restrict a, float * restrict b) {
a = __builtin_assume_aligned (a, 16);
b = __builtin_assume_aligned (b, 16);
#define SIZE 4
int i;
for (i=0; i<SIZE; i++) {
a[i] += b[i];
}
}
Výsledek nyní obsahuje instrukci MOVAPS pro čtení a zápis do „zarovnané“ oblasti paměti. Nejenom to – výsledek je kratší o jednu instrukci, protože ADDPS může jako druhý operand mít zarovnanou adresu a nikoli pouze XMM registr:
add_arrays(float* restrict, float* restrict):
movaps xmm0, XMMWORD PTR [rdi]
addps xmm0, XMMWORD PTR [rsi]
movaps XMMWORD PTR [rdi], xmm0
ret
13. Vliv změny velikosti pole na kód vytvořený překladačem
Pro delší pole (SIZE je nastaveno na větší hodnotu) s počtem prvků dělitelných šestnácti můžeme při povolení automatických vektorizací a současně i zajištění zarovnání polí atd. získat tento výsledek se smyčkou, ve které se vždy sčítají čtveřice prvků vektorů:
add_arrays(float* restrict, float* restrict):
xor eax, eax
.L2:
movaps xmm0, XMMWORD PTR [rdi+rax]
addps xmm0, XMMWORD PTR [rsi+rax]
movaps XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 64
jne .L2
ret
Pokud počet prvků není dělitelný čtyřmi, například je nastaven na 17, musí se do kódu vložit i instrukce pro dokončení výpočtů pro poslední 1 až 3 prvky. To je patrné z následujícího kódu pro SIZE=17:
add_arrays(float*, float*):
xor eax, eax
.L2:
movaps xmm0, XMMWORD PTR [rdi+rax]
addps xmm0, XMMWORD PTR [rsi+rax]
movaps XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 64
jne .L2
movss xmm0, DWORD PTR [rdi+64]
addss xmm0, DWORD PTR [rsi+64]
movss DWORD PTR [rdi+64], xmm0
ret
Jen pro ilustraci doplňme, jak dlouhý kód musí být vygenerován ve chvíli, kdy překladač nemá informace o zarovnání ani o aliasu ukazatelů:
add_arrays(float*, float*):
lea rdx, [rsi+4]
mov rax, rdi
sub rax, rdx
cmp rax, 8
mov eax, 0
jbe .L2
.L3:
movups xmm0, XMMWORD PTR [rdi+rax]
movups xmm1, XMMWORD PTR [rsi+rax]
addps xmm0, xmm1
movups XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 64
jne .L3
movss xmm0, DWORD PTR [rdi+64]
addss xmm0, DWORD PTR [rsi+64]
movss DWORD PTR [rdi+64], xmm0
ret
.L2:
movss xmm0, DWORD PTR [rdi+rax]
addss xmm0, DWORD PTR [rsi+rax]
movss DWORD PTR [rdi+rax], xmm0
add rax, 4
cmp rax, 68
jne .L2
ret
14. Výpočet druhé odmocniny všech prvků polí
Automatická vektorizace kódu může proběhnout i ve chvíli, kdy se pokusíme o překlad funkce, ve které se pro všechny prvky pole počítá například jejich převrácená hodnota nebo druhá odmocnina. Je tomu tak z toho důvodu, že v instrukční sadě SSE existují příslušné SIMD instrukce SQRTSS/SQRTPS, RSQRTSS/RSQRTPS a RCPSS/RCPPS určené pro přesný výpočet či o „odhad“ výsledku s předem známou přesností. Pokusme se tedy přeložit následující funkci se snadno pochopitelným algoritmem:
#include <math.h>
void array_sqrt(float *a) {
#define SIZE 24
int i;
for (i=0; i<SIZE; i++) {
a[i] = sqrt(a[i]);
}
}
15. Výsledek překladu se zákazem a povolením vektorizace
Výsledek překladu bude značně záviset na použitých přepínačích. Při zákazu vektorizace dostaneme tento výsledek:
array_sqrt(float*):
push rbx
pxor xmm1, xmm1
lea rbx, [rdi+96]
sub rsp, 16
.L5:
movss xmm0, DWORD PTR [rdi]
ucomiss xmm1, xmm0
ja .L8
sqrtss xmm0, xmm0
.L4:
movss DWORD PTR [rdi], xmm0
add rdi, 4
cmp rdi, rbx
jne .L5
add rsp, 16
pop rbx
ret
.L8:
mov QWORD PTR [rsp+8], rdi
call sqrtf
mov rdi, QWORD PTR [rsp+8]
pxor xmm1, xmm1
jmp .L4
Vygenerovaný kód je tak dlouhý z toho důvodu, že se v případě výpočtu druhé odmocniny (například) ze záporných čísel musí nastavit errno atd. Pokud tuto hodnotu (tedy test, jak dopadly předchozí FPU operace) stejně v našem algoritmu nepoužijeme, lze překladači přes přepínač -fno-math-errno sdělit, že případné chyby nemusí zaznamenávat. Výsledný kód bude kratší a čitelnější:
array_sqrt(float*):
lea rax, [rdi+96]
.L2:
movss xmm0, DWORD PTR [rdi]
add rdi, 4
sqrtss xmm0, xmm0
movss DWORD PTR [rdi-4], xmm0
cmp rdi, rax
jne .L2
ret
Povolením vektorizace získáme ještě optimálnější variantu, ve které se provádí výpočty druhé odmocniny vždy pro čtyři hodnoty současně:
array_sqrt(float*):
lea rax, [rdi+96]
.L2:
movups xmm1, XMMWORD PTR [rdi]
add rdi, 16
sqrtps xmm0, xmm1
movups XMMWORD PTR [rdi-16], xmm0
cmp rax, rdi
jne .L2
ret
16. Vliv přepínače -ffast-math
Mohlo by se zdát, že v předchozí kapitole uvedený vektorizovaný kód již není možné žádným způsobem vylepšit. Ovšem v případě, že jsou vyžadovány skutečně velmi rychlé výpočty, může být užitečné použít přepínač -ffast-math, který překladači umožní vygenerovat instrukce, které nebudou provádět výpočty přesně podle normy IEEE 754. Výsledkem bude v našem konkrétním případě kód, ve kterém se nejprve provede přibližný (ale o to rychlejší) výpočet převrácené hodnoty druhé odmocniny a poté se z této hodnoty (a s využitím dvou vhodně zvolených konstant) dopočítají druhé odmocniny:
array_sqrt(float*):
movss xmm5, DWORD PTR .LC1[rip]
movss xmm4, DWORD PTR .LC3[rip]
lea rax, [rdi+96]
pxor xmm3, xmm3
shufps xmm5, xmm5, 0
shufps xmm4, xmm4, 0
.L2:
movups xmm1, XMMWORD PTR [rdi]
movaps xmm2, xmm3
add rdi, 16
rsqrtps xmm0, xmm1
cmpneqps xmm2, xmm1
andps xmm0, xmm2
mulps xmm1, xmm0
mulps xmm0, xmm1
mulps xmm1, xmm4
addps xmm0, xmm5
mulps xmm0, xmm1
movups XMMWORD PTR [rdi-16], xmm0
cmp rax, rdi
jne .L2
ret
.LC1:
.long -1069547520
.LC3:
.long -1090519040
17. Automatická vektorizace u operací typu reduce
Většina demonstračních příkladů, které jsme si až doposud uvedli, vektorizovala programový kód zapsaný s využitím počítané programové smyčky, ve které se postupně prochází (a popř. modifikují) prvky jednoho pole nebo dvou polí. Ovšem současné verze překladačů dokážou vektorizovat i takové programové smyčky, v nichž se například prvky sčítají, počítá se skalární součin dvou vektorů, vyhledává se nejnižší nebo nejvyšší hodnota atd. Takové operace se nazývají operace typu reduce a se způsobem jejich vektorizace se seznámíme příště.
18. Repositář s demonstračními příklady
Demonstrační příklady naprogramované v jazyku, které jsou určené pro překlad s využitím assembleru gcc, byly uloženy do Git repositáře, který je dostupný na adrese https://github.com/tisnik/8bit-fame. Kromě zdrojových kódů příkladů jsou do repositáře přidány i výsledky překladu do assembleru v syntaxi kompatibilní s Intelem. Jednotlivé demonstrační příklady si můžete v případě potřeby stáhnout i jednotlivě bez nutnosti klonovat celý (dnes již poměrně rozsáhlý) repositář:
19. Seznam všech předchozích částí tohoto seriálu a článků o SIMD instrukcích
Podporou SIMD instrukcí na úrovni intrinsic jsme se už na Rootu zabývali, stejně jako samotnými SIMD instrukcemi na úrovni assembleru. Pro úplnost jsou v této příloze uvedeny odkazy na příslušné články:
- Užitečné rozšíření GCC: podpora SIMD (vektorových) instrukcí
https://www.root.cz/clanky/uzitecne-rozsireni-gcc-podpora-simd-vektorovych-instrukci/ - Užitečné rozšíření GCC – podpora SIMD (vektorových) instrukcí: nedostatky technologie
https://www.root.cz/clanky/uzitecne-rozsireni-gcc-podpora-simd-vektorovych-instrukci-nedostatky-technologie/ - Podpora SIMD (vektorových) instrukcí na RISCových procesorech
https://www.root.cz/clanky/podpora-simd-vektorovych-instrukci-na-riscovych-procesorech/ - Podpora SIMD operací v GCC s využitím intrinsic pro nízkoúrovňové optimalizace
https://www.root.cz/clanky/podpora-simd-operaci-v-gcc-s-vyuzitim-intrinsic-pro-nizkourovnove-optimalizace/ - Podpora SIMD operací v GCC s využitím intrinsic: technologie SSE
https://www.root.cz/clanky/podpora-simd-operaci-v-gcc-s-vyuzitim-intrinsic-technologie-sse/ - Rozšíření instrukční sady „Advanced Vector Extensions“ na platformě x86–64
https://www.root.cz/clanky/rozsireni-instrukcni-sady-advanced-vector-extensions-na-platforme-x86–64/ - Rozšíření instrukční sady F16C, FMA a AVX-512 na platformě x86–64
https://www.root.cz/clanky/rozsireni-instrukcni-sady-f16c-fma-a-avx-512-na-platforme-x86–64/ - Rozšíření instrukční sady AVX-512 na platformě x86–64 (dokončení)
https://www.root.cz/clanky/rozsireni-instrukcni-sady-avx-512-na-platforme-x86–64-dokonceni/ - SIMD instrukce na platformě 80×86: instrukční sada MMX
https://www.root.cz/clanky/simd-instrukce-na-platforme-80×86-instrukcni-sada-mmx/ - SIMD instrukce na 80×86: dokončení popisu MMX, instrukce 3DNow!
https://www.root.cz/clanky/simd-instrukce-na-80–86-dokonceni-popisu-mmx-instrukce-3dnow/ - SIMD instrukce v rozšíření SSE
https://www.root.cz/clanky/simd-instrukce-v-rozsireni-sse/ - SIMD instrukce v rozšíření SSE (2. část)
https://www.root.cz/clanky/simd-instrukce-v-rozsireni-sse-2-cast/ - Pokročilejší SSE operace: přeskupení, promíchání a rozbalování prvků vektorů
https://www.root.cz/clanky/pokrocilejsi-sse-operace-preskupeni-promichani-a-rozbalovani-prvku-vektoru/ - Od instrukční sady SSE k sadě SSE2
https://www.root.cz/clanky/od-instrukcni-sady-sse-k-sade-sse2/
20. Odkazy na Internetu
- Auto-vectorization in GCC
https://gcc.gnu.org/projects/tree-ssa/vectorization.html - GCC documentation: Extensions to the C Language Family
https://gcc.gnu.org/onlinedocs/gcc/C-Extensions.html#C-Extensions - GCC documentation: Using Vector Instructions through Built-in Functions
https://gcc.gnu.org/onlinedocs/gcc/Vector-Extensions.html - SSE (Streaming SIMD Extentions)
http://www.songho.ca/misc/sse/sse.html - Timothy A. Chagnon: SSE and SSE2
http://www.cs.drexel.edu/~tc365/mpi-wht/sse.pdf - Intel corporation: Extending the Worldr's Most Popular Processor Architecture
http://download.intel.com/technology/architecture/new-instructions-paper.pdf - SIMD architectures:
http://arstechnica.com/old/content/2000/03/simd.ars/ - Tour of the Black Holes of Computing!: Floating Point
http://www.cs.hmc.edu/~geoff/classes/hmc.cs105…/slides/class02_floats.ppt - 3Dnow! Technology Manual
AMD Inc., 2000 - Intel MMXTM Technology Overview
Intel corporation, 1996 - MultiMedia eXtensions
http://softpixel.com/~cwright/programming/simd/mmx.phpi - AMD K5 („K5“ / „5k86“)
http://www.pcguide.com/ref/cpu/fam/g5K5-c.html - Sixth Generation Processors
http://www.pcguide.com/ref/cpu/fam/g6.htm - Great Microprocessors of the Past and Present
http://www.cpushack.com/CPU/cpu1.html - Very long instruction word (Wikipedia)
http://en.wikipedia.org/wiki/Very_long_instruction_word - CPU design (Wikipedia)
http://en.wikipedia.org/wiki/CPU_design - Bulldozer (microarchitecture)
https://en.wikipedia.org/wiki/Bulldozer_(microarchitecture) - SIMD Instructions Considered Harmful
https://www.sigarch.org/simd-instructions-considered-harmful/ - GCC Compiler Intrinsics
https://iq.opengenus.org/gcc-compiler-intrinsics/ - Scalable_Vector_Extension_(SVE)
https://en.wikipedia.org/wiki/AArch64#Scalable_Vector_Extension_(SVE) - Improve the Multimedia User Experience
https://www.arm.com/technologies/neon - NEON Technology (stránky ARM)
https://developer.arm.com/technologies/neon - SIMD Assembly Tutorial: ARM NEON – Xiph.org
https://people.xiph.org/~tterribe/daala/neon_tutorial.pdf - Ne10
http://projectne10.github.io/Ne10/ - NEON and Floating-Point architecture
http://infocenter.arm.com/help/index.jsp?topic=/com.arm.doc.den0024a/BABIGHEB.html - An Introduction to ARM NEON
http://peterdn.com/post/an-introduction-to-ARM-NEON.aspx - ARM NEON Intrinsics Reference
http://infocenter.arm.com/help/topic/com.arm.doc.ihi0073a/IHI0073A_arm_neon_intrinsics_ref.pdf - Arm Neon Intrinsics vs hand assembly
https://stackoverflow.com/questions/9828567/arm-neon-intrinsics-vs-hand-assembly - ARM NEON Optimization. An Example
http://hilbert-space.de/?p=22 - AArch64 NEON instruction format
https://developer.arm.com/docs/den0024/latest/7-aarch64-floating-point-and-neon/73-aarch64-neon-instruction-format - ARM SIMD instructions
https://developer.arm.com/documentation/dht0002/a/Introducing-NEON/What-is-SIMD-/ARM-SIMD-instructions - Learn the architecture – Migrate Neon to SVE Version 1.0
https://developer.arm.com/documentation/102131/0100/?lang=en - 1.2.2. Comparison between NEON technology and other SIMD solutions
https://developer.arm.com/documentation/den0018/a/Introduction/Comparison-between-ARM-NEON-technology-and-other-implementations/Comparison-between-NEON-technology-and-other-SIMD-solutions?lang=en - NEON Programmer’s Guide
https://documentation-service.arm.com/static/63299276e68c6809a6b41308 - Brain Floating Point – nový formát uložení čísel pro strojové učení a chytrá čidla
https://www.root.cz/clanky/brain-floating-point-ndash-novy-format-ulozeni-cisel-pro-strojove-uceni-a-chytra-cidla/ - Other Built-in Functions Provided by GCC
https://gcc.gnu.org/onlinedocs/gcc/Other-Builtins.html - GCC: 6.60 Built-in Functions Specific to Particular Target Machines
https://gcc.gnu.org/onlinedocs/gcc/Target-Builtins.html#Target-Builtins - Advanced Vector Extensions
https://en.wikipedia.org/wiki/Advanced_Vector_Extensions - Top 10 Craziest Assembly Language Instructions
https://www.youtube.com/watch?v=Wz_xJPN7lAY - Intel x86: let's take a look at one of the most complex instruction set!
https://www.youtube.com/watch?v=KBLy23B38-c - x64 Assembly Tutorial 58: Intro to AVX
https://www.youtube.com/watch?v=yAvuHd8cBJY - AVX512 (1 of 3): Introduction and Overview
https://www.youtube.com/watch?v=D-mM6X5×nTY - AVX512 (2 of 3): Programming AVX512 in 3 Different Ways
https://www.youtube.com/watch?v=I3efQKLgsjM - AVX512 (3 of 3): Deep Dive into AVX512 Mechanisms
https://www.youtube.com/watch?v=543a1b-cPmU - AVX-512
https://en.wikipedia.org/wiki/AVX-512 - AVX-512
https://iq.opengenus.org/avx512/ - SIMD Algorithms Youtube course
https://denisyaroshevskiy.github.io/presentations/ - Compiler explorer
https://godbolt.org/ - Restricting pointers
https://gcc.gnu.org/onlinedocs/gcc/Restricted-Pointers.html - Does the restrict keyword provide significant benefits in gcc/g++
https://stackoverflow.com/questions/1965487/does-the-restrict-keyword-provide-significant-benefits-in-gcc-g - Demystifying The Restrict Keyword
https://cellperformance.beyond3d.com/articles/2006/05/demystifying-the-restrict-keyword.html - Basics of Vectorization for Fortran Applications
https://inria.hal.science/hal-01688488/document
