Effektive Zahlengenerierung in einem bestimmten Intervall

Bild

Die überwiegende Mehrheit meiner Beiträge zur Zufallszahlengenerierung befasste sich hauptsächlich mit den Eigenschaften verschiedener Generierungsschemata. Dies kann sich als unerwartet herausstellen, aber die Leistung des Randomisierungsalgorithmus hängt möglicherweise nicht vom gewählten Generierungsschema ab, sondern von anderen Faktoren. In diesem Beitrag (der mich von einem ausgezeichneten Artikel von Daniel Lemyr inspiriert hat ) werden wir die Hauptgründe für den Rückgang der Leistung bei der Erzeugung von Zufallszahlen untersuchen, die häufig die Leistung der PRN-Engine überwiegen.

Stellen Sie sich diese Situation vor:

Als Hausaufgabe implementieren Juan und Sasha denselben zufälligen Algorithmus in C ++, der auf demselben Universitätscomputer und mit einem Datensatz ausgeführt wird. Ihr Code ist nahezu identisch und unterscheidet sich nur in der Erzeugung von Zufallszahlen. Juan hat es eilig für seinen Musikunterricht, also hat er sich einfach für den Wirbelwind von Mersenne entschieden. Sasha hingegen verbrachte einige zusätzliche Stunden mit Nachforschungen. Sasha führte Benchmarks für mehrere der schnellsten PRNGs durch, die er kürzlich von sozialen Netzwerken gelernt hatte, und wählte die schnellsten aus. Während des Treffens war Sasha ungeduldig zu prahlen und fragte Juan: "Welches PRNG-System haben Sie verwendet?"

"Ich persönlich habe gerade den Mersenne-Wirbel genommen - er ist in die Sprache eingebaut und scheint ziemlich gut zu funktionieren."

"Ha!" Antwortete Sasha. „Ich habe jsf32 . Es ist viel schneller als der alte und langsame Mersenne-Wirbelwind! Mein Programm läuft in 3 Minuten 15 Sekunden! ”

"Hmm, nicht schlecht, aber meiner schafft es in weniger als einer Minute", sagt Juan und zuckt die Achseln. „Na dann muss ich zum Konzert gehen. Kommst du mit mir? "

"Nein", antwortet Sasha. "Ich ... äh ... muss meinen Code noch einmal ansehen."

Diese unangenehme fiktive Situation ist nicht besonders fiktiv; es basiert auf realen Ergebnissen. Wenn Ihr zufälliger Algorithmus nicht so schnell läuft, wie wir es möchten, und der Engpass in der Erzeugung von Zufallszahlen zu liegen scheint, liegt das Problem seltsamerweise möglicherweise nicht im Zufallszahlengenerator!

Einführung: Zufallszahlen in der Praxis


Die meisten modernen Zufallszahlengeneratoren hoher Qualität erzeugen Maschinenwörter, die mit Zufallsbits gefüllt sind, dh sie erzeugen normalerweise Zahlen im Intervall [0..2 32 ) oder [0..2 64 ). In vielen Fällen benötigen Benutzer jedoch Zahlen in einem bestimmten Intervall. Um beispielsweise einen Würfel zu werfen oder eine zufällige Spielkarte auszuwählen, werden Zahlen in kleinen konstanten Intervallen benötigt. Viele Algorithmen, von Mischen und Reservoir-Sampling bis hin zu randomisierten binären Suchbäumen, erfordern jedoch Zahlen aus anderen Intervallen.

Methoden


Wir werden uns viele verschiedene Methoden ansehen. Um die Diskussion zu vereinfachen, werden anstelle von Zahlen im Intervall [ i .. j ) oder [ i .. j ] Zahlen im Intervall [0 .. k ) generiert. Mit einem solchen Schema können wir beispielsweise Zahlen im Intervall [ i .. j ) erzeugen, indem wir k = j - i setzen , eine Zahl im Intervall [0 .. k ) erzeugen und dann i hinzufügen.

Integrierte C ++ - Tools


Viele Sprachen verfügen über integrierte Tools zum Abrufen einer Zufallszahl in einem bestimmten Intervall. Um beispielsweise eine Karte mit 52 Karten in Skriptsprachen wie Perl und Python aus einem Deck zu entfernen, können Sie int(rand(52)) bzw. random.randint(0,52) schreiben. [Anmerkung CryptoPirate- Benutzer: Es scheint mir hier ein Fehler zu sein, in Python generiert Randint (a, b) Zahlen von a nach b, einschließlich b. Und da sich 52 Karten im Deck befinden und die erste "0" ist, sollte sie zufällig sein. Randint (0,51) .] In C ++ können wir uniform_int_distribution gleiche uniform_int_distribution .

