Capturar y calcular el desbordamiento durante la multiplicación de dos enteros grandes

Estoy buscando una solución eficiente (opcionalmente estándar, elegante y fácil de implementar) para multiplicar números relativamente grandes y almacenar el resultado en uno o varios enteros:

Digamos que tengo dos enteros de 64 bits declarados así:

uint64_t a = xxx, b = yyy; 

Cuando hago a * b , ¿cómo puedo detectar si la operación da como resultado un desbordamiento y en este caso almacenar el transporte en alguna parte?

Tenga en cuenta que no quiero usar ninguna biblioteca de números grandes ya que tengo restricciones en la forma en que almaceno los números.

1. Detectando el desbordamiento :

 x = a * b; if (a != 0 && x / a != b) { // overflow handling } 

Editar: división fija por 0 (gracias Mark!)

2. Computar el acarreo es bastante complicado. Un enfoque es dividir ambos operandos en medias palabras, luego aplicar una multiplicación larga a las medias palabras:

 uint64_t hi(uint64_t x) { return x >> 32; } uint64_t lo(uint64_t x) { return ((1L << 32) - 1) & x; } void multiply(uint64_t a, uint64_t b) { // actually uint32_t would do, but the casting is annoying uint64_t s0, s1, s2, s3; uint64_t x = lo(a) * lo(b); s0 = lo(x); x = hi(a) * lo(b) + hi(x); s1 = lo(x); s2 = hi(x); x = s1 + lo(a) * hi(b); s1 = lo(x); x = s2 + hi(a) * hi(b) + hi(x); s2 = lo(x); s3 = hi(x); uint64_t result = s1 << 32 | s0; uint64_t carry = s3 << 32 | s2; } 

Para ver que ninguna de las sums parciales en sí pueden desbordarse, consideramos el peor de los casos:

  x = s2 + hi(a) * hi(b) + hi(x) 

Deje B = 1 << 32 . Entonces tenemos

  x <= (B - 1) + (B - 1)(B - 1) + (B - 1) <= B*B - 1 < B*B 

Creo que esto funcionará, al menos maneja el caso de prueba de Sjlver. Aparte de eso, no está probado (y puede que ni siquiera compile, ya que no tengo un comstackdor de C ++ a mano).

La idea es usar el siguiente hecho, que es cierto para la operación integral:

a*b > c si y solo si a > c/b

/ es la división integral aquí.

El pseudocódigo para verificar contra el desbordamiento de números positivos sigue:

if (a> max_int64 / b) luego “overflow” else else “ok” .

Para manejar ceros y números negativos, debe agregar más controles.

El código C para a y b no negativo sigue:

 if (b > 0 && a > 18446744073709551615 / b) { // overflow handling }; else { c = a * b; } 

Nota:

 18446744073709551615 == (1<<64)-1 

Para calcular el acarreo, podemos usar el enfoque para dividir el número en dos dígitos de 32 y multiplicarlos a medida que hacemos esto en el papel. Necesitamos dividir los números para evitar el desbordamiento.

El código sigue:

 // split input numbers into 32-bit digits uint64_t a0 = a & ((1LL<<32)-1); uint64_t a1 = a >> 32; uint64_t b0 = b & ((1LL<<32)-1); uint64_t b1 = b >> 32; // The following 3 lines of code is to calculate the carry of d1 // (d1 - 32-bit second digit of result, and it can be calculated as d1=d11+d12), // but to avoid overflow. // Actually rewriting the following 2 lines: // uint64_t d1 = (a0 * b0 >> 32) + a1 * b0 + a0 * b1; // uint64_t c1 = d1 >> 32; uint64_t d11 = a1 * b0 + (a0 * b0 >> 32); uint64_t d12 = a0 * b1; uint64_t c1 = (d11 > 18446744073709551615 - d12) ? 1 : 0; uint64_t d2 = a1 * b1 + c1; uint64_t carry = d2; // needed carry stored here 

Aunque ha habido varias otras respuestas a esta pregunta, varias de ellas tienen un código que no se ha probado completamente, y hasta el momento nadie ha comparado adecuadamente las diferentes opciones posibles.

Por esa razón, escribí y probé varias implementaciones posibles (la última se basa en este código de OpenBSD, discutido en Reddit aquí ). Aquí está el código:

 /* Multiply with overflow checking, emulating clang's builtin function * * __builtin_umull_overflow * * This code benchmarks five possible schemes for doing so. */ #include  #include  #include  #include  #include  #ifndef BOOL #define BOOL int #endif // Option 1, check for overflow a wider type // - Often fastest and the least code, especially on modern compilers // - When long is a 64-bit int, requires compiler support for 128-bits // ints (requires GCC >= 3.0 or Clang) #if LONG_BIT > 32 typedef __uint128_t long_overflow_t ; #else typedef uint64_t long_overflow_t; #endif BOOL umull_overflow1(unsigned long lhs, unsigned long rhs, unsigned long* result) { long_overflow_t prod = (long_overflow_t)lhs * (long_overflow_t)rhs; *result = (unsigned long) prod; return (prod >> LONG_BIT) != 0; } // Option 2, perform long multiplication using a smaller type // - Sometimes the fastest (eg, when mulitply on longs is a library // call). // - Performs at most three multiplies, and sometimes only performs one. // - Highly portable code; works no matter how many bits unsigned long is BOOL umull_overflow2(unsigned long lhs, unsigned long rhs, unsigned long* result) { const unsigned long HALFSIZE_MAX = (1ul << LONG_BIT/2) - 1ul; unsigned long lhs_high = lhs >> LONG_BIT/2; unsigned long lhs_low = lhs & HALFSIZE_MAX; unsigned long rhs_high = rhs >> LONG_BIT/2; unsigned long rhs_low = rhs & HALFSIZE_MAX; unsigned long bot_bits = lhs_low * rhs_low; if (!(lhs_high || rhs_high)) { *result = bot_bits; return 0; } BOOL overflowed = lhs_high && rhs_high; unsigned long mid_bits1 = lhs_low * rhs_high; unsigned long mid_bits2 = lhs_high * rhs_low; *result = bot_bits + ((mid_bits1+mid_bits2) << LONG_BIT/2); return overflowed || *result < bot_bits || (mid_bits1 >> LONG_BIT/2) != 0 || (mid_bits2 >> LONG_BIT/2) != 0; } // Option 3, perform long multiplication using a smaller type (this code is // very similar to option 2, but calculates overflow using a different but // equivalent method). // - Sometimes the fastest (eg, when mulitply on longs is a library // call; clang likes this code). // - Performs at most three multiplies, and sometimes only performs one. // - Highly portable code; works no matter how many bits unsigned long is BOOL umull_overflow3(unsigned long lhs, unsigned long rhs, unsigned long* result) { const unsigned long HALFSIZE_MAX = (1ul << LONG_BIT/2) - 1ul; unsigned long lhs_high = lhs >> LONG_BIT/2; unsigned long lhs_low = lhs & HALFSIZE_MAX; unsigned long rhs_high = rhs >> LONG_BIT/2; unsigned long rhs_low = rhs & HALFSIZE_MAX; unsigned long lowbits = lhs_low * rhs_low; if (!(lhs_high || rhs_high)) { *result = lowbits; return 0; } BOOL overflowed = lhs_high && rhs_high; unsigned long midbits1 = lhs_low * rhs_high; unsigned long midbits2 = lhs_high * rhs_low; unsigned long midbits = midbits1 + midbits2; overflowed = overflowed || midbits < midbits1 || midbits > HALFSIZE_MAX; unsigned long product = lowbits + (midbits << LONG_BIT/2); overflowed = overflowed || product < lowbits; *result = product; return overflowed; } // Option 4, checks for overflow using division // - Checks for overflow using division // - Division is slow, especially if it is a library call BOOL umull_overflow4(unsigned long lhs, unsigned long rhs, unsigned long* result) { *result = lhs * rhs; return rhs > 0 && (SIZE_MAX / rhs) < lhs; } // Option 5, checks for overflow using division // - Checks for overflow using division // - Avoids division when the numbers are "small enough" to trivially // rule out overflow // - Division is slow, especially if it is a library call BOOL umull_overflow5(unsigned long lhs, unsigned long rhs, unsigned long* result) { const unsigned long MUL_NO_OVERFLOW = (1ul << LONG_BIT/2) - 1ul; *result = lhs * rhs; return (lhs >= MUL_NO_OVERFLOW || rhs >= MUL_NO_OVERFLOW) && rhs > 0 && SIZE_MAX / rhs < lhs; } #ifndef umull_overflow #define umull_overflow2 #endif /* * This benchmark code performs a multiply at all bit sizes, * essentially assuming that sizes are logarithmically distributed. */ int main() { unsigned long i, j, k; int count = 0; unsigned long mult; unsigned long total = 0; for (k = 0; k < 0x40000000 / LONG_BIT / LONG_BIT; ++k) for (i = 0; i != LONG_MAX; i = i*2+1) for (j = 0; j != LONG_MAX; j = j*2+1) { count += umull_overflow(i+k, j+k, &mult); total += mult; } printf("%d overflows (total %lu)\n", count, total); } 

Aquí están los resultados, las pruebas con varios comstackdores y sistemas que tengo (en este caso, todas las pruebas se realizaron en OS X, pero los resultados deberían ser similares en los sistemas BSD o Linux):

 +------------------+----------+----------+----------+----------+----------+ | | Option 1 | Option 2 | Option 3 | Option 4 | Option 5 | | | BigInt | LngMult1 | LngMult2 | Div | OptDiv | +------------------+----------+----------+----------+----------+----------+ | Clang 3.5 i386 | 1.610 | 3.217 | 3.129 | 4.405 | 4.398 | | GCC 4.9.0 i386 | 1.488 | 3.469 | 5.853 | 4.704 | 4.712 | | GCC 4.2.1 i386 | 2.842 | 4.022 | 3.629 | 4.160 | 4.696 | | GCC 4.2.1 PPC32 | 8.227 | 7.756 | 7.242 | 20.632 | 20.481 | | GCC 3.3 PPC32 | 5.684 | 9.804 | 11.525 | 21.734 | 22.517 | +------------------+----------+----------+----------+----------+----------+ | Clang 3.5 x86_64 | 1.584 | 2.472 | 2.449 | 9.246 | 7.280 | | GCC 4.9 x86_64 | 1.414 | 2.623 | 4.327 | 9.047 | 7.538 | | GCC 4.2.1 x86_64 | 2.143 | 2.618 | 2.750 | 9.510 | 7.389 | | GCC 4.2.1 PPC64 | 13.178 | 8.994 | 8.567 | 37.504 | 29.851 | +------------------+----------+----------+----------+----------+----------+ 

En base a estos resultados, podemos sacar algunas conclusiones:

  • Claramente, el enfoque basado en divisiones, aunque simple y portátil, es lento.
  • Ninguna técnica es un claro ganador en todos los casos.
  • En los comstackdores modernos, el enfoque use-a-larger-int es el mejor, si puede usarlo
  • En comstackdores más antiguos, el enfoque de multiplicación larga es el mejor
  • Sorprendentemente, GCC 4.9.0 tiene regresiones de rendimiento sobre GCC 4.2.1, y GCC 4.2.1 tiene regresiones de rendimiento sobre GCC 3.3

Una versión que también funciona cuando a == 0:

  x = a * b; if (a != 0 && x / a != b) { // overflow handling } 

Si no solo necesita detectar el desbordamiento, sino también capturar el transporte, es mejor que desglose sus números en partes de 32 bits. El código es una pesadilla; lo que sigue es solo un boceto:

 #include  uint64_t mul(uint64_t a, uint64_t b) { uint32_t ah = a >> 32; uint32_t al = a; // truncates: now a = al + 2**32 * ah uint32_t bh = b >> 32; uint32_t bl = b; // truncates: now b = bl + 2**32 * bh // a * b = 2**64 * ah * bh + 2**32 * (ah * bl + bh * al) + al * bl uint64_t partial = (uint64_t) al * (uint64_t) bl; uint64_t mid1 = (uint64_t) ah * (uint64_t) bl; uint64_t mid2 = (uint64_t) al * (uint64_t) bh; uint64_t carry = (uint64_t) ah * (uint64_t) bh; // add high parts of mid1 and mid2 to carry // add low parts of mid1 and mid2 to partial, carrying // any carry bits into carry... } 

El problema no es solo los productos parciales, sino el hecho de que cualquiera de las sums puede desbordarse.

Si tuviera que hacer esto de verdad, escribiría una rutina de multiplicación extendida en el lenguaje ensamblador local. Es decir, por ejemplo, multiplicar dos enteros de 64 bits para obtener un resultado de 128 bits, que se almacena en dos registros de 64 bits. Todo el hardware razonable proporciona esta funcionalidad en una sola instrucción de multiplicación nativa, no solo es accesible desde C.

Este es uno de esos raros casos en los que la solución más elegante y fácil de progtwigr es utilizar lenguaje ensamblador. Pero ciertamente no es portátil 🙁

He estado trabajando con este problema en estos días y debo decir que me ha impresionado la cantidad de veces que he visto a gente decir que la mejor manera de saber si ha habido un desbordamiento es dividir el resultado, eso es totalmente ineficiente y innecesario. El punto de esta función es que debe ser lo más rápido posible.

Hay dos opciones para la detección de desbordamiento:

1º- Si es posible, crea la variable de resultado dos veces más grande que los multiplicadores, por ejemplo:

 struct INT32struct {INT16 high, low;}; typedef union { struct INT32struct s; INT32 ll; } INT32union; INT16 mulFunction(INT16 a, INT16 b) { INT32union result.ll = a * b; //32Bits result if(result.s.high > 0) Overflow(); return (result.s.low) } 

Sabrá inmediatamente si ha habido un desbordamiento, y el código es el más rápido posible sin escribirlo en código de máquina. Dependiendo del comstackdor, este código se puede mejorar en código de máquina.

2º- Es imposible crear una variable de resultado dos veces más grande que la variable de multiplicadores: Entonces debería jugar con las condiciones si para determinar la mejor ruta. Continuando con el ejemplo:

 INT32 mulFunction(INT32 a, INT32 b) { INT32union s_a.ll = abs(a); INT32union s_b.ll = abs(b); //32Bits result INT32union result; if(s_a.s.hi > 0 && s_b.s.hi > 0) { Overflow(); } else if (s_a.s.hi > 0) { INT32union res1.ll = s_a.s.hi * s_b.s.lo; INT32union res2.ll = s_a.s.lo * s_b.s.lo; if (res1.hi == 0) { result.s.lo = res1.s.lo + res2.s.hi; if (result.s.hi == 0) { result.s.ll = result.s.lo << 16 + res2.s.lo; if ((ashi >> 15) ^ (bshi >> 15) == 1) { result.s.ll = -result.s.ll; } return result.s.ll }else { Overflow(); } }else { Overflow(); } }else if (s_b.s.hi > 0) { //Same code changing a with b }else { return (s_a.lo * s_b.lo); } } 

Espero que este código te ayude a tener un progtwig bastante eficiente y espero que el código sea claro, si no, pondré algunos comentarios.

atentamente.

Quizás la mejor forma de resolver este problema es tener una función que multiplique dos UInt64 y obtenga un par de UInt64, una parte superior y una parte inferior del resultado UInt128. Aquí está la solución, que incluye una función, que muestra el resultado en hexadecimal. Supongo que quizás prefiera una solución C ++, pero tengo una solución Swift que muestra cómo administrar el problema:

 func hex128 (_ hi: UInt64, _ lo: UInt64) -> String { var s: String = String(format: "%08X", hi >> 32) + String(format: "%08X", hi & 0xFFFFFFFF) + String(format: "%08X", lo >> 32) + String(format: "%08X", lo & 0xFFFFFFFF) return (s) } func mul64to128 (_ multiplier: UInt64, _ multiplicand : UInt64) -> (result_hi: UInt64, result_lo: UInt64) { let x: UInt64 = multiplier let x_lo: UInt64 = (x & 0xffffffff) let x_hi: UInt64 = x >> 32 let y: UInt64 = multiplicand let y_lo: UInt64 = (y & 0xffffffff) let y_hi: UInt64 = y >> 32 let mul_lo: UInt64 = (x_lo * y_lo) let mul_hi: UInt64 = (x_hi * y_lo) + (mul_lo >> 32) let mul_carry: UInt64 = (x_lo * y_hi) + (mul_hi & 0xffffffff) let result_hi: UInt64 = (x_hi * y_hi) + (mul_hi >> 32) + (mul_carry >> 32) let result_lo: UInt64 = (mul_carry << 32) + (mul_lo & 0xffffffff) return (result_hi, result_lo) } 

Aquí hay un ejemplo para verificar, que la función funciona:

 var c: UInt64 = 0 var d: UInt64 = 0 (c, d) = mul64to128(0x1234567890123456, 0x9876543210987654) // 0AD77D742CE3C72E45FD10D81D28D038 is the result of the above example print(hex128(c, d)) (c, d) = mul64to128(0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF) // FFFFFFFFFFFFFFFE0000000000000001 is the result of the above example print(hex128(c, d)) 

Si solo desea detectar desbordamiento, ¿qué hay de convertir a doble, hacer la multiplicación y si

| x | <2 ^ 53, convierta a int64

| x | <2 ^ 63, haz la multiplicación usando int64

de lo contrario, producirá el error que desee?

Esto parece funcionar:

 int64_t safemult(int64_t a, int64_t b) { double dx; dx = (double)a * (double)b; if ( fabs(dx) < (double)9007199254740992 ) return (int64_t)dx; if ( (double)INT64_MAX < fabs(dx) ) return INT64_MAX; return a*b; } 

Aquí hay un truco para detectar si la multiplicación de dos enteros sin signo se desborda.

Hacemos la observación de que si multiplicamos un número binario de N-bit-wide con un número binario de M-bit-wide, el producto no tiene más de N + M bits.

Por ejemplo, si se nos pide que multipliquemos un número de tres bits con un número de veintinueve bits, sabemos que esto no sobrepasa los treinta y dos bits.

 #include  #include  int might_be_mul_oflow(unsigned long a, unsigned long b) { if (!a || !b) return 0; a = a | (a >> 1) | (a >> 2) | (a >> 4) | (a >> 8) | (a >> 16) | (a >> 32); b = b | (b >> 1) | (b >> 2) | (b >> 4) | (b >> 8) | (b >> 16) | (b >> 32); for (;;) { unsigned long na = a << 1; if (na <= a) break; a = na; } return (a & b) ? 1 : 0; } int main(int argc, char **argv) { unsigned long a, b; char *endptr; if (argc < 3) { printf("supply two unsigned long integers in C form\n"); return EXIT_FAILURE; } a = strtoul(argv[1], &endptr, 0); if (*endptr != 0) { printf("%s is garbage\n", argv[1]); return EXIT_FAILURE; } b = strtoul(argv[2], &endptr, 0); if (*endptr != 0) { printf("%s is garbage\n", argv[2]); return EXIT_FAILURE; } if (might_be_mul_oflow(a, b)) printf("might be multiplication overflow\n"); { unsigned long c = a * b; printf("%lu * %lu = %lu\n", a, b, c); if (a != 0 && c / a != b) printf("confirmed multiplication overflow\n"); } return 0; } 

Un puñado de pruebas: (en el sistema de 64 bits):

 $ ./ flujo 0x3 0x3FFFFFFFFFFFFFFF
 3 * 4611686018427387903 = 13835058055282163709

 $ ./ flujo 0x7 0x3FFFFFFFFFFFFFFF
 podría ser desbordamiento de multiplicación
 7 * 4611686018427387903 = 13835058055282163705
 desbordamiento de multiplicación confirmado

 $ ./uflow 0x4 0x3FFFFFFFFFFFFFFF
 podría ser desbordamiento de multiplicación
 4 * 4611686018427387903 = 18446744073709551612

 $ ./ flujo 0x5 0x3FFFFFFFFFFFFFFF
 podría ser desbordamiento de multiplicación
 5 * 4611686018427387903 = 4611686018427387899
 desbordamiento de multiplicación confirmado 

Los pasos en might_be_mul_oflow son casi seguramente más lentos que solo hacer la prueba de división, al menos en procesadores convencionales utilizados en estaciones de trabajo de escritorio, servidores y dispositivos móviles. En chips sin buen soporte de división, podría ser útil.


Se me ocurre que hay otra forma de hacer esta prueba de rechazo temprano.

  1. Comenzamos con un par de números arng y brng que se inicializan a 0x7FFF...FFFF y 1 .

  2. Si a <= arng b <= brng podemos concluir que no hay desbordamiento.

  3. De lo contrario, desplazamos arng a la derecha, y brng a brng a la izquierda, agregando un bit a brng , de modo que sean 0x3FFF...FFFF y 3 .

  4. Si arng es cero, termina; de lo contrario, repita a las 2.

La función ahora se ve así:

 int might_be_mul_oflow(unsigned long a, unsigned long b) { if (!a || !b) return 0; { unsigned long arng = ULONG_MAX >> 1; unsigned long brng = 1; while (arng != 0) { if (a <= arng && b <= brng) return 0; arng >>= 1; brng <<= 1; brng |= 1; } return 1; } }