Wyjaśniono szybką metodę zaokrąglania liczby podwójnej do 32-bitowej liczby int

169

Czytając kod źródłowy Lua , zauważyłem, że Lua używa a, macroaby zaokrąglić a doubledo 32-bitowego int. Wyodrębniłem macroi wygląda to tak:

union i_cast {double d; int i[2]};
#define double2int(i, d, t)  \
    {volatile union i_cast u; u.d = (d) + 6755399441055744.0; \
    (i) = (t)u.i[ENDIANLOC];}

Tutaj ENDIANLOCdefiniuje się jako endianness , 0dla little endian, 1dla big endian. Lua ostrożnie radzi sobie z endianizmem. toznacza typ liczby całkowitej, na przykład intlub unsigned int.

Zrobiłem trochę badań i istnieje prostszy format, macroktóry wykorzystuje tę samą myśl:

#define double2int(i, d) \
    {double t = ((d) + 6755399441055744.0); i = *((int *)(&t));}

Lub w stylu C ++:

inline int double2int(double d)
{
    d += 6755399441055744.0;
    return reinterpret_cast<int&>(d);
}

Ta sztuczka może działać na każdej maszynie używającej IEEE 754 (co oznacza dzisiaj prawie każdą maszynę). Działa zarówno dla liczb dodatnich, jak i ujemnych, a zaokrąglanie odbywa się zgodnie z regułą bankiera . (Nie jest to zaskakujące, ponieważ jest zgodne z IEEE 754.)

Napisałem mały program do testowania:

int main()
{
    double d = -12345678.9;
    int i;
    double2int(i, d)
    printf("%d\n", i);
    return 0;
}

I wyprowadza -12345679, zgodnie z oczekiwaniami.

Chciałbym szczegółowo wyjaśnić, jak macrodziała ta sztuczka . Magiczna liczba 6755399441055744.0to w rzeczywistości 2^51 + 2^52lub 1.5 * 2^52, 1.5aw systemie dwójkowym może być reprezentowana jako 1.1. Kiedy do tej magicznej liczby zostanie dodana jakakolwiek 32-bitowa liczba całkowita, cóż, zgubiłem się stąd. Jak działa ta sztuczka?

PS: To jest w kodzie źródłowym Lua, Llimits.h .

AKTUALIZACJA :

  1. Jak wskazuje @Mysticial, ta metoda nie ogranicza się do wersji 32-bitowej int, można ją również rozszerzyć do wersji 64-bitowej, into ile liczba mieści się w zakresie 2 ^ 52. ( macroPotrzebuje pewnych modyfikacji.)
  2. Niektóre materiały mówią, że ta metoda nie może być używana w Direct3D .
  3. Podczas pracy z asemblerem Microsoft dla x86, jest jeszcze szybszy macrozapis assembly(jest to również wyodrębnione ze źródła Lua):

    #define double2int(i,n)  __asm {__asm fld n   __asm fistp i}
  4. Istnieje podobna liczba magiczna dla liczby pojedynczej precyzji: 1.5 * 2 ^23

Yu Hao
źródło
3
„szybko” w porównaniu z czym?
Cory Nelson
3
@CoryNelson Szybko w porównaniu do prostej obsady. Ta metoda, jeśli jest prawidłowo wdrożona (z elementami wewnętrznymi SSE), jest dosłownie sto razy szybsza niż rzutowanie. (co wywołuje nieprzyjemne wywołanie funkcji do dość kosztownego kodu konwersji)
Mysticial
2
Racja - widzę, że jest szybszy niż ftoi. Ale jeśli mówisz o SSE, dlaczego nie skorzystać po prostu z jednej instrukcji CVTTSD2SI?
Cory Nelson
3
@tmyklebu Wiele z tych przypadków użycia double -> int64rzeczywiście mieści się w 2^52zakresie. Są one szczególnie powszechne podczas wykonywania zwojów całkowitych przy użyciu zmiennoprzecinkowych FFT.
Mysticial
7
@MSalters Niekoniecznie prawda. Rzutowanie musi być zgodne ze specyfikacją języka - w tym z odpowiednią obsługą przypadków przepełnienia i NAN. (lub cokolwiek określi kompilator w przypadku IB lub UB) Te kontrole są zwykle bardzo kosztowne. Sztuczka wspomniana w tym pytaniu całkowicie ignoruje takie przypadki narożne. Więc jeśli chcesz szybkości, a twoja aplikacja nie dba o takie przypadki narożne (lub nigdy nie napotyka), ten hack jest idealnie odpowiedni.
Mysticial

