Определение машинного эпсилон для типов double и float в С++

Question

Поставлена стандартная для новичка задача - определить машинный эпсилон для типов данных float и double на C++. Код для float выглядит так:

    i = 0;
    epsilon_f = 1.0;
    

    while(1.0 + epsilon_f > 1.0)
    {
        epsilon_f = epsilon_f / 2.0;
        i++;
    }
    cout << "Machine epsilon for float type: " << epsilon_f << endl;
    cout << "The amount of iterations for float type: " << i << endl;

Вывод программы следующий:

Machine epsilon for float type: 1.11022e-16

The amount of iterations for float type: 53

Но насколько я понимаю, такое значение машинного эпсилон у типа double (только в два раза больше, так как цикл while сработал лишний раз). В чём ошибка? Как исправить код? Может, неправильная работа с типами данных?

Answer 1

Просто вы работаете с double — так как все литералы именно double по умолчанию.

Сравните вывод в этом случае:

int i = 0;
float epsilon_f = 1.0;


while(1.0f + epsilon_f > 1.0f)
{
    epsilon_f = epsilon_f / 2.0f;
    i++;
}
cout << "Machine epsilon for float type: " << epsilon_f << endl;
cout << "The amount of iterations for float type: " << i << endl;

Здесь мы явно указываем тип float.

См. https://ideone.com/luOC4F

Answer 2

Не совсем ответ, а скорее заметка по вопросу точности.

Объявляя вещественные переменные как float (32 бит) или double (64 бит), многие почему-то забывают (или не все знают), что процессор (сопроцессор) всегда преобразует эти переменные в 80-битный тип long double - это происходит автоматически при их загрузке в 80-битные регистры стека сопроцессора (при выгрузке значений в память, происходит обратное преобразование). И все вычисления производятся над 80-битными вещественными значениями. Реальная же точность операций задаётся определёнными битами управляющего слова в регистре CWR. По умолчанию эта точность соответствует типу double. А значит, даже если исходные значения хранятся во float, процессор всё равно обрабатывает их в регистрах с точностью double. По-этому, чтобы по-настоящему ограничить точность вычислений в 32 бита, необходимо соответствующим образом изменить значение регистра CWR, либо в ассемблере (с помощью инструкций fstcw, fldcw), либо вызвав библиотечную функцию _control87().

Так как же в этом случае отработает следующий код? Ведь вычисления всё равно имеют точность 64 бит!

    int i = 0;
    float epsilon_f = 1.0;

    while(1.0f + epsilon_f > 1.0f)
    {
        epsilon_f = epsilon_f / 2.0f;
        i++;
    }

Я откомпилировал его, и отработал он, как ни странно, корректно:

Machine epsilon for float type: 5.96046e-008
The amount of iterations for float type: 24

В чём же дело? Ответ - в дизассемблере:

.text:004A62E6                 xor     eax, eax        ; int i = 0
.text:004A62E8                 fld1                    ; epsilon_f = 1.0
.text:004A62EA                 lea     esi, [esi+0]
.text:004A62F0
.text:004A62F0 Do_While:                               ; CODE XREF: _main+2A↓j
.text:004A62F0                 inc     eax             ; do {i++
.text:004A62F1                 fmul    ds:const_0_5f   ; epsilon_f = epsilon_f*0.5
.text:004A62F7                 cmp     eax, 24         ; } while (i!=24)
.text:004A62FA                 jnz     short Do_While  ; do {i++

Короче, компилятор оказался настолько умён, что понял, что от него требуется и просто заменил цикл while(1.0f + epsilon_f > 1.0f) {...} на цикл:

    do
    {
        i++;
        epsilon_f = epsilon_f*0.5;
    } while (i!=24);

Вот такие чудеса оптимизации.

БЛОГ НА HUSL

Определение машинного эпсилон для типов double и float в С++

Ответы (2 шт):