So vermeiden Sie einen Überlauf in expr. A B C D

Ich muss einen Ausdruck berechnen, der so aussieht:
A*B - C*Dwo ihre Typen sind: signed long long int A, B, C, D;
Jede Zahl kann wirklich groß sein (ohne ihren Typ zu überschreiten). Während A*B könnte Überlauf verursachen, gleichzeitig Ausdruck A*B - C*D kann wirklich klein sein. Wie kann ich es richtig berechnen?

Zum Beispiel: MAX * MAX - (MAX - 1) * (MAX + 1) == 1wo MAX = LLONG_MAX - n und n – irgendeine natürliche Zahl.

Das erscheint mir zu banal. Aber A*B ist derjenige, der überlaufen könnte.

Sie könnten Folgendes tun, ohne an Genauigkeit zu verlieren

A*B - C*D = A(D+E) - (A+F)D
          = AD + AE - AD - DF
          = AE - DF
             ^smaller quantities E & F

E = B - D (hence, far smaller than B)
F = C - A (hence, far smaller than C)

Diese Zerlegung kann sein weiter gemacht.
Wie @Gian betonte, muss während des Subtraktionsvorgangs möglicherweise darauf geachtet werden, wenn der Typ unsigned long long ist.

Bei dem Fall, den Sie in der Frage haben, dauert es beispielsweise nur eine Iteration,

 MAX * MAX - (MAX - 1) * (MAX + 1)
  A     B       C           D

E = B - D = -1
F = C - A = -1

AE - DF = {MAX * -1} - {(MAX + 1) * -1} = -MAX + MAX + 1 = 1

Die einfachste und allgemeinste Lösung besteht darin, eine Darstellung zu verwenden, die nicht überlaufen kann, entweder durch Verwendung einer Bibliothek mit langen Ganzzahlen (z http://gmplib.org/) oder die Verwendung einer Struktur oder eines Arrays darstellen und eine Art lange Multiplikation implementieren (dh jede Zahl in zwei 32-Bit-Hälften trennen und die Multiplikation wie folgt durchführen:

(R1 + R2 * 2^32 + R3 * 2^64 + R4 * 2^96) = R = A*B = (A1 + A2 * 2^32) * (B1 + B2 * 2^32) 
R1 = (A1*B1) % 2^32
R2 = ((A1*B1) / 2^32 + (A1*B2) % 2^32 + (A2*B1) % 2^32) % 2^32
R3 = (((A1*B1) / 2^32 + (A1*B2) % 2^32 + (A2*B1) % 2^32) / 2^32 + (A1*B2) / 2^32 + (A2*B1) / 2^32 + (A2*B2) % 2^32) %2^32
R4 = ((((A1*B1) / 2^32 + (A1*B2) % 2^32 + (A2*B1) % 2^32) / 2^32 + (A1*B2) / 2^32 + (A2*B1) / 2^32 + (A2*B2) % 2^32) / 2^32) + (A2*B2) / 2^32

Unter der Annahme, dass das Endergebnis in 64 Bit passt, benötigen Sie eigentlich nicht wirklich die meisten Bits von R3 und keines von R4

Beachten Sie, dass dies kein Standard ist, da es auf Wrap-Around-Signed-Overflow angewiesen ist. (GCC hat Compiler-Flags, die dies ermöglichen.)

Aber wenn Sie nur alle Berechnungen durchführen long longdas Ergebnis der direkten Anwendung der Formel:
(A * B - C * D) ist genau, solange das richtige Ergebnis in a passt long long.

Hier ist eine Problemumgehung, die sich nur auf das implementierungsdefinierte Verhalten des Umwandelns einer Ganzzahl ohne Vorzeichen in eine Ganzzahl mit Vorzeichen stützt. Es ist jedoch zu erwarten, dass dies heute auf fast jedem System funktioniert.

(long long)((unsigned long long)A * B - (unsigned long long)C * D)

Dies wirft die Eingaben auf unsigned long long wobei das Überlaufverhalten durch den Standard garantiert umlaufend ist. Das Zurücksetzen auf eine vorzeichenbehaftete Ganzzahl am Ende ist der von der Implementierung definierte Teil, funktioniert aber heute in fast allen Umgebungen.

Wenn Sie eine pedantischere Lösung benötigen, müssen Sie meiner Meinung nach “lange Arithmetik” verwenden.

Das sollte funktionieren (glaube ich):

signed long long int a = 0x7ffffffffffffffd;
signed long long int b = 0x7ffffffffffffffd;
signed long long int c = 0x7ffffffffffffffc;
signed long long int d = 0x7ffffffffffffffe;
signed long long int bd = b / d;
signed long long int bdmod = b % d;
signed long long int ca = c / a;
signed long long int camod = c % a;
signed long long int x = (bd - ca) * a * d - (camod * d - bdmod * a);

Hier meine Ableitung:

x = a * b - c * d
x / (a * d) = (a * b - c * d) / (a * d)
x / (a * d) = b / d - c / a

now, the integer/mod stuff:
x / (a * d) = (b / d + ( b % d ) / d) - (c / a + ( c % a ) / a )
x / (a * d) = (b / d - c / a) - ( ( c % a ) / a - ( b % d ) / d)
x = (b / d - c / a) * a * d - ( ( c % a ) * d - ( b % d ) * a)

E = max(A,B,C,D)
A1 = A -E;
B1 = B -E;
C1 = C -E;
D1 = D -E;

dann

A*B - C*D = (A1+E)*(B1+E)-(C1+E)(D1+E) = (A1+B1-C1-D1)*E + A1*B1 -C1*D1

Sie könnten erwägen, einen größten gemeinsamen Faktor für alle Ihre Werte zu berechnen und sie dann durch diesen Faktor zu dividieren, bevor Sie Ihre arithmetischen Operationen durchführen, und dann erneut multiplizieren. Dies setzt jedoch voraus, dass ein solcher Faktor existiert (z. B. wenn A, B, C und D zufällig relativ teilerfremd sind, haben sie keinen gemeinsamen Teiler).

In ähnlicher Weise könnten Sie in Betracht ziehen, auf logarithmischen Skalen zu arbeiten, aber dies wird ein wenig beängstigend sein, abhängig von der numerischen Genauigkeit.

Wenn das Ergebnis in ein long long int passt, ist der Ausdruck A*BC*D in Ordnung, da er die Arithmetik mod 2^64 ausführt und das richtige Ergebnis liefert. Das Problem besteht darin, zu wissen, ob das Ergebnis in ein long long int passt. Um dies zu erkennen, können Sie den folgenden Trick mit Doubles anwenden:

if( abs( (double)A*B - (double)C*D ) > MAX_LLONG ) 
    Overflow
else 
    return A*B-C*D;

Das Problem bei diesem Ansatz besteht darin, dass Sie durch die Genauigkeit der Mantisse der Doubles (54 Bit?) Begrenzt sind, sodass Sie die Produkte A * B und C * D auf 63 + 54 Bit (oder wahrscheinlich etwas weniger) begrenzen müssen.