Générer un entier aléatoire à partir d'une plage

J'ai besoin d'une fonction qui générerait un entier aléatoire dans une plage donnée (y compris les valeurs de bordure). Je n'ai pas d'exigences de qualité / aléatoire déraisonnables, j'ai quatre exigences:

• J'ai besoin que ce soit rapide. Mon projet doit générer des millions (voire des dizaines de millions) de nombres aléatoires et ma fonction de générateur actuelle s'est avérée être un goulot d'étranglement.
• J'ai besoin qu'il soit raisonnablement uniforme (l'utilisation de rand () est parfaitement bien).
• les plages min-max peuvent être n'importe quoi de <0, 1> à <-32727, 32727>.
• il doit pouvoir être semé.

J'ai actuellement le code C ++ suivant:

output = min + (rand() * (int)(max - min) / RAND_MAX)

Le problème est qu'il n'est pas vraiment uniforme - max est renvoyé uniquement lorsque rand () = RAND_MAX (pour Visual C ++, il est 1/32727). C'est un problème majeur pour les petites plages telles que <-1, 1>, où la dernière valeur n'est presque jamais retournée.

J'ai donc pris un stylo et du papier et suis venu avec la formule suivante (qui s'appuie sur l'astuce d'arrondi entier (int) (n + 0.5)):

Mais cela ne me donne toujours pas une distribution uniforme. Des analyses répétées avec 10000 échantillons me donnent un rapport de 37:50:13 pour les valeurs -1, 0. 1.

Pouvez-vous suggérer une meilleure formule? (ou même fonction de générateur de nombres pseudo-aléatoires entiers)

Une solution distribuée rapide, un peu meilleure que la vôtre, mais toujours pas correctement uniforme est

output = min + (rand() % static_cast<int>(max - min + 1))

Sauf lorsque la taille de la plage est une puissance de 2, cette méthode produit des nombres distribués non uniformes biaisés quelle que soit la qualité de rand(). Pour un test complet de la qualité de cette méthode, veuillez lire ceci .

La réponse C ++ la plus simple (et donc la meilleure) (en utilisant la norme 2011) est

#include <random>

std::random_device rd;     // only used once to initialise (seed) engine
std::mt19937 rng(rd());    // random-number engine used (Mersenne-Twister in this case)
std::uniform_int_distribution<int> uni(min,max); // guaranteed unbiased

auto random_integer = uni(rng);

Pas besoin de réinventer la roue. Pas besoin de s'inquiéter des préjugés. Pas besoin de s'inquiéter d'utiliser le temps comme une graine aléatoire.

Si votre compilateur prend en charge C ++ 0x et que son utilisation est une option pour vous, le nouvel en- <random>tête standard répondra probablement à vos besoins. Il a une haute qualité uniform_int_distributionqui acceptera les limites minimales et maximales (incluses selon vos besoins), et vous pouvez choisir parmi divers générateurs de nombres aléatoires pour vous connecter à cette distribution.

Voici le code qui génère un million de ints aléatoires uniformément répartis dans [-57, 365]. J'ai utilisé les nouvelles <chrono>installations std pour chronométrer, car vous avez mentionné que les performances sont une préoccupation majeure pour vous.

#include <iostream>
#include <random>
#include <chrono>

int main()
{
    typedef std::chrono::high_resolution_clock Clock;
    typedef std::chrono::duration<double> sec;
    Clock::time_point t0 = Clock::now();
    const int N = 10000000;
    typedef std::minstd_rand G;
    G g;
    typedef std::uniform_int_distribution<> D;
    D d(-57, 365);
    int c = 0;
    for (int i = 0; i < N; ++i) 
        c += d(g);
    Clock::time_point t1 = Clock::now();
    std::cout << N/sec(t1-t0).count() << " random numbers per second.\n";
    return c;
}

Pour moi (Intel Core i5 2,8 GHz), cela imprime:

2.10268e + 07 nombres aléatoires par seconde.

Vous pouvez amorcer le générateur en passant un int à son constructeur:

    G g(seed);

Si vous trouvez plus tard que int cela ne couvre pas la plage dont vous avez besoin pour votre distribution, cela peut être corrigé en changeant la uniform_int_distributionmême chose (par exemple en long long):

    typedef std::uniform_int_distribution<long long> D;

Si vous constatez plus tard que le minstd_randgénérateur n'est pas de qualité suffisante, il peut également être facilement remplacé. Par exemple:

    typedef std::mt19937 G;  // Now using mersenne_twister_engine

Avoir un contrôle séparé sur le générateur de nombres aléatoires et la distribution aléatoire peut être assez libérateur.

