Génération efficace de nombres dans un intervalle donné

image

La grande majorité de mes articles sur la génération de nombres aléatoires traitaient principalement des propriétés de divers schémas de génération. Cela peut s'avérer inattendu, mais les performances de l'algorithme de randomisation peuvent dépendre non pas du schéma de génération choisi, mais d'autres facteurs. Dans cet article (qui m'a inspiré un excellent article de Daniel Lemyr ), nous étudierons les principales raisons de la baisse des performances de génération de nombres aléatoires, qui l'emportent souvent sur les performances du moteur PRN.

Imaginez cette situation:

Comme devoirs, Juan et Sasha implémentent le même algorithme randomisé en C ++, qui s'exécutera sur le même ordinateur universitaire et avec un seul ensemble de données. Leur code est presque identique et ne diffère que par la génération de nombres aléatoires. Juan est pressé pour ses cours de musique, alors il a simplement choisi le tourbillon de Mersenne. Sasha, d'autre part, a passé quelques heures de recherche supplémentaires. Sasha a effectué des comparaisons de plusieurs des PRNG les plus rapides, qu'il a récemment appris des réseaux sociaux, et a choisi le plus rapide. Lors de la réunion, Sasha était impatient de se vanter et il a demandé à Juan: "Quel système PRNG avez-vous utilisé?"

"Personnellement, je viens de prendre le vortex de Mersenne - il est intégré dans la langue et semble fonctionner assez bien."

"Ha!", Répondit Sasha. «J'ai utilisé jsf32 . Il est bien plus rapide que l'ancien et lent tourbillon de Mersenne! Mon programme se déroule en 3 minutes 15 secondes! »

"Hmm, pas mal, mais le mien peut le faire en moins d'une minute", dit Juan en haussant les épaules. «Eh bien, je dois aller au concert. Voulez-vous venir avec moi? "

"Non", répond Sasha. "J'ai ... euh ... besoin de revoir mon code."

Cette situation fictive gênante n'est pas particulièrement fictive; il est basé sur des résultats réels. Si votre algorithme randomisé ne fonctionne pas aussi vite que nous le souhaiterions et que le goulot d'étranglement semble être la génération de nombres aléatoires, alors, curieusement, le problème peut ne pas être dans le générateur de nombres aléatoires!

Introduction: les nombres aléatoires dans la pratique


La plupart des générateurs de nombres aléatoires modernes de haute qualité créent des mots machine remplis de bits aléatoires, c'est-à-dire qu'ils génèrent généralement des nombres dans l'intervalle [0..2 32 ) ou [0..2 64 ). Mais dans de nombreux cas d'utilisation, les utilisateurs ont besoin de nombres dans un certain intervalle - par exemple, pour lancer un dé ou choisir une carte à jouer au hasard, les nombres sont nécessaires en petits intervalles constants. Cependant, de nombreux algorithmes, du mélange et de l' échantillonnage du réservoir aux arbres de recherche binaire randomisés, nécessitent des nombres provenant d'autres intervalles.

Les méthodes


Nous examinerons de nombreuses méthodes différentes. Pour simplifier la discussion, au lieu de générer des nombres dans l'intervalle [ i .. j ) ou [ i .. j ], nous allons générer des nombres dans l'intervalle [0 .. k ). Avec un tel schéma, nous pouvons, par exemple, générer des nombres dans l'intervalle [ i .. j ) en fixant k = j - i , en générant un nombre dans l'intervalle [0 .. k ), puis en y ajoutant i .

Outils C ++ intégrés


De nombreuses langues ont des outils intégrés pour obtenir un nombre aléatoire dans un intervalle spécifié. Par exemple, pour supprimer une carte d'un jeu de 52 cartes dans des langages de script tels que Perl et Python, nous pouvons écrire respectivement int(rand(52)) et random.randint(0,52) . [Remarque Utilisateur de CryptoPirate : Il me semble qu'une erreur ici, en Python randint (a, b) génère des nombres de a à b dont b. Et comme il y a 52 cartes dans le jeu et que la première est «0», ce devrait être random.randint (0,51) .] En C ++, nous pouvons utiliser uniform_int_distribution même uniform_int_distribution .

