Стандартные приёмы программирования

54

Стадарные приемы программирования

Урок 5. Сколько будет 2×2 на ЭВМ, или Машинная арфметика и как с него бороться

Знание ответа на вопрос «Сколько будет 2×2?» является современно необходимым для человека, который претендует на значение математики хотя бы в объеме начальной шкалы. А как ответает на этот вопрос ЭВМ? Неожиданность заключается в том, что на ЭВМ

$$2 \times 2 = 4$$
но $2.0 \times 2.0 \neq 4.0$, точнее $2.0 \times 2.0 \approx 4.0^*$

Дело в том, что целые числа — точные, а действительные — приближенные. Мало того, все числа еще имеют ограниченный диапазон изменения. Какие же выводы из этого должен сделать программист?

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

Действительные числа в ЭВМ представляют в так называемой форме с плавающей точкой, то есть число a представляется в виде $a \cdot n^6$, где $1/n < |a| < 1$, n — основание той системы численения, в которой работает данная ЭВМ (как правило, это двойная или шестнадцатеричная система, что, впрочем, для понимания сутки дела совершенно неважно). Предположим, что наша ЭВМ работает в десятичной системе численения и имеет 4 разряд под mantису и 1 разряд под порядок. Годда (см. рис. 1) в ячейке памяти запоминаются только цифры и знак mantисы и порядка, точка подразумевается перед первой цифрой mantиссы, а порядок показывает в какую сторону и на сколько цифр «уплывает» точка на самом деле (отсюда и термин «плавающая точка»). Таким образом, максимальное представление в нашей ЭВМ число — это $0.9999 \times 10^9$, то есть около миллиарда.

Тепер рассмотрим задачу 2 и предположим, что $n = 1000$, а все $a_i$ имеют значения около $2 \times 10^6$. Тогда попытка получить сумму элементов массива не дает результат, так как эту сумму нельзя записать в памяти ЭВМ. Чтобы избежать этого, формулу среднего арфметического надо переписать в виде

$$m = \sum_{i=1}^{n} \frac{a_i}{n}$$

и программировать деление на $n$ каждого $a_i$. Разумеется, это заменяет программу, но наша более короткая программа, приведенная в таблице, может просто не дать результат.

Тепер при том же $n = 1000$ допустим, что

$$\sum_{i=1}^{900} a_i = 9555, a_{901} = a_{902} = \dots = a_{1000} = 0.2.$$

При сложении $a_{901}$ с накопленной суммой получим $0.95552 \times 10^4$, но ЭВМ отбросит последний знак «2» (теряются разряды,ышедшие за границы разрядной сетки). Таким образом, последнее 100 компонент не внесут свой вклад в результат и получится грубая ошибка (равная $0.2 \times 100 = 20.0$).

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

Можно привести примеры ассоциативных и дистрибутных законов*).

Реальные ЭВМ работыются с граzdо большим диапазоном и большей точностью представления чисел, поэтому в большихстве случае разница между машинной и обычной арфметикой неощущается. Тем не менее известны случаи, когда программа на языке фортран при переноссе с ЭВМ БЭСМ-6 (которая работа с точностью около 12 десятичных знаков на ЭВМ серии EC (точность около 7 десятичных знаков) записывается: низкая точность EC ЭВМ приводила к ошибкам из-за накопления погрешностей.

Как правило, диапазона и точности ЭВМ хватает с избежатьом. Так, однократное сложение чисел с резко различшимся порядками нам не мешает: когда счет идет на миллионы, десяты доли не важны; а для потери точности при многократном сложении чисел близких порядков надо сделать около $10^8$ сложений на EC ЭВМ и около $10^{13}$ на БЭСМ-6. Особенности машинной арфметики начинают сказаться, если нужно точность, близкая к точности ЭВМ, или если объем вычислений подходит к критическому для данной ЭВМ.

*) Для допущенных читателей сообщим, что $2.0 \times 2.0$ все же же в точности равно 4.0, так как и 2.0 и 4.0 представляют в памяти точное. Но если один из компонентов получается вычисленем, то результат уже может оказаться не точным: $(1.0/3.0) \times 6.0 \neq 2.0$, так как промежуткий результат 1/3 не может быть представлен точно на ЭВМ, работающих в двойной или десятичной системе численения.

55

Решение подготовительных задач

Задание 1 (таблица значений $x^2$)

$$x := a;$$
while $x < = b$ do
begin writeln $(x, '**2 =', x**2);$
$$x := x + h$$
end;

Задание 2 (среднее арифметическое).

