Computación rápida de log2 para enteros de 64 bits

Un gran recurso de progtwigción, Bit Twiddling Hacks, propone ( aquí ) el siguiente método para calcular log2 de un entero de 32 bits:

#define LT(n) n, n, n, n, n, n, n, n, n, n, n, n, n, n, n, n static const char LogTable256[256] = { -1, 0, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, LT(4), LT(5), LT(5), LT(6), LT(6), LT(6), LT(6), LT(7), LT(7), LT(7), LT(7), LT(7), LT(7), LT(7), LT(7) }; unsigned int v; // 32-bit word to find the log of unsigned r; // r will be lg(v) register unsigned int t, tt; // temporaries if (tt = v >> 16) { r = (t = tt >> 8) ? 24 + LogTable256[t] : 16 + LogTable256[tt]; } else { r = (t = v >> 8) ? 8 + LogTable256[t] : LogTable256[v]; } 

y menciona que

El método de la tabla de búsqueda solo requiere unas 7 operaciones para encontrar el registro de un valor de 32 bits. Si se extiende para cantidades de 64 bits, tomaría aproximadamente 9 operaciones.

pero, por desgracia, no proporciona ninguna información adicional sobre la forma en que uno realmente debería ir para extender el algoritmo a enteros de 64 bits.

¿Alguna pista sobre cómo sería un algoritmo de 64 bits de este tipo?

Las funciones intrínsecas son realmente rápidas, pero aún son insuficientes para una implementación verdaderamente independiente de Log2 independiente del comstackdor. Entonces, en caso de que alguien esté interesado, aquí está el algoritmo de DeBruijn más rápido, libre de ramificaciones y abstracto de CPU al que he llegado mientras investigaba el tema yo solo.

 const int tab64[64] = { 63, 0, 58, 1, 59, 47, 53, 2, 60, 39, 48, 27, 54, 33, 42, 3, 61, 51, 37, 40, 49, 18, 28, 20, 55, 30, 34, 11, 43, 14, 22, 4, 62, 57, 46, 52, 38, 26, 32, 41, 50, 36, 17, 19, 29, 10, 13, 21, 56, 45, 25, 31, 35, 16, 9, 12, 44, 24, 15, 8, 23, 7, 6, 5}; int log2_64 (uint64_t value) { value |= value >> 1; value |= value >> 2; value |= value >> 4; value |= value >> 8; value |= value >> 16; value |= value >> 32; return tab64[((uint64_t)((value - (value >> 1))*0x07EDD5E59A4E28C2)) >> 58]; } 

La parte de redondear a la siguiente potencia inferior de 2 se tomó de los límites de Poder-de-2 y la parte de obtener el número de ceros finales se tomó de BitScan (el código (bb & -bb) que hay para seleccionar el único el bit más a la derecha que se establece en 1, que no es necesario después de haber redondeado el valor a la siguiente potencia de 2).

Y la implementación de 32 bits, por cierto, es

 const int tab32[32] = { 0, 9, 1, 10, 13, 21, 2, 29, 11, 14, 16, 18, 22, 25, 3, 30, 8, 12, 20, 28, 15, 17, 24, 7, 19, 27, 23, 6, 26, 5, 4, 31}; int log2_32 (uint32_t value) { value |= value >> 1; value |= value >> 2; value |= value >> 4; value |= value >> 8; value |= value >> 16; return tab32[(uint32_t)(value*0x07C4ACDD) >> 27]; } 

Como con cualquier otro método computacional, log2 requiere que el valor de entrada sea mayor que cero.

Si está utilizando GCC, una tabla de búsqueda no es necesaria en este caso.

GCC proporciona una función integrada para determinar la cantidad de ceros a la izquierda:

Función incorporada: int __builtin_clz (unsigned int x)
Devuelve el número de 0 bits principales en x, comenzando en la posición de bit más significativa. Si x es 0, el resultado no está definido.

Entonces puedes definir:

 #define LOG2(X) ((unsigned) (8*sizeof (unsigned long long) - __builtin_clzll((X)) - 1)) 

y funcionará para cualquier int largo largo sin signo. El resultado se redondea hacia abajo.

Para x86 y AMD64, GCC lo comstackrá en una instrucción bsr , por lo que la solución es muy rápida (mucho más rápido que las tablas de búsqueda).

Ejemplo de trabajo :

 #include  #define LOG2(X) ((unsigned) (8*sizeof (unsigned long long) - __builtin_clzll((X)) - 1)) int main(void) { unsigned long long input; while (scanf("%llu", &input) == 1) { printf("log(%llu) = %u\n", input, LOG2(input)); } return 0; } 

Intentaba convertir Encuentre la base de registro 2 de un entero de N bits en operaciones O (lg (N)) con multiplicación y búsqueda hasta 64 bits por fuerza bruta forzando el número mágico. Huelga decir que estaba tomando un tiempo.

Luego encontré la respuesta de Desmond y decidí probar su número mágico como punto de partida. Como tengo un procesador de 6 núcleos, lo ejecuté en paralelo comenzando en 0x07EDD5E59A4E28C2 / 6 múltiplos. Me sorprendió que encontrara algo inmediatamente. Resulta que 0x07EDD5E59A4E28C2 / 2 funcionó.