J'ai également calculé (non montré) les 4 premiers "moments" de cette distribution (en utilisant minstd_rand) et les ai comparés aux valeurs théoriques pour tenter de quantifier la qualité de la distribution:

min = -57
max = 365
mean = 154.131
x_mean = 154
var = 14931.9
x_var = 14910.7
skew = -0.00197375
x_skew = 0
kurtosis = -1.20129
x_kurtosis = -1.20001

(Le x_préfixe fait référence à "attendu")

Divisons le problème en deux parties:

• Générer un nombre aléatoire n entre 0 et (max-min).
• Ajouter min à ce nombre

La première partie est évidemment la plus difficile. Supposons que la valeur de retour de rand () soit parfaitement uniforme. L'utilisation de modulo ajoutera un biais aux premiers (RAND_MAX + 1) % (max-min+1)nombres. Donc , si nous pouvions changer comme par magie RAND_MAXàRAND_MAX - (RAND_MAX + 1) % (max-min+1) , il n'y aurait plus aucun parti pris.

Il s'avère que nous pouvons utiliser cette intuition si nous sommes prêts à autoriser le pseudo-non-déterminisme dans le temps d'exécution de notre algorithme. Chaque fois que rand () renvoie un nombre trop grand, nous demandons simplement un autre nombre aléatoire jusqu'à ce que nous en obtenions un qui soit assez petit.

Le temps d'exécution est maintenant géométriquement distribué , avec la valeur attendue 1/poù pest la probabilité d'obtenir un nombre suffisamment petit au premier essai. Puisque RAND_MAX - (RAND_MAX + 1) % (max-min+1)est toujours inférieur à (RAND_MAX + 1) / 2, nous le savons p > 1/2, donc le nombre d'itérations attendu sera toujours inférieur à deux pour n'importe quelle plage. Il devrait être possible de générer des dizaines de millions de nombres aléatoires en moins d'une seconde sur un processeur standard avec cette technique.

ÉDITER:

Bien que ce qui précède soit techniquement correct, la réponse de DSimon est probablement plus utile dans la pratique. Vous ne devriez pas implémenter cela vous-même. J'ai vu beaucoup d'implémentations d'échantillonnage de rejet et il est souvent très difficile de voir si c'est correct ou non.

Et le Mersenne Twister ? L'implémentation Boost est plutôt facile à utiliser et est bien testée dans de nombreuses applications du monde réel. Je l'ai moi-même utilisé dans plusieurs projets académiques tels que l'intelligence artificielle et les algorithmes évolutifs.

Voici leur exemple où ils créent une fonction simple pour lancer un dé à six faces:

#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_int.hpp>
#include <boost/random/variate_generator.hpp>

boost::mt19937 gen;

int roll_die() {
    boost::uniform_int<> dist(1, 6);
    boost::variate_generator<boost::mt19937&, boost::uniform_int<> > die(gen, dist);
    return die();
}

Oh, et voici un peu plus de proxénétisme de ce générateur juste au cas où vous ne seriez pas convaincu que vous devriez l'utiliser sur le bien inférieur rand():

Le Mersenne Twister est un générateur de "nombres aléatoires" inventé par Makoto Matsumoto et Takuji Nishimura; leur site Web comprend de nombreuses implémentations de l'algorithme.

Essentiellement, le Mersenne Twister est un très grand registre à décalage à rétroaction linéaire. L'algorithme fonctionne sur une graine de 19 937 bits, stockée dans un tableau de 624 éléments d'entiers non signés de 32 bits. La valeur 2 ^ 19937-1 est un nombre premier de Mersenne; la technique de manipulation de la graine est basée sur un ancien algorithme de "torsion" - d'où le nom "Mersenne Twister".

Un aspect attrayant du Mersenne Twister est son utilisation d'opérations binaires - par opposition à une multiplication chronophage - pour générer des nombres. L'algorithme a également une très longue période et une bonne granularité. Il est à la fois rapide et efficace pour les applications non cryptographiques.

int RandU(int nMin, int nMax)
{
    return nMin + (int)((double)rand() / (RAND_MAX+1) * (nMax-nMin+1));
}

Il s'agit d'un mappage de 32768 entiers en (nMax-nMin + 1) entiers. Le mappage sera assez bon si (nMax-nMin + 1) est petit (comme dans votre exigence). Notez cependant que si (nMax-nMin + 1) est grand, le mappage ne fonctionnera pas (par exemple, vous ne pouvez pas mapper 32768 valeurs à 30000 valeurs avec une probabilité égale). Si de telles plages sont nécessaires, vous devez utiliser une source aléatoire 32 bits ou 64 bits, au lieu de rand () 15 bits, ou ignorer les résultats rand () qui sont hors limites.

