Overflow di cattura e calcolo durante la moltiplicazione di due interi grandi

Sto cercando una soluzione efficiente (facoltativamente standard, elegante e facile da implementare) per moltiplicare numeri relativamente grandi e memorizzare il risultato in uno o più interi:

Diciamo che ho due interi a 64 bit dichiarati in questo modo:

uint64_t a = xxx, b = yyy; 

Quando faccio a * b , come posso rilevare se l’operazione si traduce in un trabocco e in questo caso conservare il carry da qualche parte?

Tieni presente che non voglio utilizzare alcuna libreria di numeri importanti poiché ho dei vincoli sul modo in cui memorizzo i numeri.

1. Rilevazione del trabocco :

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

Modifica: riparato per 0 (grazie Marco!)

2. Calcolare il carry è abbastanza complicato. Un approccio è quello di suddividere entrambi gli operandi in mezze parole, quindi applicare la moltiplicazione lunga alle mezze parole:

 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; } 

Per vedere che nessuno degli stessi somme parziali può traboccare, consideriamo il caso peggiore:

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

Sia B = 1 << 32 . Allora abbiamo

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

Credo che funzionerà, almeno gestirà il caso di Sjlver. A parte questo, non è testato (e potrebbe anche non essere compilato, dato che non ho più un compilatore C ++ a portata di mano).

L’idea è di usare il seguente fatto che è vero per l’operazione integrale:

a*b > c se e solo se a > c/b

/ è una divisione integrale qui.

Segue lo pseudocodice per verificare l’overflow per i numeri positivi:

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

Per gestire gli zeri e i numeri negativi è necessario aggiungere ulteriori controlli.

Il codice C per non negativi a e b segue:

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

Nota:

 18446744073709551615 == (1<<64)-1 

Per calcolare il carry possiamo usare l'approccio per dividere il numero in due 32 cifre e moltiplicarle come facciamo sul foglio. Abbiamo bisogno di dividere i numeri per evitare l'overflow.

Il codice segue:

 // 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 

Sebbene ci siano state molte altre risposte a questa domanda, molti di loro hanno un codice completamente non testato, e finora nessuno ha paragonato adeguatamente le diverse opzioni possibili.

