Czy 1.0 jest prawidłowym wyjściem z std :: gene_canonical?

Zawsze myślałem, że liczby losowe leżą między zerem a jedynką, bez1 , tj. Są to liczby z półotwartego przedziału [0,1). Potwierdza to dokumentacja na cppreference.com z dnia std::generate_canonical.

Jednak gdy uruchamiam następujący program:

#include <iostream>
#include <limits>
#include <random>

int main()
{
    std::mt19937 rng;

    std::seed_seq sequence{0, 1, 2, 3, 4, 5, 6, 7, 8, 9};
    rng.seed(sequence);
    rng.discard(12 * 629143 + 6);

    float random = std::generate_canonical<float,
                   std::numeric_limits<float>::digits>(rng);

    if (random == 1.0f)
    {
        std::cout << "Bug!\n";
    }

    return 0;
}

Daje mi następujący wynik:

Bug!

tzn. generuje mi perfekcję 1, co powoduje problemy w integracji z MC. Czy to prawidłowe zachowanie, czy po mojej stronie wystąpił błąd? Daje to ten sam wynik z G ++ 4.7.3

g++ -std=c++11 test.c && ./a.out

i brzęk 3.3

clang++ -stdlib=libc++ -std=c++11 test.c && ./a.out

Jeśli to prawidłowe zachowanie, jak mogę tego uniknąć 1?

Edycja 1 : G ++ z git wydaje się cierpieć na ten sam problem. jestem na

commit baf369d7a57fb4d0d5897b02549c3517bb8800fd
Date:   Mon Sep 1 08:26:51 2014 +0000

a kompilacja z ~/temp/prefix/bin/c++ -std=c++11 -Wl,-rpath,/home/cschwan/temp/prefix/lib64 test.c && ./a.outdaje ten sam wynik, ldddaje

linux-vdso.so.1 (0x00007fff39d0d000)
libstdc++.so.6 => /home/cschwan/temp/prefix/lib64/libstdc++.so.6 (0x00007f123d785000)
libm.so.6 => /lib64/libm.so.6 (0x000000317ea00000)
libgcc_s.so.1 => /home/cschwan/temp/prefix/lib64/libgcc_s.so.1 (0x00007f123d54e000)
libc.so.6 => /lib64/libc.so.6 (0x000000317e600000)
/lib64/ld-linux-x86-64.so.2 (0x000000317e200000)

Edycja 2 : zgłosiłem to zachowanie tutaj: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=63176

Edycja 3 : Zespół clang wydaje się być świadomy problemu: http://llvm.org/bugs/show_bug.cgi?id=18767

Problem polega na mapowaniu z kodomeny std::mt19937( std::uint_fast32_t) na float; algorytm opisany przez normę daje nieprawidłowe wyniki (niezgodne z opisem wyjścia algorytmu), gdy utrata precyzji występuje, jeśli aktualny tryb zaokrąglania IEEE754 jest inny niż zaokrąglanie do ujemnej nieskończoności (należy pamiętać, że wartość domyślna to okrągła -do najbliższego).

7549723. wyjście mt19937 z twoim nasionem to 4294967257 ( 0xffffffd9u), co po zaokrągleniu do 32-bitowej 0x1p+32liczby zmiennoprzecinkowej daje , co jest równe maksymalnej wartości mt19937, 4294967295 ( 0xffffffffu), gdy jest to również zaokrąglone do 32-bitowej liczby zmiennoprzecinkowej.

Norma mogłaby zapewnić poprawne zachowanie, gdyby określała, że ​​podczas konwersji z wyniku URNG do RealTypeof generate_canonical, zaokrąglanie ma być wykonywane w kierunku ujemnej nieskończoności; dałoby to prawidłowy wynik w tym przypadku. Jako QOI, byłoby dobrze, gdyby libstdc ++ dokonał tej zmiany.

Dzięki tej zmianie 1.0nie będzie już generowany; zamiast tego wartości graniczne 0x1.fffffep-Ndla 0 < N <= 8będą generowane częściej (w przybliżeniu 2^(8 - N - 32)na N, w zależności od rzeczywistego rozkładu MT19937).

Zalecałbym nie używać bezpośrednio floatz std::generate_canonical; raczej wygeneruj liczbę w, doublea następnie zaokrąglij w kierunku ujemnej nieskończoności:

    double rd = std::generate_canonical<double,
        std::numeric_limits<float>::digits>(rng);
    float rf = rd;
    if (rf > rd) {
      rf = std::nextafter(rf, -std::numeric_limits<float>::infinity());
    }

