[P&AM Lab] Катастрофическая потеря значимых разрядов -- пример с детальным выводом

[Grigoriy A. Sitkarev]

Приветствую всех!

В прошедшую субботу мы рассматривали один пример с вычислением формулы

f(x) = (1 - cos(x)) / x^2,

для случая, когда x близок к нулю.

Я сделал один пример, который даёт хорошее представление о том, что 
происходило с вычислением и откуда взялась ошибка. Достаточно 
скомпилировать исходник и запустить. Немного по выводу:

45) y = 9.84254003e-01
     x = 2.46085750e-04
  base 10:
   1.0    = 1.000000000
   cos(x) = 0.999999940
  base 2:
   1.0    = 1.00000000000000000000000x2^0
   cos(x) = 1.11111111111111111111111x2^-1

На индексе 45 мы видим, что для x (не самого маленького, всего 
приблизительно 0,000246) при вычислении по формуле «в лоб» получаем 
результат ~ 0,984, в то время как известен правильный ответ ~ 0,5.

relerr(y) = |y - y'| / |y| = |0.5 - 0.984| / |0.5| =~ 0,968 ,

т.е. ошибка приблизительно равна 100%.

Промежуточные числа в выводе представлены в десятичной и двоичной форме; 
расшифровка соответственно base 10 и base 2. Обратите внимание, что для 
base 2 там даётся значение мантиссы (с напечатанным скрытым битом) и 
экспоненты (мантисса умножается на степень двойки). Последние вытащены 
прямо из формата IEEE single precision.

Видно, что аппроксимация для cos(x) в мантиссе содержит только единицы. 
Вычитая из (1.0)_2 число (0.111111111111111111111111)_2 мы получаем 
ответ, где вообще не содержится корректных разрядов, а только ошибка от 
округления аппроксимации.

Здесь в примере удалось довольно плавно подойти к нарастанию ошибки, 
можно хорошо рассмотреть, что происходило по мере приближения 
аппроксимации cos(x) к единице.

Когда аппроксимация cos(x) (для x = 2.05071454e-04) становится равной 
единице, то алгоритм вообще выдаёт полную ахинею. Результат вычисления 
-- ноль.

Думаю, всем должно быть интересно посмотреть на эти «внутренности».

--
Г.А.
----------- следущая часть -----------
A non-text attachment was scrubbed...
Name: test.c
Type: text/x-csrc
Size: 1035 bytes
Desc: отсутствует
URL: <http://amplab.syktsu.ru/pipermail/lab/attachments/20120910/6c6b6298/attachment.c>