Дж. Деммель - Вычислительная линейная алгебра, страница 4

Примером подобного подхода может служить книга «Шаблоны для решения линейных систем. Строительные блоки для итерационных методов» ]24]. Аналогичный набор шаблонов разрабатывается в настоящее время для задач на собственные значения. 1.4. Пример: вычисление многочлена Мы хотим проиллюстрировать понятия теории возмущений, числа обусловленности, обратной устойчивости и анализа погрешностей округлений на примере вычисления многочлена р(х) = ~~~ агх'.

Напомним правило Горнера для вычисления многочлена: р = а« Еог» = 4 — 1 4оггп Го 0 р=х*р+а; ен4 1ог Применим зту процедуру к многочлену р (х) = (х — 2) г = хг — 18ха+ 144х~— 672х +2016ха — 4032х +5376хг — 4608хг+2304х — 512. Рассматривая нижнюю часть рисунка 1.1, мы видим, что вблизи точки х = 2 значение р(х), вычисляемое по правилу Горнера, ведет себя совершенно непредсказуемым образом и по справедливости может быть названо «шумомы Действительный график р (х) показан в верхней части рис.

1.1. 16 Глава 1. Введение Попробуем разобраться, какие выводы для себя мы должны сделать из рисунка. С этой целью посмотрим, что произойдет, если попытаться найти нуль многочлена р(х) с помощью простой программы вычисления нулей, основанной на методе бисекции (деления пополам); эта программа представлена ниже алгоритмом 1.1. Бисекция применяется к начальному интервалу [хь,„,хыдь], на котором р(х) изменяет знак (т.е. Р(хьо ) р(хыдь) < 0); тем самым, многочлен р(х) должен иметь нуль на этом интервале. Алгоритм вычисляет значение р(х ы) в средней точке интервала х ы = (хьо + хыдь)/2, а затем выясняет, в какой из половин исходного интервала — левой [хь, х„,ы] или правой [хтптз, хыдь],— р (х) изменяет знак.

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

И действительно, в зависимости от выбора начальных значений хь, и хыдь, алгоритм может сойтись к любой тпочке »области шума», т.е. отрезка от 1.95 до 2.05 (см. вопрос 1.21). Чтобы получить полное объяснение этого примера, мы должны обратиться к свойствам арифметики с плавающей точкой. Алгоритм 1.1. Вычисление нулей функции р(х) методом бисекции. ргос Ь1весС (р хыктхыдь Со1) /" найти корень уравтьения р(х) = 0 на отрезке [хьоот, хььдь] в предполозтсении, ьто р(хь ) р(хы ь) < 0 */ /" оспьанов, если корень найден с тпочностью до х Со1 */ РЬото = Р(ХЬопт) ры,ь =р(хы ь) юЫе хы ь — хь > 2 Со1 Хп»Ы = (Х!оы + Хмде)/2 Ртпи = Р(ХпьЫ) ь/ры, рптд < 0 11ьеп /» имеется корень на [хьо „х„„.,ь]»/ Хдндь — Хт»Ы Рыдь = Рпты е1де ь/р»в рыдь < 0 Йеп /" имеепься корень на [х ы,хььдь] "/ ХЬопт — Хт»Ы Рьом = Ртпы е1де /" х„„н яв яетпся корнем»/ хьопт — хптЫ химь =х ы епд ь/ епд юЫе гоо1 = (хьо + хыдь)/2 17 1.4.

Пример: вычисление многочлена х10ле 1.5 0.5 -0.5 -1.5 ! .92 1.94 1.96 1.98 2 2.02 2.04 2.06 2.08 х10 0.8 -1.6 1.92 1.94 1.96 1.98 2 2.02 2.04 2.06 2.08 Рис. 1.1. График многочлена у = (х — 2) = х — 18хв + 144х~ — 672х + 2016х 4032х + 5376х — 4608х + 2304х — 512. Значения многочлена вычислены в 8000 равноотстоящих точек посредством формулы у = (х — 2) (верхняя часть рисунка) и 9 правила Горнера (нижняя часть).

