guter zufall (multiply-with-carry prng)

Blue-Tiger

Basher schrieb:

Tim schrieb:

Mal ne andere Frage wie bekomme ich von einem "Zufalls"generator x "Zufalls"zahlen im Bereich [10 100)?

meinste das: http://www.maths.abdn.ac.uk/~igc/tch/mx4002/notes/node83.html ?

Oder sowas: http://www.mathematik.uni-ulm.de/stochastik/lehre/ss03/markov/skript/node28.html ?

SeppJ

@hustbaer: In diesem Fall eine ungümstige Methode, da du unnötigerweise Zahlen wegwirfst. Blue-Tigers Link zeigt auf ein paar bessere, sehr allgemeine Methoden, die jedoch für den Fall einer Gleichverteilung [10,100) wie Kanonen auf Spatzen sind. Wenn man (wie im Computer üblich) von einer gleichmäßigen Verteilung von Zufallszahlen X im Bereich [0,RMAX) ausgeht, bekommt man diese direkt mittels (X/RMAX)*90+10 auf den fraglichen Bereich transformiert.

hustbaer

@SeppJ:
Nö, das kannst du knicken. Üblicherweise werfen Zufallszahlengeneratoren Ganzzahlen aus. Und wenn du deine Transformation auf Ganzzahlen loslässt (und das Ergebnis wieder auf Ganzzahlen "beschneidest"), dann zerstörst du damit die Gleichverteilung.
(Mal davon abgesehen dass bei Integer-Arithmetik die Division nach der Multiplikation kommen muss, damit überhaupt irgendwas annähernd sinnvolles rauskommen kann.)

Das Ergebnis ist gleich schlecht, wie wenn du Modulo (ohne Verwerfen) als Transformation verwendest. Der einzige Unterschied ist, dass die "bevorzugten" Zahlen sich nicht wie bei Modulo alle am Anfang der Range sammeln, sondern schön über die ganze Range verteilt werden. Bringt mir aber auch nicht viel, wenn ich will, dass alle Zahlen in der Output-Range genau die gleiche Wahrscheinlichkeit haben.

----

Ich persönlich glaube dass das Verwerfen von Zahlen die einfachste und daher beste Methode ist.

Natürlich gäbe es auch andere Möglichkeiten die die Gleichverteilung nicht zerstören. z.B. könnte man einen "Range-Coder" (den Decoder-Teil genauergenommen) verwenden, und damit einen Stream aus zufälligen Bits dekodieren. (Die Wahrscheinlichkeiten/Häufigkeiten würde man dazu alle gleich annehmen, um die gewünschte Gleichverteilung zu erreichen.)

http://de.wikipedia.org/wiki/Arithmetische_Kodierung

Das würde ich dann aber wirklich mit Kanonen auf Spatzen schiessen nennen.

----

BTW: die Implementierung in der Boost macht genau das was ich vorhin skizziert habe: Zahlen verwerfen.

Basher_ist_nicht_BashAr

hustbaer schrieb:

Nö, das kannst du knicken. Üblicherweise werfen Zufallszahlengeneratoren Ganzzahlen aus.

hihi, vielleicht war [a b**)** nur ein tippfehler. damit war weder das mathematische 'halboffene intervall' noch reelle zahlen gemeint (das komma fehlt ja auch).

SeppJ

hustbaer schrieb:

@SeppJ:
Nö, das kannst du knicken. Üblicherweise werfen Zufallszahlengeneratoren Ganzzahlen aus. Und wenn du deine Transformation auf Ganzzahlen loslässt (und das Ergebnis wieder auf Ganzzahlen "beschneidest"), dann zerstörst du damit die Gleichverteilung.

Das sehe ich überhaupt nicht ein, kannst du das mal erläutern?

X aus [0,RMAX) gleichverteilt, meinetwegen nur Ganzzahlen

=> X/RMAX in [0,1) gleichverteilt (und man MUSS vorher dividieren, da normalerweise RMAX auch gleichzeitig die größte darstellbare Zahl ist). Ich nehme mal an, jeder im Forum weiß, wie man mit einem Computer irchtig dividiert.

=> (X/RMAX) * 90 in [0,90) gleichverteilt, jedoch als Fließkommazahlen

=> (X/RMAX) * 90 + 10 in [10,100) gleichverteilt, jedoch als Fließkommazahlen

Jetzt kann man natürlich nur den ganzzahligen Teil der Zahlen nehmen. Da die Werte aber gleichverteilt waren, wird dabei kein bestimmter Wert bevorzugt.

volkard

SeppJ schrieb:

Jetzt kann man natürlich nur den ganzzahligen Teil der Zahlen nehmen. Da die Werte aber gleichverteilt waren, wird dabei kein bestimmter Wert bevorzugt.