Per questo motivo, ho scritto e testato diverse possibili implementazioni (l’ultima è basata su questo codice di OpenBSD, discusso su Reddit qui ). Ecco il codice:

 /* 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); } 

Ecco i risultati, test con vari compilatori e sistemi che ho (in questo caso, tutti i test sono stati eseguiti su OS X, ma i risultati dovrebbero essere simili su sistemi 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 | +------------------+----------+----------+----------+----------+----------+ 

Sulla base di questi risultati, possiamo trarre alcune conclusioni:

  • Chiaramente, l'approccio basato sulla divisione, sebbene semplice e portatile, è lento.
  • Nessuna tecnica è un chiaro vincitore in tutti i casi.
  • Sui moderni compilatori, l'approccio use-a-large-int è il migliore, se puoi usarlo
  • Nei compilatori precedenti, l'approccio di moltiplicazione lunga è il migliore
  • Sorprendentemente, GCC 4.9.0 ha regressioni di prestazioni rispetto a GCC 4.2.1 e GCC 4.2.1 ha regressioni di prestazioni rispetto a GCC 3.3

Una versione che funziona anche quando a == 0:

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

Se è necessario non solo per rilevare l’overflow ma anche per acquisire il carry, è meglio rompere i numeri in parti a 32 bit. Il codice è un incubo; ciò che segue è solo uno schizzo:

 #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... } 

Il problema non sono solo i prodotti parziali, ma il fatto che una qualsiasi delle somme può eccedere.

Se dovessi fare questo per davvero, scriverò una routine di moltiplicazione estesa nel linguaggio dell’assemblaggio locale. Ad esempio, moltiplicare due interi a 64 bit per ottenere un risultato a 128 bit, che viene memorizzato in due registri a 64 bit. Tutto l’hardware ragionevole fornisce questa funzionalità in un’unica istruzione multipla nativa, non solo accessibile da C.

Questo è uno di quei rari casi in cui la soluzione più elegante e facile da programmare è in realtà utilizzare il linguaggio assembly. Ma non è certamente portatile 🙁

Ho lavorato con questo problema in questi giorni e devo dire che mi ha colpito il numero di volte in cui ho visto persone dire che il modo migliore per sapere se c’è stato un overflow è dividere il risultato, questo è totalmente inefficiente e non necessario. Il punto di questa funzione è che deve essere il più veloce ansible.

Ci sono due opzioni per il rilevamento dell’overflow:

1º- Se ansible, crea la variabile risultato due volte più grande dei moltiplicatori, ad esempio:

 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) } 

Saprai immediatamente se c’è stato un overflow e il codice è il più veloce ansible senza scriverlo in codice macchina. A seconda del compilatore, questo codice può essere migliorato in codice macchina.

2º- È imansible creare una variabile di risultato due volte più grande della variabile moltiplicatori: Quindi dovresti giocare con se le condizioni per determinare il percorso migliore. Continuando con l’esempio:

 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); } } 

Spero che questo codice ti aiuti ad avere un programma abbastanza efficiente e spero che il codice sia chiaro, altrimenti inserirò qualche commento.

i migliori saluti.

Forse il modo migliore per risolvere questo problema è avere una funzione, che moltiplica due UInt64 e risultati una coppia di UInt64, una parte superiore e una parte inferiore del risultato di UInt128. Ecco la soluzione, inclusa una funzione, che mostra il risultato in esadecimale. Immagino che tu preferisca una soluzione C ++, ma ho una soluzione Swift funzionante che mostra come gestire il 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) } 

Ecco un esempio per verificare che la funzione funzioni:

 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)) 

Se vuoi solo rilevare un overflow, che ne dici di convertirlo in double, fare la moltiplicazione e if

| X | <2 ^ 53, convertire in int64

| X | <2 ^ 63, effettuare la moltiplicazione usando int64

altrimenti produci qualunque errore tu voglia?

Questo sembra funzionare:

 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; } 

Ecco un trucco per rilevare se la moltiplicazione di due interi non firmati trabocca.

Facciamo l’osservazione che se moltiplichiamo un numero binario di larghezza N-bit con un numero binario a livello di bit M, il prodotto non ha più di N + M bit.

Ad esempio, se ci viene chiesto di moltiplicare un numero di tre bit con un numero di ventinove bit, sappiamo che questo non trabocca trentadue bit.

 #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; } 

Una manciata di test: (sul sistema a 64 bit):

 $ ./uflow 0x3 0x3FFFFFFFFFFFFFFFFF
 3 * 4611686018427387903 = 13835058055282163709

 $ ./uflow 0x7 0x3FFFFFFFFFFFFFFFFF
 potrebbe essere un eccesso di moltiplicazione
 7 * 4611686018427387903 = 13835058055282163705
 overflow di moltiplicazione confermato

 $ ./uflow 0x4 0x3FFFFFFFFFFFFFFFFF
 potrebbe essere un eccesso di moltiplicazione
 4 * 4611686018427387903 = 18446744073709551612

 $ ./uflow 0x5 0x3FFFFFFFFFFFFFFFFF
 potrebbe essere un eccesso di moltiplicazione
 5 * 4611686018427387903 = 4611686018427387899
 overflow di moltiplicazione confermato 

I passaggi in might_be_mul_oflow sono quasi certamente più lenti del semplice test di divisione, almeno sui processori mainstream utilizzati nelle workstation desktop, server e dispositivi mobili. Su chip senza un buon supporto di divisione, potrebbe essere utile.


Mi viene in mente che esiste un altro modo per fare questo test di rigetto.

  1. Iniziamo con una coppia di numeri arng e arng che sono inizializzati su 0x7FFF...FFFF e 1 .

  2. Se a <= arng e b <= brng possiamo concludere che non vi è overflow.

  3. Altrimenti, arng verso destra, e brng verso sinistra, aggiungendo un bit a brng , in modo che siano 0x3FFF...FFFF e 3 .

  4. Se arng è zero, termina; altrimenti ripetere a 2.

La funzione ora sembra:

 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; } }