Nehmen wir an, wir haben einen binären Zufallszahlengenerator, int r(); , der eine Null oder eine Eins mit einer Wahrscheinlichkeit von 0,5 zurückgibt.
Ich habe Boost.Random angeschaut und sie erzeugen, sagen wir, 32 Bits und machen so etwas (Pseudocode):
%Vor%Ich habe ernsthafte Zweifel daran. Ein Double hat 53 Bit Mantisse, und 32 Bit können nie richtig eine völlig zufällige Mantisse erzeugen, unter anderem Rundungsfehler usw.
Was wäre ein schneller Weg, um einen gleichmäßig verteilten float oder double im halboffenen Bereich [min, max) zu erstellen, unter der Annahme von IEEE754? Der Schwerpunkt liegt hier auf Korrektheit der Verteilung , nicht auf Geschwindigkeit.
Um korrekt zu definieren, wäre die korrekte Verteilung gleich der, die wir bekommen würden, wenn wir einen unendlich genau gleichverteilten Zufallszahlengenerator nehmen würden und für jede Zahl würden wir auf die nächste IEEE754-Darstellung runden, wenn diese Darstellung wäre immer noch in [min, max) , sonst würde die Zahl nicht für die Verteilung zählen.
PS .: Ich würde mich auch für die richtigen Lösungen für offene Bereiche interessieren.
Hier ist ein korrekter Ansatz ohne Effizienzsteigerung.
Wir beginnen mit einer Bignum-Klasse und dann mit einem rationalen Wrapper dieser Bignums.
Wir erzeugen einen Bereich, der "ausreichend größer als" unser [min, max) -Bereich ist, so dass die Rundung unserer smaller_min und bigger_max Fließkommawerte außerhalb dieses Bereichs erzeugt, in unserem auf dem Bignum basierenden Rational.
Jetzt unterteilen wir den Bereich in zwei Teile perfekt in der Mitte (was wir tun können, da wir ein rationales Bignum-System haben). Wir wählen einen der beiden Teile zufällig.
Wenn nach dem Runden der obere und untere Bereich des ausgewählten Bereichs (A) außerhalb von [min, max) liegt (auf der gleichen Seite, wohlgemerkt!), lehnen Sie ab und starten von vorne.
Wenn (B) der obere und untere Bereich Ihres Bereichs auf den gleichen double (oder float , wenn Sie einen Float zurückgeben) rundet, sind Sie fertig und Sie geben diesen Wert zurück.
Ansonsten (C) rekrutieren Sie sich auf diesen neuen, kleineren Bereich (Unterteilen, Zufallsauswahl, Test).
Es gibt keine Garantie dafür, dass diese Prozedur angehalten wird, da Sie entweder kontinuierlich zwischen zwei Rundungen double s bis zur "Kante" vorbohren können, oder Sie können ständig Werte außerhalb des Bereichs [min, max) auswählen. Die Wahrscheinlichkeit dafür ist (niemals haltend), jedoch Null (unter der Annahme eines guten Zufallszahlengenerators und einer [min, max) von Nicht-Null-Größe).
Dies funktioniert auch für (min, max) , oder wählt sogar eine Zahl in der gerundeten, ausreichend dicken Cantor-Menge aus. Solange das Maß des gültigen Bereichs von Realen, der auf die korrekten Gleitkommawerte aufrundet, ungleich Null ist und der Bereich eine kompakte Unterstützung hat, kann diese Prozedur ausgeführt werden und hat eine Wahrscheinlichkeit von 100% der Beendigung, aber keine harte Obergrenze gebunden an die benötigte Zeit kann gemacht werden.
Das Problem hier ist, dass in IEEE754 die Doubles, die dargestellt werden können, nicht gleich verteilt sind. Das heißt, wenn wir einen Generator haben, der reelle Zahlen erzeugt, sagen wir in (0,1) und dann IEEE754-darstellbaren Zahlen zuordnet, wird das Ergebnis nicht gleich verteilt.
Also müssen wir "equi-distribution" definieren. Das heißt, unter der Annahme, dass jede IEEE754-Nummer nur ein Repräsentant für die Wahrscheinlichkeit ist, in dem durch die IEEE754-Rundung definierten Intervall zu liegen, erzeugt die Prozedur der ersten Erzeugung gleichverteilter "Zahlen" und der Runde nach IEEE754 (per Definition) ein " Equi-Distribution "von IEEE754-Nummern.
Daher glaube ich, dass die obige Formel beliebig nahe einer solchen Verteilung wird, wenn wir nur die Genauigkeit hoch genug wählen. Wenn wir das Problem darauf beschränken, eine Zahl in [0,1] zu finden, bedeutet dies, dass man sich auf den Satz von denominierten IEEE 754-Zahlen beschränkt, die Eins-zu-Eins zu einer 53-Bit-Ganzzahl sind. Daher sollte es schnell und korrekt sein, nur die Mantisse durch einen 53-Bit-Binär-Zufallsgenerator zu erzeugen.
Die IEEE 754-Arithmetik ist immer "Arithmetik mit unendlicher Genauigkeit gefolgt von Rundung", dh die IEEE754-Zahl, die ein b darstellt, ist diejenige, die einem b am nächsten kommt (anders ausgedrückt, man kann sich ein * vorstellen) b mit unendlicher Genauigkeit berechnet, dann auf die geschlossene IEEE754-Zahl gerundet). Daher glaube ich, dass min + (max-min) * x, wobei x eine nummerierte Zahl ist, ein praktikabler Ansatz ist.
(Anmerkung: Wie aus meinem Kommentar hervorgeht, war mir zuerst nicht bewusst, dass Sie auf den Fall mit min und max von 0,1 verweisen. Die denormalisierten Zahlen haben die Eigenschaft, dass sie gleichmäßig verteilt sind Durch die Abbildung der 53 Bits auf die Mantisse können Sie die Fließkomma-Arithmetik verwenden, da sie bis zur Maschinenpräzision korrekt ist.Wenn Sie das umgekehrte Mapping verwenden, werden Sie die Gleichverteilung wiederherstellen.
Siehe diese Frage für einen anderen Aspekt dieses Problems: Scaling Int uniform Zufallsbereich in Double One
std::uniform_real_distribution .
Es gibt ein wirklich gutes Gespräch von S.T.L. von der diesjährigen Going Native-Konferenz, die erklärt, warum Sie die Standard-Distributionen wann immer möglich verwenden sollten. Kurz gesagt, handgewalzter Code tendiert dazu, von lächerlich schlechter Qualität zu sein (think std::rand() % 100 ), oder hat subtilere Uniformitätsfehler, wie in (std::rand() * 1.0 / RAND_MAX) * 99 , das ist das Beispiel, das in der Rede gegeben wird, und ist ein Spezialfall der Code in der Frage geschrieben.
EDIT: Ich habe mir die Implementierung von std::uniform_real_distribution von libstdc ++ angeschaut, und das habe ich gefunden:
Die Implementierung erzeugt eine Zahl im Bereich [dist_min, dist_max) , indem eine einfache lineare Transformation von einer Zahl verwendet wird, die im Bereich [0, 1) erzeugt wird. Er generiert diese Quellennummer mit Hilfe von std::generate_canonical , deren Implementierung hier zu finden ist (am Ende der Datei). std::generate_canonical bestimmt die Anzahl der Male (bezeichnet als k ) Der Bereich der Verteilung , ausgedrückt als Integer und hier als r * bezeichnet, passt in die Mantisse des Zieltyps . Was es dann tut, ist im Wesentlichen eine Zahl in [0, r) für jedes r -große Segment der Mantisse zu erzeugen und unter Verwendung von Arithmetik jedes Segment entsprechend aufzufüllen. Die Formel für den resultierenden Wert kann als
Dabei ist X eine stochastische Variable in [0, r) . Jede Division durch den Bereich entspricht einer Verschiebung um die Anzahl von Bits, die verwendet werden, um sie darzustellen (d. H.% Co_de%) und füllt so das entsprechende Mantissensegment. Auf diese Weise wird die gesamte Genauigkeit des Zieltyps verwendet, und da der Bereich des Ergebnisses log2(r) ist, bleibt der Exponent [0, 1) ** (modulo bias) und Sie erhalten nicht die Gleichmäßigkeitsprobleme, die Sie haben wenn Sie mit dem Exponenten herumspielen.
Ich würde nicht glauben, dass diese Methode kryptographisch sicher ist (und ich habe Verdacht auf mögliche Fehler bei der Berechnung der Größe von 0 ), aber ich kann mir vorstellen, dass es wesentlich zuverlässiger ist Einheitlichkeit als die Boost-Implementierung, die Sie gepostet haben, und definitiv besser als herumspielen mit r .
Es ist erwähnenswert, dass der Boost-Code tatsächlich ein degenerierter Fall dieses Algorithmus ist, wobei std::rand bedeutet, dass er äquivalent ist wenn der Eingabebereich mindestens 23 Bits benötigt Größe (IEE 754 Single-Precision) oder mindestens 52 Bit (Double-Precision). Dies bedeutet eine minimale Reichweite von ~ 8,4 Millionen bzw. ~ 4,5e15. Angesichts dieser Informationen glaube ich nicht, dass die Boost-Implementierung, wenn Sie einen Binärgenerator verwenden, ziemlich dafür ausgelegt ist, ihn zu schneiden.
Nach einem kurzen Blick auf die Implementierung von libc ++ sieht es so aus, als ob sie was verwenden ist der gleiche Algorithmus, etwas anders implementiert.
(*) k = 1 ist eigentlich der Bereich der Eingabe plus eins . Dies ermöglicht die Verwendung des r -Wertes des Urng als gültige Eingabe.
(**) Streng genommen ist der codierte Exponent nicht max , da IEEE 754 eine implizite führende 1 vor der Radix des Signifikanden codiert. Konzeptionell ist dies jedoch für diesen Algorithmus irrelevant.
AFAIK, der korrekte (und wahrscheinlich auch schnellste) Weg besteht darin, zuerst eine 64-Bit-Ganzzahl ohne Vorzeichen zu erzeugen, wobei die 52 Bruchbits zufällige Bits sind und der Exponent 1023 ist, was bei Eingabe in ein (IEEE 754) -Verdoppel sei ein gleichmäßig verteilter Zufallswert im Bereich [1.0, 2.0]. Der letzte Schritt besteht darin, 1,0 davon zu subtrahieren, was zu einem gleichmäßig verteilten zufälligen Doppelwert im Bereich [0.0, 1.0] führt.
Im Pseudocode:
rndDouble = bitCastUInt64ToDouble (1023 & lt; & lt; 52 | rndUInt64 & amp; 0xfffffffffffff) - 1,0
Diese Methode wird hier erwähnt: Ссылка (Siehe "Generieren von gleichmäßigen Doubles im Einheitenintervall")
BEARBEITEN: Die empfohlene Methode wurde seitdem geändert in: (x & gt; & gt; 11) * (1. / (UINT64_C (1) & lt; & lt; 53))
Siehe obigen Link für Details.