Dlaczego GCC nie optymalizuje * a * a * a * a * a do (a * a * a) * (a * a * a)?

Robię optymalizację numeryczną w aplikacji naukowej. Zauważyłem tylko, że GCC zoptymalizuje wywołanie pow(a,2), kompilując je a*a, ale wywołanie pow(a,6)nie jest zoptymalizowane i faktycznie wywoła funkcję biblioteki pow, co znacznie spowalnia działanie. (Natomiast kompilator Intel C ++ , wykonywalny icc, wyeliminuje wywołanie biblioteki pow(a,6)).

Jestem ciekaw co o to, że kiedy otrzymuje pow(a,6)z a*a*a*a*a*aużyciem GCC 4.5.1 i opcje „ -O3 -lm -funroll-loops -msse4”, używa 5 mulsdwskazówek:

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13

podczas gdy jeśli napiszę (a*a*a)*(a*a*a), będzie produkować

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm13, %xmm13

co zmniejsza liczbę instrukcji mnożenia do 3. iccma podobne zachowanie.

Dlaczego kompilatory nie rozpoznają tej sztuczki optymalizacji?

Ponieważ matematyka zmiennoprzecinkowa nie jest asocjacyjna . Sposób grupowania operandów w mnożeniu zmiennoprzecinkowym ma wpływ na dokładność liczbową odpowiedzi.

W rezultacie większość kompilatorów jest bardzo konserwatywnych w kwestii zmiany kolejności obliczeń zmiennoprzecinkowych, chyba że mogą być pewni, że odpowiedź pozostanie taka sama lub jeśli nie powiesz im, że nie zależy ci na dokładności liczbowej. Na przykład: opcja gcc gcc, który pozwala na ponowne skojarzenie operacji zmiennoprzecinkowych, a nawet opcja, która pozwala jeszcze bardziej agresywny kompromisów dokładności wobec prędkości.-fassociative-math-ffast-math

Lambdageek słusznie zauważa, że ponieważ asocjatywność nie trzymać dla liczb zmiennoprzecinkowych, „Optymalizacja” oda*a*a*a*a*ado(a*a*a)*(a*a*a)może zmienić wartość. Dlatego jest niedozwolony przez C99 (chyba że użytkownik wyraźnie na to zezwolił, poprzez flagę kompilatora lub pragma). Ogólnie zakłada się, że programista napisał to, co zrobiła bez powodu, i kompilator powinien to uszanować. Jeśli chcesz(a*a*a)*(a*a*a), napisz to.

Może to jednak być trudny do napisania; dlaczego kompilator nie może po prostu zrobić tego, co uważasz za właściwe pow(a,6)? Ponieważ byłoby to niewłaściwe . Na platformie z dobrą biblioteką matematyczną pow(a,6)jest znacznie bardziej dokładny niż jeden a*a*a*a*a*alub (a*a*a)*(a*a*a). Aby dostarczyć trochę danych, przeprowadziłem mały eksperyment na moim Macu Pro, mierząc najgorszy błąd w ocenie ^ 6 dla wszystkich liczb zmiennoprzecinkowych pojedynczej precyzji między [1,2]:

worst relative error using    powf(a, 6.f): 5.96e-08
worst relative error using (a*a*a)*(a*a*a): 2.94e-07
worst relative error using     a*a*a*a*a*a: 2.58e-07

Użycie powzamiast drzewa mnożenia zmniejsza błąd związany czterokrotnie . Kompilatory nie powinny (i generalnie nie robią) „optymalizacji”, które zwiększają błąd, chyba że użytkownik na to zezwoli (np. Via -ffast-math).

Zauważ, że GCC stanowi __builtin_powi(x,n)alternatywę dla pow( ), która powinna wygenerować wbudowane drzewo mnożenia. Użyj tego, jeśli chcesz obniżyć dokładność pod względem wydajności, ale nie chcesz włączać szybkiej matematyki.

Inny podobny przypadek: większość kompilatory nie będzie optymalizować a + b + c + ddo (a + b) + (c + d)(jest to optymalizacja ponieważ drugi wyrażenie może być lepiej potokowym) i ocenia je jako dane (tj (((a + b) + c) + d)). Jest to również spowodowane przypadkami narożnymi:

float a = 1e35, b = 1e-5, c = -1e35, d = 1e-5;
printf("%e %e\n", a + b + c + d, (a + b) + (c + d));

To wychodzi 1.000000e-05 0.000000e+00

