Как избежать переполнения в expr. A * B - C * D

161

Мне нужно вычислить выражение, которое выглядит следующим образом:, A*B - C*Dгде их типы: signed long long int A, B, C, D; Каждое число может быть очень большим (не выходя за пределы его типа). Хотя A*Bможет вызвать переполнение, в то же время выражение A*B - C*Dможет быть очень маленьким. Как я могу вычислить это правильно?

Например:, MAX * MAX - (MAX - 1) * (MAX + 1) == 1где MAX = LLONG_MAX - nи n - некоторое натуральное число.

NGix
источник
17
Насколько важна точность?
Анируд Раманатан
1
@Cthulhu, отличный вопрос. Он мог попытаться сделать эквивалентную функцию, используя меньшее число, разделив их все на 10 или что-то еще, а затем умножив результат.
Крис
4
Vars A, B, C, D подписаны. Это подразумевает, что A - Cможет переполниться. Это вопрос, который стоит рассмотреть, или вы знаете, что это не произойдет с вашими данными?
Уильям Моррис
2
@MooingDuck, но вы можете заранее проверить, не переполнится ли операция stackoverflow.com/a/3224630/158285
bradgonesurfing
1
@ Крис: Нет, я говорю, что нет никакого портативного способа проверить, произошло ли переполнение со знаком. (Brad правильно , что вы можете переносимо обнаружить , что это будет происходить). Использование встроенной сборки - один из многих непереносимых способов проверки.
Mooing Duck

Ответы:

120

Это кажется слишком тривиальным, я думаю. Но A*Bтот, который может переполниться.

Вы могли бы сделать следующее, не теряя точности

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)

Это разложение может быть сделано в дальнейшем .
Как отметил @Gian, во время операции вычитания может потребоваться осторожность, если тип длинный без знака.


Например, в случае, если у вас есть вопрос, это займет всего одну итерацию,

 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
Анирудх Раманатан
источник
4
@Caleb, просто примените тот же алгоритм кC*D
Крис
2
Я думаю, вы должны объяснить, что Е представляет.
Калеб
7
И long long и double - это 64 бита. Поскольку double должен выделить несколько битов для показателя степени, он имеет меньший диапазон возможных значений без потери точности.
Джим Гаррисон
3
@Cthulhu - мне кажется, что это сработало бы, только если все числа очень велики ... например, у вас все равно будет переполнение {A, B, C, D} = {MAX, MAX, MAX, 2}. ОП говорит: «Каждое число может быть действительно большим», но из постановки задачи не ясно, что каждое число должно быть действительно большим.
Кевин К
4
Что если A,B,C,Dхотя бы один из них отрицательный? Не будет Eили Fбудет еще больше тогда?
Supr
68

Самое простое и наиболее общее решение состоит в том, чтобы использовать представление, которое не может переполниться, либо с помощью библиотеки длинных целых чисел (например, http://gmplib.org/ ), либо с помощью представления структуры или массива и реализации своего рода длинного умножения ( т.е. разделение каждого числа на две 32-битные половинки и выполнение умножения, как показано ниже:

(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

Предполагая, что конечный результат умещается в 64 бита, вам на самом деле не нужно большинство битов R3 и ни одного из R4

Офир
источник
8
Вышеприведенный расчет на самом деле не так сложен, как кажется, на самом деле это простое длинное умножение в базе 2 ^ 32, и код на C должен выглядеть проще. Также будет хорошей идеей создать универсальные функции для выполнения этой работы в вашей программе.
Офир
46

Обратите внимание, что это не является стандартным, поскольку оно основано на переполнении со знаком. (GCC имеет флаги компилятора, которые позволяют это.)

Но если вы просто выполните все вычисления long long, результат применения формулы напрямую:
(A * B - C * D)будет точным, если правильный результат вписывается в long long.


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

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

Это unsigned long longприводит к тому, что стандарт переполнения гарантирует, что поведение переполнения гарантируется. Возвращение к целому числу со знаком в конце - это часть, определяемая реализацией, но сегодня она будет работать практически во всех средах.


Если вам нужно более педантичное решение, я думаю, что вы должны использовать «длинную арифметику»

Риад
источник
+1 Ты единственный, кто это заметил. Единственная сложная задача - настроить компилятор для выполнения переполнения и проверки, действительно ли правильный результат вписывается в long long.
Мистика
2
Даже наивная версия без каких-либо хитростей будет делать правильные вещи в большинстве реализаций; это не гарантируется стандартом, но вам придется найти машину с дополнением 1 или какое-то другое странное устройство, чтобы оно не сработало.
Хоббс
1
Я думаю, что это важный ответ. Я согласен, что это может быть неправильное программирование, предполагающее поведение, специфичное для реализации, но каждый инженер должен понимать арифметику по модулю и то, как получить правильные флаги компилятора для обеспечения согласованного поведения, если производительность важна. Инженеры DSP полагаются на такое поведение для реализаций фильтра с фиксированной запятой, для которых принятый ответ будет иметь недопустимую производительность.
Питер М
18

Это должно работать (я думаю):

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

Вот мой вывод:

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)
paquetp
источник
1
Спасибо @bradgonesurfing - не могли бы вы предоставить такой вклад? Я обновил свой ответ, выполнил его, и он работает (bd и ca равны 0) ...
paquetp
1
Хммм. Теперь я думаю об этом, может быть, нет. Вырожденный случай с d = 1 и a = 1 и b = maxint и c = maxint все еще работает. Круто :)
bradgonesurfing
1
@paquetp: a = 1, b = 0x7fffffffffffffff, c = -0x7fffffffffffffff, d = 1 (примечание c является отрицательным). Умно, хотя, я уверен, что ваш код правильно обрабатывает все положительные числа.
Mooing Duck
3
@MooingDuck, но окончательный ответ для вашего набора также переполнен, так что это недопустимая настройка. Это работает, только если каждая сторона имеет один и тот же знак, поэтому полученное вычитание находится в пределах диапазона.
Bradgonesurfing
1
Что-то странное в StackOverflow, когда этот ответ, который является самым простым и лучшим, получил такой низкий балл по сравнению с ответом с самым высоким баллом.
Bradgonesurfing
9

Вы могли бы рассмотреть вычисление как наибольший общий фактор для всех ваших значений, а затем разделить их на этот коэффициент, прежде чем выполнять арифметические операции, а затем снова умножить. Это предполагает , что такой фактор существует, однако (например, если A, B, Cи Dслучается быть относительно простыми, они не будут иметь общий фактор).

Точно так же вы могли бы рассмотреть возможность работы с логарифмическими масштабами, но это будет немного страшно, учитывая числовую точность.

Джан
источник
1
Логарифмирование кажется хорошим, если long doubleдоступно. В этом случае может быть достигнут приемлемый уровень точности (и результат может быть округлен).
9

Если результат помещается в long long int, тогда выражение A * BC * D нормально, поскольку оно выполняет арифметический мод 2 ^ 64 и даст правильный результат. Проблема в том, чтобы узнать, подходит ли результат в long long int. Чтобы обнаружить это, вы можете использовать следующий трюк с использованием double:

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

Проблема этого подхода заключается в том, что вы ограничены точностью мантиссы двойных чисел (54 бита?), Поэтому вам нужно ограничить произведения A * B и C * D до 63 + 54 бита (или, возможно, чуть меньше).

Эстебан Креспи
источник
Это самый практичный пример. Очистить и дает правильный ответ (или выдает исключение, когда входные данные плохие).
Марк Лаката
1
Красиво и элегантно! Вы не попали в ловушку, в которую попали другие. Еще одна вещь: я бы поспорил, что есть несколько примеров, когда двойной расчет ниже MAX_LLONG только из-за ошибок округления. Мой математический инстинкт подсказывает мне, что вы должны вместо этого вычислить разницу между двойным и длинным результатом и сравнить это с MAX_LLONG / 2 или чем-то еще. Эта разница является ошибками округления двойного вычисления и плюс переполнение и обычно должна быть относительно низкой, но в случае, о котором я упоминал, она будет большой. Но сейчас мне лень выяснить наверняка. :-)
Hans-Peter Störr
9
E = max(A,B,C,D)
A1 = A -E;
B1 = B -E;
C1 = C -E;
D1 = D -E;

затем

A*B - C*D = (A1+E)*(B1+E)-(C1+E)(D1+E) = (A1+B1-C1-D1)*E + A1*B1 -C1*D1
Ландри
источник
7