Voici une version non biaisée qui génère des nombres dans [low, high]:

int r;
do {
  r = rand();
} while (r < ((unsigned int)(RAND_MAX) + 1) % (high + 1 - low));
return r % (high + 1 - low) + low;

Si votre plage est raisonnablement petite, il n'y a aucune raison de mettre en cache le côté droit de la comparaison dans la doboucle.

Je recommande la bibliothèque Boost.Random , elle est très détaillée et bien documentée, vous permet de spécifier explicitement la distribution que vous voulez, et dans les scénarios non cryptographiques, elle peut en fait surpasser une implémentation typique de la bibliothèque C.

supposons que min et max sont des valeurs int, [et] signifie inclure cette valeur, (et) signifie ne pas inclure cette valeur, en utilisant ci-dessus pour obtenir la bonne valeur en utilisant c ++ rand ()

référence: pour () [] définir, visitez:

https://en.wikipedia.org/wiki/Interval_(mathematics)

pour la fonction rand et srand ou la définition RAND_MAX, visitez:

http://en.cppreference.com/w/cpp/numeric/random/rand

[min max]

int randNum = rand() % (max - min + 1) + min

(min max]

int randNum = rand() % (max - min) + min + 1

[min max)

int randNum = rand() % (max - min) + min

(min max)

int randNum = rand() % (max - min - 1) + min + 1

Dans ce thread, l'échantillonnage de rejet a déjà été discuté, mais je voulais suggérer une optimisation basée sur le fait qu'elle rand() % 2^somethingn'introduit aucun biais comme déjà mentionné ci-dessus.

L'algorithme est vraiment simple:

• calculer la plus petite puissance de 2 supérieure à la longueur de l'intervalle
• randomiser un nombre dans ce "nouvel" intervalle
• renvoie ce nombre s'il est inférieur à la longueur de l'intervalle d'origine
  • rejeter autrement

Voici mon exemple de code:

int randInInterval(int min, int max) {
    int intervalLen = max - min + 1;
    //now calculate the smallest power of 2 that is >= than `intervalLen`
    int ceilingPowerOf2 = pow(2, ceil(log2(intervalLen)));

    int randomNumber = rand() % ceilingPowerOf2; //this is "as uniform as rand()"

    if (randomNumber < intervalLen)
        return min + randomNumber;      //ok!
    return randInInterval(min, max);    //reject sample and try again
} 

Cela fonctionne bien en particulier pour les petits intervalles, car la puissance de 2 sera "plus proche" de la longueur réelle de l'intervalle, et donc le nombre d'échecs sera plus petit.

PS
Évidemment, éviter la récursivité serait plus efficace (pas besoin de calculer au-delà du plafond du journal ..) mais j'ai pensé que c'était plus lisible pour cet exemple.

Notez que dans la plupart des suggestions, la valeur aléatoire initiale que vous avez obtenue de la fonction rand (), qui est généralement de 0 à RAND_MAX, est simplement gaspillée. Vous n'en créez qu'un seul nombre aléatoire, alors qu'il existe une procédure sonore qui peut vous en donner plus.

Supposons que vous vouliez une région [min, max] de nombres aléatoires entiers. Nous partons de [0, max-min]

Prendre la base b = max-min + 1

Commencez par représenter un nombre obtenu de rand () en base b.

De cette façon, vous avez floor (log (b, RAND_MAX)) parce que chaque chiffre de la base b, sauf peut-être le dernier, représente un nombre aléatoire dans la plage [0, max-min].

Bien sûr, le décalage final vers [min, max] est simple pour chaque nombre aléatoire r + min.

int n = NUM_DIGIT-1;
while(n >= 0)
{
    r[n] = res % b;
    res -= r[n];
    res /= b;
    n--;
}

Si NUM_DIGIT est le nombre de chiffres en base b que vous pouvez extraire et c'est

NUM_DIGIT = floor(log(b,RAND_MAX))

alors ce qui précède est une implémentation simple d'extraction de NUM_DIGIT nombres aléatoires de 0 à b-1 à partir d'un nombre aléatoire RAND_MAX fournissant b <RAND_MAX.

La formule pour cela est très simple, alors essayez cette expression,

 int num = (int) rand() % (max - min) + min;  
 //Where rand() returns a random number between 0.0 and 1.0

L'expression suivante doit être impartiale si je ne me trompe pas:

std::floor( ( max - min + 1.0 ) * rand() ) + min;

Je suppose ici que rand () vous donne une valeur aléatoire comprise entre 0,0 et 1,0 NON compris 1,0 et que max et min sont des entiers avec la condition que min <max.