Поведение пользовательских плаваний NaN в Python и Numpy

Мне нужно упаковать некоторую дополнительную информацию в значения NaN с плавающей запятой. Я использую float IEEE 754 с одной точностью (32-битные поплавки) в Python. Как Python и NumPy обрабатывают эти значения?

Теория

Стандарт IEEE 754-2008, кажется, считает, что число действительно не является числом, если установлены биты экспоненты (23..30), и по крайней мере один из значащих бит установлен. Таким образом, если мы преобразуем float в 32-битное целочисленное представление, то все, что удовлетворяет следующим условиям:

  • i & 0x7f800000 == 0x7f800000
  • i & 0x007fffff != 0

Это оставило бы меня с большим выбором. Тем не менее, стандарт, кажется, говорит, что старший бит значимости is_quiet и должен быть установлен, чтобы избежать исключений в вычислениях.

Практические тесты

Python 2.7

Чтобы убедиться, я провел несколько тестов с интересными результатами:

import math
import struct

std_nan = struct.unpack("f4", struct.pack("I", 0x7fc00000))[0]
spec_nan = struct.unpack("f4", struct.pack("I", 0x7f800001))[0]
spec2_nan = struct.unpack("f4", struct.pack("I", 0x7fc00001))[0]

print "{:08x}".format(struct.unpack("I", struct.pack("f4", std_nan))[0])
print "{:08x}".format(struct.unpack("I", struct.pack("f4", spec_nan))[0])
print "{:08x}".format(struct.unpack("I", struct.pack("f4", spec2_nan))[0])

Это дает:

7fc00000
7fc00001 <<< should be 7f800001
7fc00001

Это и некоторые дополнительные тесты, похоже, предполагают, что что-то (struct.unpack?) всегда устанавливает бит is_quiet.

NumPy

Я пробовал то же самое с NumPy, потому что там я всегда могу полагаться на конверсии, не меняя ни одного бита:

import numpy as np

intarr = np.array([0x7f800001], dtype='uint32')
f = np.fromstring(intarr.tostring(), dtype='f4')
print np.isnan(f)

Это дает:

RuntimeWarning: invalid value encountered in isnan
[True]

но если значение заменено на 0x7fc00001, ошибки нет.

Гипотезы

Оба Python и NumPy будут довольны, если я установлю is_quiet и использую остальные биты для своих целей. Python обрабатывает бит сам по себе, NumPy полагается на языковые реализации более низкого уровня и/или на реализацию аппаратного FP.

Вопрос

Является ли моя гипотеза правильной, и может ли она быть доказана или опровергнута какой-либо официальной документацией? Или это одна из тех вещей, которые зависят от платформы?

Я нашел здесь что-то очень важное: Как отличить разные типы NaN-плавающих в Python, но я не мог найти никакого официального слова о том, как переносить дополнительную информацию NaNs следует обрабатывать в Python или NumPy.

Ответы

Ответ 1

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

Поскольку NumPy и Python обрабатывают числа совершенно по-другому, этот ответ состоит из двух частей.

Что действительно происходит в Python и NumPy с NaNs

NumPy

Это может быть немного специфично для платформы, но на большинстве платформ NumPy использует gcc встроенный isnan, который, в свою очередь, делает что-то быстро. Предупреждения во время выполнения исходят от более глубоких уровней, от аппаратного обеспечения в большинстве случаев. (NumPy может использовать несколько методов определения состояния NaN, таких как x!= X, который работает на платформах AMD 64, но с gcc он меньше gcc, который, вероятно, использует довольно короткий код для этой цели.)

Итак, в теории нет способа гарантировать, как NumPy обрабатывает NaN, но на практике на более общих платформах он будет делать, как говорит стандарт, потому что это то, что делает оборудование. Сам NumPy не заботится о типах NaN. (За исключением некоторых поддерживаемых NumPy не-hw-поддерживаемых типов данных и платформ.)

Python

Здесь история становится интересной. Если платформа поддерживает поплавки IEEE (в большинстве случаев), Python использует библиотеку C для арифметики с плавающей запятой и, следовательно, почти в большинстве случаев аппаратные инструкции. Таким образом, для NumPy не должно быть никаких различий.

