Определяем машинный ноль, машинную бесконечность и машинный эпсилон

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

Поэтому повторим в этой заметке несколько азбучных истин о представлении вещественных чисел в компьютере и правилах выполнения операций с ними.

В IBM-совместимой ЭВМ для вещественных чисел используется двоичная система счисления и принята форма представления чисел с плавающей точкой вида x = m*2p, где мантисса m = ± (g1*2-1 + g2*2-2 + ... + gt*2-t), g1, ..., gt - двоичные цифры, причём, g1=1, а целое значение p называется двоичным порядком. Количество цифр t, которое отводится для записи мантиссы, называется разрядностью мантиссы. Диапазон представления чисел в ЭВМ ограничен конечной разрядностью мантиссы и значением числа p.

Все представимые на ЭВМ вещественные числа x удовлетворяют неравенствам 0 < X0 ≤ |x| < X∞, где X0 = 2-pmax+1, X∞ = 2pmax, а значение pmax соответствует разрядности вычислительной системы.

Все числа, по модулю большие X∞, не представимы на ЭВМ и рассматриваются как машинная бесконечность. Все числа, по модулю меньшие X0, для компьютера не отличаются от нуля и рассматриваются как машинный ноль. Машинным эпсилон εM называется относительная точность ЭВМ, то есть граница относительной погрешности представления вещественных чисел. Можно показать, что εM ≈ 2-t. Пусть x* = m*2p. Тогда граница абсолютной погрешности представления этого числа равна Δ(x*) ≈ 2-t-1*2p. Поскольку 1/2≤m<1, то величина относительной погрешности представления оценивается как δ(x*) ≈ Δ(x*) / |x*| ≈ (2-t-1*2p) / (m*2p) = 2-t-1 / m ≤ 2-t-1 / 2-1 = 2-t.

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

Примем следующие способы определения приближённых значений искомых величин:

• положим X∞ = 2n, где n - первое натуральное число, при котором произошло переполнение;
• положим X0 = 2-m, где m – первое натуральное число , при котором 2-m совпадает с нулем;
• положим εM = 2-k, где k – наибольшее натуральное число, при котором сумма вычисленного значения 1+2-k ещё больше 1. Фактически, εM есть граница относительной погрешности представления числа x* ≈ 1.

Дальше задаём это в нужной среде (пакете) и подбираем значения параметров, вот пример для моего MathCAD 15:

машинный ноль, машинная бесконечность и машинный эпсилон в MathCAD 15

А вот что вышло в Visual Studio 2010 при использовании проекта Windows Forms, C++/CLI, библиотеки System::Math и типа данных long double:

Inf: 1024
Zero: 1075
Eps: 53

Код:

//Функции для подсчёта
long double inf (int n) { return Math::Pow(2.,n); }
long double zero (int m) { return Math::Pow(2.,-m); }
long double eps (int k) { return 1.+Math::Pow(2.,-k); }

//...
//Расчёт, сделанный по нажатию кнопки с выводом результатов в метку label1
 label1->Text = "";

 int n=1,m=1,k=1;
 long double res;
 while (1) {
  res=inf(n);
  if (res==Double::PositiveInfinity) break;
  else n++;
 };
 label1->Text +=  "Inf: " + n + Environment::NewLine;

 while (1) {
  res=zero(m);
  if (res==0.) break;
  else m++;
 };
 label1->Text +=  "Zero: " + m + Environment::NewLine;

 while (1) {
  res=eps(k);
  if (res==1.) break;
  else k++;
 };
 label1->Text +=  "Eps: " + k + Environment::NewLine;

Ну и пара стандартных напоминаний напоследок:

К вещественным значениям в общем случае неприменима операция == ("сравнение") из-за неточного представления этих значений в памяти компьютера. Поэтому для вещественных переменных отношение вида a==b обычно заменяется на fabs(a-b)≤eps, где fabs() - функция вычисления модуля вещественного числа, а eps - малая величина, определяющая допустимую погрешность.

Допустимую погрешность можно ввести в расчёт также через стандартный метод округления round, например, левый расчёт произведения чисел в MathCAD не даст нуля, а правый - да:

учёт погрешностей через метод round (Mathcad)

27.10.2015, 17:32 [29113 просмотров]