18 Глава 1. Введение 1.5. Арифметика с плавающей точкой Число — 3.1416 можно записать в математической нотации как — 0.31416 Х 10 1 г~ ~ ~ показатель мантисса основание знак Схожее представление, называемое представлением с плавающей точкой, используется в компьютерах. Правда, основанием обычно является число 2 (с такими исключениями, как 16 для машины 1ВМ 370 и 10 для некоторых решающих таблиц и большинства калькуляторов). Например, 0.10101з х 2з = 5 25~о. Число с плавающей точкой называется нормализованным, если старший разряд его мантиссы отличен от нуля. К примеру, число 0.10101з х 2з нормализовано, а число 0.010101з х 2~ нет.

Обычно числа с плавающей точкой нормализуют, что дает выигрыш в двух отношениях: всякое ненулевое число с плавающей точкой в этом случае имеет единственное представление в виде строки битов, и старший разряд двоичной мантиссы можно не хранить в явном виде (поскольку он всегда равен 1); за счет сэкономленного бита можно удлинить мантиссу. Самыми важными среди параметров, описывающих числа с плавающей точкой, являются основание, число разрядов (битов) мантиссы, определяющее точность представления, и число разрядов (битов) покыателя, определяющее область изменения показателей и, тем самым, наибольшее и наименьшее из представимых чисел. Различные арифметики с плавающей точкой разнятся между собой также способом округления вычисленных результатов, решениями, которые принимаются в отношении чисел, слишком близких к нулю (машинные нули) или слишком больших (переполнения), наличием или отсутствием символов ~ оо и некоторых полезных нечисел.

Нечисло, называемое еще неопределенной величиной или специальным операндом, иногда обозначается Ха1ч. Все эти понятия обсуждаются ниже. Рассмотрим прежде всего вопрос о точности представления чисел. К примеру, число 0.31416 х 10~ имеет пять десятичных разрядов, поэтому в нем могла быть потеряна информация, меньшая, чем 0.5 х 10 4. Это означает, что если х — вещественное число, для которого наилучшим пятиразрядным представлением является 0.31416 х 10, то относительная ошибка представления для 0.31416 х 10~ может быть оценена как )х — 0.31416 х 10~! 0.5 х 10 4 4 0.31416 х 10' 0.31416 х 10' Максимум относительной ошибки представления, возможной для нормализованного числа, реализуе..ся в числе 0.10000 х 10', представляющем собой наиболее точное пятиразрядное представление для всех чисел в интервале от 0.999995 до 1.00005.

Относительная ошибка для этого числа оценивается величиной 0.5 10 4. Более общо, для арифметики с плавающей точкой, имеющей р-разрядную мантиссу и основание )1, максимальнав относительная ошибка представления равна 0.5 х ф г. Эта величина, к тому же, равна половине 19 1.5. Арифметика с плаватотей точкой расстояния между 1 и ближайшим ббльшим числом с плавающей точкой, т.е. числом 1+)3' ". На протяжении истории компьютеров использовались многие различные комбинации основания, длины мантиссы и области показателей. К счастью, в настоящее время почти общепринятым является 1ЕЕЕ-стандарт двоичпой арифметики.

Он реализован рабочими станциями 8пп, ПЕС, НР и 1ВМ, а также всеми персональными компьютерами. 1ЕЕЕ-арифметика предусматривает два типа чисел с плавающей точкой: числа обычной точности (с 32-битовым представлением) и числа двойной точности (64 бита). ! ЕЕЕ число обычной точности двоичная точка Пусть в представлении 1ЕЕЕ-числа обычной точности в, е и 1' ( 1 суть, соответственно, 1-битовый знак, 8-битовый показатель и 23-битовая мантисса; тогда само число равно ( — 1)' 2' 'гт (1+ 1).