Falsch.
Wenn ich 100 Gäste habe und nur 90 Zimmer, MÜSSEN gelegentliche Doppelbelegungen vorkommen und die Gleichverteilung ist kaputt.
Aber auch wenn ich 1000000 Betten habe und 90 Zimmer, kann es keine Gleichverteilung geben, weil 90 einfach kein Teiler von 1000000 ist.
Wenn ich aber nur 10 Gäste ablehne (mit einem Einnahmeverlust von (10/1000000) einem Tausendstel Prozent, kann ich die verbleibenden 999990 Gäste prima auf 90 Zimmer gleichverteilen.

volkard

SeppJ schrieb:

Das sehe ich überhaupt nicht ein, kannst du das mal erläutern?

Das Skriptum war ein wenig Weltfremd. Das ist normal.

X aus [0,RMAX) gleichverteilt, meinetwegen nur Ganzzahlen

OK. Also wie rand().

=> X/RMAX in [0,1) gleichverteilt

Nö. Nur RMAX ausgesuchte Werte kommen vor und nicht alle reellen Zahlen zwischen 0 und 1. Das ist ein riesiger Unterschied.

=> (X/RMAX) * 90 in [0,90) gleichverteilt, jedoch als Fließkommazahlen

Nicht haltbar, da die Vorbedingung schon falsch war.

=> (X/RMAX) * 90 + 10 in [10,100) gleichverteilt, jedoch als Fließkommazahlen

Nicht haltbar, da die Vorbedingung schon falsch war.

Jetzt kann man natürlich nur den ganzzahligen Teil der Zahlen nehmen. Da die Werte aber gleichverteilt waren, wird dabei kein bestimmter Wert bevorzugt.

Nicht haltbar, da die Vorbedingung schon falsch war.

SeppJ

Ah, jetzt verstehe ich, was ihr gemeint habt. Ja, da habt ihr wohl Recht. Ich habe nicht bedacht, dass ich mit der Transformation nur (quasi-)kontinuierliche Verteilungen erzeugen kann, was für Ganzzahlen zwischen 10 und 100 nicht gegeben ist.

Aber es muss doch was besseres geben als das was hustbaer vorgeschlagen hat. Wenn ich im ungünstigsten Fall Ganzzahlen zwischen 0 und RMAX/2+1 suche, müsste ich ja fast die Hälfte der Zahlen wegwerfen. Mir fällt spontan jedoch auch nichts besseres ein.

life

Wenn Java das auch so macht, wirds schon richtig sein ;):

public int nextInt(int n)

    Returns a pseudorandom, uniformly distributed int value between 0 (inclusive) and the specified value (exclusive), drawn from this random number generator's sequence. The general contract of nextInt is that one int value in the specified range is pseudorandomly generated and returned. All n possible int values are produced with (approximately) equal probability. The method nextInt(int n) is implemented by class Random as follows:

         public int nextInt(int n) {
             if (n<=0)
        		throw new IllegalArgumentException("n must be positive");

             if ((n & -n) == n)  // i.e., n is a power of 2
                 return (int)((n * (long)next(31)) >> 31);

             int bits, val;
             do {
                 bits = next(31);
                 val = bits % n;
             } while(bits - val + (n-1) < 0);
             return val;
         }

    The hedge "approximately" is used in the foregoing description only because the next method is only approximately an unbiased source of independently chosen bits. If it were a perfect source of randomly chosen bits, then the algorithm shown would choose int values from the stated range with perfect uniformity.

    The algorithm is slightly tricky. It rejects values that would result in an uneven distribution (due to the fact that 2^31 is not divisible by n). The probability of a value being rejected depends on n. The worst case is n=2^30+1, for which the probability of a reject is 1/2, and the expected number of iterations before the loop terminates is 2.

    The algorithm treats the case where n is a power of two specially: it returns the correct number of high-order bits from the underlying pseudo-random number generator. In the absence of special treatment, the correct number of low-order bits would be returned. Linear congruential pseudo-random number generators such as the one implemented by this class are known to have short periods in the sequence of values of their low-order bits. Thus, this special case greatly increases the length of the sequence of values returned by successive calls to this method if n is a small power of two.

    Parameters:
        n - the bound on the random number to be returned. Must be positive. 
    Returns:
        a pseudorandom, uniformly distributed int value between 0 (inclusive) and n (exclusive). 
    Throws:
        IllegalArgumentException - n is not positive.
    Since:
        1.2

volkard

SeppJ schrieb:

Aber es muss doch was besseres geben als das was hustbaer vorgeschlagen hat. Wenn ich im ungünstigsten Fall Ganzzahlen zwischen 0 und RMAX/2+1 suche, müsste ich ja fast die Hälfte der Zahlen wegwerfen. Mir fällt spontan jedoch auch nichts besseres ein.

Du sagst "Nur die Hälfte? :("
Sag mal "Nur die Hälfte! :)"

Außerdem braucht man das selten. Ungefähr immer geht rand()%49, außer ganz manchmal, dann nimmt man rand()/double(RAND_MAX+1)*49, und nur wenn das auch nicht geht, nimmt man die 1337-Version. Normalerweise stört an der Modulo-Version nämlich gar nicht die Ungleichverteilung, sondern die Schieflage, weil sich die volleren Zimmer links sammeln und dort Slums entstehen, da bringt die double-Version schon Abhilfe.

volkard

//Multiply-With-Carry Pseudo Random Number Generator
class Random {
private:
	static UInt64 const a=4294967118;
	UInt64 x;
public:
	Random(UInt64 seed=0) {
		x=seed+!seed;
	}
	UInt32 operator()() {
		x=a*(x&0xffffffff)+(x>>32);
		return x;
	}
};

//Numerical Receipes sagt:
// Wenn der Seed in 1 <= x <= 2^32-1 liegt 
// und sowohl a*b-1 als auch (a*b-2)/2 Primzahlen sind, 
// dann ist die Periodenlänge gleich (a*b-2)/2. 
//Volkard sagt: 
// Die größte Zahl, die diese Bedingungen erfüllt, ist 
// obiges a=2^32-178 mit beinahe 2^63 als Periodenlänge.

Basher_ist_nicht_BashAr

nimm doch den: http://burtleburtle.net/bob/rand/isaacafa.html