Различные значения в зависимости от установленных флагов исключения с плавающей запятой

Короткий вопрос:

Как установка флага исключения _EM_INVALID на FPU может привести к другим значениям?

Длинный вопрос:

В нашем проекте мы отключили исключения с плавающей запятой в нашей сборке Release, но включили ZERODIVIDE, INVALID и OVERFLOW, используя _controlfp_s () в нашей сборке Debug. Это для того, чтобы выявлять ошибки, если они есть.

Однако мы также хотели бы, чтобы результаты численных расчетов (включая алгоритмы оптимизации, инверсию матриц, Монте-Карло и все такое) были согласованы между сборкой отладки и выпуска, чтобы упростить отладку.

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

Это распространяется на разницу в 0,5% в результирующем значении.

Приведенный ниже код даст показанный вывод программы при добавлении ее в новое решение в Visual Studio 2005, Windows XP и компиляции в конфигурации отладки. (Release даст другой результат, но это потому, что оптимизатор повторно использует результат первого вызова log ().)

Я надеюсь, что кто-нибудь сможет пролить свет на это. Спасибо.

/*
Program output:

Xi, 3893f76f, 7.4555176582633598
K,  c0a682c7, 7.44466687218

Untouched
x,  da8caea1, 0.0014564635732296288

Invalid exception on
x,  da8caea2, 0.001456463573229629

Invalid exception off
x,  da8caea1, 0.0014564635732296288
*/

#include <float.h>
#include <math.h>
#include <limits>
#include <iostream>
#include <iomanip>
using namespace std;

int main()
{
    unsigned uMaskOld  = 0;
    errno_t err;

    cout << std::setprecision (numeric_limits<double>::digits10 + 2);

    double Xi = 7.4555176582633598;
    double K  = 7.44466687218;
    double x;

    cout << "Xi, " << hex << setw(8) << setfill('0') << *(unsigned*)(&Xi) << ", " << dec << Xi << endl; 
    cout << "K,  " << hex << setw(8) << setfill('0') << *(unsigned*)(&K) << ", " << dec << K << endl; 
    cout << endl;

    cout << "Untouched" << endl;
    x = log(Xi/K);
    cout << "x,  " << hex << setw(8) << setfill('0') << *(unsigned*)(&x) << ", " << dec << x << endl; 
    cout << endl;

    cout << "Invalid exception on" << endl;

    ::_clearfp();
    err = ::_controlfp_s(&uMaskOld, 0, _EM_INVALID);

    x = log(Xi/K);
    cout << "x,  " << hex << setw(8) << setfill('0') << *(unsigned*)(&x) << ", " << dec << x << endl; 
    cout << endl;

    cout << "Invalid exception off" << endl;

    ::_clearfp();
    err = ::_controlfp_s(&uMaskOld, _EM_INVALID, _EM_INVALID);

    x = log(Xi/K);
    cout << "x,  " << hex << setw(8) << setfill('0') << *(unsigned*)(&x) << ", " << dec << x << endl; 
    cout << endl;

    return 0;
}

person nielses    schedule 07.08.2013    source источник
comment
спецификация IEEE делает, или, по крайней мере, последняя из тех, что я видел, диктует разные результаты в зависимости от того, включены исключения или нет. но это было для таких вещей, как деление на ноль, переполнение и тому подобное. Точно так же, если у вас уже есть сигнальный nan и вы используете его в качестве входа, вы можете получить тихий nan как выход следующей операции. флаги исключения также могут определять, какой тип наночастиц вы получите, но я просто предполагаю, что у меня больше нет доступа к спецификации.   -  person old_timer    schedule 07.08.2013
comment
Да, я мог понять, что могу получить другие результаты, если сделаю что-то, что вызовет недопустимое исключение. Как операции с сигнальными значениями NaN, бесконечность - бесконечность, 0,0 * бесконечность, 0,0 / 0,0, квадратный корень из отрицательного числа и т.д. Мои цифры относительно благоприятны.   -  person nielses    schedule 07.08.2013
comment
Возможно ли, что при вычислении журнала используются разные алгоритмы в зависимости от настроек флага с плавающей запятой? Он может использовать более быстрый алгоритм, если он может иметь бесконечное количество промежуточных результатов, не вызывая исключения.   -  person Patricia Shanahan    schedule 07.08.2013
comment
Да, это возможно. Я надеюсь, что кто-нибудь, у кого есть документация о внутренней работе IEEE-754 и / или x87, сможет объяснить, является ли это причиной.   -  person nielses    schedule 08.08.2013


Ответы (1)


Это не полный ответ, но для комментария он слишком длинный.

Я предлагаю вам изолировать код, который выполняет сомнительные вычисления, и поместить его в подпрограмму, предпочтительно в исходный модуль, который компилируется отдельно. Что-то типа:

void foo(void)
{
    double Xi = 7.4555176582633598;
    double K  = 7.44466687218;
    double x;
    x = log(Xi/K);
    …Insert output statements here…
}

Затем вы вызываете процедуру с другими настройками:

cout << "Untouched:\n";
foo();

cout << "Invalid exception on:\n";
…Change FP state…
foo();

Это гарантирует, что одни и те же инструкции выполняются в каждом случае, исключая возможность того, что компилятор по какой-то причине сгенерировал отдельный код для каждой последовательности. То, как вы скомпилировали код, я подозреваю, возможно, компилятор мог использовать 80-битную арифметику в одном случае и 64-битную арифметику в другом, или, возможно, использовал 80-битную арифметику в целом, но преобразовал некоторый результат в 64-битный в одном случае, но не в другом

Как только это будет сделано, вы можете дополнительно разделить и изолировать код. Например, попробуйте оценить Xi/K один раз перед любым из тестов, сохранить это в double и передать его foo в качестве параметра. Проверяет, отличается ли вызов log в зависимости от состояния с плавающей запятой. Я подозреваю, что это так, потому что вряд ли операция разделения будет отличаться.

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

Случайные примечания:

Если вы знаете, что Xi и K близки друг к другу, лучше вычислить log(Xi/K) как log1p((Xi-K)/K). Когда Xi и K близки друг к другу, вычитание Xi-K является точным (без ошибок), а в частном больше полезных битов (1, о которой мы уже знали, и некоторые нулевые биты, следующие за ним, пропали).

Тот факт, что небольшие изменения в вашей среде с плавающей запятой вызывают изменение вашего результата на 0,5%, означает, что ваши вычисления очень чувствительны к ошибкам. Это говорит о том, что даже если вы сделаете свои результаты воспроизводимыми, ошибки, которые обязательно существуют в арифметике с плавающей запятой, приведут к тому, что ваш результат будет неточным. То есть последняя ошибка все равно будет существовать, просто она не будет обращена на ваше внимание разницей между двумя разными способами расчета.

В вашей реализации на C ++ видно, что unsigned - это четыре байта, а double - восемь байтов. Таким образом, при печати кодировки double путем присвоения ей псевдонима unsigned половина битов опускается. Вместо этого вы должны преобразовать указатель на double в указатель на const char и вывести sizeof(double) байт.

person Eric Postpischil    schedule 07.08.2013