Le code C ++ pour implémenter cette approche est simple:

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

Habituellement, l'une des techniques décrites ci-dessous est utilisée dans les outils intégrés, mais la plupart des utilisateurs utilisent simplement ces outils, sans penser à ce qui se passe «sous le capot», estimant que ces outils sont correctement conçus et assez efficaces. En C ++, les outils intégrés sont plus complexes car ils devraient pouvoir fonctionner avec des moteurs de génération assez arbitraires - un générateur qui produit des valeurs dans la plage de -3 à 17 peut être tout à fait valide et peut être utilisé avec std::uniform_int_distribution pour créer des nombres dans n'importe quel intervalle, par exemple [0..1000). Autrement dit, les outils C ++ intégrés sont trop compliqués pour la plupart des cas où ils sont utilisés.

Le reste classique de la division (asymétrique)


Passons d'une approche trop simplifiée à une approche trop simpliste.

Lorsque j'ai étudié la programmation, nous avons généré des nombres dans l'intervalle (par exemple, pour sélectionner une carte dans un jeu de 52 cartes) en utilisant l'opérateur restant. Pour obtenir le nombre dans l'intervalle [0..52), nous avons écrit rand() % 52 .

En C ++, cette approche peut être implémentée comme suit:

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

Malgré la simplicité de cette approche, elle démontre la raison pour laquelle obtenir des nombres dans le bon intervalle est généralement une tâche lente - elle nécessite une division (pour calculer le reste obtenu par l'opérateur % ). La division est généralement au moins un ordre de grandeur plus lente que les autres opérations arithmétiques, de sorte qu'une seule opération arithmétique prend plus de temps que tout le travail effectué par un PRNG rapide.

Mais en plus de la faible vitesse, elle est également asymétrique . Pour comprendre pourquoi rand() % 52 renvoie des nombres asymétriques, supposons que rand() crée des nombres dans l'intervalle [0..2 32 ), et notez que 52 ne divise pas complètement 2 32 , il le divise 82 595 524 fois avec le reste 48. Autrement dit, si nous utilisons rand() % 52 , nous aurons 82 595 525 façons de sélectionner les 48 premières cartes du jeu et seulement 82 595 524 façons de sélectionner les quatre dernières cartes. En d'autres termes, il y a un biais de 0,00000121% par rapport à ces quatre dernières cartes (peut-être que ce sont des rois!). Quand j'étais étudiant et que j'écrivais des devoirs sur le lancer de dés ou le tirage de cartes, personne ne se souciait généralement de ces petites distorsions, mais avec une augmentation de l'intervalle, la distorsion croît de façon linéaire. Pour un PRNG 32 bits, un intervalle limité de moins de 2 24 a un biais de moins de 0,5%, mais au-dessus de 2 31 un décalage de 50% - certains nombres reviendront deux fois plus souvent que d'autres.

Dans cet article, nous considérerons principalement les techniques qui utilisent des stratégies pour éliminer une erreur systématique, mais il vaut probablement la peine de dire que pour un PRNG 64 bits, la valeur de biais dans les applications normales est probablement négligeable.

Un autre problème peut être que certains générateurs ont des bits faibles faibles. Par exemple, les familles GPRS Xoroshiro + et Xoshiro + ont des bits bas qui ne passent pas les tests statistiques. Lorsque nous exécutons % 52 (car 52 est pair), nous passons le bit de poids faible directement à la sortie.

Multiplier les nombres à virgule flottante (asymétrique)


Une autre technique courante consiste à utiliser un PRNG qui génère des nombres à virgule flottante dans l'intervalle [0..1) avec la conversion ultérieure de ces nombres à l'intervalle souhaité. Cette approche est utilisée en Perl, il est recommandé d' utiliser int(rand(10)) pour générer un entier dans l'intervalle [0..10) en générant un nombre à virgule flottante suivi d'un arrondi.

En C ++, cette approche est écrite comme ceci:

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

(Notez que 0x1.0p-32 est une constante binaire à virgule flottante pour 2 -32 , que nous utilisons pour convertir un entier aléatoire dans l'intervalle [0..2 32 ) pour doubler dans l'intervalle unitaire; à la place, nous pouvons effectuer une telle conversion en utilisant ldexp(rng(), -32) , mais quand j'ai évalué cette approche, cela s'est avéré beaucoup plus lent.)

Cette approche est aussi biaisée que le reste classique de la division, mais le biais apparaît différemment. Par exemple, si nous sélectionnions des nombres dans l'intervalle [0..52), alors les nombres 0, 13, 26 et 39 se produiraient une fois moins souvent que les autres.

Cette version, lors de la généralisation à 64 bits, est encore plus désagréable, car elle nécessite un type à virgule flottante dont la mantisse est d'au moins 64 bits. Sur les machines x86 avec Linux et macOS, nous pouvons utiliser le long double pour tirer parti des nombres à virgule flottante x86 de précision accrue qui ont une mantisse 64 bits, mais le long double pas universellement porté sur tous les systèmes - dans certains systèmes, le long double équivaut au double .

Il y a un bon côté - cette approche est plus rapide que les solutions résiduelles pour les PRNG avec des bits faibles faibles.

Multiplication entière (asymétrique)


La méthode de multiplication peut être adaptée à l'arithmétique à virgule fixe plutôt qu'à flottante. En fait, nous multiplions simplement constamment par 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; } 

