Соответствие стандарту IEEE-754 от половины до четности

Стандартная библиотека C обеспечивает round, lround, а также llround семейство функций в C99. Однако эти функции не соответствуют стандарту IEEE-754, потому что они не реализуют «округление банкира» от полупериодов до уровня, предписанного IEEE …

Нет смысла говорить о том, является ли отдельная функция «совместимой с IEEE-754». Соответствие IEEE-754 требует, чтобы был доступен набор операций типов данных с определенной семантикой. Он не требует, чтобы эти типы или операции имели конкретные имена, а также не требует, чтобы только эти операции будут доступны. Реализация может предоставлять любые дополнительные функции, которые она хочет, и при этом быть совместимой. Если реализация хочет обеспечить округление до нечетности, округление в случайном порядке, округление от нуля и прерывание, если неточное, она может сделать это.

Что IEEE-754 фактически требует для округления, так это то, что следующие шесть операций являются предоставлена:

convertToIntegerTiesToEven(Икс)

convertToIntegerTowardZero(Икс)

convertToIntegerTowardPositive(Икс)

convertToIntegerTowardNegative(Икс)

convertToIntegerTiesToAway(Икс)

convertToIntegerExact(Икс)

В C и C ++ последние пять из этих операций связаны с trunc, ceil, floor, round, а также rint функции соответственно. C11 и C ++ 14 не имеют привязки для первого, но будущие ревизии будут использовать roundeven, Как вы видете, round на самом деле является одной из обязательных операций.

Тем не мение, roundeven недоступен в текущих реализациях, что приводит нас к следующей части вашего вопроса:

Обычным специальным способом реализации округления в C является выражение (int)(x + 0.5f) которая, несмотря на неправильность в строгой математике IEEE-754, обычно переводится компиляторами в правильную cvtss2si инструкция. Однако это, безусловно, не переносимое предположение.

Проблемы с этим выражением выходят далеко за рамки «строгой математики IEEE-754». Это совершенно неверно для отрицательного x, дает неправильный ответ для nextDown(0.5)и превращает все нечетные целые числа в бинаде 2 ** 23 в четные целые числа. Любой компилятор, который переводит его в cvtss2si ужасно, ужасно сломан. Если у вас есть пример этого, я бы хотел это увидеть.

Как я могу реализовать функцию, которая будет округлять любое значение с плавающей запятой с семантикой половин к четности?

Как njuffa отмеченный в комментарии, вы можете убедиться, что режим округления по умолчанию установлен и используется rint (или же lrint, так как это звучит так, как будто вы действительно хотите получить целочисленный результат), или вы можете реализовать свою собственную функцию округления, вызвав round а затем исправить на полпути, как gnasher729 предлагает. Как только привязки n1778 для C будут приняты, вы сможете использовать roundeven или же fromfp функции для выполнения этой операции без необходимости управления режимом округления.

использование remainder(double x, 1.0) из стандартной библиотеки C. Это не зависит от текущего режима округления.

Остальные функции вычисляют остаток x REM y, требуемый IEC 60559

remainder() здесь полезно, так как отвечает требованиям OP к четным требованиям.

double round_to_nearest_ties_to_even(double x) {
x -= remainder(x, 1.0);
return x;
}

Тестовый код

void rtest(double x) {
double round_half_to_even = round_to_nearest_ties_to_even(x);
printf("x:%25.17le   z:%25.17le \n", x, round_half_to_even);
}

void rtest3(double x) {
rtest(nextafter(x, -1.0/0.0));
rtest(x);
rtest(nextafter(x, +1.0/0.0));
}

int main(void) {
rtest3(-DBL_MAX);
rtest3(-2.0);
rtest3(-1.5);
rtest3(-1.0);
rtest3(-0.5);
rtest(nextafter(-0.0, -DBL_MAX));
rtest(-0.0);
rtest(0.0);
rtest(nextafter(0.0, +DBL_MAX));
rtest3(0.5);
rtest3(1.0);
rtest3(1.5);
rtest3(2.0);
rtest3(DBL_MAX);
rtest3(0.0/0.0);
return 0;
}

Выход

