История страницы

Срок выполнения: 13–20 ноября

Задание

1. Составьте алгоритм вычисления формулы Тейлора для своего варианта.
2. Напишите функцию вычисления "машинного эпсилон".
3. Напишите функцию, которая эффективно вычисляет ряд, заданный вариантом, с точностью машинного эпсилон (оно должно вычисляться с помощью функции из пункта 2).
4. Напишите программу, которая вычисляет заданную функцию дважды: напрямую по формуле из варианта и с помощью вашей подпрограммы из пункта 3. Выведите полученные значения и разность между ними.

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

Математические функции (<cmath>): http://en.cppreference.com/w/cpp/numeric/math
(обращайте внимание на пометки:
since C++11 (или C++14) – средства из нового стандарта, пока присутствуют не во всех реализациях библиотеки;
until C++11 – средства из старого стандарта, которые изменены или заменены в новом)

Замечания о точности

"Машинное эпсилон" несёт в себе информацию о точности вычислений в конкретном формате с плавающей запятой. Но его вычисление (как и другие подобные расчёты) зависит от машинного кода, создаваемого компилятором. Так на Intel (и совместимых машинах) регистры арифметического устройства для вещественных чисел хранят по 80 битов (всего), но модуль умеет дополнительно сохранять результат в память в 32- (float) и 64-разрядном (double) форматах (и, конечно, загружать такие числа). Компилятор (пример – 32-битный GCC) может создать код, который будет считать всё на регистрах (включая сравнения). В итоге мы получим "эпсилон" для 80-битового формата. Сохранив его в float или double, мы лишь приведём его к ближайшему представимому в запрошенном формате числу (обрежем значимую часть). Полученное значение будет неверным.

float eps(float startVal) {
    float sv = startVal;
    while (true) {
       float val = sv * 0.5f;
       if (startVal + val == startVal) return sv;
       sv = val;
   }
}

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

float pass(float x) { return x; }
float eps(float startVal) {
    float sv = startVal;
    while (true) {
       float val = sv * 0.5f;
       if (pass(startVal + val) == startVal) return sv;
       sv = val;
   }
}

Отметим, что при оптимизации компиляторы могут встраивать тело функции в вызывающий код вместо обычного вызова, тогда от функции pass не останется и следа. Более надёжный способ – поместить функцию pass в другой исходный файл.

Что интересно, 64-разрядная версия GCC даст верный результат и в первом варианте, потому что она создаст код не для вещественного сопроцессора, а для векторного модуля (SSE). SSE работает только с классическими типами float и double, поэтому результаты будут верно обрезаться в любом случае.

Отсюда мораль:

1. Для подстраховки значение нужно обрезать в любом случае – это более универсальный способ.
2. Вычисления желательно вести в одном типе данных. Это относится как к типам переменных, так и констант (!!!).
3. В настоящих программах лучше использовать библиотечные средства. Например: std::numeric_limits<float>::epsilon().

"Эпсилон" даёт нам двоичную точность (2^-x). Нередко хочется выяснить десятичную, например для вывода. Выяснить, скольким знакам мы можем доверять. Для этого мы должны узнать, какой десятичной степени соответствует это число.

Для этого нужно вычислить выражение: n = -log₁₀(eps).

Значение будет дробным. Чтобы гарантированно получить все необходимые знаки, обычно его округляют в большую сторону – ceil(n). При этом нужно помнить, что достоверным является количество знаков, на единицу меньшее – floor(n).

О типах данных

Обычно для более быстрых вычислений выбирают меньший тип (например, float), для более точных – больший (например, double). Обычно можно работать с float, а в критичных местах использовать double для промежуточных данных расчёта.

Дальнейшее зависит от машины и компилятора. Как сказано выше, для Intel компилятор GCC выбирает FPU для 32-разрядных программ, SSE – для 64-разрядных.

Если вычисления идут на FPU, надо помнить, что при работе с регистрами выбор типа теряет остроту – все вычисления ведутся в 80-битном формате. Число обрезается только при сохранении в память. Обычно это происходит при передаче параметров в функцию, при записи не в локальные переменные (например, в глобальные) и при нехватке регистров.

SSE работает с float. Начиная со второй версии (SSE2), и с double. При этом размер данных и скорость работы для двух типов действительно различаются. SSE быстрее и удобнее, чем FPU (а также может работать с векторами), но выполняет только арифметические и логические операции, сравнения и вычисляет квадратный корень. Прочие встроенные функции (например, синус, арктангенс, логарифм) и тип long double есть только в FPU.

В настройках компилятора Си/Си++ обычно можно выбрать, какой модуль использовать. Приведём ключи для GCC.
Выбор FPU: -mfpmath=387
Выбор SSE: -msse -msse2 -mfpmath=sse

Расширение SSE появилось в процессоре Pentium III, а SSE2 – в Pentium 4.

Примеры

float.cpp

Варианты

Расчёт номера своего варианта см. на странице Успеваемость студентов.

float.pdf