Szablony wyrażeń zamiast kodu odręcznego

Question

Obecnie piszę bibliotekę wyrażeń szablonów C++ i porównuję niektóre wystąpienia z ręcznie pisanym kodem na poziomie zespołu. Funkcja ręcznie napisany jest następujący:

spinor multiply(vector const& a, vector const& b) 
{ 
     spinor result = { 
       a.at<1>() * b.at<1>() - a.at<2>() * b.at<2>() 
          - a.at<4>() * b.at<4>() - a.at<8>() * b.at<8>(), 
       a.at<1>() * b.at<2>() - a.at<2>() * b.at<1>(), 
       a.at<1>() * b.at<4>() - a.at<4>() * b.at<1>(), 
       a.at<1>() * b.at<8>() - a.at<8>() * b.at<1>(), 
       a.at<2>() * b.at<4>() - a.at<4>() * b.at<2>(), 
       a.at<2>() * b.at<8>() - a.at<8>() * b.at<2>(), 
       a.at<4>() * b.at<8>() - a.at<8>() * b.at<4>() 
     }; 

     return result; 
} 

Klasa vector tylko wrapper ponad cztery dwuosobowe, które można odczytać za pomocą funkcji at<index>() użytkownika. Ze względu na decyzje projektowe indeksy dla czterech składników to 1, 2, 4, 8, które są dostępne pod numerem at<index>() zamiast zwykłego 0, 1, 2, 3.

Celem tej funkcji jest zwrócenie wyniku mnożenia dwóch wektorów (w przestrzeni Minkowskiego). Jeśli znasz geometrię algebry, zobaczysz produkt kropkowy (pierwszy komponent result, symetryczny pod wymianą a i b) i produkt klinowy (reszta składników, antysymetryczna w wymianie a i b). Jeśli nie jesteś zaznajomiony z Geometryczną Algebrą, po prostu weź tę funkcję jako receptę na pomnożenie wektorów.

Gdybym skompilować funkcję powyżej GCC 4.7 i spojrzeć na demontażu danego przez objdump -SC a.out to daje mi następujący wynik:

400bc0: movsd 0x8(%rsi),%xmm6 
400bc5: mov %rdi,%rax 
400bc8: movsd (%rsi),%xmm8 
400bcd: movsd 0x8(%rdx),%xmm5 
400bd2: movapd %xmm6,%xmm9 
400bd7: movsd (%rdx),%xmm7 
400bdb: movapd %xmm8,%xmm0 
400be0: mulsd %xmm5,%xmm9 
400be5: movsd 0x10(%rsi),%xmm4 
400bea: mulsd %xmm7,%xmm0 
400bee: movsd 0x10(%rdx),%xmm1 
400bf3: movsd 0x18(%rdx),%xmm3 
400bf8: movsd 0x18(%rsi),%xmm2 
400bfd: subsd %xmm9,%xmm0 
400c02: movapd %xmm4,%xmm9 
400c07: mulsd %xmm1,%xmm9 
400c0c: subsd %xmm9,%xmm0 
400c11: movapd %xmm3,%xmm9 
400c16: mulsd %xmm2,%xmm9 
400c1b: subsd %xmm9,%xmm0 
400c20: movapd %xmm6,%xmm9 
400c25: mulsd %xmm7,%xmm9 
400c2a: movsd %xmm0,(%rdi) 
400c2e: movapd %xmm5,%xmm0 
400c32: mulsd %xmm8,%xmm0 
400c37: subsd %xmm9,%xmm0 
400c3c: movapd %xmm4,%xmm9 
400c41: mulsd %xmm7,%xmm9 
400c46: mulsd %xmm2,%xmm7 
400c4a: movsd %xmm0,0x8(%rdi) 
400c4f: movapd %xmm1,%xmm0 
400c53: mulsd %xmm8,%xmm0 
400c58: mulsd %xmm3,%xmm8 
400c5d: subsd %xmm9,%xmm0 
400c62: subsd %xmm7,%xmm8 
400c67: movapd %xmm4,%xmm7 
400c6b: mulsd %xmm5,%xmm7 
400c6f: movsd %xmm0,0x10(%rdi) 
400c74: mulsd %xmm2,%xmm5 
400c78: movapd %xmm1,%xmm0 
400c7c: mulsd %xmm6,%xmm0 
400c80: movsd %xmm8,0x18(%rdi) 
400c86: mulsd %xmm3,%xmm6 
400c8a: mulsd %xmm2,%xmm1 
400c8e: mulsd %xmm4,%xmm3 
400c92: subsd %xmm7,%xmm0 
400c96: subsd %xmm5,%xmm6 
400c9a: subsd %xmm1,%xmm3 
400c9e: movsd %xmm0,0x20(%rdi) 
400ca3: movsd %xmm6,0x28(%rdi) 
400ca8: movsd %xmm3,0x30(%rdi) 
400cad: retq 
400cae: nop 
400caf: nop 

To wygląda całkiem nieźle mi - składników pierwszy (%rsi) a drugie (%rdx) wektory są dostępne tylko raz, a rzeczywiste obliczenia są wykonywane tylko w rejestrach. Na końcu wynik jest zapisywany pod adresem w rejestrze %rdi. Ponieważ jest to pierwszy rejestr argumentów, myślę, że tutaj zastosowana jest optymalizacja wartości zwracanej.

Porównaj to z następującej listy w wersji szablonu wyrażenie funkcji powyżej:

