Estoy tratando de encontrar la forma más eficiente de calcular el módulo 255 de un entero sin signo de 32 bits. Mi enfoque principal es encontrar un algoritmo que funcione bien en las plataformas x86 y ARM con miras a la aplicabilidad más allá de eso. En primer lugar, estoy tratando de evitar las operaciones de memoria (que podrían ser costosas), por lo que estoy buscando enfoques poco complicados mientras evito las tablas. También estoy tratando de evitar operaciones potencialmente costosas, como bifurcaciones y multiplicaciones, y minimizar la cantidad de operaciones y registros utilizados.
El siguiente código ISO-C99 captura las ocho variantes que probé hasta ahora. Incluye un marco para pruebas exhaustivas. Me aferré a esta medida de tiempo de ejecución cruda que parece funcionar lo suficientemente bien como para obtener una primera impresión de rendimiento. En las pocas plataformas que probé (todas con multiplicaciones rápidas de enteros), las variantes WARREN_MUL_SHR_2 , WARREN_MUL_SHR_1 y DIGIT_SUM_CARRY_OUT_1 parecen ser las más eficaces. Mis experimentos muestran que los compiladores x86, ARM, PowerPC y MIPS que probé en Compiler Explorer hacen muy buen uso de las funciones específicas de la plataforma, como LEA de tres entradas, instrucciones de expansión de bytes, multiplicación y acumulación y predicación de instrucciones.
La variante NAIVE_USING_DIV utiliza una división de enteros, multiplica hacia atrás con el divisor seguido de la resta. Este es el caso base. Los compiladores modernos saben cómo implementar eficientemente la división de enteros sin signo por 255 (a través de la multiplicación) y usarán un reemplazo discreto para backmultiply cuando corresponda. Para calcular el módulo base-1 uno puede sumar los dígitos de la base y luego doblar el resultado. Por ejemplo, 3334 mod 9: suma 3+3+3+4 = 13, doblar 1+3 = 4. Si el resultado después de doblar es base-1 , necesitamos generar 0 en su lugar. DIGIT_SUM_THEN_FOLD usa este método.
A. Cockburn, "Implementación eficiente del algoritmo de suma de comprobación del protocolo de transporte OSI utilizando aritmética de 8/16 bits", ACM SIGCOMM Computer Communication Review , vol. 17, No. 3, julio/agosto. 1987, págs. 13-20
mostró una forma diferente de sumar dígitos módulo base-1 manera eficiente en el contexto de un cálculo de suma de comprobación módulo 255. Calcule una suma de los dígitos por bytes y, después de cada suma, agregue también cualquier resultado de la suma. Así que esto sería una secuencia ADD a, b , ADC a, 0 . Escribiendo la cadena de suma para esto usando base 256 dígitos queda claro que el cálculo es básicamente una multiplicación con 0x0101 ... 0101 . El resultado estará en la posición del dígito más significativo, excepto que se necesita capturar el acarreo de la suma en esa posición por separado. Este método solo funciona cuando un dígito base comprende 2k bits. Aquí tenemos k=3 . Probé tres formas diferentes de reasignar un resultado de base-1 a 0, lo que resultó en variantes DIGIT_SUM_CARRY_OUT_1 , DIGIT_SUM_CARRY_OUT_2 , DIGIT_SUM_CARRY_OUT_3 .
Joe Keane demostró un enfoque intrigante para calcular módulo-63 de manera eficiente en el grupo de noticias comp.lang.c el 09/07/1995. Si bien el participante del hilo Peter L. Montgomery demostró que el algoritmo era correcto, desafortunadamente el Sr. Keane no respondió a las solicitudes para explicar su derivación. Este algoritmo también se reproduce en Hacker's Delight 2nd ed de H. Warren. Pude extenderlo, de forma puramente mecánica , a módulo-127 y módulo-255. Esta es la variante KEANE_MAGIC (apropiadamente nombrada). Actualización: desde que publiqué originalmente esta pregunta, descubrí que el enfoque de Keane es básicamente una implementación inteligente de punto fijo de lo siguiente: return (uint32_t)(fmod (x * 256.0 / 255.0 + 0.5, 256.0) * (255.0 / 256.0)); . Esto lo convierte en un pariente cercano de la siguiente variante.
Henry S. Warren, Hacker's Delight 2ª ed. , pags. 272 muestra un algoritmo de "multiplicar desplazamiento a la derecha", presumiblemente ideado por el propio autor, que se basa en la propiedad matemática de que n mod 2 k-1 = piso (2 k / 2 k-1 * n) mod 2 k . El cálculo de punto fijo se usa para multiplicar con el factor 2 k / 2 k-1 . Construí dos variantes de esto que difieren en cómo manejan el mapeo de un resultado preliminar de base-1 a 0. Estas son variantes WARREN_MUL_SHR_1 y WARREN_MUL_SHR_2 .
¿Existen algoritmos para el cálculo del módulo 255 que sean aún más eficientes que los tres principales contendientes que he identificado hasta ahora, en particular para plataformas con multiplicaciones de enteros lentas? Una modificación eficiente del algoritmo sin multiplicación de Keane para la suma de cuatro dígitos de base 256 parece ser de particular interés en este contexto.
#include <stdio.h> #include <stdlib.h> #include <stdint.h> #define NAIVE_USING_DIV (1) #define DIGIT_SUM_THEN_FOLD (2) #define DIGIT_SUM_CARRY_OUT_1 (3) #define DIGIT_SUM_CARRY_OUT_2 (4) #define DIGIT_SUM_CARRY_OUT_3 (5) #define KEANE_MAGIC (6) // Joe Keane, comp.lang.c, 1995/07/09 #define WARREN_MUL_SHR_1 (7) // Hacker's Delight, 2nd ed., p. 272 #define WARREN_MUL_SHR_2 (8) // Hacker's Delight, 2nd ed., p. 272 #define VARIANT (WARREN_MUL_SHR_2) uint32_t mod255 (uint32_t x) { #if VARIANT == NAIVE_USING_DIV return x - 255 * (x / 255); #elif VARIANT == DIGIT_SUM_THEN_FOLD x = (x & 0xffff) + (x >> 16); x = (x & 0xff) + (x >> 8); x = (x & 0xff) + (x >> 8) + 1; x = (x & 0xff) + (x >> 8) - 1; return x; #elif VARIANT == DIGIT_SUM_CARRY_OUT_1 uint32_t t; t = 0x01010101 * x; t = (t >> 24) + (t < x); if (t == 255) t = 0; return t; #elif VARIANT == DIGIT_SUM_CARRY_OUT_2 uint32_t t; t = 0x01010101 * x; t = (t >> 24) + (t < x) + 1; t = (t & 0xff) + (t >> 8) - 1; return t; #elif VARIANT == DIGIT_SUM_CARRY_OUT_3 uint32_t t; t = 0x01010101 * x; t = (t >> 24) + (t < x); t = t & ((t - 255) >> 8); return t; #elif VARIANT == KEANE_MAGIC x = (((x >> 16) + x) >> 14) + (x << 2); x = ((x >> 8) + x + 2) & 0x3ff; x = (x - (x >> 8)) >> 2; return x; #elif VARIANT == WARREN_MUL_SHR_1 x = (0x01010101 * x + (x >> 8)) >> 24; x = x & ((x - 255) >> 8); return x; #elif VARIANT == WARREN_MUL_SHR_2 x = (0x01010101 * x + (x >> 8)) >> 24; if (x == 255) x = 0; return x; #else #error unknown VARIANT #endif } uint32_t ref_mod255 (uint32_t x) { volatile uint32_t t = x; t = t % 255; return t; } // timing with microsecond resolution #if defined(_WIN32) #if !defined(WIN32_LEAN_AND_MEAN) #define WIN32_LEAN_AND_MEAN #endif #include <windows.h> double second (void) { LARGE_INTEGER t; static double oofreq; static int checkedForHighResTimer; static BOOL hasHighResTimer; if (!checkedForHighResTimer) { hasHighResTimer = QueryPerformanceFrequency (&t); oofreq = 1.0 / (double)t.QuadPart; checkedForHighResTimer = 1; } if (hasHighResTimer) { QueryPerformanceCounter (&t); return (double)t.QuadPart * oofreq; } else { return (double)GetTickCount() * 1.0e-3; } } #elif defined(__linux__) || defined(__APPLE__) #include <stddef.h> #include <sys/time.h> double second (void) { struct timeval tv; gettimeofday(&tv, NULL); return (double)tv.tv_sec + (double)tv.tv_usec * 1.0e-6; } #else #error unsupported platform #endif int main (void) { double start, stop; uint32_t res, ref, x = 0; printf ("Testing VARIANT = %d\n", VARIANT); start = second(); do { res = mod255 (x); ref = ref_mod255 (x); if (res != ref) { printf ("error @ %08x: res=%08x ref=%08x\n", x, res, ref); return EXIT_FAILURE; } x++; } while (x); stop = second(); printf ("test passed\n"); printf ("elapsed = %.6f seconds\n", stop - start); return EXIT_SUCCESS; }Esta es mi idea de cómo funcionan las respuestas más rápidas. Todavía no sé si Keane se puede mejorar o generalizar fácilmente.
Dado un entero x ≥ 0, sea q = ⌊x/255⌋ (en C, q = x / 255; ) y r = x − 255 q (en C, r = x % 255; ) de modo que q ≥ 0 y 0 ≤ r < 255 son números enteros y x = 255 q + r.
Este método evalúa (x + ⌊x/255⌋) mod 2 8 (en C, (x + x / 255) & 0xff ), lo que equivale a (255 q + r + q) mod 2 8 = (2 8 q + r ) mod 2 8 = r.
Tenga en cuenta que x + ⌊x/255⌋ = ⌊x + x/255⌋ = ⌊(2 8 /255) x⌋, donde el primer paso se deriva de que x es un número entero. Este método usa el multiplicador (2 0 + 2 −8 + 2 −16 + 2 −24 + 2 −32 ) en lugar de 2 8 /255, que es la suma de la serie infinita 2 0 + 2 −8 + 2 −16 + 2 −24 + 2 −32 + …. Dado que la aproximación es ligeramente inferior, este método debe detectar el residuo 2 8 − 1 = 255.
La intuición de este método es calcular y = (2 8 /255) x mod 2 8 , que es igual a (2 8 /255) (255 q + r) mod 2 8 = (2 8 q + (2 8 /255) r) mod 2 8 = (2 8 /255) r, y devuelve y − y/2 8 , que es igual a r.
Dado que estas fórmulas no usan el hecho de que ⌊(2 8/255 ) r⌋ = r, Keane puede cambiar de 2 8 a 2 10 para dos bits de protección. Idealmente, estos siempre serían cero, pero debido al truncamiento de punto fijo y una aproximación de 2 10/255 , no lo son. Keane suma 2 para pasar del truncamiento al redondeo, lo que también evita el caso especial de Warren.
Este método usa el multiplicador 2 2 (2 0 + 2 −8 + 2 −16 + 2 −24 + 2 −32 + 2 −40 ) = 2 2 (2 0 + 2 −16 + 2 −32 ) (2 0 + 2 −8 ). El enunciado C x = (((x >> 16) + x) >> 14) + (x << 2); calcula x′ = ⌊2 2 (2 0 + 2 −16 + 2 −32 ) x⌋ modificación 2 32 . Entonces ((x >> 8) + x) & 0x3ff es x′′ = ⌊(2 0 + 2 −8 ) x′⌋ mod 2 10 .
No tengo tiempo en este momento para hacer el análisis de errores formalmente. De manera informal, el intervalo de error del primer cálculo tiene un ancho < 1; el segundo, ancho < 2 + 2 −8 ; el tercero, ancho < ((2 − 2 −8 ) + 1)/2 2 < 1, lo que permite un redondeo correcto.
Con respecto a las mejoras, el término 2 −40 de la aproximación no parece necesario (?), pero bien podríamos tenerlo a menos que podamos eliminar el término 2 −32 . Dejar caer 2 −32 empuja la calidad de la aproximación fuera de especificación.
Supongo que probablemente no esté buscando soluciones que requieran una multiplicación rápida de 64 bits, pero para que conste:
return (x * 0x101010101010102ULL) >> 56;Para enteros arbitrarios sin signo, x y n , evaluar la expresión módulo x % n implica (conceptualmente, al menos), tres operaciones: división, multiplicación y resta:
quotient = x / n; product = quotient * n; modulus = x - product;Sin embargo, cuando n es una potencia de 2 ( n = 2 p ), el módulo se puede determinar mucho más rápidamente, simplemente ocultando todos los bits excepto los p inferiores.
En la mayoría de las CPU, la suma, la resta y el enmascaramiento de bits son operaciones muy 'baratas' (rápidas), la multiplicación es más 'cara' y la división es muy costosa, pero tenga en cuenta que la mayoría de los compiladores de optimización convertirán la división por una constante de tiempo de compilación en una multiplicación (por una constante diferente) y un cambio de bit ( vide infra ).
Por lo tanto, si podemos convertir nuestro módulo 255 en un módulo 256, sin demasiada sobrecarga, probablemente podamos acelerar el proceso. Podemos hacer esto al notar que x % n es equivalente a (x + x / n) % (n + 1) † . Así, nuestras operaciones conceptuales ahora son: división, adición y enmascaramiento.
En el caso específico de enmascarar los 8 bits inferiores, las CPU basadas en x86/x64 (¿y otras?) probablemente podrán realizar una mayor optimización, ya que pueden acceder a versiones de 8 bits de (la mayoría) de los registros.
Esto es lo que genera el compilador clang-cl para una función ingenua de módulo 255 (argumento pasado en ecx y devuelto en eax ):
unsigned Naive255(unsigned x) { return x % 255; } mov edx, ecx mov eax, 2155905153 ; imul rax, rdx ; Replacing the IDIV with IMUL and SHR shr rax, 39 ; mov edx, eax shl edx, 8 sub eax, edx add eax, ecxY aquí está el código (claramente más rápido) generado usando el 'truco' descrito anteriormente:
unsigned Trick255(unsigned x) { return (x + x / 255) & 0xFF; } mov eax, ecx mov edx, 2155905153 imul rdx, rax shr rdx, 39 add edx, ecx movzx eax, dl ; Faster than an explicit AND mask?La prueba de este código en una plataforma Windows-10 (64 bits) (CPU Intel® Core™ i7-8550U) muestra que supera significativamente (pero no enormemente) a los otros algoritmos presentados en la pregunta.
† La respuesta dada por David Eisenstat explica cómo/por qué esta equivalencia es válida.
Este método (mejorado ligeramente desde la edición anterior) mezcla a Warren y Keane. En mi computadora portátil, es más rápido que Keane pero no tan rápido como un multiplicador y cambio de 64 bits. Evita la multiplicación pero se beneficia de una sola instrucción de rotación. A diferencia de la versión original, probablemente esté bien en RISC-V.
Al igual que Warren, este método aproxima ⌊(256/255) x mod 256⌋ en 8,24 punto fijo. Mod 256, cada byte b contribuye con un término (256/255) b, que es aproximadamente b.bbb base 256. La versión original de este método simplemente suma las cuatro rotaciones de bytes. (Iré a la versión revisada en un momento.) Esta suma siempre subestima el valor real, pero por menos de 4 unidades en el último lugar. Al sumar 4/2 −24 antes de truncar, garantizamos la respuesta correcta como en Keane.
La versión revisada ahorra trabajo al relajar la calidad de la aproximación. Escribimos (256/255) x = (257/256) (65536/65535) x, evaluamos (65536/65535) x en punto fijo 16.16 (es decir, sumamos x a su rotación de 16 bits), y luego multiplicamos por 257 /256 y mod por 256 en punto fijo 8.24. La primera multiplicación tiene un error de menos de 2 unidades en el último lugar de 16,16, y la segunda es exacta (!). La suma subestima en menos de (2/2 16 ) (257/256), por lo que un término constante de 514/2 24 es suficiente para corregir el truncamiento. También es posible usar un valor mayor en caso de que un operando inmediato diferente sea más eficiente.
uint32_t mod255(uint32_t x) { x += (x << 16) | (x >> 16); return ((x << 8) + x + 514) >> 24; }Si tuviéramos un método incorporado, intrínseco o optimizado para una sola instrucción addc, se podría usar la aritmética de 32 bits de la siguiente manera:
uint32_t carry = 0; // sum up top and bottom 16 bits while generating carry out x = __builtin_addc(x, x<<16, carry, &carry); x &= 0xffff0000; // store the previous carry to bit 0 while adding // bits 16:23 over bits 24:31, and producing one more carry x = __builtin_addc(x, x << 8, carry, &carry); x = __builtin_addc(x, x >> 24, carry, &carry); x &= 0x0000ffff; // actually 0x1ff is enough // final correction for 0<=x<=257, ie min(x,x-255) x = x < x-255 ? x : x - 255; En Arm64, al menos la instrucción regular de suma puede tomar la forma de add r0, r1, r2 LSL 16 ; el enmascaramiento con bits consecutivos inmediatos o borradores es una sola instrucción bfi r0, wzr, #start_bit, #length .
Para el cálculo paralelo, no se puede usar esa multiplicación de ampliación eficiente. En su lugar, se puede dividir y vencer mientras se calculan los acarreos: comenzando con 16 elementos uint32_t interpretados como 16+16 elementos uint16_t, luego pasando a la aritmética uint8_t, uno puede calcular un resultado en un poco menos de una instrucción.
a0 = vld2q_u16(ptr); // split input to top16+bot16 bits a1 = vld2q_u16(ptr + 8); // load more inputs auto b0 = vaddq_u16(a0.val[0], a0.val[1]); auto b1 = vaddq_u16(a1.val[0], a1.val[1]); auto c0 = vcltq_u16(b0, a0.val[1]); // 8 carries auto c1 = vcltq_u16(b1, a1.val[1]); // 8 more carries b0 = vsubq_u16(b0, c0); b1 = vsubq_u16(b1, c1); auto d = vuzpq_u8(b0, b1); auto result = vaddq_u8(d.val[0], d.val[1]); auto carry = vcltq_u8(result, d.val[1]); result = vsubq_u8(result, carry); auto is_255 = vceqq_u8(result, vdupq_n_u8(255)); result = vbicq_u8(result, is_255);