Fortran (zaprojektowany do obliczeń naukowych) ma wbudowany operator mocy i, o ile wiem, kompilatory Fortran zwykle optymalizują podnoszenie do mocy całkowitych w podobny sposób, jak to opisujesz. C / C ++ niestety nie ma operatora mocy, tylko funkcję biblioteki pow(). Nie uniemożliwia to inteligentnym kompilatorom powspecjalnego traktowania i obliczania go szybciej w szczególnych przypadkach, ale wygląda na to, że robią to rzadziej ...

Kilka lat temu próbowałem uczynić wygodniejszym obliczanie mocy całkowitych w optymalny sposób i wymyśliłem następujące. Jest to C ++, a nie C, i nadal zależy od tego, czy kompilator jest dość inteligentny w zakresie optymalizacji / wstawiania rzeczy. W każdym razie, mam nadzieję, że okaże się przydatny w praktyce:

template<unsigned N> struct power_impl;

template<unsigned N> struct power_impl {
    template<typename T>
    static T calc(const T &x) {
        if (N%2 == 0)
            return power_impl<N/2>::calc(x*x);
        else if (N%3 == 0)
            return power_impl<N/3>::calc(x*x*x);
        return power_impl<N-1>::calc(x)*x;
    }
};

template<> struct power_impl<0> {
    template<typename T>
    static T calc(const T &) { return 1; }
};

template<unsigned N, typename T>
inline T power(const T &x) {
    return power_impl<N>::calc(x);
}

_{Wyjaśnienie dla ciekawskich: nie znajduje to optymalnego sposobu obliczania mocy, ale ponieważ znalezienie optymalnego rozwiązania jest problemem NP-zupełnym i warto to robić tylko dla małych mocy (w przeciwieństwie do używania pow), nie ma powodu do zamieszania ze szczegółami.}

Następnie użyj go jako power<6>(a).

Ułatwia to wpisywanie mocy (nie trzeba przeliterować 6 as za pomocą parens) i pozwala na tego rodzaju optymalizację bez -ffast-mathprzypadku, gdy masz coś zależnego od precyzji, np. Sumowanie skompensowane (przykład, w którym niezbędna jest kolejność operacji) .

Prawdopodobnie możesz również zapomnieć, że jest to C ++ i po prostu użyć go w programie C (jeśli kompiluje się z kompilatorem C ++).

Mam nadzieję, że to może być przydatne.

EDYTOWAĆ:

Oto, co otrzymuję od mojego kompilatora:

dla a*a*a*a*a*a,

    movapd  %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0

dla (a*a*a)*(a*a*a),

    movapd  %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm1, %xmm0
    mulsd   %xmm0, %xmm0

dla power<6>(a),

    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm0, %xmm1

GCC jest rzeczywiście zoptymalizować a*a*a*a*a*ado (a*a*a)*(a*a*a)kiedy jest liczbą całkowitą. Próbowałem z tym poleceniem:

$ echo 'int f(int x) { return x*x*x*x*x*x; }' | gcc -o - -O2 -S -masm=intel -x c -

Istnieje wiele flag gcc, ale nic szczególnego. Oni mają na myśli: Czytaj ze standardowego; użyj poziomu optymalizacji O2; wypisuje listę języka asemblera zamiast pliku binarnego; lista powinna używać składni języka asemblera Intel; wejście jest w języku C (zwykle język wywodzi się z rozszerzenia pliku wejściowego, ale nie ma rozszerzenia pliku podczas odczytu ze standardowego wejścia); i napisz na standardowe wyjście.

Oto ważna część wyników. Zanotowałem to z kilkoma komentarzami wskazującymi, co się dzieje w języku asemblera:

; x is in edi to begin with.  eax will be used as a temporary register.
mov  eax, edi  ; temp = x
imul eax, edi  ; temp = x * temp
imul eax, edi  ; temp = x * temp
imul eax, eax  ; temp = temp * temp

Korzystam z systemu GCC na Linux Mint 16 Petra, pochodnej Ubuntu. Oto wersja gcc:

$ gcc --version
gcc (Ubuntu/Linaro 4.8.1-10ubuntu9) 4.8.1

Jak zauważyli inni plakaty, ta opcja nie jest możliwa w liczbach zmiennoprzecinkowych, ponieważ arytmetyka zmiennoprzecinkowa nie jest asocjacyjna.

Ponieważ 32-bitowa liczba zmiennoprzecinkowa - taka jak 1.024 - nie jest równa 1.024. W komputerze 1.024 to przedział: od (1.024-e) do (1.024 + e), gdzie „e” oznacza błąd. Niektórzy ludzie nie zdają sobie z tego sprawy i wierzą również, że * w * a oznacza pomnożenie liczb o dowolnej dokładności bez żadnych błędów związanych z tymi liczbami. Powodem, dla którego niektórzy nie zdają sobie z tego sprawy, są być może obliczenia matematyczne, które wykonywali w szkołach podstawowych: praca tylko z idealnymi liczbami bez błędów i przekonanie, że po prostu pomiń „e” podczas mnożenia. Nie widzą „e” domyślnie w „float a = 1.2”, „a * a * a” i podobnych kodach C.