$$s := 0;$$
for $i := 1$ to $n$ do $s := s + a[i]$:
$$s := s/n$$;

Задание 3 (попск данного числа).

$$k := 0; i := 1;$$
while $(i < = n)$ AND $(k = 0)$ do
begin if $a[i] = x$ then $k := i;$
$$i := i + 1$$
end;

В упомянутом примере с программой на фортеран нужна была точность $10^{-5}$ при точности ЭВМ равной $10^{-7}$, а точность $10^{-3}$ программы на EC ЭВМ обеспечивала.

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

a) Доказать, что для данного алгоритма погрешности машиной арифметики не выходят за допустимые пределы.

b) Доказать то же при дополнительных условиях на исходные данные (установить, что при ваших данных все будет в порядке, хотя этого гарантировать в общем случае нельзя).

c) Разработать алгоритм, компенсирующий погрешности.

e) Использовать вычисления с увеличенной точностью (реализованные аппаратно или программно), в частности, просчитать программу на другой модели ЭВМ с более точной арифметикой.

На практике чаще всего используется варианты e) или e). В некоторых языках программирования (фортран и ПЛ/1, например) есть возможность указать, что вычисления надо вести с увеличенной точностью. Это на EC ЭВМ реализуется аппаратными средствами (с неокоторым замедленным программным), а на БЭСМ-6 — программно (so значительным замедленным программным).

Действительные числа лучше не сравнивать на равенство,

так как с почти стоппроцетной вероятностью получим результат «не равно». Поэтому оператор в задание 3

if $a[i] = x$ then…
надо заменить на
if abs $(a[i] - x) < eps$ then…

где eps — константа, описанная, например, как eps = 0.0001 (в зависимости от точности ЭВМ и нужной вам точности). Тогда действие числа считаются равными, если они различаются не более чем на некоторое ε.

В задание 1, если $a = 0, b = 1, a h = 0,1$, то у нас нет гарантии, что мы получим все 11 значений функций: после 10-кратного сложения $x$ с $h$ может получить значение чуть меньше $b$, или чуть больше $b$, и в последнем случае для него уже цикл работает не будет. Для достижения нужного результата нужно либо организовать цикл с условиям while $(x < b + h/2)$, либо заранее вычислить число повторений цикла $n := \text{round}((b - a)h)$ и организовать цикл for $i := o$ to $n$ do… Однако и этот путь может не привести к желаемому результату: если значение $n$ (необходимое число шагов) получится достаточно большим, то при многократном выполненном оператора $x := x + h$ может накопиться погрешность, превосходящая значение $h$. Тогда цикл, организованный через while, выполнится на один (или даже несколько) раз меньше или больше, а цикл, организованный через for, выполнится нужное количество раз, но последнее значение будет весы далеко от $b$. Для компенсации этого нам придется организовать цикл следующим образом:

$$n := \text{round}((b - a)/h);$$
for $i := 0$ to $n$ do
begin $x := a + i*h;$
writeln $(x, '**2 =', x**2)$
end;

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

На трех простых примерах мы рассмотрели выполнение варианта e) — разработку специального алгоритма. Обратите внимание, что варианты e) и e) требуют выполнения в качестве составной части пунктов $a$ и $b$, чтобы убедиться, что разработанный алгоритм дает нужный эффект. Иного это просто. Так, например, оченьindo, что

$$\sum_{i=1}^{k} \frac{a_i}{n} < \max\{a_i\} (k < n),$$

позному при считете по формуле $\sum \frac{a_i}{n}$ переполнение нам нерозит — доказан варианта $a$. Но если $n \cdot \max\{a_i\}$ меньше максимального представимого числа, то можно работать по формуле $\frac{1}{n} \sum_{i=1}^{k} a_i$ — это варианта e):

при наших данных работает можно, хотя, вообще говоря, алгоритм тант опасность. Доказательства здесь были выполнены только экзизно, более подробно их предлагается выполнить самостояльно (для этих простых случаяев вам это доступно). В других случаях нужные доказательства $a$ или $b$ могут оказаться сложными.

КонтROLное задание

5.1. Дана ЭВМ, работающая в $p$-чной системе численения с длиной мантиссы в $n$ разрядов. Для задания 1, разработанной выше, определяте, при каких соответствии между $a, b, h, p$ и $n$ можно строить таблицу циклом while с последовательным увеличением $x$ на величины $n$ с тем, чтобы гарантировать правильное повторение цикла.

(Подготовительные задач к следующему (последнему) уроку не дается!)
Л. Ф. Штернберг