Así que aquí está el código para 0x07EDD5E59A4E28C2 que le ahorra un cambio y resta:

 int LogBase2(uint64_t n) { static const int table[64] = { 0, 58, 1, 59, 47, 53, 2, 60, 39, 48, 27, 54, 33, 42, 3, 61, 51, 37, 40, 49, 18, 28, 20, 55, 30, 34, 11, 43, 14, 22, 4, 62, 57, 46, 52, 38, 26, 32, 41, 50, 36, 17, 19, 29, 10, 13, 21, 56, 45, 25, 31, 35, 16, 9, 12, 44, 24, 15, 8, 23, 7, 6, 5, 63 }; n |= n >> 1; n |= n >> 2; n |= n >> 4; n |= n >> 8; n |= n >> 16; n |= n >> 32; return table[(n * 0x03f6eaf2cd271461) >> 58]; } 

Logaritmo de enteros de base 2

Esto es lo que hago para enteros sin signo de 64 bits. Esto calcula el piso del logaritmo de base 2, que es equivalente al índice del bit más significativo. Este método es increíblemente rápido para números grandes porque usa un ciclo desenrollado que se ejecuta siempre en log₂64 = 6 pasos.

Esencialmente, lo que hace es restar cuadrados progresivamente más pequeños en la secuencia {0 ≤ k ≤ 5: 2 ^ (2 ^ k)} = {2³², 2¹⁶, 2⁸, 2⁴, 2², 2¹} = {4294967296, 65536, 256 , 16, 4, 2, 1} y sum los exponentes k de los valores restados.

 int uint64_log2(uint64_t n) { #define S(k) if (n >= (UINT64_C(1) < < k)) { i += k; n >>= k; } int i = -(n == 0); S(32); S(16); S(8); S(4); S(2); S(1); return i; #undef S } 

Tenga en cuenta que esto devuelve -1 si se le da una entrada inválida de 0 (que es lo que la inicial -(n == 0) está buscando). Si nunca espera invocarlo con n == 0 , puede sustituir int i = 0; para el inicializador y agregar assert(n != 0); en la entrada a la función.

Logaritmo de enteros de base 10

Los logaritmos de enteros de base 10 se pueden calcular utilizando de manera similar, con el cuadrado más grande para probar siendo 10¹⁶ porque log₁₀2⁶⁴ .2 19.2659 … (Nota: Esta no es la forma más rápida de lograr un logaritmo de enteros de base 10, porque usa división de enteros, que es intrínsecamente lento. Una implementación más rápida sería usar un acumulador con valores que crecen exponencialmente, y compararlos con el acumulador, realizando de hecho una especie de búsqueda binaria).

 int uint64_log10(uint64_t n) { #define S(k, m) if (n >= UINT64_C(m)) { i += k; n /= UINT64_C(m); } int i = -(n == 0); S(16,10000000000000000); S(8,100000000); S(4,10000); S(2,100); S(1,10); return i; #undef S } 

Aquí hay una extensión bastante compacta y rápida, sin temporales adicionales:

 r = 0; /* If its wider than 32 bits, then we already know that log >= 32. So store it in R. */ if (v >> 32) { r = 32; v >>= 32; } /* Now do the exact same thing as the 32 bit algorithm, except we ADD to R this time. */ if (tt = v >> 16) { r += (t = tt >> 8) ? 24 + LogTable256[t] : 16 + LogTable256[tt]; } else { r += (t = v >> 8) ? 8 + LogTable256[t] : LogTable256[v]; } 

Aquí hay uno construido con una cadena de if s, una vez más sin temporales adicionales. Puede que no sea el más rápido sin embargo.

  if (tt = v >> 48) { r = (t = tt >> 8) ? 56 + LogTable256[t] : 48 + LogTable256[tt]; } else if (tt = v >> 32) { r = (t = tt >> 8) ? 40 + LogTable256[t] : 32 + LogTable256[tt]; } else if (tt = v >> 16) { r = (t = tt >> 8) ? 24 + LogTable256[t] : 16 + LogTable256[tt]; } else { r = (t = v >> 8) ? 8 + LogTable256[t] : LogTable256[v]; } 

El algoritmo básicamente descubre qué byte contiene el bit 1 más significativo, y luego busca ese byte en la búsqueda para encontrar el registro del byte y luego lo agrega a la posición del byte.

Aquí hay una versión algo simplificada del algoritmo de 32 bits:

 if (tt = v >> 16) { if (t = tt >> 8) { r = 24 + LogTable256[t]; } else { r = 16 + LogTable256[tt]; } } else { if (t = v >> 8) { r = 8 + LogTable256[t]; } else { r = LogTable256[v]; } } 

Este es el algoritmo equivalente de 64 bits:

 if (ttt = v >> 32) { if (tt = ttt >> 16) { if (t = tt >> 8) { r = 56 + LogTable256[t]; } else { r = 48 + LogTable256[tt]; } } else { if (t = ttt >> 8) { r = 40 + LogTable256[t]; } else { r = 32 + LogTable256[ttt]; } } } else { if (tt = v >> 16) { if (t = tt >> 8) { r = 24 + LogTable256[t]; } else { r = 16 + LogTable256[tt]; } } else { if (t = v >> 8) { r = 8 + LogTable256[t]; } else { r = LogTable256[v]; } } } 

Se me ocurrió un algoritmo para cualquier tipo de tamaño que creo que es mejor que el original.

 unsigned int v = 42; unsigned int r = 0; unsigned int b; for (b = sizeof(v) < < 2; b; b = b >> 1) { if (v >> b) { v = v >> b; r += b; } } 

Nota: b = sizeof(v) < < 2 establece b en la mitad del número de bits en v. Usé el desplazamiento en lugar de la multiplicación aquí (solo porque me dio la gana).

Podría agregar una tabla de búsqueda a ese algoritmo para acelerarlo posiblemente, pero es más una prueba de concepto.