Il peut sembler que cette version nécessite une arithmétique 64 bits, sur les processeurs x86, un bon compilateur compilera ce code dans une instruction mult 32 bits (ce qui nous donne deux valeurs de sortie 32 bits, dont l'une est la valeur de retour). On peut s'attendre à ce que cette version soit rapide, mais elle est asymétrique exactement comme la méthode de multiplication des nombres à virgule flottante.

Division de baisse (pas de biais)


Nous pouvons modifier le schéma de multiplication à virgule flottante en schéma basé sur la division. Au lieu de multiplier x * range / 2**32 nous calculons x / (2**32 / range) . Comme nous travaillons avec l'arithmétique des nombres entiers, l'arrondi dans cette version sera effectué différemment et générera parfois des valeurs en dehors de l'intervalle souhaité. Si nous rejetons ces valeurs (par exemple, nous en débarrassons et générons de nouvelles valeurs), nous obtenons ainsi une technique sans distorsion.

Par exemple, dans le cas du retrait d'une carte à l'aide d'un PRNG 32 bits, nous pouvons générer un nombre 32 bits et le diviser par 2 32/52 = 82 595 524 pour sélectionner une carte. Cette technique fonctionne si la valeur aléatoire du PRNG 32 bits est inférieure à 52 × 82595524 = 2 32/32 - 48. Si la valeur aléatoire du PRNR est l'une des 48 dernières valeurs de la partie supérieure de l'intervalle du générateur, vous devez alors la supprimer et en chercher une autre.

Notre code pour cette version utilise une astuce avec la division de 2 32 par range sans utiliser de mathématiques 64 bits. Pour un calcul direct de 2**32 / range nous devons représenter le nombre 2 32 , qui est trop grand (par un!) Pour représenter un entier 32 bits. Au lieu de cela, nous prenons en compte que pour les entiers non signés, la range opération de négation unaire calcule une valeur positive de 2 32 - range ; en divisant cette valeur par range , nous obtenons une réponse inférieure à 2**32 / range .

Par conséquent, le code C ++ pour générer des nombres à l'aide de la division et de la suppression ressemble à ceci:

 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; } } 

Bien sûr, cette approche nécessite deux opérations lentes basées sur la division, qui sont généralement plus lentes que les autres opérations arithmétiques, vous ne devez donc pas vous attendre à ce qu'elle soit rapide.