Odpowiedzi:

161

A doublejest reprezentowane w następujący sposób:

podwójna reprezentacja

i można go postrzegać jako dwie 32-bitowe liczby całkowite; teraz, intwe wszystkich wersjach twojego kodu (zakładając, że jest to 32-bitowy int), jest ten po prawej stronie rysunku, więc to, co robisz na końcu, to po prostu pobranie najniższych 32 bitów mantysy.


A teraz do magicznej liczby; jak poprawnie powiedziałeś, 6755399441055744 to 2 ^ 51 + 2 ^ 52; dodanie takiej liczby wymusza doubleprzejście do „słodkiego zakresu” między 2 ^ 52 a 2 ^ 53, co, jak wyjaśnia Wikipedia tutaj , ma interesującą właściwość:

Pomiędzy 2 52 = 4,503,599,627,370,496 a 2 53 = 9,007,199,254,740,992 liczbami, które można przedstawić, są dokładnie liczby całkowite

Wynika to z faktu, że mantysa ma 52 bity szerokości.

Innym interesującym faktem dotyczącym dodawania 2 51 +2 52 jest to, że wpływa na mantysę tylko w dwóch najwyższych bitach - które i tak są odrzucane, ponieważ bierzemy tylko jej najniższe 32 bity.


Wreszcie: znak.

Zmiennoprzecinkowe IEEE 754 używają reprezentacji wielkości i znaku, podczas gdy liczby całkowite na "normalnych" maszynach używają arytmetyki dopełnienia do 2; jak to jest tutaj traktowane?

Rozmawialiśmy tylko o dodatnich liczbach całkowitych; teraz przypuśćmy, że mamy do czynienia z liczbą ujemną w zakresie reprezentowanym przez 32-bit int, a więc mniejszą (w wartości bezwzględnej) niż (-2 ^ 31 + 1); nazwij to -a. Taka liczba jest oczywiście dodatnia przez dodanie magicznej liczby, a wynikowa wartość to 2 52 +2 51 + (- a).

Co otrzymamy, jeśli zinterpretujemy mantysę w reprezentacji dopełnienia 2? Musi być wynikiem sumy uzupełnień do 2 (2 52 +2 51 ) i (-a). Ponownie, pierwszy człon wpływa tylko na dwa górne bity, a to, co pozostaje w bitach 0 ~ 50, jest uzupełnieniem do dwójki (-a) (znowu, bez dwóch górnych bitów).

Ponieważ redukcja liczby uzupełnienia do 2 do mniejszej szerokości odbywa się po prostu przez odcięcie dodatkowych bitów po lewej stronie, biorąc dolne 32 bity daje nam poprawnie (-a) w 32-bitowej arytmetyce uzupełnienia do 2.

Matteo Italia
źródło
"" "Innym interesującym faktem związanym z dodaniem 2 ^ 51 + 2 ^ 52 jest to, że wpływa na mantysę tylko w dwóch najwyższych bitach - które i tak są odrzucane, ponieważ bierzemy tylko najniższe 32 bity." "" Co to jest? Dodanie tego może przesunąć całą mantysę!
YvesgereY
@John: oczywiście cały sens ich dodawania polega na wymuszeniu, aby wartość znajdowała się w tym zakresie, co oczywiście może skutkować przesunięciem mantysy (między innymi rzeczami) w stosunku do wartości pierwotnej. Mówiłem tutaj, że kiedy już znajdziesz się w tym zakresie, jedynymi bitami różniącymi się od odpowiadających 53-bitowych liczb całkowitych są bity 51 i 52, które i tak są odrzucane.
Matteo Italia
2
Dla tych, którzy chcieliby przejść na konwersję, int64_tmożesz to zrobić, przesuwając mantysę w lewo, a następnie w prawo o 13 bitów. Spowoduje to usunięcie wykładnika i dwóch bitów z „magicznej” liczby, ale zachowa i przeniesie znak na całą 64-bitową liczbę całkowitą ze znakiem. union { double d; int64_t l; } magic; magic.d = input + 6755399441055744.0; magic.l <<= 13; magic.l >>= 13;
Wojciech Migda