Максимальная относительная ошибка представления равна 2 г4 яг 6 10 г, а область положительных нормализованных чисел простирается от 2 тгг (,порог машинного нуля) до 2~~~ (2 — 2 гг) 2тгг (порог переполнения) или, приблизительно, от 10 гг до 10гг. Расположение этих чисел с плавающей точкой на вещественной числовой прямой показано на рис. 1.2 (где, в целях иллюстрации, для мантисс принято 3-битовое представление). ! ЕЕЕ число двойной точности двоичная точка Пусть в представлении 1ЕЕЕ-числа двойной точности г, е и 1 < 1 суть, соответственно, 1-битовый знак, 11-битовый показатель и 52-битовая мантисса; тогда само число равно ( — 1)' 2' твгг .

(1+ 1"). Максимальная относительная ошибка представления равна 2 гг — 10 тч, а границами области показателей являются 2-твгг ~парве натаиннвгв нуля) и 2тогг . (2 2 — гг) 2твг4 (порог переполнения), т.е., приблизительно, 10 гвг и 10гег. Пусть символ О обозначает одну из четырех бинарных операций +, —, г и /. Если точный результат вычисления а О Ь не может быть представлен как число с плавающей точкой, то, прежде чем записать его в память или регистр, нужно аппроксимировать его каким-либо близко расположенным числом с плавающей точкой. Это приближение будем обозначать через Я(а О Ь).

Разность (а О Ь) — Ца О Ь) называется погрешностью округления. Если Я(а О Ь) всегда является ближайшим числом с плавающей точкой к числу а О Ь, то будем говорить, что арифметика округляет правильно (или, просто, округляет). 1ЕЕЕ- арифметика обладает этим привлекательным свойством. (Если число аОЬ находится точно посредине между двумя соседними числами с плавающей точкой, то из двух возможных значений для Я(а О Ь) 1ЕЕЕ-арифметика выбирает число с нулевым последним разрядом мантиссы; такой выбор называется округлением до блилсайшего четттого.) Предположим, что используемая арифметика округляет правильно и число а О Ь не выходит за пределы области допусти- 20 Глава 1. Введение мых показателей (в противном случае, мы получили бы машинный нуль или переполнение).

Тогда можно записать Я(а О Ь) как Я(а О Ь) = (а О Ь)(1+ б), (1 1) где (б) не превышает числа г, называемого машинным эпсилон или машинной точностью и обозначаемого тасЬерэ. Поскольку округление производится с наивысшей возможной точностью, то г равно максимальной относительной ошибке представления 0.5 Д' г. 1ЕЕЕ-арифметика обеспечивает также, что Я( „/а) = ч/а(1+ б), где ~б~ < г. Эти соотношения определяют наиболее распространенную модель анализа погрешностей округлений, и именно такая модель используется в данной книге. Почти аналогичные соотношения справедливы для комплексной арифметики с плавающей точкой (см. вопрос 1.12). Отметим, однако, что формула (1.1) игнорирует некоторые заслуживающие интереса детали.

1.5.1. 0 некоторых важных деталях 1ЕЕЕ-арифметика включает в себя также субнормальные числа, т. е. ненормализованные числа с плавающей точкой, имеющие наименьший возможный показатель. Субнормальные числа суть очень малые числа, находящиеся между нулем и наименьшим нормализованным числом с плавающей точкой (см. рис. 1.2). Наличие таких чисел в арифметике означает, что разность Я(х — у) никогда не может быть машинным нулем. Результатом является замечательное свойство, состоящее в том, что равенство Я(х — у) = 0 возможно тогда и только тогда, когда х = у. Если бы это свойство отсутствовало, то формулу (1.1) пришлось бы изменить так, чтобы она учитывала эффект машинных нулей. Модифицированная формула имела бы вид Я(а О Ь) = (а О Ь) (1 + б) + и, где по-прежнему ~б( < г, а величина (у! ограничена очень малым числом, равным наибольшей ошибке, которую может повлечь за собой появление машинного нуля (для 1ЕЕЕ-арифметикн обычной точности это 2 шо и 10 4э, а для 1ЕЕЕ-арифметики двойной точности — 2 штэ и 10 гг4).