Warum ist die neue Zufallsbibliothek besser als std :: rand ()?

82

Also sah ich einen Vortrag namens rand (), der als schädlich angesehen wurde, und er befürwortete die Verwendung des Motorverteilungsparadigmas der Zufallszahlengenerierung gegenüber dem einfachen std::rand()Plusmodul-Paradigma.

Ich wollte jedoch die Fehler aus std::rand()erster Hand sehen, also machte ich ein kurzes Experiment:

  1. Im Grunde genommen habe ich 2 Funktionen getRandNum_Old()und getRandNum_New()daß eine Zufallszahl zwischen 0 und 5 einschließlich erzeugt wurde , std::rand()und std::mt19937+ std::uniform_int_distributionsind.
  2. Dann habe ich 960.000 (durch 6 teilbare) Zufallszahlen auf die "alte" Weise generiert und die Häufigkeiten der Zahlen 0-5 aufgezeichnet. Dann habe ich die Standardabweichung dieser Frequenzen berechnet. Was ich suche, ist eine möglichst geringe Standardabweichung, da dies passieren würde, wenn die Verteilung wirklich gleichmäßig wäre.
  3. Ich habe diese Simulation 1000 Mal ausgeführt und die Standardabweichung für jede Simulation aufgezeichnet. Ich habe auch die Zeit in Millisekunden aufgezeichnet.
  4. Danach habe ich genau das gleiche noch einmal gemacht, aber diesmal habe ich Zufallszahlen auf "neue" Weise generiert.
  5. Schließlich berechnete ich den Mittelwert und die Standardabweichung der Liste der Standardabweichungen für den alten und den neuen Weg sowie den Mittelwert und die Standardabweichung für die Liste der Zeiten, die sowohl für den alten als auch für den neuen Weg genommen wurden.

Hier waren die Ergebnisse:

[OLD WAY]
Spread
       mean:  346.554406
    std dev:  110.318361
Time Taken (ms)
       mean:  6.662910
    std dev:  0.366301

[NEW WAY]
Spread
       mean:  350.346792
    std dev:  110.449190
Time Taken (ms)
       mean:  28.053907
    std dev:  0.654964

Überraschenderweise war die Gesamtverteilung der Rollen für beide Methoden gleich. Dh std::mt19937+ std::uniform_int_distributionwar nicht "einheitlicher" als einfach std::rand()+ %. Eine andere Beobachtung, die ich machte, war, dass der neue ungefähr 4x langsamer war als der alte Weg. Insgesamt schien es, als würde ich enorme Geschwindigkeitskosten für fast keinen Qualitätsgewinn zahlen.

Ist mein Experiment irgendwie fehlerhaft? Oder ist das std::rand()wirklich nicht so schlimm und vielleicht sogar noch besser?

Als Referenz ist hier der Code, den ich in seiner Gesamtheit verwendet habe:

#include <cstdio>
#include <random>
#include <algorithm>
#include <chrono>

int getRandNum_Old() {
    static bool init = false;
    if (!init) {
        std::srand(time(nullptr)); // Seed std::rand
        init = true;
    }

    return std::rand() % 6;
}

int getRandNum_New() {
    static bool init = false;
    static std::random_device rd;
    static std::mt19937 eng;
    static std::uniform_int_distribution<int> dist(0,5);
    if (!init) {
        eng.seed(rd()); // Seed random engine
        init = true;
    }

    return dist(eng);
}

template <typename T>
double mean(T* data, int n) {
    double m = 0;
    std::for_each(data, data+n, [&](T x){ m += x; });
    m /= n;
    return m;
}

template <typename T>
double stdDev(T* data, int n) {
    double m = mean(data, n);
    double sd = 0.0;
    std::for_each(data, data+n, [&](T x){ sd += ((x-m) * (x-m)); });
    sd /= n;
    sd = sqrt(sd);
    return sd;
}

