Quel est le moyen le plus rapide de renvoyer les positions de tous les bits définis dans un entier de 64 bits?

Question

J'ai besoin d'un moyen rapide pour obtenir la position de tous les bits d'un entier de 64 bits. Par exemple, x = 123703, je voudrais remplir un tableau à l' idx[] = {0, 1, 2, 4, 5, 8, 9, 13, 14, 15, 16}. Nous pouvons supposer que nous savons que le nombre de bits à priori. Ce sera appelée à 10^12 - 10^15 fois, donc, la vitesse est de l'essence. Le plus rapide de la réponse que j'ai trouvé jusqu'à présent est la suivante monstruosité, qui utilise chaque octet de l'entier de 64 bits comme un index dans les tables qui donnent le nombre de bits de cet octet et les positions de ceux:

int64_t x;            // this is the input
unsigned char idx[K]; // this is the array of K bits that are set
unsigned char *dst=idx, *src;
unsigned char zero, one, two, three, four, five;  // these hold the 0th-5th bytes
zero  =  x & 0x0000000000FFUL;
one   = (x & 0x00000000FF00UL) >> 8;
two   = (x & 0x000000FF0000UL) >> 16;
three = (x & 0x0000FF000000UL) >> 24;
four  = (x & 0x00FF00000000UL) >> 32;
five  = (x & 0xFF0000000000UL) >> 40;
src=tab0+tabofs[zero ]; COPY(dst, src, n[zero ]);
src=tab1+tabofs[one  ]; COPY(dst, src, n[one  ]);
src=tab2+tabofs[two  ]; COPY(dst, src, n[two  ]);
src=tab3+tabofs[three]; COPY(dst, src, n[three]);
src=tab4+tabofs[four ]; COPY(dst, src, n[four ]);
src=tab5+tabofs[five ]; COPY(dst, src, n[five ]);

où COPY est une instruction switch pour copier jusqu'à 8 octets, n tableau est le nombre de bits dans un octet, et tabofs donne le décalage en tabX, qui détient les positions des bits du X-ème octet. C'est environ 3x plus rapide que le déroulé de la boucle des méthodes basées sur avec __builtin_ctz() sur mon Xeon E5-2609. (Voir ci-dessous). Je suis actuellement en itérant x dans le vocabulaire de la commande pour un nombre donné de bits définis.

Est-il un meilleur moyen?

EDIT: Ajout d'un exemple (que j'ai ensuite fixé). Le code complet est disponible ici: http://pastebin.com/79X8XL2P . Remarque: avec GCC -O2 semble pour l'optimiser à l'écart, mais d'Intel compilateur (que j'ai utilisé pour composer il) n'est pas...

Aussi, permettez-moi de vous donner quelques informations supplémentaires pour répondre à certains des commentaires ci-dessous. L'objectif est d'effectuer un test statistique sur chaque sous-ensemble de K variables d'un univers de N variables explicatives; l'objectif spécifique est maintenant de N=41, mais je peux voir certains projets nécessitant N jusqu'à 45-50. Le test consiste à factoriser les données correspondantes submatrix. En pseudo-code, quelque chose comme ceci:

double doTest(double *data, int64_t model) {
  int nidx, idx[];
  double submatrix[][];
  nidx = getIndices(model, idx);  // get the locations of ones in model
  // copy data into submatrix
  for(int i=0; i<nidx; i++) {
    for(int j=0; j<nidx; j++) {
      submatrix[i][j] = data[idx[i]][idx[j]];
    }
  }
  factorize(submatrix, nidx);
  return the_answer;
}