400cb0: mov (%rsi),%rdx 
400cb3: mov 0x8(%rsi),%rax 
400cb7: movsd 0x1f1(%rip),%xmm4  # 400eb0 <_IO_stdin_used+0x10> 
400cbe: 
400cbf: movsd 0x10(%rdx),%xmm3 
400cc4: movsd 0x18(%rdx),%xmm0 
400cc9: mulsd 0x10(%rax),%xmm3 
400cce: xorpd %xmm4,%xmm0 
400cd2: mulsd 0x18(%rax),%xmm0 
400cd7: movsd 0x8(%rdx),%xmm2 
400cdc: movsd (%rdx),%xmm1 
400ce0: mulsd 0x8(%rax),%xmm2 
400ce5: mulsd (%rax),%xmm1 
400ce9: subsd %xmm3,%xmm0 
400ced: subsd %xmm2,%xmm0 
400cf1: addsd %xmm0,%xmm1 
400cf5: movsd %xmm1,(%rdi) 
400cf9: movsd (%rdx),%xmm0 
400cfd: movsd 0x8(%rdx),%xmm1 
400d02: mulsd 0x8(%rax),%xmm0 
400d07: mulsd (%rax),%xmm1 
400d0b: subsd %xmm1,%xmm0 
400d0f: movsd %xmm0,0x8(%rdi) 
400d14: movsd (%rdx),%xmm0 
400d18: movsd 0x10(%rdx),%xmm1 
400d1d: mulsd 0x10(%rax),%xmm0 
400d22: mulsd (%rax),%xmm1 
400d26: subsd %xmm1,%xmm0 
400d2a: movsd %xmm0,0x10(%rdi) 
400d2f: movsd 0x8(%rdx),%xmm0 
400d34: movsd 0x10(%rdx),%xmm1 
400d39: mulsd 0x10(%rax),%xmm0 
400d3e: mulsd 0x8(%rax),%xmm1 
400d43: subsd %xmm1,%xmm0 
400d47: movsd %xmm0,0x18(%rdi) 
400d4c: movsd (%rdx),%xmm0 
400d50: movsd 0x18(%rdx),%xmm1 
400d55: mulsd 0x18(%rax),%xmm0 
400d5a: mulsd (%rax),%xmm1 
400d5e: subsd %xmm1,%xmm0 
400d62: movsd %xmm0,0x20(%rdi) 
400d67: movsd 0x8(%rdx),%xmm0 
400d6c: movsd 0x18(%rdx),%xmm1 
400d71: mulsd 0x18(%rax),%xmm0 
400d76: mulsd 0x8(%rax),%xmm1 
400d7b: subsd %xmm1,%xmm0 
400d7f: movsd %xmm0,0x28(%rdi) 
400d84: movsd 0x10(%rdx),%xmm0 
400d89: movsd 0x18(%rdx),%xmm1 
400d8e: mulsd 0x18(%rax),%xmm0 
400d93: mulsd 0x10(%rax),%xmm1 
400d98: subsd %xmm1,%xmm0 
400d9c: movsd %xmm0,0x30(%rdi) 
400da1: retq 

Podpis tej funkcji jest

spinor<product<vector, vector>>(product<vector, vector> const&) 

Mam nadzieję, że mi zaufać, że zarówno wersja dać takie same wynik. Pierwsze dwie linie wyodrębniają pierwszy i drugi wektor, które są przechowywane jako odniesienia w product. Zastanawiałem się nad następującymi rzeczami:

• Co robi movsd 0x1f1(%rip),%xmm4 w połączeniu z xorpd %xmm4,%xmm0? Już się dowiedziałem, że nazywa się to "względnym adresowaniem RIP", zobacz http://www.x86-64.org/documentation/assembly.html
• Dlaczego GCC nie wykorzystuje więcej rejestrów, np. buforować 0x10(%rax), który jest czytany cztery razy?

ja również oceniać w obie funkcje generując 100000000 losowych wektorów i poświęcenie czasu obie funkcje potrzebne:

ET: 7.5 sec 
HW: 6.8 sec 

Funkcja ręcznie napisane jest o 10% szybciej. Czy ktokolwiek ma doświadczenie z szablonami ekspresji i wie, jak sprawić, by były bardziej zbliżone do pisanego odręcznie odpowiednika?

Answer 1

Byłoby jasne, czy my wiemy na pewno zawartość adres 0x400eb0, ale podejrzewam, że jest 0x8000 0000 0000 0000 8000 0000 0000 0000 lub podobne (ewentualnie z wiodącym 0, ponieważ kod nie jest wektorowy), napisany w 128 bitową int.

W takim przypadku xorpd zmienia znak drugiego operandu.

Przyczyna, dla której odczytany rejestr nie jest buforowany - lepiej zapytać o to na liście mailingowej gcc-help. Prawdopodobnie kompilator nie może udowodnić, że dwa wektory lub wynik pośredni nie są aliasami.

Ale wbrew ogólnej opinii, kompilatory nie optymalizują się zawsze idealnie, ale tylko lepiej niż 90% (lub 99%?) Wszystkich programistów (jeśli próbują napisać zespół), a czasami (rzadko) produkują bardzo wolno kod.

Ale twoje podejście jest bardzo dobre - analiza porównawcza i szukanie wygenerowanego kodu obiektu jest właściwą rzeczą, jeśli chcesz zoptymalizować.

PS: Kod może zostać przyspieszony za pomocą instrukcji wektorowych (mulpd zamiast mulsd), co powoduje pomnożenie dwóch lub czterech kopii w jednym przejściu), czyli SSE lub AVX. Ale potrzebne są pewne instrukcje, aby przetasować wartości w odpowiednie miejsca w rejestrach, więc wzmocnienie jest zawsze wolniejsze niż dwa lub cztery razy.