C ++ - Code zur Implementierung dieses Ansatzes ist einfach:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { std::uniform_int_distribution<uint32_t> dist(0, range-1); return dist(rng); } 

Normalerweise wird eine der unten beschriebenen Techniken in den integrierten Tools verwendet, aber die meisten Benutzer verwenden diese Tools einfach, ohne darüber nachzudenken, was „unter der Haube“ passiert, und glauben, dass diese Tools korrekt entworfen und recht effektiv sind. In C ++ sind integrierte Tools komplexer, da sie mit ziemlich willkürlichen Generierungs-Engines arbeiten können sollten - ein Generator, der Werte im Bereich von -3 bis 17 erzeugt, kann durchaus gültig sein und mit std::uniform_int_distribution , um Zahlen in jedem Intervall zu erstellen. zum Beispiel [0..1000]. Das heißt, die integrierten C ++ - Tools sind für die meisten Fälle, in denen sie verwendet werden, zu kompliziert.

Der klassische Rest der Division (schief)


Gehen wir von einem stark vereinfachten zu einem zu vereinfachten Ansatz über.

Als ich Programmieren studierte, generierten wir mit dem Restoperator Zahlen in dem Intervall (zum Beispiel, um eine Karte in einem Kartenspiel mit 52 Karten auszuwählen). Um die Zahl im Intervall [0..52) zu erhalten, haben wir rand() % 52 .

In C ++ kann dieser Ansatz wie folgt implementiert werden:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { return rng() % range; } 

Trotz der Einfachheit dieses Ansatzes zeigt es den Grund, warum das Abrufen von Zahlen im richtigen Intervall normalerweise eine langsame Aufgabe ist - es erfordert eine Division (um den vom % -Operator erhaltenen Rest zu berechnen). Das Teilen ist normalerweise mindestens eine Größenordnung langsamer als andere arithmetische Operationen, so dass eine einzelne arithmetische Operation länger dauert als alle Arbeiten, die von einem schnellen PRNG ausgeführt werden.

Neben der geringen Geschwindigkeit ist es aber auch schief . Um zu verstehen, warum rand() % 52 verzerrte Zahlen zurückgibt, nehmen wir an, dass rand() Zahlen im Intervall [0..2 32 ) erstellt, und beachten Sie, dass 52 2 32 nicht vollständig teilt, sondern 82 595 524 mal mit dem Rest teilt 48. Das heißt, wenn wir rand() % 52 , haben wir 82 595 525 Möglichkeiten, die ersten 48 Karten aus dem Stapel auszuwählen, und nur 82 595 524 Möglichkeiten, die letzten vier Karten auszuwählen. Mit anderen Worten, es gibt einen Versatz von 0,00000121% gegenüber diesen letzten vier Karten (vielleicht sind dies Könige!). Als ich Student war und Hausaufgaben über das Werfen von Würfeln oder das Zeichnen von Karten schrieb, kümmerte sich normalerweise niemand um so kleine Verzerrungen, aber mit zunehmendem Intervall wächst die Verzerrung linear. Bei einem 32-Bit-PRNG hat ein begrenztes Intervall von weniger als 2 24 einen Versatz von weniger als 0,5%, aber über 2 31 einen Versatz von 50% - einige Zahlen werden doppelt so oft zurückgegeben wie andere.

In diesem Artikel werden wir hauptsächlich Techniken betrachten, die Strategien verwenden, um einen systematischen Fehler zu beseitigen. Es ist jedoch wahrscheinlich erwähnenswert, dass bei einem 64-Bit-PRNG der Versatzwert in normalen Anwendungen wahrscheinlich vernachlässigbar ist.

Ein weiteres Problem kann sein, dass einige Generatoren schwache niedrige Bits haben. Beispielsweise haben die GPRS-Familien Xoroshiro + und Xoshiro + niedrige Bits, die statistische Tests nicht bestehen. Wenn wir % 52 ausführen (weil 52 gerade ist), übergeben wir das niederwertige Bit direkt an den Ausgang.

Gleitkommazahlen multiplizieren (schief)