int main() {
    const int N = 960000; // Number of trials
    const int M = 1000;   // Number of simulations
    const int D = 6;      // Num sides on die

    /* Do the things the "old" way (blech) */

    int freqList_Old[D];
    double stdDevList_Old[M];
    double timeTakenList_Old[M];

    for (int j = 0; j < M; j++) {
        auto start = std::chrono::high_resolution_clock::now();
        std::fill_n(freqList_Old, D, 0);
        for (int i = 0; i < N; i++) {
            int roll = getRandNum_Old();
            freqList_Old[roll] += 1;
        }
        stdDevList_Old[j] = stdDev(freqList_Old, D);
        auto end = std::chrono::high_resolution_clock::now();
        auto dur = std::chrono::duration_cast<std::chrono::microseconds>(end-start);
        double timeTaken = dur.count() / 1000.0;
        timeTakenList_Old[j] = timeTaken;
    }

    /* Do the things the cool new way! */

    int freqList_New[D];
    double stdDevList_New[M];
    double timeTakenList_New[M];

    for (int j = 0; j < M; j++) {
        auto start = std::chrono::high_resolution_clock::now();
        std::fill_n(freqList_New, D, 0);
        for (int i = 0; i < N; i++) {
            int roll = getRandNum_New();
            freqList_New[roll] += 1;
        }
        stdDevList_New[j] = stdDev(freqList_New, D);
        auto end = std::chrono::high_resolution_clock::now();
        auto dur = std::chrono::duration_cast<std::chrono::microseconds>(end-start);
        double timeTaken = dur.count() / 1000.0;
        timeTakenList_New[j] = timeTaken;
    }

    /* Display Results */

    printf("[OLD WAY]\n");
    printf("Spread\n");
    printf("       mean:  %.6f\n", mean(stdDevList_Old, M));
    printf("    std dev:  %.6f\n", stdDev(stdDevList_Old, M));
    printf("Time Taken (ms)\n");
    printf("       mean:  %.6f\n", mean(timeTakenList_Old, M));
    printf("    std dev:  %.6f\n", stdDev(timeTakenList_Old, M));
    printf("\n");
    printf("[NEW WAY]\n");
    printf("Spread\n");
    printf("       mean:  %.6f\n", mean(stdDevList_New, M));
    printf("    std dev:  %.6f\n", stdDev(stdDevList_New, M));
    printf("Time Taken (ms)\n");
    printf("       mean:  %.6f\n", mean(timeTakenList_New, M));
    printf("    std dev:  %.6f\n", stdDev(timeTakenList_New, M));
}
rcplusplus
quelle
32
Das ist so ziemlich der Grund, warum dieser Rat existiert. Wenn Sie nicht wissen, wie Sie das RNG auf ausreichende Entropie testen sollen oder ob es für Ihr Programm wichtig ist oder nicht, sollten Sie davon ausgehen, dass std :: rand () nicht gut genug ist. en.wikipedia.org/wiki/Entropy_(computing)
Hans Passant
4
Das Endergebnis davon, ob rand()es gut genug ist, hängt weitgehend davon ab, wofür Sie die Sammlung von Zufallszahlen verwenden. Wenn Sie eine bestimmte Art der Zufallsverteilung benötigen, ist die Bibliotheksimplementierung natürlich besser. Wenn Sie einfach Zufallszahlen benötigen und sich nicht um die "Zufälligkeit" oder die Art der Verteilung kümmern, rand()ist dies in Ordnung. Passen Sie das richtige Werkzeug an den jeweiligen Auftrag an.
David C. Rankin
2
möglicher Betrüger: stackoverflow.com/questions/52869166/… Ich möchte diesen einfach nicht hämmern, also verzichte ich darauf, tatsächlich abzustimmen.
Bolov
18
for (i=0; i<k*n; i++) a[i]=i%n;erzeugt den gleichen exakten Mittelwert und die gleiche Standardabweichung wie das beste RNG da draußen. Wenn dies für Ihre Anwendung gut genug ist, verwenden Sie einfach diese Sequenz.
n. 'Pronomen' m.
3
"Standardabweichung so gering wie möglich" - Nr. Das ist falsch. Sie erwarten, dass die Frequenzen etwas anders sind - bei sqrt (Frequenz) handelt es sich um die Standardabweichung, die Sie erwarten. Der "inkrementierende Zähler", den nm erzeugt, hat einen viel niedrigeren SD (und ist ein sehr schlechter Rng).
Martin Bonner unterstützt Monica