Вы можете записать каждое число в массив, каждый элемент является цифрой и выполнять вычисления в виде полиномов . Возьмите получившийся полином, который является массивом, и вычислите результат, умножив каждый элемент массива на 10 на степень позиции в массиве (первая позиция является наибольшей, а последняя - нулем).

Число 123может быть выражено как:

123 = 100 * 1 + 10 * 2 + 3

для которого вы просто создаете массив [1 2 3].

Вы делаете это для всех чисел A, B, C и D, а затем умножаете их на полиномы. Получив полученный многочлен, вы просто восстанавливаете число по нему.

Михай
источник
2
не знаю что это, но я должен найти. ставить :) . это решение головы, пока я делаю покупки с моей девушкой :)
Михай
вы реализуете bignums в массиве base10. GMP - это качественная библиотека bignum, использующая базу 4294967296. НАМНОГО быстрее. Хотя нет никакого отрицательного ответа, потому что ответ правильный и полезный.
Mooing Duck
Спасибо :) . Полезно знать, что это один из способов сделать это, но есть лучшие способы, поэтому не делайте так. по крайней мере, не в этой ситуации :)
Михай
в любом случае ... используя это решение, вы можете вычислить намного большее число, чем любой примитивный тип, выделенный жирным шрифтом (например, 100-значное число s), и сохранить результат в виде массива. это заслуживает голосования: p
Михай
Я не уверен, что это произойдет, так как этот метод (хотя и эффективный и относительно простой для понимания) требует много памяти и медленно работает.
Mooing Duck
6

Пока signed long long intне проведут A*B, двое из них будут. Таким образом, A*Bможно разложить на три разных показателя степени, подходящих для любого из них signed long long int.

A1=A>>32;
A0=A & 0xffffffff;
B1=B>>32;
B0=B & 0xffffffff;

AB_0=A0*B0;
AB_1=A0*B1+A1*B0;
AB_2=A1*B1;

То же самое для C*D .

Следуя прямому пути, вычитание может быть выполнено для каждой пары, AB_iа CD_iтакже с использованием дополнительного бита переноса (точно 1-битное целое число) для каждой. Поэтому, если мы скажем E = A * BC * D, вы получите что-то вроде:

E_00=AB_0-CD_0 
E_01=(AB_0 > CD_0) == (AB_0 - CD_0 < 0) ? 0 : 1  // carry bit if overflow
E_10=AB_1-CD_1 
...

Мы продолжаем, перенося верхнюю половину E_10в E_20(сдвинуть на 32 и добавить, затем стереть верхнюю половину E_10).

Теперь вы можете избавиться от бита для переноса E_11, добавив его с нужным знаком (полученным из части без переноса) в E_20. Если это вызывает переполнение, результат не будет соответствовать.

E_10теперь достаточно места, чтобы взять верхнюю половину E_00 (сдвиг, добавление, стирание) и бит переносаE_01 .

E_10 может быть больше теперь снова, поэтому мы повторяем передачу E_20 .

В этот момент E_20должен стать ноль, иначе результат не будет соответствовать. В E_10результате переноса верхняя половина также пуста.

Последний шаг - перенести нижнюю половину E_20 в E_10снова.

Если ожидание, которое E=A*B+C*Dбудет соответствовать ожиданиям, у signed long long intнас теперь есть

E_20=0
E_10=0
E_00=E
dronus
источник
1
На самом деле это упрощенная формула, которую можно получить, если использовать формулу умножения Офира и удалить все бесполезные временные результаты.
Дронус
3

Если вы знаете, что конечный результат представлен в целочисленном типе, вы можете быстро выполнить этот расчет, используя код ниже. Поскольку в стандарте C указано, что арифметика без знака является арифметикой по модулю и не переполняется, вы можете использовать тип без знака для выполнения вычислений.

В следующем коде предполагается, что существует неподписанный тип такой же ширины, и что подписанный тип использует все битовые шаблоны для представления значений (нет представлений ловушек, минимум подписанного типа является отрицательным от половины модуля беззнакового типа). Если это не выполняется в реализации C, для этого можно внести простые корректировки в процедуру ConvertToSigned.

