Tirage rapide entre bitsets.

Bonjour.

Je fais des simulations individus centrées pour le boulot (comprenne qui pourra). Lors de ces simulations, les individus se reproduisent, donnent des enfants et transmettent leur gènes.
J'ai modélisé une série de gènes via des bitset, et lors de la reproduction, un enfant vas donc "hériter" (à ne pas confondre avec ce que l'on appelle héritage en c++), d'une des versions (A ou B) que chacun de ses parents porte (et qui en fait provienne de ses grands parents si on remonte...).
Grosso modo, j'ai donc besoin de remplir un bitset, en choisissant une fois sur deux une valeur issu d'un autre bitset (A ou B).
J'ai finalement opté pour un remplissage complet de mon bitset de destination par l'un des deux possibles (la A), puis par une boucle, avec un tirage aléatoire qui décide si je remplace cette valeur par l'autre possible( la B).
Cependant, si cette méthode est légèrement plus rapide que mon "if else" précédent, cela fait plus que doubler le temps de simulation (comparé à si je remplis simplement avec le A).
Ce morceau de code me parait très simple est difficilement optimisable, mais comme je n'ai pas une formation d'informaticien, je passe peut être (surement ?) à côté de quelque chose et je vous soumet donc le problème.
Pour précision, j'utilise la bibliothèque GSL (peut être améliorable aussi de ce côté là).

Code:

---

1 2 3 4 5 6 7 8 9 10 11

neutres_A = maman.neutres_A; neutres_B = papa.neutres_A;   for (int i=0;i<32;i++) { if (gsl_rng_uniform_int(gene,2)==0) { neutres_A[i] = maman.neutres_B[i]; } if (gsl_rng_uniform_int(gene,2)==0) { neutres_B[i] = papa.neutres_B[i]; } }

---

bonjour
tout d'abord concernant l'indentation, sache que pour présenter du code (en gardant son indentation) tu peux utiliser la balise code (sans l'espace juste avant le crochet fermant qui m'a servie à empêcher l'interprétation de la balise pour l'afficher).

Code:

---

1 2 3 4 5 6

[CODE ] ton code avec indentation respectée [/code ]

---

C'est normal que ton exécution soit doublée, à la place d'utiliser une instruction pour remplir ta variable, tu vas utiliser plusieurs instructions pour définir dans quel cas tu te trouves et quel valeur affecter.

Si tu n'as pas besoin de ce facteur aléatoire, tu peux utiliser une ternaire qui sera optimisée par le compilo :)

Hello,

Je ne connais pas la lib gsl, mais si cette instruction : gsl_rng_uniform_int prend beaucoup de temps, bah tu peux deja en enlever une des 2 : tu fais deux fois le MEME test sur gsl_rng_uniform_int(gene,2)

Donc remplace comme ceci :

Code:

---

1 2 3 4 5 6 7 8

for (int i = 0; i < 32; i++) { if (gsl_rng_uniform_int(gene,2) ==0 ) { neutres_A[i] = maman.neutres_B[i]; neutres_B[i] = papa.neutres_B[i]; } }

---

NB: je vois pas en quoi 'gene' depend de "i", donc s'il ne depend pas de 'i', tu peux meme gagner 63x le temps d execution de gsl_rng_uniform_int en faisant comme suivant :

Code:

---

1 2 3 4 5 6 7 8 9 10 11 12 13

// si gsl_uniform depend pas de i, tu peux meme faire le test en dehors // de ta boucle for, ca t'evite de tester 32 fois au lieu d'une ... if (gsl_mg_uniform_int(gene,2) == 0) { for (int i = 0; i < 32; i++) { // if (gsl_rng_uniform_int(gene,2)==0) // { neutres_A[i] = maman.neutres_B[i]; neutres_B[i] = papa.neutres_B[i]; // } } }

---

Ptet j'ai mal compris ton probleme, mais j'pense que ca devrait etre mieux.. Aussi, si tu te sers beaucoup de gsl_mg_uniform_int(gene,2), mets le dans une variable annexe pour pas avoir besoin de le recalculer a chaque fois, et recalcule le juste si tu sais qu il a ete modifie ;)

Ciao

Merci pour vos réponse. Mais je me suis sans doute mal expliqué.

J'ai besoin de ce double tirage, puisque chaque "gene" présent sur le bitset (ici au nombre de 32) doit être obtenu du père ou de la mère, mais indépendamment les uns des autres.
Il y a peut être une confusion avec le fait que dans : gsl_mg_uniform_int(gene,2) il y a gene, mais c'est simplement l'abréviation de "générateur" car la fonction demande à chaque fois son générateur.