Eine andere übliche Technik ist die Verwendung eines PRNG, das Gleitkommazahlen im Intervall [0..1) mit der anschließenden Umwandlung dieser Zahlen in das gewünschte Intervall erzeugt. Dieser Ansatz wird in Perl verwendet. Es wird empfohlen, int(rand(10)) zu verwenden, um eine Ganzzahl im Intervall [0..10) zu generieren, indem eine Gleitkommazahl generiert und anschließend abgerundet wird.

In C ++ ist dieser Ansatz folgendermaßen geschrieben:

 static uint32_t bounded_rand(rng_t& rng, uint32_t range) { double zeroone = 0x1.0p-32 * rng(); return range * zeroone; } 

(Beachten Sie, dass 0x1.0p-32 eine binäre Gleitkommakonstante für 2 -32 ist , mit der wir eine zufällige Ganzzahl im Intervall [0..2 32 konvertieren), um sie im Einheitsintervall zu verdoppeln. Stattdessen können wir eine solche Konvertierung mit ldexp(rng(), -32) , aber als ich diesen Ansatz ldexp(rng(), -32) , stellte sich heraus, dass er viel langsamer ist.)

Dieser Ansatz ist genauso verzerrt wie der klassische Rest der Teilung, aber die Abweichung erscheint anders. Wenn wir beispielsweise Zahlen im Intervall [0..52] auswählen würden, würden die Zahlen 0, 13, 26 und 39 einmal seltener vorkommen als andere.

Diese Version ist bei Verallgemeinerung auf 64 Bit noch unangenehmer, da sie einen Gleitkommatyp erfordert, dessen Mantisse mindestens 64 Bit beträgt. Auf x86-Computern mit Linux und macOS können wir long double , um die präziseren x86-Gleitkommazahlen mit einer 64-Bit-Mantisse zu nutzen. long double nicht universell auf alle Systeme portiert. In einigen Systemen entspricht long double double .

Es gibt eine gute Seite - dieser Ansatz ist schneller als Restlösungen für PRNGs mit schwachen niedrigen Bits.

Ganzzahlige Multiplikation (schief)


Die Multiplikationsmethode kann eher an feste als an Gleitkomma-Arithmetik angepasst werden. Tatsächlich multiplizieren wir ständig mit 2 32 ,

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { uint32_t x = rng(); uint64_t m = uint64_t(x) * uint64_t(range); return m >> 32; } 

Es scheint, dass diese Version 64-Bit-Arithmetik erfordert. Auf x86-Prozessoren kompiliert ein guter Compiler diesen Code in einen 32-Bit- mult Befehl (der uns zwei 32-Bit-Ausgabewerte gibt, von denen einer der Rückgabewert ist). Es ist zu erwarten, dass diese Version schnell ist, sie ist jedoch genau wie die Methode zum Multiplizieren von Gleitkommazahlen verzerrt.

Drop Division (kein Versatz)


Wir können das Gleitkomma-Multiplikationsschema in ein dividationsbasiertes Schema umwandeln. Anstatt x * range / 2**32 multiplizieren x * range / 2**32 berechnen wir x / (2**32 / range) . Da wir mit Ganzzahlarithmetik arbeiten, wird die Rundung in dieser Version anders ausgeführt und generiert manchmal Werte außerhalb des gewünschten Intervalls. Wenn wir diese Werte verwerfen (z. B. entfernen und neue Werte generieren), erhalten wir als Ergebnis eine Technik ohne Verzerrungen.

Wenn Sie beispielsweise eine Karte mit einem 32-Bit-PRNG herausziehen, können Sie eine 32-Bit-Zahl generieren und durch 2 32/52 = 82 595 524 teilen, um eine Karte auszuwählen. Diese Technik funktioniert, wenn der Zufallswert aus dem 32-Bit-PRNG kleiner als 52 × 82595524 = 2 32/32 - 48 ist. Wenn der Zufallswert aus dem PRNR einer der letzten 48 Werte des oberen Teils des Generatorintervalls ist, müssen Sie ihn verwerfen und nach einem anderen suchen.

Unser Code für diese Version verwendet einen Trick, bei dem 2 32 durch den range ohne 64-Bit-Mathematik zu verwenden. Für die direkte Berechnung von 2**32 / range wir die Zahl 2 32 darstellen , die zu groß ist (um eins!), Um als 32-Bit-Ganzzahl dargestellt zu werden. Stattdessen berücksichtigen wir, dass der unäre Negationsoperationsbereich für vorzeichenlose Ganzzahlen einen positiven Wert von 2 32 berechnet - range ; Wenn wir diesen Wert durch den range teilen, erhalten wir eine Antwort von weniger als 2**32 / range .