Следующее использует signed charи unsigned charдля демонстрации кода. Для вашей реализации измените определение Signedна typedef signed long long int Signed;и определение Unsignedна typedef unsigned long long int Unsigned;.

#include <limits.h>
#include <stdio.h>
#include <stdlib.h>


//  Define the signed and unsigned types we wish to use.
typedef signed char   Signed;
typedef unsigned char Unsigned;

//  uHalfModulus is half the modulus of the unsigned type.
static const Unsigned uHalfModulus = UCHAR_MAX/2+1;

//  sHalfModulus is the negation of half the modulus of the unsigned type.
static const Signed   sHalfModulus = -1 - (Signed) (UCHAR_MAX/2);


/*  Map the unsigned value to the signed value that is the same modulo the
    modulus of the unsigned type.  If the input x maps to a positive value, we
    simply return x.  If it maps to a negative value, we return x minus the
    modulus of the unsigned type.

    In most C implementations, this routine could simply be "return x;".
    However, this version uses several steps to convert x to a negative value
    so that overflow is avoided.
*/
static Signed ConvertToSigned(Unsigned x)
{
    /*  If x is representable in the signed type, return it.  (In some
        implementations, 
    */
    if (x < uHalfModulus)
        return x;

    /*  Otherwise, return x minus the modulus of the unsigned type, taking
        care not to overflow the signed type.
    */
    return (Signed) (x - uHalfModulus) - sHalfModulus;
}


/*  Calculate A*B - C*D given that the result is representable as a Signed
    value.
*/
static signed char Calculate(Signed A, Signed B, Signed C, Signed D)
{
    /*  Map signed values to unsigned values.  Positive values are unaltered.
        Negative values have the modulus of the unsigned type added.  Because
        we do modulo arithmetic below, adding the modulus does not change the
        final result.
    */
    Unsigned a = A;
    Unsigned b = B;
    Unsigned c = C;
    Unsigned d = D;

    //  Calculate with modulo arithmetic.
    Unsigned t = a*b - c*d;

    //  Map the unsigned value to the corresponding signed value.
    return ConvertToSigned(t);
}


int main()
{
    //  Test every combination of inputs for signed char.
    for (int A = SCHAR_MIN; A <= SCHAR_MAX; ++A)
    for (int B = SCHAR_MIN; B <= SCHAR_MAX; ++B)
    for (int C = SCHAR_MIN; C <= SCHAR_MAX; ++C)
    for (int D = SCHAR_MIN; D <= SCHAR_MAX; ++D)
    {
        //  Use int to calculate the expected result.
        int t0 = A*B - C*D;

        //  If the result is not representable in signed char, skip this case.
        if (t0 < SCHAR_MIN || SCHAR_MAX < t0)
            continue;

        //  Calculate the result with the sample code.
        int t1 = Calculate(A, B, C, D);

        //  Test the result for errors.
        if (t0 != t1)
        {
            printf("%d*%d - %d*%d = %d, but %d was returned.\n",
                A, B, C, D, t0, t1);
            exit(EXIT_FAILURE);
        }
    }
    return 0;
}
Эрик Постпищил
источник
2

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

AB - CD
= [ A(B - N) - C( D - M )] + [AN - CM]

= ( AK - CJ ) + ( AN - CM)

    where K = B - N
          J = D - M

Если компоненты все еще переполнены, вы можете рекурсивно разбить их на более мелкие компоненты, а затем рекомбинировать.

bradgonesurfing
источник
Это может или не может быть правильно, но определенно сбивает с толку. Вы определяете, Kа Jпочему бы Nи нет M. Кроме того, я думаю, что вы разбиваете уравнение на большие части. Так как ваш шаг 3 такой же, как вопрос ОП, за исключением более сложного (AK-CJ)->(AB-CD)
Mooing Duck
N не упрощается ни от чего. Это просто число, вычтенное из A, чтобы сделать его меньше. На самом деле это похожее, но худшее решение для Paquetp. Здесь я использую вычитание вместо целочисленного деления, чтобы сделать его меньше.
Bradgonesurfing
2

Возможно, я не охватил все граничные случаи и не провел тщательного тестирования, но это реализует технику, которую я помню, использовал в 80-х годах при попытке выполнить 32-разрядное целочисленное вычисление на 16-разрядном процессоре. По сути, вы разбиваете 32 бита на два 16-битных блока и работаете с ними отдельно.