Gdyby większość programistów rozpoznała (i była w stanie wykonać) koncepcję, że wyrażenie C a * a * a * a * a * a tak naprawdę nie działa z liczbami idealnymi, kompilator GCC byłby WOLNY w celu optymalizacji „a * a” * a * a * a * a „powiedzieć” t = (a * a); t * t * t ”, co wymaga mniejszej liczby mnożenia. Niestety kompilator GCC nie wie, czy programista piszący kod uważa, że ​​„a” jest liczbą z błędem, czy bez. I tak GCC zrobi tylko to, jak wygląda kod źródłowy - ponieważ tak widzi GCC „gołym okiem”.

... kiedy już wiesz, jakim jesteś programistą , możesz użyć przełącznika „-ffast-matematyki”, aby powiedzieć GCC, że „Hej, GCC, wiem, co robię!”. Umożliwi to GCC konwersję * a * a * a * a * a na inny fragment tekstu - wygląda inaczej niż * a * a * a * a * a - ale nadal oblicza liczbę w przedziale błędu a * a * a * a * a * a. Jest to w porządku, ponieważ już wiesz, że pracujesz w interwałach, a nie w liczbach idealnych.

Żaden plakat nie wspomniał jeszcze o skracaniu wyrażeń zmiennoprzecinkowych (norma ISO C, 6.5p8 i 7.12.2). Jeśli FP_CONTRACTpragma jest ustawiona na ON, kompilator może traktować wyrażenie takie jak a*a*a*a*a*apojedyncza operacja, tak jakby zostało ocenione dokładnie za pomocą pojedynczego zaokrąglenia. Na przykład kompilator może zastąpić go wewnętrzną funkcją zasilania, która jest zarówno szybsza, jak i dokładniejsza. Jest to szczególnie interesujące, ponieważ zachowanie jest częściowo kontrolowane przez programistę bezpośrednio w kodzie źródłowym, podczas gdy opcje kompilatora dostarczone przez użytkownika końcowego mogą czasami być używane nieprawidłowo.

Domyślny stan FP_CONTRACTpragmy jest zdefiniowany w implementacji, dzięki czemu kompilator może domyślnie wykonywać takie optymalizacje. W związku z tym kod przenośny, który musi ściśle przestrzegać reguł IEEE 754, powinien wyraźnie to ustawić OFF.

Jeśli kompilator nie obsługuje tej pragmy, musi być konserwatywny, unikając takiej optymalizacji, na wypadek, gdyby programista zdecydował się ją ustawić OFF.

GCC nie obsługuje tej pragmy, ale przy domyślnych opcjach zakłada, że ​​tak jest ON; tak więc dla celów ze sprzętową FMA, jeśli chce się zapobiec transformacji a*b+cdo fma (a, b, c), należy zapewnić opcję, taką jak -ffp-contract=off(jawnie ustawić pragmę OFF) lub -std=c99(poinformować GCC, aby dostosowała się do niektórych Wersja standardowa C, tutaj C99, dlatego należy postępować zgodnie z powyższym akapitem). W przeszłości ta ostatnia opcja nie uniemożliwiała transformacji, co oznacza, że ​​GCC nie było zgodne w tym punkcie: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=37845

Jak zauważył Lambdageek, mnożenie zmiennoprzecinkowe nie jest asocjacyjne i można uzyskać mniejszą dokładność, ale także, gdy można uzyskać lepszą dokładność, można argumentować przeciwko optymalizacji, ponieważ chcesz deterministycznej aplikacji. Na przykład w grze klient / serwer do symulacji, w której każdy klient musi symulować ten sam świat, w którym obliczenia zmiennoprzecinkowe są deterministyczne.

Funkcje biblioteki, takie jak „pow”, są zwykle starannie tworzone, aby uzyskać minimalny możliwy błąd (w ogólnym przypadku). Zwykle osiąga się to w przybliżeniu funkcji z splajnami (zgodnie z komentarzem Pascala wydaje się, że najczęstszą implementacją jest algorytm Remeza )

zasadniczo następująca operacja:

pow(x,y);

ma nieodłączny błąd o wielkości w przybliżeniu tej samej wielkości co błąd w dowolnym pojedynczym pomnożeniu lub dzieleniu .

Podczas następujących operacji:

float a=someValue;
float b=a*a*a*a*a*a;