x:                     -inf   z:                     -inf
x:-1.79769313486231571e+308   z:-1.79769313486231571e+308
x:-1.79769313486231551e+308   z:-1.79769313486231551e+308
x: -2.00000000000000044e+00   z: -2.00000000000000000e+00
x: -2.00000000000000000e+00   z: -2.00000000000000000e+00
x: -1.99999999999999978e+00   z: -2.00000000000000000e+00
x: -1.50000000000000022e+00   z: -2.00000000000000000e+00
x: -1.50000000000000000e+00   z: -2.00000000000000000e+00 tie to even
x: -1.49999999999999978e+00   z: -1.00000000000000000e+00
x: -1.00000000000000022e+00   z: -1.00000000000000000e+00
x: -1.00000000000000000e+00   z: -1.00000000000000000e+00
x: -9.99999999999999889e-01   z: -1.00000000000000000e+00
x: -5.00000000000000111e-01   z: -1.00000000000000000e+00
x: -5.00000000000000000e-01   z:  0.00000000000000000e+00 tie to even
x: -4.99999999999999944e-01   z:  0.00000000000000000e+00
x:-4.94065645841246544e-324   z:  0.00000000000000000e+00
x: -0.00000000000000000e+00   z:  0.00000000000000000e+00
x:  0.00000000000000000e+00   z:  0.00000000000000000e+00
x: 4.94065645841246544e-324   z:  0.00000000000000000e+00
x:  4.99999999999999944e-01   z:  0.00000000000000000e+00
x:  5.00000000000000000e-01   z:  0.00000000000000000e+00 tie to even
x:  5.00000000000000111e-01   z:  1.00000000000000000e+00
x:  9.99999999999999889e-01   z:  1.00000000000000000e+00
x:  1.00000000000000000e+00   z:  1.00000000000000000e+00
x:  1.00000000000000022e+00   z:  1.00000000000000000e+00
x:  1.49999999999999978e+00   z:  1.00000000000000000e+00
x:  1.50000000000000000e+00   z:  2.00000000000000000e+00 tie to even
x:  1.50000000000000022e+00   z:  2.00000000000000000e+00
x:  1.99999999999999978e+00   z:  2.00000000000000000e+00
x:  2.00000000000000000e+00   z:  2.00000000000000000e+00
x:  2.00000000000000044e+00   z:  2.00000000000000000e+00
x: 1.79769313486231551e+308   z: 1.79769313486231551e+308
x: 1.79769313486231571e+308   z: 1.79769313486231571e+308
x:                      inf   z:                      inf
x:                      nan   z:                      nan
x:                      nan   z:                      nan
x:                      nan   z:                      nan

float Тип данных может представлять все целые числа, но не дроби, в диапазоне от 8388608.0f до 16777216.0f. любой float числа, которые больше 8388607.5f, являются целыми числами, и округление не требуется. Добавление 8388608.0f к любому неотрицательному float который меньше этого, даст целое число, которое будет округлено в соответствии с текущим режимом округления (обычно с округлением от половины до четности). Вычитание 8388608.0f приведет к получению правильно округленной версии оригинала (при условии, что он находится в подходящем диапазоне).

Таким образом, должно быть возможно сделать что-то вроде:

float round(float f)
{
if (!(f > -8388608.0f && f < 8388608.0f)) // Return true for NaN
return f;
else if (f > 0)
return (float)(f+8388608.0f)-8388608.0f;
else
return (float)(f-8388608.0f)+8388608.0f;
}

и использовать преимущества естественного округления при сложении, не прибегая к какому-либо другому «округлению до целого».

Ниже приводится простая реализация от половины до вечера Программа, которая соответствует стандарту IEEE округления.

Логика:
ошибка = 0,00001

1. число = 2,5
2. темп = этаж (2,5)% 2 = 2% 2 = 0
3. Икс = -1 + темп = -1
4. х * ошибка + номер = 2,40009
5. раунд (2.40009) = 2

Замечания: Ошибка здесь равна 0,00001, т. Е. Если происходит 2,500001, то она округляется до 2 вместо 3

Реализация Python 2.7:

temp = (number)
rounded_number = int( round(-1+ temp%2)*0.00001 + temp )

Реализация C ++: (Использование math.h для напольной функции)

float temp = (number)
int rounded_number = (int)( (-1+ temp%2)*0.00001 + temp + 0.5)

Результат, который это дало бы, был бы следующим согласно. стандартам:

(3.5) -> 4

(2.5) -> 2

Изменить 1: Как отмечает @Mark Dickinson в комментариях. Ошибка может быть изменена в соответствии с требованиями вашего кода для ее стандартизации. Для python, чтобы превратить его в наименьшее возможное значение с плавающей запятой, вы можете сделать следующее.

import sys
error = sys.float_info.min