Daher sieht der C ++ - Code zum Generieren von Zahlen mithilfe von Division und Drop folgendermaßen aus:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { // calculates divisor = 2**32 / range uint32_t divisor = ((-range) / range) + 1; if (divisor == 0) // overflow, it's really 2**32 return 0; for (;;) { uint32_t val = rng() / divisor; if (val < range) return val; } } 

Natürlich erfordert dieser Ansatz zwei langsame Operationen basierend auf der Division, die normalerweise langsamer sind als andere arithmetische Operationen, daher sollten Sie nicht erwarten, dass sie schnell sind.

Der Rest der Division (doppelt) ohne Verzerrungen - OpenBSD-Technik


Wir können auch den Drop-Ansatz verwenden, um den Versatz bei der klassischen Teilungsrestmethode zu beseitigen. Im Beispiel mit Spielkarten müssen wir wieder 48 Werte fallen lassen. In dieser Version werden die ersten 48 Werte (äquivalent) verworfen, anstatt die letzten 48 Werte zu verwerfen.

Hier ist die Implementierung dieses Ansatzes in C ++:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { // calculates 2**32 % range uint32_t t = (-range) % range; for (;;) { uint32_t r = rng(); if (r >= t) return r % range; } } 

Diese Technik entfernt den Versatz, erfordert jedoch zwei zeitaufwändige Teilungsvorgänge mit dem Rest jedes Ausgabewerts (und möglicherweise benötigen Sie einen internen Generator, um mehrere Zahlen zu erstellen). Daher ist zu erwarten, dass die Methode ungefähr zweimal langsamer ist als der klassische Skew-Ansatz.

arc4random_uniform OpenBSD- arc4random_uniform (der auch unter OS X und iOS verwendet wird) verwendet diese Strategie.

Rest der Teilung (einzeln) ohne Versatz - Java-Methodik


Java verwendet einen anderen Ansatz zum Generieren einer Zahl in einem Intervall, in dem nur eine Restdivisionsoperation verwendet wird, mit Ausnahme relativ seltener Fälle, in denen das Ergebnis verworfen wird. Code:

 static uint32_t bounded_rand(rng_t& rng, uint32_t range) { uint32_t x, r; do { x = rng(); r = x % range; } while (x - r > (-range)); return r; } 

Um zu verstehen, warum diese Option funktioniert, müssen Sie ein wenig nachdenken. Im Gegensatz zur vorherigen Version, die auf Residuen basiert und die Verzerrung beseitigt, indem ein Teil der niedrigsten Werte aus der internen Generierungsmaschine entfernt wird, filtert diese Version die Werte aus dem oberen Teil des Motorintervalls.

Skew Integer Multiplication - Lemira-Methode


Ähnlich wie wir die Vorspannung aus der Restdivisionsmethode entfernt haben, können wir die Vorspannung aus der Ganzzahlmultiplikationstechnik entfernen. Diese Technik wurde von Lemyr erfunden.

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { uint32_t t = (-range) % range; do { uint32_t x = rng(); uint64_t m = uint64_t(x) * uint64_t(range); uint32_t l = uint32_t(m); } while (l < t); return m >> 32; } 

Drop Bitmaske (kein Versatz) - Apple-Technik


Bei unserem letzten Ansatz werden Divisions- und Restvorgänge vollständig eliminiert. Stattdessen wird eine einfache Maskierungsoperation verwendet, um eine Zufallszahl im Intervall [0..2 k ) zu erhalten, wobei k der kleinste Wert ist, so dass 2 k größer als das Intervall ist. Wenn der Wert für unser Intervall zu groß ist, verwerfen wir ihn und versuchen, einen anderen zu erhalten. Der Code wird unten angezeigt:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { uint32_t mask = ~uint32_t(0); --range; mask >>= __builtin_clz(range|1); uint32_t x; do { x = rng() & mask; } while (x > range); return x; } 

Dieser Ansatz wurde von Apple übernommen, als (in der Version macOS Sierra) eine eigene Überarbeitung des arc4random_uniform Codes durchgeführt wurde.

Benchmarking grundlegender Techniken


Jetzt haben wir mehrere Ansätze, die bewertet werden können. Wenn wir uns Sorgen über die Kosten eines einzelnen Geschäftsbereichs machen, wird das Benchmarking leider nicht trivial. Kein Benchmark kann alle Faktoren berücksichtigen, die den Anwendungsbereich beeinflussen, und es gibt keine Garantie dafür, dass die beste Option für Ihre Anwendung mit Sicherheit die beste für meine ist.

Wir verwenden drei Benchmarks und testen die Techniken mit vielen verschiedenen PRNGs.

Benchmark Large-Shuffle


Der wahrscheinlich offensichtlichste Maßstab ist das Mischen. In diesem Benchmark simulieren wir das Mischen in großem Maßstab. Um ein Array der Größe N zu sortieren , müssen wir Zahlen in den Intervallen [0 .. N ), [0 .. ( N -1)), ..., [0..1) generieren. In diesem Benchmark nehmen wir an, dass N die maximal mögliche Zahl ist (für uint32_t es 2 32 -1). Code:

 for (uint32_t i = 0xffffffff; i > 0; --i) { uint32_t bval = bounded_rand(rng, i); assert(bval < i); sum += bval; } 

Beachten Sie, dass wir jede Zahl „verwenden“ und zur sum addieren (damit sie nicht durch Optimierung weggeworfen wird), aber keine Mischung durchführen, um uns auf die Erzeugung von Zahlen zu konzentrieren.

Zum Testen der 64-Bit-Generierung haben wir einen ähnlichen Test, aber es ist unpraktisch, einen Test durchzuführen, der dem Mischen eines Arrays mit einer Größe von 2 64 - 1 entspricht (da es viele tausend Jahre dauern wird, bis dieser größere Benchmark abgeschlossen ist). Stattdessen kreuzen wir das gesamte 64-Bit-Intervall, generieren jedoch die gleiche Anzahl von Ausgabewerten wie im 32-Bit-Test. Code:

 for (uint32_t i = 0xffffffff; i > 0; --i) { uint64_t bound = (uint64_t(i)<<32) | i; uint64_t bval = bounded_rand(rng, bound ); assert(bval < bound); sum += bval; } 

Mersenne-Wirbelergebnisse

Die unten gezeigten Ergebnisse zeigen die Leistung dieses Benchmarks für jede der Methoden, die wir bei Verwendung des Mersenne-Wirbels und beim Testen des 32-Bit-Codes (unter Verwendung von std::mt19937 aus libstdc++ ) und eines ähnlichen 64-Bit-Codes (unter Verwendung von std:mt19937_64 aus libstdc++ ) Das Ergebnis ist das geometrische Mittel von 15 Läufen mit unterschiedlichen Startwerten, das dann normalisiert wird, so dass die klassische Teilungsrestmethode eine einzige Laufzeit hat.


Es scheint, dass wir klare Antworten zur Leistung haben - es scheint, dass Sie Techniken für ihre Perfektion entwickeln und sich fragen können, woran libstdc++ Entwickler gedacht haben, als sie eine so schreckliche Implementierung für 32-Bit-Zahlen geschrieben haben. Wie so oft beim Benchmarking ist die Situation jedoch komplizierter, als es aus diesen Ergebnissen hervorgeht. Erstens besteht das Risiko, dass die Ergebnisse spezifisch für den Mersenne-Wirbel sind, sodass wir die vielen getesteten PRNGs erweitern werden. Zweitens kann es ein subtiles Problem mit dem Benchmark selbst geben. Lassen Sie uns zuerst die erste Frage behandeln.

Ergebnisse verschiedener PRNGs

Wir werden 32-Bit- arc4_rand32 mit arc4_rand32 , chacha8r , gjrand32 , jsf32 , mt19937 , pcg32 , pcg32_fast , sfc32 , splitmix32 , xoroshiro64+ , xorshift*64/32 xoshiro128+ , xoshiro128+ und xoshiro128** und 64-Bit- gjrand64 xoshiro128** jsf64 , mcg128 , mcg128_fast , mt19937_64 , pcg64 , pcg64_fast , sfc64 , splitmix64 , xoroshiro128+ , xorshift*128/64 xoshiro256+ , xoshiro256+ und xoshiro256* . Diese Kits geben uns einige langsame PRNs und viele sehr schnelle.

Hier sind die Ergebnisse:


Wir können die Hauptunterschiede zu den Ergebnissen mit dem Mersenne-Wirbel sehen. Schnellere PRNGs verschieben das Gleichgewicht in Richtung des Begrenzungscodes, und daher wird der Unterschied zwischen den verschiedenen Ansätzen stärker, insbesondere im Fall von 64-Bit-PRNRs. Mit einer größeren libstc++ von libstc++ scheint libstc++ Implementierung nicht mehr so ​​schrecklich zu sein.

Schlussfolgerungen

In dieser Benchmark gewinnt der Ansatz, der auf Multiplikation mit Bias basiert, deutlich an Geschwindigkeit. Es gibt viele Situationen, in denen die Grenzen im Verhältnis zur Größe des PRNG klein sind und die Leistung absolut kritisch ist. In solchen Situationen ist es unwahrscheinlich, dass eine leichte Abweichung spürbare Auswirkungen hat, die PRNG-Geschwindigkeit jedoch. Ein solches Beispiel ist Quicksort mit einem zufälligen Referenzpunkt. Von den verzerrten Methoden sieht die Bitmaskentechnik vielversprechend aus.

Bevor wir jedoch ernsthafte Schlussfolgerungen ziehen, müssen wir auf das große Problem dieses Benchmarks hinweisen - die meiste Zeit wird an sehr hohen Grenzen verbracht, was höchstwahrscheinlich großen Intervallen eine übermäßige Bedeutung beimisst. Deshalb müssen wir zum zweiten Benchmark gehen.

Benchmark Small-Shuffle


Dieser Benchmark ähnelt dem vorherigen, führt jedoch viel weniger „Array-Mixing“ (mehrere) durch. Code:

 for (uint32_t j = 0; j < 0xffff; ++j) { for (uint32_t i = 0xffff; i > 0; --i) { uint32_t bval = bounded_rand(rng, i); assert(bval < i); sum += bval; } } 

Mersenne-Wirbelergebnisse


Ergebnisse verschiedener PRNGs


Schlussfolgerungen

Dieser Benchmark vermeidet eine zu starke Betonung großer Grenzen und spiegelt die realen Anwendungsfälle genauer wider, verwirft jedoch jetzt große Grenzen vollständig.

Benchmark für alle Intervalle


Diese Benchmark zielt darauf ab, die Nachteile der beiden vorherigen zu vermeiden. Er führt Tests bei jeder Größe der Zweierpotenz durch, so dass jede Größe vorhanden ist, aber sein Einfluss wird nicht überschätzt.

 for (uint32_t bit = 1; bit != 0; bit <<= 1) { for (uint32_t i = 0; i < 0x1000000; ++i) { uint32_t bound = bit | (i & (bit - 1)); uint32_t bval = bounded_rand(rng, bound); assert(bval < bound); sum += bval; } } 

Mersenne-Wirbelergebnisse


Ergebnisse verschiedener PRNGs


Schlussfolgerungen

Viele unserer Ergebnisse bleiben unverändert. Die Skew-Methode ist schnell, wenn wir den Fehler ertragen können, und das Bitmaskenschema scheint eine gute gemittelte Wahl zu sein.

Wir könnten dies beenden, wenn wir nicht zurückkehren, unseren Code kritisch betrachten und Änderungen daran vornehmen möchten.

Verbesserungen vornehmen


Bis zu diesem Punkt erforderten alle Methoden zur Beseitigung von Versatz die Verwendung einer zusätzlichen Teilungsrestoperation, weshalb sie viel langsamer als Versatzmethoden durchgeführt werden. Es wäre hilfreich, wenn wir diesen Vorteil reduzieren könnten.

Schnellerer Schwellenwert-Drop


Einige unserer Algorithmen haben Code, der einen Schwellenwert verwendet, zum Beispiel:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { // calculates 2**32 % range uint32_t t = (-range) % range; for (;;) { uint32_t r = rng(); if (r >= t) return r % range; } } 

Wenn rangedie Anzahl im Vergleich zum PRNG-Ausgabeintervall klein ist, ist sie meistens viel größer als der Schwellenwert. Das heißt, wenn wir eine vorläufige Schätzung des Schwellenwerts hinzufügen können, die etwas höher sein kann, sparen wir den teuren Vorgang, den Rest der Division zu übernehmen.

Der folgende Code behandelt diese Aufgabe:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { uint32_t r = rng(); if (r < range) { uint32_t t = (-range) % range; while (r < t) r = rng(); } return r % range; } 

Diese Änderung kann sowohl auf "Double Mod ohne Verzerrungen" (siehe oben) als auch auf "Integer Multiplication ohne Verzerrungen" angewendet werden. Die Idee wurde von Lemir erfunden, der sie auf die zweite Methode anwendete (aber nicht auf die erste).

Large-Shuffle-Benchmark-Ergebnisse

Diese Optimierung führt zu einer signifikanten Verbesserung der Ergebnisse des 64-Bit-Benchmarks (bei dem der Mod noch langsamer ist), verschlechtert jedoch die Leistung im 32-Bit-Benchmark geringfügig. Trotz der Verbesserungen gewinnt die Bitmaskenmethode immer noch.


Small-Shuffle-Benchmark-Ergebnisse

Andererseits beschleunigt diese Änderung den Small-Shuffle-Benchmark sowohl für die Methode der ganzzahligen Multiplikation als auch für die Methode des doppelten Restes der Division erheblich. In beiden Fällen verschiebt sich ihre Leistung ohne Verzerrungen näher an die Ergebnisse der Optionen. Die Leistung der Doppelrestmethode (OpenBSD) entspricht jetzt fast der Leistung der Einzelrestmethode (Java).


Benchmark-Ergebnisse für alle Intervalle

Wir sehen eine ähnliche Verbesserung des Benchmarks für alle Intervalle.


Es sieht so aus, als könnten wir einen neuen universellen Gewinner bekannt geben: eine optimierte Methode zum Multiplizieren von Lemire-Ganzzahlen ohne Verzerrung.

Division Restoptimierung


Normalerweise a % berfordert eine Berechnung eine Division, aber in Situationen, in denen das a < bErgebnis einfach ist aund keine Division erforderlich ist. Und wenn a/2 < b, ist das Ergebnis einfach a - b. Daher statt zu rechnen

 a %= b; 

wir können erfüllen

 if (a >= b) { a -= b; if (a >= b) a %= b; } 

Die Kosten für die Teilung sind so hoch, dass eine Erhöhung der Kosten für diesen komplexeren Code sich durch Zeitersparnis aufgrund fehlender Teilung rechtfertigen kann.

Large-Shuffle-Benchmark-Ergebnisse

Durch Hinzufügen dieser Optimierung werden die Ergebnisse des Large-Shuffle-Benchmarks erheblich verbessert. Dies macht sich wiederum im 64-Bit-Code bemerkbar, bei dem die Übernahme des Restbetrags teurer ist. Die Doppelrestmethode (OpenBSD-Stil) zeigt Versionen mit Optimierungen für nur eine Restoperation und für beide.


In diesem Benchmark ist die Bitmaske immer noch der Gewinner, aber die Grenze zwischen ihr und Lemiras Ansatz hat sich erheblich verengt.

Small-Shuffle-Benchmark-Ergebnisse

Das Hinzufügen dieser Optimierung erhöht nicht die Leistung des Small-Shuffle-Benchmarks, sodass die Frage nur bleibt, ob dadurch erhebliche Kosten entstehen. In einigen Fällen nein, in anderen steigen die Kosten leicht an.


Benchmark-Ergebnisse für alle Intervalle

Im Benchmark für alle Intervalle sind die Änderungen ebenfalls gering.


Bonus: PRSP-Vergleichsergebnisse


Der Hauptgrund für die Verwendung vieler PRNGs zum Testen von Nummernschemata in Intervallen bestand darin, eine versehentliche Verzerrung der Ergebnisse aufgrund der Besonderheiten des Betriebs einzelner PRNG-Schemata zu vermeiden. Wir können jedoch dieselben Ergebnisse interner Tests verwenden, um die Generierungsschemata selbst zu vergleichen.

PRNG mit 32-Bit-Ausgabe

Die folgende Grafik zeigt die Leistung verschiedener 32-Bit-Generierungsschemata, gemittelt für alle Methoden und fünfzehn Läufe, normalisiert auf die Leistung des 32-Bit-Mersenne-Wirbels:


Einerseits bin ich froh zu sehen, dass es pcg32_fastwirklich schnell ist - es wurde nur von einer kleinen Version von Xoroshiro besiegt (die keine statistischen Tests besteht). Dies zeigt aber auch, warum ich mich wegen der Leistung moderner Hochleistungs-Universal-PRSPs selten aufrege - der Unterschied zwischen den verschiedenen Methoden ist sehr unbedeutend. Insbesondere die Leistung der schnellsten vier Schaltkreise unterscheidet sich um weniger als 5%, und ich glaube, dass dies einfach durch „Rauschen“ verursacht wird.

PRNG mit der Ausgabe von 64-Bit-Zahlen

Die Grafik zeigt die Leistung verschiedener 64-Bit-Generierungsschemata, gemittelt unter allen Techniken und fünfzehn Läufen, normalisiert auf die Leistung des 32-Bit-Mersenne-Wirbels. Es mag seltsam erscheinen, dass die Normalisierung mit dem 32-Bit-Mersenne-Wirbel durchgeführt wird, aber dies ermöglicht es uns, die zusätzlichen Kosten für die Verwendung der 64-Bit-Generierung in Fällen zu sehen, in denen die 32-Bit-Generierung ausreichend ist.


Diese Ergebnisse bestätigen, dass es mcg128_fastunglaublich schnell ist, aber die letzten vier Techniken unterscheiden sich wiederum nur um etwa 5%, sodass es schwierig ist, aus den schnellsten Methoden zu wählen. pcg64und pcg64_fastmüssen langsamer sein mcg128_fast, da ihre Basisgeneratoren lineare kongruente 128-Bit-Generatoren (LCG) und multiplikative kongruente 128-Bit-Generatoren (MCG, MCG) sind. Trotz der Tatsache, dass sie nicht die schnellsten Techniken in diesem Satz sind, sind sie immer pcg64noch 20% schneller als der 64-Bit-Mersenne-Wirbel.

Aber was vielleicht noch wichtiger ist, diese Ergebnisse zeigen auch, dass ein 64-Bit-PRNG normalerweise langsamer ist als ein 32-Bit-PRNG, wenn Sie keine 64-Bit-Ausgabe benötigen.

Schlussfolgerungen


An unseren Benchmarks können wir erkennen, dass der Übergang von standardmäßig verwendeten PRNGs (z. B. dem 32-Bit-Mersenne-Wirbel) zu schnelleren PRNPs die Ausführungszeit von Benchmarks um 45% reduzierte. Der Übergang von der Standardmethode zum Ermitteln der Zahl im Intervall zu unserer schnellsten Methode ermöglichte es uns jedoch, die Benchmark-Zeit um etwa 66% zu reduzieren. mit anderen Worten, bis zu einem Drittel der ursprünglichen Zeit.

Die schnellste Methode (ohne Verzerrungen) ist die Lemira-Methode (mit meiner zusätzlichen Optimierung). Da ist er:

 uint32_t bounded_rand(rng_t& rng, uint32_t range) { uint32_t x = rng(); uint64_t m = uint64_t(x) * uint64_t(range); uint32_t l = uint32_t(m); if (l < range) { uint32_t t = -range; if (t >= range) { t -= range; if (t >= range) t %= range; } while (l < t) { x = rng(); m = uint64_t(x) * uint64_t(range); l = uint32_t(m); } } return m >> 32; } 

Die Verwendung der Lemira-Methode verbessert die Leistung der meisten randomisierten Algorithmen mehr als der Wechsel von einer Engine mit schneller Generation zu einer schnelleren.

Anhänge: Testnotizen


Der Code aller Tests wird auf GitHub veröffentlicht . Insgesamt habe ich 23 Methoden zur bounded_randVerwendung von 26 verschiedenen PRNs (13 32-Bit-PRNs und 13 64-Bit-PRNs) in zwei Compilern (GCC 8 und LLVM 6) getestet , wodurch ich jeweils 26 * 23 * 2 = 1196 ausführbare Dateien erhielt Davon wurde es mit denselben 15 Samen durchgeführt, was 1196 * 15 = 17.940 eindeutige Testläufe ergibt, in denen jeweils drei Benchmarks kombiniert werden. Grundsätzlich habe ich Tests auf einem 48-Kern-Computer mit vier 2,1-GHz-Xeon E7-4830v3-Prozessoren durchgeführt. Das Durchführen eines vollständigen Satzes von Tests dauerte etwas weniger als einen Monat Prozessorzeit.

Am Ende kehren wir von der Einleitung des Artikels zur Situation zurück. Stellen Sie sich vor, Sasha jsf32.STD-libc++und Juan -mt19937.BIASED_FP_MULT_SCALE. In Benchmark 3 benötigt letzteres 69,6% weniger Zeit. Das heißt, die Zeit von dieser fiktiven Situation basiert auf Daten aus der Realität.

Source: https://habr.com/ru/post/de455702/


All Articles