Le reste de la division (double) sans distorsions - technique OpenBSD


Nous pouvons également adopter l'approche de baisse pour éliminer le biais dans la méthode du reste de division classique. Dans l'exemple avec les cartes à jouer, nous devons à nouveau supprimer 48 valeurs. Dans cette version, au lieu de supprimer les 48 dernières valeurs, nous rejetons (de manière équivalente) les 48 premières valeurs.

Voici l'implémentation de cette approche en 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; } } 

Cette technique supprime le biais, mais elle nécessite deux opérations de division chronophages avec le reste de chaque valeur de sortie (et vous pouvez avoir besoin d'un générateur interne pour créer plusieurs nombres). Par conséquent, il faut s'attendre à ce que la méthode soit environ deux fois plus lente que l'approche asymétrique classique.

arc4random_uniform OpenBSD arc4random_uniform (qui est également utilisé sur OS X et iOS) utilise cette stratégie.

Reste de division (simple) sans biais - méthodologie Java


Java utilise une approche différente pour générer un nombre dans un intervalle qui utilise une seule opération de division restante, à l'exception de cas assez rares de rejet du résultat. 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; } 

Pour comprendre pourquoi cette option fonctionne, vous devez réfléchir un peu. Contrairement à la version précédente, basée sur les résidus, qui élimine les biais en supprimant une partie des valeurs les plus basses du moteur de génération interne, cette version filtre les valeurs de la partie supérieure de l'intervalle du moteur.

Multiplication d'entiers asymétriques - Méthode Lemira


De la même manière que nous avons supprimé le biais de la méthode du reste de la division, nous pouvons éliminer le biais de la technique de multiplication des nombres entiers. Cette technique a été inventée par Lemyr.

 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 bitmask (no skew) - Apple technique


Dans notre dernière approche, la division et les opérations restantes sont complètement éliminées. Au lieu de cela, il utilise une opération de masquage simple pour obtenir un nombre aléatoire dans l'intervalle [0..2 k ), où k est la plus petite valeur, telle que 2 k est supérieur à l'intervalle. Si la valeur est trop grande pour notre intervalle, nous la rejetons et essayons d'en obtenir une autre. Le code est illustré ci-dessous:

 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; } 

Cette approche a été adoptée par Apple lorsque (dans la version macOS Sierra) elle a effectué sa propre révision du code arc4random_uniform .

Analyse comparative des techniques de base


Nous avons maintenant plusieurs approches qui peuvent être évaluées. Malheureusement, lorsque nous sommes préoccupés par les coûts d'une opération à division unique, l'analyse comparative devient une chose non triviale. Aucune référence ne peut prendre en compte tous les facteurs affectant le domaine d'application, et rien ne garantit que la meilleure option pour votre application sera certainement la meilleure pour la mienne.

Nous utilisons trois repères et testons les techniques avec de nombreux PRNG différents.

Benchmark Large-Shuffle


La référence la plus évidente est probablement le mixage. Dans cette référence, nous simulons un mixage à grande échelle. Pour trier un tableau de taille N, nous devons générer des nombres dans les intervalles [0 .. N ), [0 .. ( N -1)), ..., [0..1). Dans ce cas-ci, nous supposerons que N est le nombre maximum possible (pour uint32_t il est 2 32 -1). Code:

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

Notez que nous «utilisons» chaque nombre en l'ajoutant à la sum (afin qu'il ne soit pas jeté par optimisation), mais nous n'effectuons aucun mixage pour nous concentrer sur la génération de nombres.