J'ai codé une version de ce pour un processeur Intel Phi conseil que doivent remplir les N=41 cas dans environ 15 jours, dont ~5-10% du temps est consacré à un naïf getIndices() donc dès le départ une version plus rapide pourrait sauver une journée ou plus. Je suis en train de travailler sur une mise en œuvre pour NVidia Kepler trop, mais malheureusement, le problème que j'ai (ridicule nombre de petites opérations matricielles) n'est pas idéal pour le matériel (ridiculement grandes opérations matricielles). Cela dit, ce papier présente une solution qui semble atteindre des centaines de GFLOPS/s sur les matrices de ma taille, de façon dynamique, en dérouler les boucles et la réalisation de l'ensemble de la factorisation dans les registres, avec la réserve que les dimensions de la matrice être définies au moment de la compilation. (Ce déroulement de la boucle devrait aider à réduire les frais généraux et d'améliorer la vectorisation du Phi version trop, getIndices() deviendra de plus en plus important!) Donc maintenant, je pense à mon noyau devrait ressembler à:

double *data;  // move data to GPU/Phi once into shared memory
template<unsigned int K> double doTestUnrolled(int *idx) {
  double submatrix[K][K];
  // copy data into submatrix
  #pragma unroll
  for(int i=0; i<K; i++) {
    #pragma unroll
    for(int j=0; j<K; j++) {
      submatrix[i][j] = data[idx[i]][idx[j]];
    }
  }
  factorizeUnrolled<K>(submatrix);
  return the_answer;
}

Le Phi version résout chaque modèle dans un "cilk_for" en boucle à partir du modèle=0 à 2^N (ou, plutôt, un sous-ensemble de test), mais maintenant, dans l'ordre des lots de travail pour le GPU et d'amortir le noyau de lancement de frais généraux que j'ai pour itérer les numéros de modèle lexicographique de commande pour chaque K=1 à 41 bits (comme doynax noté).

EDIT 2: Maintenant que les vacances, voici quelques résultats sur mon Xeon E5-2602 à l'aide de la cpi version 15. Le code que j'ai utilisé pour indice de référence est ici: http://pastebin.com/XvrGQUat. J'ai effectuer les bits d'extraction sur des entiers qui ont exactement K bits, donc il y a une surcharge pour le lexicographique itération mesurée dans la "Base" de la colonne dans le tableau ci-dessous. Celles-ci sont effectuées 2^30 fois avec N=48 (à répéter au besoin).

"CTZ" est une boucle qui utilise la gcc intrinsèque __builtin_ctzll d'obtenir le plus petit bit afin de définir:

for(int i=0; i<K; i++) {
    idx[i] = __builtin_ctzll(tmp);
    lb = tmp & -tmp;    // get lowest bit
    tmp ^= lb;      // remove lowest bit from tmp
} 

Marque est Marque est dépourvu de branches sur une boucle:

for(int i=0; i<K; i++) {
    *dst = i;
    dst += x & 1;
    x >>= 1;
} 

Tab1 est mon origine basée sur la table de code avec la copie de la macro:

#define COPY(d, s, n) \
switch(n) { \
case 8: *(d++) = *(s++); \
case 7: *(d++) = *(s++); \
case 6: *(d++) = *(s++); \
case 5: *(d++) = *(s++); \
case 4: *(d++) = *(s++); \
case 3: *(d++) = *(s++); \
case 2: *(d++) = *(s++); \
case 1: *(d++) = *(s++); \
case 0: break;        \
}

Tab2 est le même code que Tab1, mais la copie de la macro passe juste à 8 octets en un seul exemplaire (prendre des idées de doynax et Lưu Vĩnh Phúc... mais cette remarque n'est pas de s'assurer de l'alignement):

#define COPY2(d, s, n) { *((uint64_t *)d) = *((uint64_t *)s); d+=n; }

Voici les résultats. Je suppose que ma demande initiale, Tab1 est 3x plus rapide que CTZ ne tient que pour les grands K (où j'ai été le tester). Marque boucle est plus rapide que mon code d'origine, mais se débarrasser de la branche dans l' COPY2 macro prend le gâteau pour K > 8.

 K    Base    CTZ   Mark   Tab1   Tab2
001  4.97s  6.42s  6.66s 18.23s 12.77s
002  4.95s  8.49s  7.28s 19.50s 12.33s
004  4.95s  9.83s  8.68s 19.74s 11.92s
006  4.95s 16.86s  9.53s 20.48s 11.66s
008  4.95s 19.21s 13.87s 20.77s 11.92s
010  4.95s 21.53s 13.09s 21.02s 11.28s
015  4.95s 32.64s 17.75s 23.30s 10.98s
020  4.99s 42.00s 21.75s 27.15s 10.96s
030  5.00s 100.64s 35.48s 35.84s 11.07s
040  5.01s 131.96s 44.55s 44.51s 11.58s

Answer 1

Je crois que la clé de la performance ici est de se concentrer sur le problème plutôt que sur la micro-optimisation de l'extraction de positions de bits d'un entier aléatoire.

À en juger par votre exemple de code précédente et DONC, la question que vous êtes l'énumération de tous les mots avec K bits dans l'ordre, et d'en extraire les bits les indices de ces. Cela simplifie grandement les choses.

Si oui, alors au lieu de reconstruire la position du bit à chaque itération essayer directement l'incrémentation de la position dans le tableau de bits. La moitié du temps, il s'agira d'une seule itération de boucle et d'un incrément.

Quelque chose le long de ces lignes:

// Walk through all len-bit words with num-bits set in order
void enumerate(size_t num, size_t len) {
    size_t i;
    unsigned int bitpos[64 + 1];

    // Seed with the lowest word plus a sentinel
    for(i = 0; i < num; ++i)
        bitpos[i] = i;
    bitpos[i] = 0;

    // Here goes the main loop
    do {
        // Do something with the resulting data
        process(bitpos, num);

        // Increment the least-significant series of consecutive bits
        for(i = 0; bitpos[i + 1] == bitpos[i] + 1; ++i)
            bitpos[i] = i;
    // Stop on reaching the top
    } while(++bitpos[i] != len);
}

// Test function
void process(const unsigned int *bits, size_t num) {
    do
        printf("%d ", bits[--num]);
    while(num);
    putchar('\n');
}

Pas particulièrement optimisé, mais vous avez l'idée générale.

Answer 2

Voici quelque chose de très simple qui pourrait être plus rapide - aucun moyen de savoir sans test. Cela dépendra beaucoup du nombre de bits défini par rapport au nombre non défini. Vous pouvez dérouler cette opération pour supprimer complètement les branches, mais avec les processeurs actuels, je ne sais pas si cela accélérerait du tout.

 unsigned char idx[K+1]; // need one extra for overwrite protection
unsigned char *dst=idx;
for (unsigned char i = 0; i < 50; i++)
{
    *dst = i;
    dst += x & 1;
    x >>= 1;
}
 

PS votre échantillon de sortie dans la question est faux, voir http://ideone.com/2o032E

Answer 3

Comme un minimum de modifications:

int64_t x;            
char idx[K+1];
char *dst=idx;
const int BITS = 8;
for (int i = 0 ; i < 64+BITS; i += BITS) {
  int y = (x & ((1<<BITS)-1));
  char* end = strcat(dst, tab[y]); // tab[y] is a _string_
  for (; dst != end; ++dst)
  {
    *dst += (i - 1); // tab[] is null-terminated so bit positions are 1 to BITS.
  }
  x >>= BITS;
}

Le choix de l' BITS détermine la taille de la table. 8, 13 et 16 sont le choix logique. Chaque entrée est une chaîne de caractères, sans fin et contenant des positions de bits à 1 de décalage. I. e. onglet[5] est - "\x03\x01". La boucle interne corrige ce décalage.

Légèrement plus efficace: remplacer l' strcat et l'intérieur de la boucle par

char const* ptr = tab[y];
while (*ptr)
{
   *dst++ = *ptr++ + (i-1);
}

Déroulement de la boucle peut être un peu de douleur si la boucle contient des branches, parce que la copie de ceux de la branche des déclarations qui n'aide pas la branche prédicteur. Je vais heureux de laisser le choix à l'compilateur.

Une chose que j'envisage, c'est que tab[y] est un tableau de pointeurs vers des chaînes de caractères. Ceux-ci sont très similaires: "\x1" est un suffixe d' "\x3\x1". En fait, chaque chaîne qui ne commence pas par" "\x8" est un suffixe d'une chaîne qui ne. Je me demande combien de chaînes uniques dont vous avez besoin, et à quel degré, tab[y] est en fait nécessaire. E. g. par la logique ci-dessus, tab[128+x] == tab[x]-1.

[modifier]

Tant pis, vous avez certainement besoin 128 onglet entrées commençant par "\x8" depuis, ils ne sont jamais le suffixe d'une autre chaîne. Encore, l' tab[128+x] == tab[x]-1 règle signifie que vous pouvez économiser de la moitié des entrées, mais au prix de deux instructions supplémentaires: char const* ptr = tab[x & 0x7F] - ((x>>7) & 1). (Configurer tab[] de point après l' \x8)

Answer 4

À l'aide de char ne serait pas vous aider à augmenter la vitesse, mais en fait souvent besoin de plus de ANDing et signe/zero extension pendant le calcul. Seulement dans le cas de très grands tableaux qui doit tenir dans le cache, les petites int types doivent être utilisés

Une autre chose que vous pouvez améliorer est la COPIE de la macro. Au lieu de copier l'octet-par-octet, copiez le mot en entier si possible

inline COPY(unsigned char *dst, unsigned char *src, int n)
{
switch(n) { // remember to align dst and src when declaring
case 8:
    *((int64_t*)dst) = *((int64_t*)src);
    break;
case 7:
    *((int32_t*)dst) = *((int32_t*)src);
    *((int16_t*)(dst + 4)) = *((int32_t*)(src + 4));
    dst[6] = src[6];
    break;
case 6:
    *((int32_t*)dst) = *((int32_t*)src);
    *((int16_t*)(dst + 4)) = *((int32_t*)(src + 4));
    break;
case 5:
    *((int32_t*)dst) = *((int32_t*)src);
    dst[4] = src[4];
    break;
case 4:
    *((int32_t*)dst) = *((int32_t*)src);
    break;
case 3:
    *((int16_t*)dst) = *((int16_t*)src);
    dst[2] = src[2];
    break;
case 2:
    *((int16_t*)dst) = *((int16_t*)src);
    break;
case 1:
    dst[0] = src[0];
    break;
case 0:
    break;
}

Aussi, depuis tabofs[x] et n[x] est souvent accès à des proches les uns des autres, essayez de le mettre à proximité de la mémoire pour s'assurer qu'ils sont toujours dans le cache en même temps

typedef struct TAB_N
{
    int16_t n, tabofs;
} tab_n[256];

src=tab0+tab_n[b0].tabofs; COPY(dst, src, tab_n[b0].n);
src=tab0+tab_n[b1].tabofs; COPY(dst, src, tab_n[b1].n);
src=tab0+tab_n[b2].tabofs; COPY(dst, src, tab_n[b2].n);
src=tab0+tab_n[b3].tabofs; COPY(dst, src, tab_n[b3].n);
src=tab0+tab_n[b4].tabofs; COPY(dst, src, tab_n[b4].n);
src=tab0+tab_n[b5].tabofs; COPY(dst, src, tab_n[b5].n);

Dernier mais non le moins, gettimeofday n'est pas pour la performance de comptage. Utilisation QueryPerformanceCounter au lieu de cela, il est beaucoup plus précis

Answer 5

Votre code utilise une table d'index à 1 octet (256 entrées). Vous pouvez l'accélérer d'un facteur 2 si vous utilisez une table d'index à 2 octets (65 536 entrées).

Malheureusement, vous ne pouvez probablement pas étendre cela davantage - pour une taille de table de 3 octets, 16 Mo, ne rentrant probablement pas dans le cache local de la CPU, ce qui ralentirait les choses.