За исключением... В Python обычно нет такой вещи, как 32-битный float. В плавающих объектах Python используется C double, который представляет собой 64-разрядный формат. Как преобразовать специальные NaN между этими форматами? Чтобы увидеть, что происходит на практике, следующий небольшой код C помогает:

/* nantest.c - Test floating point nan behaviour with type casts */

#include <stdio.h>
#include <stdint.h>

static uint32_t u1 = 0x7fc00000;
static uint32_t u2 = 0x7f800001;
static uint32_t u3 = 0x7fc00001;

int main(void)
    {
    float f1, f2, f3;
    float f1p, f2p, f3p;
    double d1, d2, d3;
    uint32_t u1p, u2p, u3p;
    uint64_t l1, l2, l3;

    // Convert uint32 -> float
    f1 = *(float *)&u1; f2 = *(float *)&u2; f3 = *(float *)&u3;

    // Convert float -> double (type cast, real conversion)
    d1 = (double)f1; d2 = (double)f2; d3 = (double)f3;

    // Convert the doubles into long ints
    l1 = *(uint64_t *)&d1; l2 = *(uint64_t *)&d2; l3 = *(uint64_t *)&d3;

    // Convert the doubles back to floats
    f1p = (float)d1; f2p = (float)d2; f3p = (float)d3;

    // Convert the floats back to uints
    u1p = *(uint32_t *)&f1p; u2p = *(uint32_t *)&f2p; u3p = *(uint32_t *)&f3p;

    printf("%f (%08x) -> %lf (%016llx) -> %f (%08x)\n", f1, u1, d1, l1, f1p, u1p);
    printf("%f (%08x) -> %lf (%016llx) -> %f (%08x)\n", f2, u2, d2, l2, f2p, u2p);
    printf("%f (%08x) -> %lf (%016llx) -> %f (%08x)\n", f3, u3, d3, l3, f3p, u3p);

    return 0;
    }

Отпечатки:

nan (7fc00000) -> nan (7ff8000000000000) -> nan (7fc00000)
nan (7f800001) -> nan (7ff8000020000000) -> nan (7fc00001)
nan (7fc00001) -> nan (7ff8000020000000) -> nan (7fc00001)

Посмотрев на строку 2, очевидно, что мы имеем то же явление, что и у Python. Таким образом, это преобразование в double, которое вводит дополнительный бит is_quiet сразу после экспоненты в 64-битной версии.

Это звучит немного странно, но на самом деле стандарт говорит (IEEE 754-2008, раздел 6.2.3):

Преобразование тихого NaN из более узкого формата в более широкий формат в том же самом основании, а затем обратно в тот же более узкий формат, не должно изменять тихую полезную нагрузку NaN каким-либо образом, кроме как сделать ее канонической.

Это ничего не говорит о распространении сигнальных NaN. Однако это объясняется разделом 6.2.1.:

Для двоичных форматов полезная нагрузка кодируется в p - 2 младших значащих битах конечного значения поля.

P выше - точность, 24 бит для 32-битного поплавка. Итак, моя ошибка заключалась в использовании сигнальных NaN для полезной нагрузки.

Резюме

Я получил следующие домашние очки:

  • использование qNaNs (тихие NaN) поддерживается и поощряется IEEE 754-2008
  • Нечетные результаты были связаны с тем, что я пытался использовать sNaN и преобразовывать типы, в результате которых был установлен бит is_quiet
  • Оба NumPy и Python действуют в соответствии с IEEE 754 на наиболее распространенных платформах.
  • реализация сильно опирается на базовую реализацию C и, таким образом, гарантирует очень мало (в Python есть даже некоторый код, который признает, что NaN не обрабатываются, поскольку они должны быть на некоторых платформах).
  • единственный безопасный способ справиться с этим - сделать немного DIY с полезной нагрузкой

Однако есть одна вещь, которая реализована ни в Python, ни в NumPy (ни на каком-либо другом языке, с которым я столкнулся). Раздел 5.12.1:

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