Pour simplifié, j'ai un bitset A, un bitset B, et je veux remplir un bitset C de même longueur avec pour chaque bit soit ce qu'il y a dans A, soit ce qu'il y a dans B. Et je cherche quelque chose de plus rapide (quitte à ce que ce soit une fourberie). Mais qui ne dénature pas le fait que ce qu'il va y avoir en position X de mon bitset est "indépendant" de ce qu'il y a en position X+1.
(relativement indépendant puisque ça viens potentiellement du même endroit...)

Je ne sais pas si c'est plus clair comme cela. Je ne sais pas ce qu'est un ternaire par contre.

Si tu veux utiliser de l'aléatoire, regarde du coter de la fonction rand().

Elle permet d'avoir un nombre random, tu appliques un modulo dessus, et hop c'est plier :):

Code:

---

1 2 3 4 5 6 7 8 9

int i = rand()%2; if (i == 0) { //un petit traitement } else { //deuxieme traitement }

---

Bonsoir,

tu parles de bits de partout mais je ne vois pas le moindre opérateur de bits ou bitset ? :calim2:

J'imagine que neutre_A et neutre_B sont des bitsets.

L'avantage d'utiliser un bitset est la compression de tes données en mémoire : la plus petite taille pour une variable en C/C++ est en général de 8 bits, ce qui veut dire que si tu utilisais à la place un tableau de booléens, celui-ci occuperait 8 fois plus de place en mémoire que le bitset. C'est donc très avantageux si tu es limité par la mémoire disponible. En revanche, cette limitation de 8 bits minimum par variable n'est pas là pour rien : c'est l'unité avec laquelle travaille le CPU (en réalité ça sera même un multiple de 8, selon l'architecture : 32 ou 64 bits). Manipuler des variables de 8 bits est ainsi beaucoup plus rapide que la magouille utilisée par le bitset, et l'optimum en terme de performance est de manipuler des variables de taille machine, soit 32 ou 64 bits.

Dans ton cas, il existe une possibilité pour avoir le meilleur des deux mondes. L'astuce va consister à essayer de ne jamais manipuler individuellement les éléments du bitset. Je m'explique.

Prenons un cas simple où tu n'as que 8 gènes. Tu as donc deux std::bitset nommés A et B de 8 éléments, qui vont occuper chacun 8 bits en mémoire, soit 1 octet (la taille d'un char). Tu veux ensuite construire un nouveau bitset de la même taille, qui, pour chaque gène, va piocher ses valeurs aléatoirement dans A et B.

Un cas concret simple :

Code:

---

1 2 3

A = {1,0,0,1,1,1,0,1} B = {0,1,0,0,1,0,1,0} => C = {0,1,0,0,1,1,0,1}

---

Ici le fils a récupéré la première moitié de ses gènes de B, et la seconde moitié de A. De manière générale, tu peux représenter cette opération par l'intermédiaire d'un masque M de sorte que :

Code:

---

C = (A & M) | (B & ~M);

---

Ce qui est équivalent à :

Code:

---

1 2 3 4 5

for (int i = 0; i < A.size(); ++i) { C[i] = (A[i] && M[i]) || (B[i] && !M[i]); // ou encore, plus simplement : C[i] = M[i] ? A[i] : B[i]; }

---

Pour notre cas particulier, M = {0,0,0,0,1,1,1,1}. Utiliser la première écriture sera bien plus rapide car, au lieu de faire explicitement une opération sur chaque bit, le bitset peut utiliser d'autres opérations qui vont agir sur un groupe de N bits en même temps, et diviser donc par N le nombre d'instructions CPU pour un même calcul. Si tu utilises std::bitset, tu peux déjà écrire le code ci-dessus sans problème : c'est implémenté de base.

Le dernier point technique va être de générer un bitset aléatoire rapidement. Ici, au lieu de générer un unsigned long aléatoire entre 0 et 1 pour chaque bit individuellement (64 bits de perdus!), on peut utiliser tous les bits émis par un appel du générateur aléatoire entre 0 et (unsigned long)(-1) et les ranger au fur et a mesure dans le bitset jusqu'à ce qu'il soit rempli. Un exemple de code :

Code:

---

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23

template<std::size_t N, typename Gen> std::bitset<N> random_bitset(Gen& g) { std::uniform_int_distribution<unsigned long> dist(0, -1);   // 'b' est initialisé à zéro std::bitset<N> b;   for (std::size_t i = 0; i < ceil(N/(8.0*sizeof(unsigned long))); ++i) { // On génère les sizeof(unsigned long) premiers bits de 'tmp' aléatoirement std::bitset<N> tmp(dist(g));   // On décale ces bits pour qu'ils se trouvent à un endroit où on n'a encore rien généré tmp <<= 8*i*sizeof(unsigned long);   // Et on range le tout dans le bitset final b |= tmp;   // Autrement dit : // b |= std::bitset<N>(dist(g)) << 8*i*sizeof(unsigned long); }   return b; }

---

Si tu fais tout ça, tu as gagné du temps pour deux raison :

1. tu laisse std::bitset optimiser pour toi les opérations bits par bits (il va convertir ça en une poignée d'instructions CPU pour la totalité du calcul)
2. tu appelles le générateur de nombres aléatoires moins souvent

J'ai fait un benchmark pour tester. Avec ces astuces, j'obtiens un code 5 fois plus rapide (en optimisant aussi la fonction de génération pour le cas particulier où ton bitset contient moins de bits qu'un unsigned long, ce qui est le cas dans ton exemple avec 32 gènes). Si tu veux tester chez toi (ça compile avec gcc 4.7 et -std=c++0x) :

Code:

---

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106

#include <iostream> #include <fstream> #include <bitset> #include <random> #include <chrono> #include <array>   // ------------------------------ // Générateur de bitset aléatoire // ------------------------------   // Version générique, pour un gros bitset template<std::size_t N, typename Gen> std::bitset<N> random_bitset_impl(Gen& g, std::true_type) { std::uniform_int_distribution<unsigned long> dist(0, -1); std::bitset<N> b;   for (std::size_t i = 0; i < ceil(N/(8.0*sizeof(unsigned long))); ++i) { b |= std::bitset<N>(dist(g)) << 8*i*sizeof(unsigned long); }   return b; }   // Version optimisée, pour les petits bitsets (taille < sizeof(unsigned long)) template<std::size_t N, typename Gen> std::bitset<N> random_bitset_impl(Gen& g, std::false_type) { std::uniform_int_distribution<unsigned long> dist(0, -1); return std::bitset<N>(dist(g)); }   // Fonction interface, qui va choisir entre la version optimisée ou la version générique template<std::size_t N, typename Gen> std::bitset<N> random_bitset(Gen& g) { return random_bitset_impl<N>(g, std::integral_constant<bool,(N > 8*sizeof(unsigned long))>()); }   // Génère juste 'true' ou 'false' template<typename Gen> bool random_bool(Gen& g) { std::uniform_int_distribution<unsigned long> dist(0, 1); return dist(g); }   // -------------------------------------------- // Un petit raccourcis pour mesurer le temps :) // --------------------------------------------   double now() { return std::chrono::duration_cast<std::chrono::microseconds>( std::chrono::high_resolution_clock::now().time_since_epoch() ).count()*1e-6; }   // -------------- // Programme test // --------------   int main() { std::mt19937 seed(time(0));   { // Vérifie que tous les bits ont bien une chance sur deux d'être 0 ou 1 const std::size_t N = 130; const std::size_t M = 10000; std::array<std::size_t,N> counts; for (std::size_t i = 0; i < M; ++i) { std::bitset<N> b = random_bitset<N>(seed); for (std::size_t j = 0; j < N; ++j) { counts[j] += b[j]; } }   for (std::size_t j = 0; j < N; ++j) { std::cout << counts[j]/double(M) << std::endl; // doit être proche de 0.5 } }   // Papa et maman std::bitset<8> a("10011101"); std::bitset<8> b("01001010");   // Version optimisée double start = now(); for (std::size_t i = 0; i < 100000; ++i) { std::bitset<8> m = random_bitset<8>(seed); std::bitset<8> c = (a & m) | (b & ~m); }   std::cout << now() - start << std::endl; // 0.017097 chez moi   // Version de base start = now(); for (std::size_t i = 0; i < 100000; ++i) { std::bitset<8> c = a; for (std::size_t j = 0; j < 8; ++j) { if (random_bool(seed)) { c[j] = b[j]; } } }   std::cout << now() - start << std::endl; // 0.089474 chez moi   return 0; }

---

Sans passer par std::bitset (en gérant toi même les bits dans un tableau de unsigned long par exemple), je pense tu peux coder quelque chose d'encore plus rapide, mais ça demande plus d'investissements (à la fois pour écrire le truc, mais aussi pour s'assurer que tout fonctionne bien).

Dernière remarque : je ne sais pas quel est ton but au final, mais en règle générale si tu veux des résultats qui soient statistiquement propres, il faut éviter d'utiliser des générateurs aléatoires trop simples (comme rand(), qui n'est bon que pour des usages basiques (et encore, même là on peut avoir des surprises !)). Ici j'ai choisi un Mersenne twister : c'est ce qu'on recommande souvent comme bon compromis entre rapidité et qualité. Il est disponible dans la STL avec le standard C+11, mais tu dois avoir ça aussi dans la GSL.

Merci (beaucoup) pour cette réponse très complète que je vais regarder de prêt et essayer de transposer dans mon petit code. Et merci d'avoir compris qu'en effet les "neutres" étaient les bitsets.
J'essayerais de penser à venir raconter le résultat final et à passer en "résolu".