public class DoubleMaths {
  private static class SplitLong {
    // High half (or integral part).
    private final long h;
    // Low half.
    private final long l;
    // Split.
    private static final int SPLIT = (Long.SIZE / 2);

    // Make from an existing pair.
    private SplitLong(long h, long l) {
      // Let l overflow into h.
      this.h = h + (l >> SPLIT);
      this.l = l % (1l << SPLIT);
    }

    public SplitLong(long v) {
      h = v >> SPLIT;
      l = v % (1l << SPLIT);
    }

    public long longValue() {
      return (h << SPLIT) + l;
    }

    public SplitLong add ( SplitLong b ) {
      // TODO: Check for overflow.
      return new SplitLong ( longValue() + b.longValue() );
    }

    public SplitLong sub ( SplitLong b ) {
      // TODO: Check for overflow.
      return new SplitLong ( longValue() - b.longValue() );
    }

    public SplitLong mul ( SplitLong b ) {
      /*
       * e.g. 10 * 15 = 150
       * 
       * Divide 10 and 15 by 5
       * 
       * 2 * 3 = 5
       * 
       * Must therefore multiply up by 5 * 5 = 25
       * 
       * 5 * 25 = 150
       */
      long lbl = l * b.l;
      long hbh = h * b.h;
      long lbh = l * b.h;
      long hbl = h * b.l;
      return new SplitLong ( lbh + hbl, lbl + hbh );
    }

    @Override
    public String toString () {
      return Long.toHexString(h)+"|"+Long.toHexString(l);
    }
  }

  // I'll use long and int but this can apply just as easily to long-long and long.
  // The aim is to calculate A*B - C*D without overflow.
  static final long A = Long.MAX_VALUE;
  static final long B = Long.MAX_VALUE - 1;
  static final long C = Long.MAX_VALUE;
  static final long D = Long.MAX_VALUE - 2;

  public static void main(String[] args) throws InterruptedException {
    // First do it with BigIntegers to get what the result should be.
    BigInteger a = BigInteger.valueOf(A);
    BigInteger b = BigInteger.valueOf(B);
    BigInteger c = BigInteger.valueOf(C);
    BigInteger d = BigInteger.valueOf(D);
    BigInteger answer = a.multiply(b).subtract(c.multiply(d));
    System.out.println("A*B - C*D = "+answer+" = "+answer.toString(16));

    // Make one and test its integrity.
    SplitLong sla = new SplitLong(A);
    System.out.println("A="+Long.toHexString(A)+" ("+sla.toString()+") = "+Long.toHexString(sla.longValue()));

    // Start small.
    SplitLong sl10 = new SplitLong(10);
    SplitLong sl15 = new SplitLong(15);
    SplitLong sl150 = sl10.mul(sl15);
    System.out.println("10="+sl10.longValue()+"("+sl10.toString()+") * 15="+sl15.longValue()+"("+sl15.toString()+") = "+sl150.longValue() + " ("+sl150.toString()+")");

    // The real thing.
    SplitLong slb = new SplitLong(B);
    SplitLong slc = new SplitLong(C);
    SplitLong sld = new SplitLong(D);
    System.out.println("B="+Long.toHexString(B)+" ("+slb.toString()+") = "+Long.toHexString(slb.longValue()));
    System.out.println("C="+Long.toHexString(C)+" ("+slc.toString()+") = "+Long.toHexString(slc.longValue()));
    System.out.println("D="+Long.toHexString(D)+" ("+sld.toString()+") = "+Long.toHexString(sld.longValue()));
    SplitLong sanswer = sla.mul(slb).sub(slc.mul(sld));
    System.out.println("A*B - C*D = "+sanswer+" = "+sanswer.longValue());

  }

}

Печать:

A*B - C*D = 9223372036854775807 = 7fffffffffffffff
A=7fffffffffffffff (7fffffff|ffffffff) = 7fffffffffffffff
10=10(0|a) * 15=15(0|f) = 150 (0|96)
B=7ffffffffffffffe (7fffffff|fffffffe) = 7ffffffffffffffe
C=7fffffffffffffff (7fffffff|ffffffff) = 7fffffffffffffff
D=7ffffffffffffffd (7fffffff|fffffffd) = 7ffffffffffffffd
A*B - C*D = 7fffffff|ffffffff = 9223372036854775807