Pour tester la génération 64 bits, nous avons un test similaire, mais il ne sera pas pratique d'effectuer un test correspondant au mélange d'un tableau de taille 2 64 - 1 (car il faudra plusieurs milliers d'années pour compléter ce plus grand benchmark). Au lieu de cela, nous traversons l'intégralité de l'intervalle 64 bits, mais générons le même nombre de valeurs de sortie que dans le test 32 bits. 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; } 

Résultats du tourbillon de Mersenne

Les résultats ci-dessous montrent les performances de cette référence pour chacune des méthodes que nous avons examinées lors de l'utilisation du vortex de Mersenne et du test sur 32 bits (en utilisant std::mt19937 de libstdc++ ) et du code 64 bits similaire (en utilisant std:mt19937_64 de libstdc++ ) Les résultats sont la moyenne géométrique de 15 séries avec différentes valeurs de graine, qui est ensuite normalisée de sorte que la méthode du reste de la division classique ait un seul temps d'exécution.


Il peut sembler que nous ayons des réponses claires sur les performances - il semble que vous pouvez construire des techniques pour leur perfection et vous demander à quoi pensaient les développeurs de libstdc++ lorsqu'ils ont écrit une implémentation aussi terrible pour les nombres 32 bits. Mais, comme c'est souvent le cas avec l'analyse comparative, la situation est plus compliquée qu'il n'y paraît d'après ces résultats. Premièrement, il existe un risque que les résultats soient spécifiques au vortex de Mersenne, nous allons donc étendre les nombreux PRNG testés. Deuxièmement, il peut y avoir un problème subtil avec l'indice de référence lui-même. Voyons d'abord la première question.

Résultats des différents PRNG

Nous testerons des arc4_rand32 32 bits avec les arc4_rand32 , chacha8r , gjrand32 , jsf32 , mt19937 , pcg32 , pcg32_fast , sfc32 , splitmix32 , xoroshiro64+ , xorshift*64/32 xoshiro128+ , xoshiro128+ et xoshiro128** , et xoshiro128** gjrand64 jsf64 , mcg128 , mcg128_fast , mt19937_64 , pcg64 , pcg64_fast , sfc64 , splitmix64 , xoroshiro128+ , xorshift*128/64 xoshiro256+ , xoshiro256+ et xoshiro256* . Ces kits nous donneront des PRN lents et beaucoup de très rapides.

Voici les résultats:


Nous pouvons voir les principales différences par rapport aux résultats avec le vortex de Mersenne. Des PRNG plus rapides déplacent l'équilibre vers le code de délimitation, et donc la différence entre les différentes approches devient plus prononcée, en particulier dans le cas des PRNR 64 bits. Avec un ensemble plus large de libstc++ implémentation de libstc++ cesse de paraître si terrible.

Conclusions

Dans ce benchmark d'une marge importante, l'approche basée sur la multiplication avec biais gagne en vitesse. Il existe de nombreuses situations dans lesquelles les limites seront petites par rapport à la taille du PRNG, et les performances sont absolument critiques. Dans de telles situations, un léger biais est peu susceptible d'avoir un effet notable, mais la vitesse PRNG aura. Un tel exemple est Quicksort avec un point de référence aléatoire. Parmi les méthodes biaisées, la technique du masque bitmap semble prometteuse.

Mais avant de tirer des conclusions sérieuses, nous devons souligner l'énorme problème de cette référence - la plupart du temps est consacré à des limites très élevées, ce qui donne très probablement une importance excessive à de grands intervalles. Par conséquent, nous devons passer au deuxième repère.

Benchmark Small-Shuffle


Ce benchmark est similaire au précédent, mais effectue beaucoup moins de «mixage matriciel» (multiple). 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; } } 

Résultats du tourbillon de Mersenne


Résultats des différents PRNG


Conclusions

Cette référence évite de trop insister sur les grandes frontières et reflète plus précisément les cas d'utilisation réels, mais rejette désormais complètement les grandes frontières.

Benchmark pour tous les intervalles


Cette référence vise à éviter les inconvénients des deux précédents; il effectue des tests à chaque taille de la puissance de deux pour que chaque taille soit présente, mais son influence n'est pas surestimée.

 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; } } 

Résultats du tourbillon de Mersenne


Résultats des différents PRNG


Conclusions

Bon nombre de nos constatations demeurent inchangées. La méthode asymétrique est rapide si nous pouvons supporter l'erreur, et le schéma de masque de bits semble être un bon choix moyen.