Ten problem może również wystąpić w przypadku std::uniform_real_distribution<float>; rozwiązanie jest takie samo, aby wyspecjalizować rozkład na doublei zaokrąglić wynik w kierunku ujemnej nieskończoności w float.

Zgodnie ze standardem 1.0nie obowiązuje.

C ++ 11 §26.5.7.2 Szablon funkcji generation_canonical

Każda funkcja, której wystąpienie pochodzi z szablonu opisanego w tej sekcji 26.5.7.2, odwzorowuje wynik jednego lub więcej wywołań dostarczonego jednolitego generatora liczb losowych gna jeden element członkowski określonego typu RealType tak, że jeśli wartości g _i utworzone przez gsą równomiernie rozłożone, wyniki instancji t _j , 0 ≤ t _j <1 , są rozłożone tak równomiernie, jak to możliwe, jak określono poniżej.

Właśnie natknąłem się na podobne pytanie uniform_real_distribution , a oto jak interpretuję oszczędne sformułowanie normy na ten temat:

Standard zawsze definiuje funkcje matematyczne w kategoriach matematycznych , nigdy w kategoriach zmiennoprzecinkowych IEEE (ponieważ Standard nadal udaje, że zmiennoprzecinkowe mogą nie oznaczać zmiennoprzecinkowych IEEE). Tak więc za każdym razem, gdy zobaczysz sformułowanie matematyczne w standardzie, oznacza to prawdziwą matematykę , a nie o IEEE.

Standard mówi, że oba uniform_real_distribution<T>(0,1)(g)i generate_canonical<T,1000>(g)powinny zwracać wartości w półotwartym zakresie [0,1). Ale to są wartości matematyczne . Kiedy weźmiesz liczbę rzeczywistą z półotwartego zakresu [0,1) i przedstawisz ją jako zmiennoprzecinkową IEEE, cóż, znaczny ułamek czasu będzie zaokrąglany w górę T(1.0).

Kiedy Tjest float(24 bity mantysy), spodziewamy się zobaczyć uniform_real_distribution<float>(0,1)(g) == 1.0fokoło 1 na 2 ^ 25 razy. Moje eksperymenty siłowe z libc ++ potwierdzają to oczekiwanie.

template<class F>
void test(long long N, const F& get_a_float) {
    int count = 0;
    for (long long i = 0; i < N; ++i) {
        float f = get_a_float();
        if (f == 1.0f) {
            ++count;
        }
    }
    printf("Expected %d '1.0' results; got %d in practice\n", (int)(N >> 25), count);
}

int main() {
    std::mt19937 g(std::random_device{}());
    auto N = (1uLL << 29);
    test(N, [&g]() { return std::uniform_real_distribution<float>(0,1)(g); });
    test(N, [&g]() { return std::generate_canonical<float, 32>(g); });
}

Przykładowe dane wyjściowe:

Expected 16 '1.0' results; got 19 in practice
Expected 16 '1.0' results; got 11 in practice

Kiedy Tjest double(53 bity mantysy), spodziewamy się zobaczyć uniform_real_distribution<double>(0,1)(g) == 1.0około 1 na 2 ^ 54 razy. Nie mam cierpliwości, aby sprawdzić to oczekiwanie. :)

Rozumiem, że takie zachowanie jest w porządku. To może obrazić nasze poczucie „half-open-rangeness” twierdząc, że dystrybucja powrót numerów „mniej niż 1,0” można w liczbach obie strony fakt, że są równe do 1.0; ale to są dwa różne znaczenia "1.0", rozumiesz? Pierwsza to matematyczna 1.0; druga to liczba zmiennoprzecinkowa o pojedynczej precyzji IEEE 1.0. Od dziesięcioleci uczono nas, aby nie porównywać liczb zmiennoprzecinkowych w celu uzyskania dokładnej równości.

Bez względu na algorytm, do którego wprowadzisz liczby losowe, nie będzie dbać, czy czasami otrzyma dokładnie 1.0. Z liczbą zmiennoprzecinkową nie można nic zrobić poza operacjami matematycznymi, a gdy tylko wykonasz jakąś operację matematyczną, twój kod będzie musiał radzić sobie z zaokrąglaniem. Nawet gdybyś mógł to zasadnie założyć generate_canonical<float,1000>(g) != 1.0f, nadal nie byłbyś w stanie tego założyć generate_canonical<float,1000>(g) + 1.0f != 2.0f- z powodu zaokrąglania. Po prostu nie możesz od tego uciec; więc dlaczego mielibyśmy udawać w tym jednym przypadku, że możesz?