который выглядит для меня, как будто он работает.

Бьюсь об заклад, я пропустил некоторые тонкости, такие как наблюдение за переполнением знака и т. Д., Но я думаю, что суть есть.

OldCurmudgeon
источник
1
Я думаю, что это реализация того, что предложил @Ofir.
OldCurmudgeon
2

Для полноты картины, поскольку никто не упомянул об этом, некоторые компиляторы (например, GCC) фактически предоставляют вам 128-битное целое число в настоящее время.

Таким образом, простое решение может быть:

(long long)((__int128)A * B - (__int128)C * D)
я код 4 еды
источник
1

AB-CD = (AB-CD) * AC / AC = (B/C-D/A)*A*C, Ни B/Cни D/Aможет переполниться, поэтому не рассчитать (B/C-D/A)первый. Поскольку окончательный результат не будет переполнен в соответствии с вашим определением, вы можете безопасно выполнить оставшиеся умножения и рассчитать, (B/C-D/A)*A*Cкакой результат является требуемым.

Обратите внимание: если ваш ввод также может быть очень маленьким , B/Cили D/Aможет переполниться. Если это возможно, могут потребоваться более сложные манипуляции в соответствии с входной проверкой.

SomeWittyUsername
источник
2
Это не сработает, так как целочисленное деление теряет информацию (часть результата)
Офир
@ Да, верно, но вы не можете съесть торт и не трогать его. Вы должны платить либо с точностью, либо с помощью дополнительных ресурсов (как вы предложили в своем ответе). Мой ответ имеет математическую природу, а ваш - компьютерный. Каждый может быть правильным в зависимости от обстоятельств.
SomeWittyUsername
2
Вы правы - я должен был сформулировать это так - не даст точного результата, а не сработает, поскольку математика верна. Тем не менее, обратите внимание, что в случаях, которые могут заинтересовать отправителя вопроса (например, в примере в вопросе), ошибка, вероятно, будет удивительно большой - гораздо большей, чем может быть приемлемо для любого практического применения. В любом случае - это был проницательный ответ, и я не должен был использовать этот язык.
Офир
@Ofir Я не думаю, что ваш язык неуместен. ОП явно запрашивал «правильный» расчет, а не тот, который потерял бы точность из-за того, что выполнялся в условиях ограниченных ресурсов.
user4815162342
1

Выберите K = a big number(например K = A - sqrt(A))

A*B - C*D = (A-K)*(B-K) - (C-K)*(D-K) + K*(A-C+B-D); // Avoid overflow.

Зачем?

(A-K)*(B-K) = A*B - K*(A+B) + K^2
(C-K)*(D-K) = C*D - K*(C+D) + K^2

=>
(A-K)*(B-K) - (C-K)*(D-K) = A*B - K*(A+B) + K^2 - {C*D - K*(C+D) + K^2}
(A-K)*(B-K) - (C-K)*(D-K) = A*B - C*D - K*(A+B) + K*(C+D) + K^2 - K^2
(A-K)*(B-K) - (C-K)*(D-K) = A*B - C*D - K*(A+B-C-D)

=>
A*B - C*D = (A-K)*(B-K) - (C-K)*(D-K) + K*(A+B-C-D)

=>
A*B - C*D = (A-K)*(B-K) - (C-K)*(D-K) + K*(A-C+B-D)

Обратите внимание, что, поскольку A, B, C и D являются большими числами, следовательно, A-Cи B-Dявляются маленькими числами.

Амир Саниян
источник
Как вы выбираете K на практике? Кроме того, K * (A-C + BD) может все еще переполниться.
YLC
@ylc: Выберите K = sqrt (A), не так уж A-C+B-Dмало. Поскольку A, B, C и D - большие числа, таким образом, AC - небольшое число.
Амир Саниян
Если вы выберете K = sqrt (A) , то (AK) * (BK) может снова переполниться.
2012 г.
@ylc: ОК! Я изменяю это на A - sqrt(A):)
Амир Саниян
Тогда K * (A-C + BD) может переполниться.
YLC