Nous pourrions mettre fin à cela si nous ne voulions pas revenir en arrière, jeter un regard critique sur notre code et y apporter des modifications.

Apportez des améliorations


Jusqu'à ce point, toutes les méthodes d'élimination de l'inclinaison nécessitaient l'utilisation d'une opération de division résiduelle supplémentaire, c'est pourquoi elles sont exécutées beaucoup plus lentement que les méthodes d'inclinaison. Il serait utile de réduire cet avantage.

Chute basée sur un seuil plus rapide


Certains de nos algorithmes ont du code utilisant une valeur de seuil, par exemple:

 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; } } 

Lorsqu'il est rangepetit par rapport à l'intervalle de sortie PRNG, le plus souvent le nombre sera beaucoup plus grand que le seuil. Autrement dit, si nous pouvons ajouter une estimation préliminaire du seuil, qui peut être un peu plus, nous économiserons sur l'opération coûteuse de prendre le reste de la division.

Le code suivant gère cette tâche:

 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; } 

Cette modification peut être appliquée à la fois au «double Mod sans distorsions» (voir ci-dessus) et à la «multiplication entière sans distorsions». L'idée a été inventée par Lemir, qui l'a appliquée à la deuxième méthode (mais pas à la première).

Résultats de référence à grande échelle

Cette optimisation conduit à une amélioration significative des résultats du benchmark 64 bits (dans lequel le mod est encore plus lent), mais dégrade en fait légèrement les performances du benchmark 32 bits. Malgré les améliorations, la méthode bitmask gagne toujours.


Résultats des tests de référence à petite échelle

D'un autre côté, ce changement accélère considérablement le repère de petit shuffle pour la méthode de multiplication entière et la méthode du double reste de division. Dans les deux cas, leurs performances se rapprochent davantage des résultats des options sans distorsions. Les performances de la méthode des résidus doubles (OpenBSD) sont désormais quasiment égales aux performances de la méthode des résidus simples (Java).


Résultats de référence pour tous les intervalles

Nous constatons une amélioration similaire de l'indice de référence pour tous les intervalles.


Il semble que nous pouvons annoncer un nouveau gagnant universel: une méthode optimisée pour multiplier les entiers Lemire sans biais.

Optimisation du reste de la division


En règle générale, un calcul a % bnécessite une division, mais dans les situations où le a < brésultat est simple aet la division n'est pas requise. Et quand a/2 < b, le résultat est simple a - b. Par conséquent, au lieu de calculer

 a %= b; 

nous pouvons accomplir

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

Le coût de la division est si important que l'augmentation du coût de ce code plus complexe peut se justifier par un gain de temps dû au manque de division.

Résultats de référence à grande échelle

L'ajout de cette optimisation améliore considérablement les résultats du test de référence à grand shuffle. Ceci est encore plus visible dans le code 64 bits, où l'opération de prise du reste est plus coûteuse. La méthode du double reste (style OpenBSD) affiche les versions avec optimisations pour une seule opération restante et pour les deux.


Dans cette référence, le masque de bits est toujours le gagnant, mais la frontière entre celui-ci et l'approche de Lemira s'est considérablement rétrécie.

Résultats des tests de référence à petite échelle

L'ajout de cette optimisation n'augmente pas les performances de la référence de petit shuffle, donc la question reste de savoir si elle ajoute des coûts importants. Dans certains cas, non; dans d'autres, les coûts augmentent légèrement.


Résultats de référence pour tous les intervalles

Dans la référence pour tous les intervalles, les changements sont également faibles.


Bonus: résultats de la comparaison du DSRP


La raison principale de l'utilisation d'un grand nombre de PRNG pour tester les schémas de numérotation à intervalles était d'éviter la distorsion involontaire des résultats en raison des particularités du fonctionnement des schémas PRNG individuels. Mais nous pouvons utiliser les mêmes résultats de tests internes pour comparer les schémas de génération eux-mêmes.