Antworten:

106

So ziemlich jede Implementierung von "alt" rand()verwendet eine LCG ; Während sie im Allgemeinen nicht die besten Generatoren sind, werden sie bei einem so grundlegenden Test normalerweise nicht versagen - Mittelwert und Standardabweichung werden im Allgemeinen selbst von den schlechtesten PRNGs richtig eingestellt.

Häufige Fehler bei "schlechten" - aber häufig genug - rand()Implementierungen sind:

  • geringe Zufälligkeit von Bits niedriger Ordnung;
  • kurze Zeit;
  • niedrig RAND_MAX;
  • Eine gewisse Korrelation zwischen aufeinanderfolgenden Extraktionen (im Allgemeinen erzeugen LCGs Zahlen, die sich auf einer begrenzten Anzahl von Hyperebenen befinden, obwohl dies irgendwie gemildert werden kann).

Dennoch ist keines davon spezifisch für die API von rand(). Eine bestimmte Implementierung könnte einen Generator der Xorshift-Familie hinter srand/ platzieren randund algoritmisch gesehen ein PRNG auf dem neuesten Stand der Technik ohne Änderungen der Schnittstelle erhalten, sodass kein Test wie der von Ihnen durchgeführte eine Schwäche in der Ausgabe zeigen würde.

Bearbeiten: @R. stellt korrekt fest, dass die rand/ srand-Schnittstelle durch die Tatsache begrenzt ist, dass srandein unsigned intGenerator erforderlich ist, sodass jeder Generator, den eine Implementierung hinter sich lässt, an sich auf UINT_MAXmögliche Start-Seeds (und damit generierte Sequenzen) beschränkt ist. Dies ist in der Tat wahr, obwohl die API trivial erweitert werden könnte, srandum eine unsigned long longoder eine separate srand(unsigned char *, size_t)Überladung hinzuzufügen .


In der Tat ist das eigentliche Problem rand()nicht viel von der Umsetzung im Prinzip, aber:

  • Abwärtskompatibilität; Viele aktuelle Implementierungen verwenden suboptimale Generatoren, normalerweise mit schlecht gewählten Parametern. Ein berüchtigtes Beispiel ist Visual C ++ mit einem RAND_MAXWert von nur 32767. Dies kann jedoch nicht einfach geändert werden, da dies die Kompatibilität mit der Vergangenheit beeinträchtigen würde - Personen, die srandeinen festen Startwert für reproduzierbare Simulationen verwenden, wären nicht allzu glücklich (tatsächlich IIRC) Die oben genannte Implementierung geht auf frühere Versionen von Microsoft C (oder sogar Lattice C) aus der Mitte der achtziger Jahre zurück.
  • vereinfachte Schnittstelle; rand()stellt einen einzelnen Generator mit dem globalen Status für das gesamte Programm bereit. Dies ist zwar für viele einfache Anwendungsfälle vollkommen in Ordnung (und tatsächlich recht praktisch), wirft jedoch Probleme auf:

    • mit Multithread-Code: Um dies zu beheben, benötigen Sie entweder einen globalen Mutex, der alles ohne Grund verlangsamt und die Wahrscheinlichkeit einer Wiederholbarkeit zunichte macht, da die Reihenfolge der Aufrufe selbst zufällig wird, oder einen Thread-lokalen Status. Letzteres wurde von mehreren Implementierungen übernommen (insbesondere Visual C ++).
    • Wenn Sie eine "private", reproduzierbare Sequenz in einem bestimmten Modul Ihres Programms wünschen, die den globalen Status nicht beeinflusst.

Schließlich der randStand der Dinge:

  • spezifiziert keine tatsächliche Implementierung (der C-Standard bietet nur eine Beispielimplementierung), daher muss jedes Programm, das reproduzierbare Ausgabe (oder ein PRNG von bekannter Qualität) für verschiedene Compiler erzeugen soll, seinen eigenen Generator rollen;
  • bietet keine plattformübergreifende Methode, um einen anständigen Startwert zu erhalten ( time(NULL)ist nicht, da es nicht granular genug ist und oft - denken Sie an eingebettete Geräte ohne RTC - nicht einmal zufällig genug).

