Поймайте и вычислите переполнение во время умножения двух больших целых чисел

Идея состоит в том, чтобы использовать следующий факт, который верен для интегральной операции:

a * b> c тогда и только тогда, когда a> c / b

/ здесь целочисленное деление.

Псевдокод для проверки переполнения для положительных чисел следующий:

if (a> max_int64 / b) then "overflow" else "ok" .

Для обработки нулей и отрицательных чисел вам следует добавить дополнительные проверки.

Код C для неотрицательных a и b следует:

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

Примечание:

18446744073709551615 == (1<<64)-1

Чтобы вычислить перенос, мы можем использовать подход, чтобы разделить число на две 32-значные цифры и умножить их, как мы делаем это на бумаге. Нам нужно разделить числа, чтобы избежать переполнения.

Код следующий:

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

Версия, которая также работает при a == 0:

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

Если вам нужно не только обнаружить переполнение, но и зафиксировать перенос, лучше всего разбить числа на 32-битные части. Код - кошмар; далее следует лишь набросок:

#include <stdint.h>

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

Проблема заключается не только в частичных произведениях, но и в том, что любая из сумм может переполниться.

Если бы мне пришлось делать это по-настоящему, я бы написал процедуру расширенного умножения на местный язык ассемблера. То есть, например, умножьте два 64-битных целых числа, чтобы получить 128-битный результат, который сохраняется в двух 64-битных регистрах. Все разумное оборудование обеспечивает эту функциональность в единственной собственной инструкции умножения - она ​​доступна не только из C.

Это один из тех редких случаев, когда наиболее элегантным и простым в программировании решением является фактически использование языка ассемблера. Но это уж точно не переносно: