При использовании удвоений, почему (x/(y * z)) не совпадает (x/y/z)?

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

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

@Test
public void shouldEqual() {
  double expected = 450.00d / (7d * 60);  // 1.0714285714285714
  double actual = 450.00d / 7d / 60;      // 1.0714285714285716

  assertThat(actual).isEqualTo(expected);
}

Но с этим выходом не получается:

org.junit.ComparisonFailure: 
Expected :1.0714285714285714
Actual   :1.0714285714285716

Может кто-нибудь объяснить подробно, что происходит под капотом, чтобы привести к тому, что значение на 1.000000000000000 X будет другим?

Некоторые из вопросов, которые я ищу в ответе, следующие: Где теряется точность? Какой метод является предпочтительным и почему? Что на самом деле правильно? (В чистой математике, оба не могут быть правы. Возможно, оба ошибаются?) Есть ли лучшее решение или метод для этих арифметических операций?

Ответ 1

Я вижу кучу вопросов, которые рассказывают вам, как обойти эту проблему, но не тот, который действительно объясняет, что происходит, кроме "ошибки округления с плавающей точкой - это плохо, m'kay?" Поэтому позвольте мне сделать снимок. Позвольте мне прежде всего отметить, что ничто в этом ответе не является специфическим для Java. Ошибка округления - это проблема, присущая любому фиксированному представлению чисел, поэтому вы получаете те же проблемы, например, в C.

Ошибка округления в десятичном типе данных

В качестве упрощенного примера предположим, что у нас есть какой-то компьютер, который изначально использует беззнаковый десятичный тип данных, пусть назовет его float6d. Длина типа данных - 6 цифр: 4, посвященная мантиссе, и 2, посвященная экспоненте. Например, число 3.142 может быть выражено как

3.142 x 10^0

который будет храниться в 6 цифрах как

Первые две цифры - показатель плюс 50, а последние четыре - мантисса. Этот тип данных может представлять любое число от 0.001 x 10^-50 до 9.999 x 10^+49.

Собственно, это не так. Он не может хранить номер. Что делать, если вы хотите представить 3.141592? Или 3.1412034? Или 3.141488906? Жесткая удача, тип данных не может хранить более четырех цифр точности, поэтому компилятор должен округлить все с большим количеством цифр, чтобы вписаться в ограничения типа данных. Если вы пишете

float6d x = 3.141592;
float6d y = 3.1412034;
float6d z = 3.141488906;

тогда компилятор преобразует каждое из этих трех значений в одно и то же внутреннее представление, 3.142 x 10^0 (которое, помните, сохраняется как 503142), так что x == y == z будет иметь значение true.

Дело в том, что существует целый ряд действительных чисел, которые все сопоставляются с одной и той же базовой последовательностью цифр (или бит, на реальном компьютере). В частности, любой x, удовлетворяющий 3.1415 <= x <= 3.1425 (предполагающий получетное округление), преобразуется в представление 503142 для хранения в памяти.

Это округление происходит каждый раз, когда ваша программа сохраняет значение с плавающей запятой в памяти. В первый раз это происходит, когда вы пишете константу в своем исходном коде, как я сделал выше, с x, y и z. Это происходит снова, когда вы выполняете арифметическую операцию, которая увеличивает количество цифр точности, превышающих то, что может представлять тип данных. Любой из этих эффектов называется ошибка округления. Это может произойти несколькими способами:

Сложение и вычитание: если одно из значений, которые вы добавляете, имеет другой показатель от другого, вы получите дополнительные цифры точности, и если их будет достаточно, наименее значимые будут необходимо отбросить. Например, 2.718 и 121.0 - оба значения, которые могут быть точно представлены в типе данных float6d. Но если вы попытаетесь добавить их вместе:
```
   1.210     x 10^2
+  0.02718   x 10^2
-------------------
   1.23718   x 10^2
```
который округляется до 1.237 x 10^2 или 123.7, отбрасывая две цифры точности.
Умножение: количество цифр в результате приблизительно равно количеству цифр в двух операндах. Это приведет к некоторой ошибке округления, если ваши операнды уже имеют много значащих цифр. Например, 121 x 2.718 дает вам
```
   1.210     x 10^2
x  0.02718   x 10^2
-------------------
   3.28878   x 10^2
```
который округляется до 3.289 x 10^2 или 328.9, снова опустив две цифры точности.

Однако полезно помнить, что если ваши операнды являются "хорошими" числами без значительных цифр, формат с плавающей запятой может, вероятно, точно представлять результат, поэтому вам не нужно иметь дело с ошибкой округления, Например, 2.3 x 140 дает
```
   1.40      x 10^2
x  0.23      x 10^2
-------------------
   3.22      x 10^2
```
который не имеет проблем округления.
Раздел: здесь все становится беспорядочным. Подразделение в значительной степени всегда приводит к некоторой сумме ошибок округления, если только количество, которое вы делите, не является базой (в этом случае деление представляет собой просто сдвиг цифр или сдвиг бит в двоичном формате). В качестве примера возьмите два очень простых числа, 3 и 7, разделите их, и вы получите
```
   3.                x 10^0
/  7.                x 10^0
----------------------------
   0.428571428571... x 10^0
```
Ближайшим значением для этого числа, которое может быть представлено как float6d, является 4.286 x 10^-1 или 0.4286, что явно отличается от точного результата.

Как мы увидим в следующем разделе, ошибка, возникающая при округлении, растет с каждой выполняемой вами операцией. Итак, , если вы работаете с "хорошими" номерами, как в вашем примере, как правило, лучше всего выполнять операции деления как можно позже, поскольку это операции, которые, скорее всего, вводят ошибку округления в вашу программу где раньше не существовало.

Анализ ошибки округления

В целом, если вы не можете предположить, что ваши цифры "хороши", ошибка округления может быть положительной или отрицательной, и очень сложно предсказать, в каком направлении она будет идти только на основе операции. Это зависит от конкретных значений. Посмотрите на этот график ошибки округления для 2.718 z как функцию z (все еще используя тип данных float6d):

На практике, когда вы работаете со значениями, использующими полную точность вашего типа данных, часто проще обрабатывать ошибку округления как случайную ошибку. Рассматривая сюжет, вы можете угадать, что величина ошибки зависит от порядка величины результата операции. В этом частном случае, когда z имеет порядок 10 ^-12.718 z также имеет порядок 10 ^-1 поэтому это будет число формы 0.XXXX. Максимальная ошибка округления - это половина последней цифры точности; в этом случае "последней цифрой точности" я имею в виду 0.0001, поэтому ошибка округления колеблется между -0.00005 и +0.00005. В точке, где 2.718 z переходит в следующий порядок величины, который равен 1/2.718 = 0.3679, вы можете видеть, что ошибка округления также скачет на порядок.

Вы можете использовать известные методы анализа ошибок , чтобы проанализировать, как случайная (или непредсказуемая) ошибка определенной величины влияет на ваш результат. В частности, для умножения или деления "средняя" относительная ошибка в вашем результате может быть аппроксимирована добавлением относительной ошибки в каждом из операндов в квадратуре, т.е. Их квадратом, их добавлением и выводом квадратного корня. С нашим типом данных float6d относительная ошибка изменяется от 0.0005 (для значения, такого как 0.101) и 0.00005 (для значения, равного 0.995).

Возьмите 0,0001 в качестве грубой средней относительной погрешности в значениях x и y. Относительная ошибка в x * y или x / y затем задается

sqrt(0.0001^2 + 0.0001^2) = 0.0001414

что является фактором sqrt(2) больше относительной ошибки в каждом из отдельных значений.

Когда дело доходит до объединения операций, вы можете применить эту формулу несколько раз, один раз для каждой операции с плавающей запятой. Так, например, для z / (x * y) относительная ошибка в x * y составляет в среднем 0,0001414 (в этом десятичном примере), а затем относительная ошибка в z / (x * y) равна

sqrt(0.0001^2 + 0.0001414^2) = 0.0001732

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

Аналогично, для z / x * y средняя относительная ошибка в z / x равна 0,0001414, а относительная ошибка в z / x * y равна

sqrt(0.0001414^2 + 0.0001^2) = 0.0001732

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

Сведения о горах

Вам может быть интересно узнать, какой конкретный расчет вы представили в вопросе, а не только в среднем. Для этого анализа перейдите в реальный мир двоичной арифметики. Номера плавающей запятой в большинстве систем и языков представлены с помощью стандарта IEEE 754. Для 64-битных чисел format указывает 52 бита, посвященных мантиссе, 11 - экспоненте, а один - знаку. Другими словами, при записи в базе 2 число с плавающей запятой является значением формы

1.1100000000000000000000000000000000000000000000000000 x 2^00000000010
                       52 bits                             11 bits

Ведущий 1 явно не хранится и составляет 53-й бит. Кроме того, вы должны отметить, что 11 бит, хранящихся для представления экспоненты, фактически являются реальным показателем плюс 1023. Например, это конкретное значение равно 7, что равно 1.75 x 2 ². Мантисса равна 1,75 в двоичном выражении или 1.11, а показатель составляет 1023 + 2 = 1025 в двоичном выражении или 10000000001, поэтому содержимое, хранящееся в памяти, составляет

01000000000111100000000000000000000000000000000000000000000000000
 ^          ^
 exponent   mantissa

но это не имеет большого значения.

В вашем примере также 450,

1.1100001000000000000000000000000000000000000000000000 x 2^00000001000

и 60,

1.1110000000000000000000000000000000000000000000000000 x 2^00000000101

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

Когда вы вычисляете первое выражение, 450/(7*60), процессор сначала выполняет умножение, получение 420 или

1.1010010000000000000000000000000000000000000000000000 x 2^00000001000

Затем он делит 450 на 420. Это дает 15/14, что составляет

1.0001001001001001001001001001001001001001001001001001001001001001001001...

в двоичном формате. Теперь спецификация языка Java говорит, что

Неточные результаты должны округляться до представляемого значения, ближайшего к бесконечно точному результату; если два ближайших представимых значения одинаково близки, выбирается единица с наименьшим значащим битом нуль. Это стандартный режим округления по стандарту IEEE 754, известный как от округленного до ближайшего.

а ближайшее представимое значение до 15/14 в 64-битном формате IEEE 754 -

1.0001001001001001001001001001001001001001001001001001 x 2^00000000000

что приблизительно равно 1.0714285714285714 в десятичной форме. (Точнее, это наименьшее точное десятичное значение, которое однозначно определяет это конкретное двоичное представление.)

С другой стороны, если вы сначала вычислите 450/7, результат будет 64.2857142857... или в двоичном формате,

1000000.01001001001001001001001001001001001001001001001001001001001001001...

для которого ближайшее представимое значение

1.0000000100100100100100100100100100100100100100100101 x 2^00000000110

который равен 64.28571428571429180465... Обратите внимание на изменение последней цифры бинарной мантиссы (по сравнению с точным значением) из-за ошибки округления. Разделение этого на 60 дает вам

1.000100100100100100100100100100100100100100100100100110011001100110011...

Посмотрите на конец: картина отличается! Это 0011, которое повторяется вместо 001, как в другом случае. Ближайшим представимым значением является

1.0001001001001001001001001001001001001001001001001010 x 2^00000000000

который отличается от другого порядка операций в последних двух битах: они 10 вместо 01. Десятичный эквивалент равен 1.0714285714285716.

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

1.0001001001001001001001001001001001001001001001001001001001001001001001...
1.0001001001001001001001001001001001001001001001001001100110011001100110...
                                                     ^ last bit of mantissa

В этом случае получается, что первый результат, численно 15/14, оказывается наиболее точным представлением точного значения. Это пример того, как уход из отдела до конца вам выгоден. Но опять же, это правило выполняется только до тех пор, пока значения, с которыми вы работаете, не используют полную точность типа данных. После того, как вы начнете работать с неточными (округленными) значениями, вы больше не защитите себя от дальнейших ошибок округления, выполнив сначала умножения.

Ответ 2

Это связано с тем, как реализуется тип double, и тот факт, что типы с плавающей точкой не обеспечивают такую же точность, как другие простые численные типы. Хотя следующий ответ более конкретно касается сумм, он также отвечает на ваш вопрос, объясняя, как нет гарантии бесконечной точности в математических операциях с плавающей запятой: Почему изменение порядка сумм возвращает другой результат?. По сути, вы никогда не должны пытаться определить равенство значений с плавающей точкой без указания допустимого предела погрешности. Библиотека Google Guava включает DoubleMath.fuzzyEquals(double, double, double), чтобы определить равенство двух значений double в пределах определенной точности. Если вы хотите ознакомиться с особенностями равенства с плавающей запятой этот сайт весьма полезен; тот же сайт также объясняет ошибки округления с плавающей запятой. В суммировании: ожидаемые и фактические значения ваших вычислений различаются из-за округления, различающегося между расчетами из-за порядка операций.

Ответ 3

Немного упростите вещи. То, что вы хотите знать, - это то, почему 450d / 420 и 450d / 7 / 60 (в частности) дают разные результаты.

Посмотрим, как выполняется деление в формате с плавающей запятой с двойной точностью IEE. Не углубляясь в детали реализации, это в основном XOR -значение бита знака, вычитание показателя делителя из показателя дивиденда, деление мантис и нормализация результата.

Сначала мы должны представить наши числа в правильном формате для double:

450    is  0 10000000111 1100001000000000000000000000000000000000000000000000

420    is  0 10000000111 1010010000000000000000000000000000000000000000000000

7      is  0 10000000001 1100000000000000000000000000000000000000000000000000

60     is  0 10000000100 1110000000000000000000000000000000000000000000000000

Пусть сначала разделите 450 на 420

Сначала появляется знаковый бит, он 0 (0 xor 0 == 0).

Затем появляется показатель. 10000000111b - 10000000111b + 1023 == 10000000111b - 10000000111b + 01111111111b == 01111111111b

Хорошо выглядящий, теперь мантисса:

1.1100001000000000000000000000000000000000000000000000 / 1.1010010000000000000000000000000000000000000000000000 == 1.1100001 / 1.101001. Есть несколько способов сделать это, я немного поговорю о них позже. Результат 1.0(001) (вы можете проверить его here).

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

0001001001001001001001001001001001001001001001001001 0 0 1

Охрана бит 0, мы не делаем округления. Результат: в двоичном формате:

0 01111111111 0001001001001001001001001001001001001001001001001001

который представляется как 1.0714285714285714 в десятичном формате.

Теперь разделите 450 на 7 по аналогии.

Знак бит = 0

Экспонент = 10000000111b - 10000000001b + 01111111111b == -01111111001b + 01111111111b + 01111111111b == 10000000101b

Мантисса = 1.1100001 / 1.11 == 1.00000(001)

Округление:

0000000100100100100100100100100100100100100100100100 1 0 0

Бит защиты установлен, круглые и липкие бит - нет. Мы округляем до ближайшего (режим по умолчанию для IEEE), и мы застреваем прямо между двумя возможными значениями, которые мы могли бы округлить до. Так как lsb 0, добавим 1. Это дает нам округлую мантису:

0000000100100100100100100100100100100100100100100101

Результат

0 10000000101 0000000100100100100100100100100100100100100100100101

который представляется как 64.28571428571429 в десятичном формате.

Теперь нам нужно разделить его на 60... Но вы уже знаете, что мы потеряли определенную точность. Разделение 450 на 420 вообще не требовало округления, но здесь нам уже приходилось округлять результат хотя бы один раз. Но, для полноты, позвольте закончить работу:

Разделение 64.28571428571429 на 60

Знак бит = 0

Экспонент = 10000000101b - 10000000100b + 01111111111b == 01111111110b

Мантисса = 1.0000000100100100100100100100100100100100100100100101 / 1.111 == 0.10001001001001001001001001001001001001001001001001001100110011

Круглый и сдвиг:

0.1000100100100100100100100100100100100100100100100100 1 1 0 0

1.0001001001001001001001001001001001001001001001001001 1 0 0

Округление, как и в предыдущем случае, получим мантиссу: 0001001001001001001001001001001001001001001001001010.

По мере того, как мы сдвигались на 1, добавляем к экспоненте, получая

Экспонент = 01111111111b

Итак, результат:

0 01111111111 0001001001001001001001001001001001001001001001001010

который представляется как 1.0714285714285716 в десятичном формате.

Тл; д-р:

Первое подразделение дало нам:

0 01111111111 0001001001001001001001001001001001001001001001001001

И последнее подразделение дало нам:

0 01111111111 0001001001001001001001001001001001001001001001001010

Разница заключается только в последних двух битах, но мы могли бы потерять больше - в конце концов, чтобы получить второй результат, нам пришлось округлить два раза вместо них!

Теперь о разделении мантиссы. Разделение с плавающей запятой осуществляется двумя основными способами.

То, что предусмотрено длинным делением IEEE (здесь, является хорошими примерами: это в основном регулярное длинное деление, но с двоичным, а не десятичным), и это довольно медленно. Это то, что сделал ваш компьютер.

Существует также более быстрый, но меньший вариант accrate, умножение на обратное. Во-первых, найдено обратное от делителя, а затем выполняется умножение.

Ответ 4

Это потому, что двойное разделение часто приводит к потере точности. Указанная потеря может варьироваться в зависимости от порядка отделов.

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

Когда вы делите на 7d * 60, вам нужно использовать только одно разделение, тем самым теряя точность только один раз.

Обратите внимание, что двойное умножение может иногда терпеть неудачу, но это гораздо менее распространено.

Ответ 5

Конечно, порядок операций, смешанный с тем, что удваивается, не является точным:

450.00d / (7d * 60) --> a = 7d * 60 --> result = 450.00d / a

против

450.00d / 7d / 60 --> a = 450.00d /7d --> result = a / 60