Daher der neue <random>Header, der versucht, dieses Durcheinander zu beheben, indem er folgende Algorithmen bereitstellt:

  • vollständig spezifiziert (damit Sie eine über den Compiler reproduzierbare Ausgabe und garantierte Eigenschaften haben können - beispielsweise die Reichweite des Generators);
  • im Allgemeinen auf dem neuesten Stand der Technik ( ab dem Zeitpunkt der Gestaltung der Bibliothek ; siehe unten);
  • Eingekapselt in Klassen (so wird Ihnen kein globaler Status aufgezwungen, wodurch Threading- und Nichtlokalitätsprobleme vollständig vermieden werden);

... und eine Standardeinstellung, random_deviceum sie zu säen.

Wenn Sie mich fragen, hätte ich mir auch eine einfache API gewünscht , die darauf aufbaut, und zwar für die Fälle "einfach", "erraten Sie eine Anzahl" (ähnlich wie Python die "komplizierte" API bereitstellt, aber auch die triviale random.randint& Co. (mit einem globalen, vorab gesetzten PRNG für uns unkomplizierte Leute, die nicht jedes Mal, wenn wir eine Nummer für die Bingokarten extrahieren möchten, in zufälligen Geräten / Motoren / Adaptern ertrinken möchten), aber es ist wahr, dass Sie es leicht können Erstellen Sie es selbst über die aktuellen Einrichtungen (während das Erstellen der "vollständigen" API über eine vereinfachte API nicht möglich wäre).


Um schließlich zu Ihrem Leistungsvergleich zurückzukehren: Wie andere angegeben haben, vergleichen Sie ein schnelles LCG mit einem langsameren (aber allgemein als besser angesehenen) Mersenne Twister. Wenn Sie mit der Qualität eines LCG einverstanden sind, können Sie std::minstd_randstattdessen verwenden std::mt19937.

In der Tat, nachdem Sie Ihre Funktion optimiert haben, std::minstd_randum nutzlose statische Variablen für die Initialisierung zu verwenden und zu vermeiden

int getRandNum_New() {
    static std::minstd_rand eng{std::random_device{}()};
    static std::uniform_int_distribution<int> dist{0, 5};
    return dist(eng);
}

Ich bekomme 9 ms (alt) vs 21 ms (neu); Schließlich, wenn ich los werde dist(was im Vergleich zum klassischen Modulo-Operator den Verteilungsversatz für den Ausgabebereich behandelt, der nicht ein Vielfaches des Eingabebereichs beträgt) und zurück zu dem komme, was Sie gerade tungetRandNum_Old()

int getRandNum_New() {
    static std::minstd_rand eng{std::random_device{}()};
    return eng() % 6;
}

Ich bekomme es auf 6 ms nach unten (also 30% schneller), wahrscheinlich , weil, im Gegensatz zu dem Aufruf rand(), std::minstd_randeinfacher zu inline ist.


Übrigens habe ich den gleichen Test mit einem handgerollten Test durchgeführt (entspricht jedoch weitgehend der Standardbibliotheksschnittstelle) XorShift64*und ist 2,3-mal schneller als rand()(3,68 ms gegenüber 8,61 ms). Angesichts der Tatsache, dass es im Gegensatz zum Mersenne Twister und den verschiedenen mitgelieferten LCGs die aktuellen Zufälligkeitstestsuiten mit Bravour besteht und unglaublich schnell ist, wundert man sich, warum es noch nicht in der Standardbibliothek enthalten ist.

Matteo Italia
quelle
3
Es ist genau die Kombination von srandund einem nicht spezifizierten Algorithmus, der std::rand in Schwierigkeiten gerät. Siehe auch meine Antwort auf eine andere Frage .
Peter O.
2
randist auf API-Ebene grundsätzlich dadurch begrenzt, dass der Startwert (und damit die Anzahl der möglichen Sequenzen, die produziert werden können) durch begrenzt ist UINT_MAX+1.
R .. GitHub STOP HELPING ICE
2
Nur eine Anmerkung: minstd ist ein schlechtes PRNG, mt19973 ist besser, aber nicht viel: pcg-random.org/… (in dieser Tabelle minstd == LCG32 / 64). Es ist ziemlich schade, dass C ++ keine hochwertigen, schnellen PRNGs wie PCG oder xoroshiro128 + bietet.
user60561
2
@MatteoItalia Wir sind uns nicht einig. Dies war auch Bjarnes Punkt. Wir wollen wirklich <random>in den Standard, aber wir möchten auch die Option "Gib mir nur eine anständige Implementierung, die ich jetzt verwenden kann". Für PRNGs und andere Dinge.
Ravnsgaard
2
Ein paar Anmerkungen: 1. Durch Ersetzen std::uniform_int_distribution<int> dist{0, 5}(eng);durch wird eng() % 6der Versatzfaktor, unter dem der std::randCode leidet, wieder eingeführt (zugegebenermaßen geringfügiger Versatz in diesem Fall, wenn die Engine 2**31 - 1Ausgänge hat und Sie diese auf 6 Buckets verteilen). 2. In Ihrem Hinweis zu " srandTakes a unsigned int", das die möglichen Ausgaben begrenzt, wie geschrieben, hat das Aussäen Ihres Motors genau std::random_device{}()das gleiche Problem. Sie benötigen eine seed_seq, um die meisten PRNGs ordnungsgemäß zu initialisieren .
ShadowRanger
6

Wenn Sie Ihr Experiment mit einem Bereich von mehr als 5 wiederholen, werden Sie wahrscheinlich unterschiedliche Ergebnisse sehen. Wenn Ihre Reichweite erheblich kleiner ist als RAND_MAXfür die meisten Anwendungen, gibt es kein Problem.

Wenn wir zum Beispiel eine RAND_MAXvon 25 rand() % 5haben, werden Zahlen mit den folgenden Frequenzen erzeugt:

0: 6
1: 5
2: 5
3: 5
4: 5

Da RAND_MAXgarantiert mehr als 32767 beträgt und der Frequenzunterschied zwischen dem unwahrscheinlichsten und dem wahrscheinlichsten nur 1 beträgt, ist die Verteilung für kleine Zahlen für die meisten Anwendungsfälle nahezu zufällig genug.

Alan Birtles
quelle
3
Dies wird in STLs zweiter Folie erklärt
Alan Birtles
4
Ok, aber ... wer ist STL? Und welche Folien? (ernste Frage)
Kebs
@kebs, Stephan Lavavej, siehe die Youtube-Referenz in der Frage.
Evg
3

Erstens ändert sich die Antwort überraschenderweise je nachdem, wofür Sie die Zufallszahl verwenden. Wenn Sie beispielsweise einen zufälligen Hintergrundfarbwechsler verwenden möchten, ist die Verwendung von rand () vollkommen in Ordnung. Wenn Sie eine Zufallszahl verwenden, um eine zufällige Pokerhand oder einen kryptografisch sicheren Schlüssel zu erstellen, ist dies nicht in Ordnung.

Vorhersagbarkeit: Die Sequenz 012345012345012345012345 ... würde eine gleichmäßige Verteilung jeder Zahl in Ihrer Stichprobe liefern, ist aber offensichtlich nicht zufällig. Damit eine Sequenz zufällig ist, kann der Wert von n + 1 nicht leicht durch den Wert von n (oder sogar durch die Werte von n, n-1, n-2, n-3 usw.) vorhergesagt werden. Dies ist eindeutig eine sich wiederholende Sequenz von den gleichen Ziffern ist ein entarteter Fall, aber eine Sequenz, die mit einem linearen Kongruenzgenerator erzeugt wird, kann einer Analyse unterzogen werden; Wenn Sie die Standardeinstellungen einer allgemeinen LCG aus einer gemeinsamen Bibliothek verwenden, kann eine böswillige Person die Sequenz ohne großen Aufwand "unterbrechen". In der Vergangenheit wurden mehrere Online-Casinos (und einige stationäre) von Maschinen mit schlechten Zufallsgeneratoren wegen Verlusten getroffen. Sogar Leute, die es besser wissen sollten, wurden eingeholt;

Verteilung: Wie im Video erwähnt, garantiert die Verwendung eines Modulos von 100 (oder eines Werts, der nicht gleichmäßig in die Länge der Sequenz teilbar ist), dass einige Ergebnisse zumindest geringfügig wahrscheinlicher werden als andere. Im Universum von 32767 möglichen Startwerten modulo 100 erscheinen die Zahlen 0 bis 66 328/327 (0,3%) häufiger als die Werte 67 bis 99; Ein Faktor, der einem Angreifer einen Vorteil verschaffen kann.

JackLThornton
quelle
1
"Vorhersagbarkeit: Die Sequenz 012345012345012345012345 ... würde Ihren Test auf" Zufälligkeit "bestehen, da es eine gleichmäßige Verteilung jeder Zahl in Ihrer Stichprobe geben würde" tatsächlich, nicht wirklich; Was er misst, ist der Standardwert der Standardwerte zwischen den Läufen, dh im Wesentlichen, wie das Histogramm der verschiedenen Läufe verteilt ist. Mit einem 012345012345012345 ... Generator wäre es immer Null.
Matteo Italia
Guter Punkt; Ich fürchte, ich habe den Code von OP etwas zu schnell durchgelesen. Ich habe meine Antwort bearbeitet, um nachzudenken.
JackLThornton
Hehe ich weiß, weil ich diesen Test auch machen wollte, und ich bemerkte, dass ich unterschiedliche Ergebnisse erhielt 😄
Matteo Italia
1

Die richtige Antwort lautet: Es kommt darauf an, was Sie unter "besser" verstehen.

Die "neuen" <random>Engines wurden vor über 13 Jahren in C ++ eingeführt, sind also nicht wirklich neu. Die C-Bibliothek rand()wurde vor Jahrzehnten eingeführt und war in dieser Zeit für eine Reihe von Dingen sehr nützlich.

Die C ++ - Standardbibliothek bietet drei Klassen von Zufallszahlengenerator-Engines: die lineare Kongruenz (von denen rand() ein Beispiel ist), die verzögerte Fibonacci und die Mersenne Twister. Es gibt Kompromisse für jede Klasse, und jede Klasse ist in gewisser Weise "am besten". Zum Beispiel haben die LCGs einen sehr kleinen Zustand und wenn die richtigen Parameter ausgewählt werden, ist dies auf modernen Desktop-Prozessoren ziemlich schnell. Die LFGs haben einen größeren Status und verwenden nur Speicherabrufe und Additionsoperationen. Sie sind daher auf eingebetteten Systemen und Mikrocontrollern, denen spezielle mathematische Hardware fehlt, sehr schnell. Das MTG hat einen riesigen Zustand und ist langsam, kann aber eine sehr große, sich nicht wiederholende Sequenz mit ausgezeichneten spektralen Eigenschaften aufweisen.

Wenn keiner der mitgelieferten Generatoren für Ihre spezifische Verwendung ausreichend ist, bietet die C ++ - Standardbibliothek auch eine Schnittstelle für einen Hardwaregenerator oder Ihre eigene benutzerdefinierte Engine. Keiner der Generatoren ist für die eigenständige Verwendung vorgesehen: Die beabsichtigte Verwendung erfolgt über ein Verteilungsobjekt, das eine zufällige Sequenz mit einer bestimmten Wahrscheinlichkeitsverteilungsfunktion bereitstellt.

Ein weiterer Vorteil von <random>over rand()besteht darin, dass rand()der globale Status verwendet wird, nicht wiedereintrittsfähig oder threadsicher ist und eine einzelne Instanz pro Prozess zulässig ist. Wenn Sie eine feinkörnige Kontrolle oder Vorhersagbarkeit benötigen (dh in der Lage sind, einen Fehler angesichts des RNG-Startzustands zu reproduzieren), rand()ist dies nutzlos. Die <random>Generatoren sind lokal instanziiert und haben einen serialisierbaren (und wiederherstellbaren) Zustand.

Stephen M. Webb
quelle