ma nieodłączny błąd większy niż 5-krotność błędu pojedynczego pomnożenia lub podziału (ponieważ łączysz 5 pomnożenia).

Kompilator powinien bardzo uważać na optymalizację:

1. jeśli optymalizacja pow(a,6)do a*a*a*a*a*atego może poprawić wydajność, ale drastycznie zmniejszyć dokładność dla liczb zmiennoprzecinkowych.
2. jeśli optymalizacja a*a*a*a*a*a dopow(a,6) niego może faktycznie zmniejszyć dokładność, ponieważ „a” było jakąś specjalną wartością, która pozwala na pomnożenie bez błędu (potęga 2 lub niewielka liczba całkowita)
3. jeśli optymalizacja pow(a,6)do (a*a*a)*(a*a*a)lub (a*a)*(a*a)*(a*a)nadal może wystąpić utrata dokładności w porównaniu do powfunkcji.

Ogólnie wiesz, że dla dowolnych wartości zmiennoprzecinkowych „pow” ma lepszą dokładność niż jakakolwiek funkcja, którą ostatecznie możesz napisać, ale w niektórych szczególnych przypadkach wielokrotne mnożenie może mieć lepszą dokładność i wydajność, to programista wybiera to, co jest bardziej odpowiednie, ostatecznie komentując kod, aby nikt inny nie „zoptymalizował” tego kodu.

Jedyną rzeczą, która ma sens (osobista opinia i najwyraźniej wybór w GCC bez konkretnej optymalizacji lub flagi kompilatora) do optymalizacji, to zastąpienie „pow (a, 2)” przez „a * a”. To byłaby jedyna rozsądna rzecz, którą powinien zrobić dostawca kompilatora.

W ogóle nie spodziewałbym się, że ta sprawa zostanie zoptymalizowana. Często zdarza się, że wyrażenie zawiera podwyrażenia, które można zgrupować w celu usunięcia całych operacji. Spodziewałbym się, że autorzy kompilatorów zainwestują swój czas w obszary, w których bardziej prawdopodobne jest zauważalne ulepszenie, niż omawianie rzadko spotykanego przypadku.

Byłem zaskoczony, gdy dowiedziałem się z innych odpowiedzi, że to wyrażenie można rzeczywiście zoptymalizować za pomocą odpowiednich przełączników kompilatora. Albo optymalizacja jest trywialna, albo jest skrajnym przypadkiem znacznie powszechniejszej optymalizacji, lub autorzy kompilatora byli bardzo dokładni.

Nie ma nic złego w udzielaniu wskazówek kompilatorowi, tak jak tutaj zrobiłeś. Jest to normalna i oczekiwana część procesu mikrooptymalizacji polegająca na zmianie kolejności instrukcji i wyrażeń w celu sprawdzenia, jakie różnice przyniosą.

Chociaż kompilator może być uzasadniony, biorąc pod uwagę dwa wyrażenia w celu dostarczenia niespójnych wyników (bez odpowiednich przełączników), nie musisz być związany tym ograniczeniem. Różnica będzie niewiarygodnie mała - do tego stopnia, że ​​jeśli różnica jest dla Ciebie ważna, nie powinieneś używać standardowej arytmetyki zmiennoprzecinkowej.

Istnieje już kilka dobrych odpowiedzi na to pytanie, ale ze względu na kompletność chciałem zauważyć, że odpowiednią sekcją normy C jest 5.1.2.2.3 / 15 (która jest taka sama jak sekcja 1.9 / 9 w C ++ 11). W tej sekcji stwierdzono, że operatory można przegrupować tylko wtedy, gdy naprawdę są asocjacyjne lub przemienne.

gcc faktycznie może przeprowadzić tę optymalizację, nawet dla liczb zmiennoprzecinkowych. Na przykład,

double foo(double a) {
  return a*a*a*a*a*a;
}

staje się

foo(double):
    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm1, %xmm0
    ret

z -O -funsafe-math-optimizations. Ta zmiana kolejności narusza IEEE-754, więc wymaga flagi.

Podpisane liczby całkowite, jak zauważył Peter Cordes w komentarzu, mogą przeprowadzić tę optymalizację bez, -funsafe-math-optimizationsponieważ zachowuje się dokładnie wtedy, gdy nie ma przepełnienia, a jeśli występuje przepełnienie, zachowanie jest niezdefiniowane. Więc dostajesz

foo(long):
    movq    %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rax, %rax
    ret

z właśnie -O. W przypadku liczb całkowitych bez znaku jest to jeszcze łatwiejsze, ponieważ działają one na modach o wartości 2, dzięki czemu można je dowolnie zmieniać nawet w przypadku przepełnienia.