PRNG avec sortie 32 bits

Le graphique ci-dessous montre les performances de divers schémas de génération 32 bits, moyennées pour toutes les méthodes et quinze exécutions, normalisées aux performances du vortex Mersenne 32 bits:


D'une part, je suis heureux de voir que c'est pcg32_fastvraiment rapide - il n'a été vaincu que par une petite version de Xoroshiro (qui ne passe pas les tests statistiques). Mais cela montre aussi pourquoi je me fâche rarement à cause des performances des DSRP polyvalents modernes à hautes performances - la différence entre les différentes méthodes est très insignifiante. En particulier, les quatre circuits les plus rapides présentent des performances inférieures à 5%, et je pense que cela est simplement dû au «bruit».

PRNG avec sortie de nombres 64 bits

Le graphique montre les performances de divers schémas de génération 64 bits en moyenne parmi toutes les techniques et quinze cycles normalisés aux performances du vortex Mersenne 32 bits. Il peut sembler étrange que la normalisation soit effectuée à l'aide du vortex Mersenne 32 bits, mais cela nous permet de voir les coûts supplémentaires liés à l'utilisation de la génération 64 bits dans les cas où la génération 32 bits est suffisante.


Ces résultats confirment qu'elle est mcg128_fastincroyablement rapide, mais les quatre dernières techniques ne diffèrent à nouveau que d'environ 5%, il est donc difficile de choisir parmi les méthodes les plus rapides. pcg64et pcg64_fastdoit être plus lent mcg128_fast, car leurs générateurs de base sont des générateurs congruents linéaires 128 bits (LCG) et des générateurs congruents multiplicatifs 128 bits (MCG, MCG). Bien qu'ils ne soient pas les techniques les plus rapides de cet ensemble, ils sont pcg64toujours 20% plus rapides que le vortex Mersenne 64 bits.

Mais peut-être plus important encore, ces résultats montrent également que si vous n'avez pas besoin d'une sortie 64 bits, un PRNG 64 bits est généralement plus lent qu'un 32 bits.

Conclusions


À partir de nos benchmarks, nous pouvons voir que la transition des PRNG standard (par exemple, le vortex Mersenne 32 bits) vers des PRNP plus rapides a réduit le temps d'exécution des benchmarks de 45%. Mais le passage de la méthode standard de recherche du nombre dans l'intervalle à notre méthode la plus rapide nous a permis de réduire le temps de référence d'environ 66%; en d'autres termes, jusqu'à un tiers de l'heure d'origine.

La méthode la plus rapide (sans distorsions) est la méthode Lemira (avec mon optimisation supplémentaire). Le voici:

 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; } 

L'utilisation de la méthode Lemira améliorera davantage les performances de la plupart des algorithmes randomisés que le passage d'un moteur de génération rapide à un moteur plus rapide.

Annexes: notes de test


Le code de tous les tests est publié sur GitHub . Au total, j'ai testé 23 méthodes pour bounded_randutiliser 26 PRN différents (13 PRN 32 bits et 13 64 bits), dans deux compilateurs (GCC 8 et LLVM 6), ce qui m'a donné 26 * 23 * 2 = 1196 fichiers exécutables, chacun dont il a été effectué avec les mêmes 15 semences, ce qui donne 1196 * 15 = 17 940 séries de tests uniques, dans chacune desquelles trois repères sont combinés. Fondamentalement, j'ai effectué des tests sur une machine à 48 cœurs avec quatre processeurs Xeon E7-4830v3 à 2,1 GHz. L'exécution d'un ensemble complet de tests a pris un peu moins d'un mois de temps processeur.

Au final, nous revenons à la situation depuis l'introduction de l'article. Imaginez que Sasha a utilisé jsf32.STD-libc++, et Juan -mt19937.BIASED_FP_MULT_SCALE. Dans le benchmark 3, ce dernier prend 69,6% de temps en moins. Autrement dit, le temps de cette situation fictive est basé sur des données de la réalité.

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


All Articles