Представление чисел в двоичном коде с плавающей запятой
Часто приходится обрабатывать очень большие числа (например, расстояние между звёздами) или наоборот очень маленькие числа (например, размеры атомов или электронов). При таких вычислениях пришлось бы использовать числа с очень большой разрядностью. В то же время нам не нужно знать расстояние между звёздами с точностью до миллиметра. Для вычислений с такими величинами числа с фиксированной запятой неэффективны.
В десятичной арифметике для записи таких чисел используется алгебраическая форма. При этом число записывается в виде мантиссы, умноженной на 10 в степени, отображающей порядок числа, Например:
Для записи двоичных чисел тоже используется такая форма записи. Она позволяет работать с числами с большим диапазоном значений Эта форма записи называется запись числа с плавающей точкой. Напомним, что мантисса не может быть больше единицы и после запятой в мантиссе не можетзаписываться ноль.
Для записи числа в формате с плавающей запятой одинарной точности требуется тридцатидвухбитовое слово. Для записи чисел с двойной точностью требуется шестидесятичетырёхбитовое слово. Чаще всего числа хранятся в нескольких соседних ячейках памяти процессора. Форматы числа в формате с плавающей запятой одинарной точности и числа в формате с плавающей запятой удвоенной точности приведены на рисунке
Рисунок 1. Форматы числа с плавающей запятой
На рисунке буквой S обозначен знак числа, 0 — это положительное число, 1 — отрицательное число.
Группа бит, обозначенная e предназначена для записи смещённого порядка числа. Смещение потребовалось, чтобы не вводить в двоичный код числа с плавающей запятой еще один знак. Смещённый порядок всегда является положительным числом. В двоичном коде одинарной точности float для записи порядка числа выделено восемь бит. Для него смещение порядка числа принято 127. Для смещённого порядка в двоичном коде числа с плавающей запятой двойной точности double отводится 11 бит. В нем смещение порядка числа составляет — 1023.
В десятичной мантиссе после запятой могут присутствовать цифры 1. 9, а в двоичной — только 1. Поэтому для хранения единицы после двоичной запятой не выделяется отдельный бит в числе с плавающей запятой. Единица подразумевается, как и двоичная запятая. Кроме того, в формате чисел с плавающей запятой принято, что мантисса всегда больше 1. То есть диапазон значений мантиссы лежит в диапазоне от 1 до 2.
Рассмотрим несколько примеров:
1) Определить число с плавающей запятой, лежащее в четырёх соседних байтах:
11000001 01001000 00000000 00000000
— Знаковый бит, равный 1 показывает, что число отрицательное.
— Экспонента 10000010 в десятичном виде соответствует числу 130. Вычтя число 127 из 130, получим число 3.
— Теперь запишем мантиссу: 1,100 1000 0000 0000 0000 0000
— И, наконец, определим десятичное число: 1100,1b = 12,5d
2) Определить число с плавающей запятой, лежащее в четырёх соседних байтах:
11000011 00110100 00000000 00000000
— Знаковый бит, равный 1 показывает, что число отрицательное.
— Экспонента 10000110 в десятичном виде соответствует числу 134. Вычтя число 127 из 134, получим число 7.
— Теперь запишем мантиссу: 1,011 0100 0000 0000 0000 0000
— И, наконец, определим десятичное число: 10110100b=180d
Для того чтобы записать ноль, в двоичном представлении числа с плавающей запятой достаточно записать в смещенный порядок число 00000000b. Значение мантиссы при этом не имеет значения. Число, в котором все байты равны 0, тоже попадает в этот диапазон значений.
Бесконечность в числе с плавающей запятой соответствует смещенному порядку 11111111b и мантиссе, равной 1,0. При этом существует минус бесконечность и плюс бесконечность (переполнение и антипереполнение), которые часто отображаются на экран монитора компьютера или дисплей микропроцессорного устройства как +INF и -INF.
Все остальные комбинации битов мантиссы числа с плавающей запятой (в том числе и все единицы) при смещенном порядке 11111111b воспринимаются языками программирования как не числа и отображаются на экран: NaN.
Понравился материал? Поделись с друзьями!
Другие виды двоичных кодов:
Целочисленные двоичные коды Представление двоичных чисел в памяти компьютера или микроконтроллера https://digteh.ru/proc/IntCod.php
Двоично-десятичный код Иногда бывает удобно хранить числа в памяти процессора в десятичном виде https://digteh.ru/proc/DecCod.php
Запись текстов двоичным кодом Представление текстов в памяти компьютеров и микроконтроллеров https://digteh.ru/proc/text.php
Системы счисления В настоящее время и в технике и в быту широко используются как позиционные, так и непозиционные системы счисления. https://digteh.ru/digital/SysSchis.php
Предыдущие версии сайта: http://neic.nsk.su/
Об авторе: к.т.н., доц., Александр Владимирович Микушин
Кандидат технических наук, доцент кафедры САПР СибГУТИ. Выпускник факультета радиосвязи и радиовещания (1982) Новосибирского электротехнического института связи (НЭИС).
А.В.Микушин длительное время проработал ведущим инженером в научно исследовательском секторе НЭИС, конструкторско технологическом центре «Сигнал», Научно производственной фирме «Булат». В процессе этой деятельности он внёс вклад в разработку систем радионавигации, радиосвязи и транкинговой связи.
Научные исследования внедрены в аппаратуре радинавигационной системы Loran-C, комплексов мобильной и транкинговой связи «Сигнал-201», авиационной системы передачи данных «Орлан-СТД», отечественном развитии системы SmarTrunkII и радиостанций специального назначения.
Представление чисел в позиционной системе исчисления
Любое n-разрядное число K в позиционной системе может быть представлено в виде следующего полинома где ai – весовые коэффициенты при соответствующих разрядах, X – основание системы.
Основание системы определяет количество чисел, которое используется в системе исчисления. Коэффициенты ai могут быть представлены только числами из этого набора. Например, в двоичной системе (основание X=2), используются числа 0 и 1, в десятичной (X=10) – 0…9, в шестнадцатеричной (X=16) – 0…9 и символы A,B,C,D,E,F.
Двоичная система исчисления является доминирующей в вычислительной технике. Она наиболее удобна для обработки данных, представленных на уровне электрических сигналов. Полином (1) в этом случае имеет вид K = an-1*2 n-1 + an-2*2 n-2 + … + a1*2 1 + a0.
В табл.1 приведен пример записи чисел в трех различных системах исчисления.
Табл.1. Числа в различных системах исчисления.
Представление чисел в дополнительном коде
Представление неотрицательных чисел в дополнительном коде (см. табл.2) совпадает с соответствующим представлением чисел в позиционной системе. Отрицательные же числа в дополнительном коде получаются в результате следующего действия — X = 2 n – X (2).
Табл.2. Представление чисел в дополнительном коде.
Десятичная система
Двоично-десятичное кодирование
Для хранения десятичных чисел часто используется двоично-десятичное кодирование (ДДК). В этом случае байты разбиваются на тетрады, в каждой из которых размещаются двоично-кодированные десятичные числа. Например, числу 100101102 в ДДК соответствует десятичное число 9610 (10012 = 910 – старший полубайт, 01102 = 610 – младший полубайт).
ДДК применяется в основном для представления пользовательской информации в удобной и понятной для человека десятичной системе исчисления. Однако формат ДДК позволяет хранить в одном байте только 100 десятичных чисел (числа 0…9910). И, кроме того, далеко не каждый микропроцессор имеет в своей системе команды для обработки чисел в ДДК (в частности, у микроконтроллеров AVR такие команды отсутствуют).
Представление чисел в форме с фиксированной запятой
Числа такого типа содержат две части: целую и дробную. Под целую часть отводиться n-m разрядов, а под дробную m разрядов. Естественно, что в случае m=0 дробная часть отсутствует и (3) переходит формулу (1).
Разряды a7a6a5a4 представляют собой целую часть числа, а a3a2a1a0 – дробную.
Формат чисел с фиксированной запятой имеет обозначение (N.Q), где N и Q число целых и дробных разрядов соответственно. Так в предыдущем примере число K имеет формат (4.4).
При программировании 8-разрядных микроконтроллеров на ассемблере в подобной форме представления, как правило, нет необходимости. Однако числа с фиксированной запятой имеют определенные преимущества при цифровой обработке сигналов в специализированных DSP (Digital Signal Processor) процессорах, где их поддержка введена аппаратном уровне.
Представление чисел в форме с плавающей запятой
где M – мантисса, E – порядок.
Таким образом, числа с плавающей запятой должны отображаться в ОЗУ микропроцессора в виде двух частей, каждая из которых может иметь соответствующее число разрядов (чаще всего кратное 8).
Как M, так и E, могут быть знаковыми числами. Мантисса должна быть представлена в виде числа с фиксированной запятой и находиться в пределах 0≤|m|≤1, а порядок M – виде целого числа.
Для повышения точности представления чисел разрядность M требуется всегда выше, чем разрядность E. Диапазон же представления чисел практически полностью определяется разрядностью E.
Программирование математических операций над числами с плавающей запятой связано со значительными трудностями и требует большого расхода памяти программ и данных. К счастью такой формат является избыточным для подавляющего большинства приложений, где используются 8-разрядные микроконтроллеры, и в основном встречается только в некоторых задачах управления и фильтрации. Числа с плавающей запятой иногда также применяют для создания и пересылки архивов данных в компьютер.
Представление символьной информации в кодировке ASCII
ASCII (American Standard Code for Information Interchange) – общепринятый стандарт кодирования информации. В соответствии с ним каждому символу из табл.3 сопоставляется индивидуальный 7-разрядный числовой код. Таблица разделена на 128 частей, но реально, только диапазону адресов 0x20…0x7E соответствует информация, которая воспринимается, как символьная. Остальные кода, размещенные по адресам 0x00…0x1F и 0x7F, отведены командам управления. Большинство из них поддерживает любая терминальная программа, совместимая с протоколом обмена ASCII. Набор из 21 команды очень универсален и предоставляет широкие возможности по управлению процессом передачи данных. Наиболее распространенные из них приведены в табл.4.
Символы в таблице ASCII расположены не случайным образом, а имеют строгий порядок. Кода символов ‘0’…‘9’ получаются добавлением смещения 0x30 к своему числовому представлению. Буквы ‘A’…‘Z’ имеют порядок следования аналогичный латинскому алфавиту и получаются добавлением смещения 0x41 к порядковому номеру буквы в алфавите. То же самое относится к набору из 26 символов нижнего регистра ‘a’…‘z’ для получения, которых необходимо смещение 0x61.
Преобразуем в строку. Часть 2. Числа с фиксированной и плавающей точкой.
Продолжаем преобразовывать всё, что можно в строки. В предыдущей статье были целые числа, теперь очередь чисел с фиксированной и плавающей точкой. Все рассмотренные примеры с фиксированной точкой используют формат с 16-ю битами для дробной части и 16-ю битами для целой части, так называемый формат Q16, однако легко могут быть адаптированы для других форматов. В качестве чисел с плавающей точкой использован 32-х разрядный float.
Числа с фиксированной точкой.
Стандартная библиотека языка Си не предоставляет никаких средств для вывода чисел с фиксированной точкой. Зато есть вывод чисел с плавающей точкой — функции *printf. Можно преобразовать число с фиксированной в число с плавающей точкой и вывести его с помощью, например, sprintf. Первый наивный способ преобразовать фиксированную точку в плавающую, состоит в том, что целочисленное представление числа с фиксированной точкой непосредственно приводится к типу с плавающей точкой, после чего делится на целочисленное представление на единицы в формате с фиксированной точкой:
Тип double, а не float, в качестве типа с плавающей точкой использован потому, что sprintf принимает именно его. Для AVR это не имеет значения, так как float и double там одно и тоже. Этот метод очень прост и быстр в реализации, к тому-же достаточно надёжен: меньше кода — меньше ошибок. На этом его преимущества заканчиваются и остаются недостатки: — медленный: AVR 4200 тактов, STM32 — 4234 тактов в среднем; — использует sprintf и тянет за собой всё её обвязку; — использует целых две операции с плавающей точкой.
Первый недостаток усугубляется третьим. Попробуем избавится от последнего недостатка, вручную преобразовав число с фиксированной точкой во float. Для этого надо найти в числе с фиксированной точкой самый старший единичный бит и сдвинуть его на место неявного старшего единичного бита числа с плавающей точкой — 23-его бита. После чего этот бит надо обнулить, так как там хранится младший бит экспоненты. Таким образом получится мантисса искомого float-а. Экспонента вычисляется по формуле: exponent = 127 + «число бит в дробной части» — 9 + «насколько бит сдвинули исходное число» — 9 здесь нужно потому, что целевой бит, куда нужно сдвинуть старший единичный бит — девятый с конца. «насколько бит сдвинули исходное число» — положительное, если сдвигали вправо, отрицательное, если влево. При сдвиге вправо могут потеряться до 8 бит младших дробной части.
Число с фиксированной точкой в этом примере без знаковое, поэтому знак не трогаем. В результате на AVR такой метод работает за 2298 тактов в среднем — почти в два раза быстрее. На STM32 всё по прежнему «плохо» — 4194 такта — операции с плавающей точкой там не являются узким местом — это sprintf так медленно работает.
Следующий вариант — это преобразовывать целую часть отдельно по одному из методов из предыдущей статьи, а дробную — методом умножения на 10. Для преобразования целой части целесообразно взять самый быстрый метод — деление на 10 сдвигами и сложениями.
Получилось неплохо никаких внешних зависимостей и кода не много. На AVR этот метод работает в среднем за 490 тактов, на STM32 — за 148. Результат отличный — быстрее, чем 32-х битное целое преобразуется по методу быстрого деления.
Также для чисел с фиксированной точкой применим метод умножения на 10. Он почти ничем не отличается от аналогичного для целых чисел. Только на этот раз для умножения исходного числа на AVR-ках я применил не тупой цикл со сдвигами и сложениями, честное 64-битное умножение по алгоритму Д. Кнута. Дело в том, что встроенное 64-битное умножение в avr-gcc не оптимизировано и работает неприлично долго порядка 800 тактов. Реализация, приведённая ниже работает за чуть более 100 тактов. Функция mul32hu возвращает старшие 32 бита от произведения, которые получаются в переменной u32[1]. Если нужно полное 64-х битное произведение, то можно вытащить и младшую часть — она находится в переменной u32[0].
На STM32 имеется инструкция умножения 32×32=>64, по этому используем встроенное умножение.
Кстати, если из метода с умножением на 10 из предыдущей статьи (utoa_fract) выкинуть цикл умножения и заменить его на mul32hu, то он будет выполняться на AVR в среднем за 648 тактов, что делает его быстрее, чем метод вычитания степеней 10 (utoa_cycle_sub).
В самом методе, как и для целых чисел, сначала умножаем исходное число на множитель, который представляет из себя степень двойки, соответствующую той битовой позиции, где будут извлекаться цифры, делённую на степень 10, такую, чтоб не было переполнения. Затем в цикле извлекаются десятичные цифры с определённой битовой позиции, в данном случае с 29 бита, и число умножается на 10. После извлечения всех цифр производится коррекция коррекция старшей цифры, так как она может оказаться больше 10. Для этой коррекции в строковом буфере предусмотрен 1 символ, который инициализируется в ‘0’. После удаляются ведущие и завершающие нули.
На AVR этот метод работает в среднем за 604 такта, на STM32 — 236. Не плохо, но всё-же хуже чем у предыдущего метода. Однако, если дробная часть будет не 16 разрядов, как в примере, а занимать все 32 бита, то умножение будет не нужно и этот метод станет проще и быстрее предыдущего. Внимательный читатель заметит эту конструкцию.
Почему нельзя оставить только закомментированный фрагмент? Потому, что в avr-gcc сдвиг 32-х разрядного числа плохо оптимизирован. Вместо того, чтоб взять старший байт и сдвинуть его на 4 разряда, компилятор генерирует цикл, в котором честно сдвигает все 4 байта на 28 разрядов. Это у него занимает 28*(4+1) + 2 = 142 цикла вместо 3 (mov, swap, andi), которые должны быть. А ведь это выполняется в цикле 9 раз, итого 1278 такта — дофига для такой мелочи. Приходится компилятору немного помогать.
Плавающая точка.
Для преобразования чисел с плавающей точкой в строку тоже можно придумать много вариантов, но я этого делать не буду. Приведу только один вариант и сравню его с sprintf. Этот вариант будет, естественно, умножение на 10, потому, что чтоб вытащить цифры из мантиссы, её всё-равно придётся умножать. В методах умножения на 10 для целых и чисел с фиксированной точкой, начальное 64 битное умножение занимало заметную долю времени. Поэтому они оказывались несколько медленнее, чем некоторые другие, несмотря на то, что умножение на 10 само по себе быстрее, чем деление. В случае с плавающей точкой он должен оказаться безоговорочным лидером. Для начала посмотрим на что способна sprintf, при выводе float-ов. Для AVR sprintf работает в среднем за 2000 такта, на STM32 — 3950. Не стоит ругать sprintf на STM32 — помним, что там преобразуется double, а не float. При этом sprintf с ключиком %g округляет результат до заданной точности, удаляет незначащие нули и выводит результат в наиболее подходящем виде: в обычном 123.45, или научном 1.2345e+2. Задача состоит в том, чтоб получить вывод идентичный sprintf, но гораздо быстрее. Разделим задачу на две части: 1 — извлечение из float-а значащих цифр и десятичной экспоненты; 2 — форматирование результата в нужном виде. Для первой части нужно: — распаковать float, вытащить из него мантиссу, двоичную экспоненту и знак; — обработать особые случаи: ±0, ±inf, nan. — получить десятичную экспоненту; — получить множитель 2 в степени экспоненты, делённая на степень 10, такую, чтоб не было переполнения; — умножить мантиссу на этот множитель; — извлечь требуемое количество значащих цифр + одну для округления; — округлить полученные цифры; — удалить незначащие нули.
Извлечение цифр и десятичной экспоненты выделим в отдельную функцию:
presc — это необходимое количество значащих цифр. Теперь распакуем float:
Финт ушами со сгвигим экспоненты на 16 бит, а потом еще на 7, опять-же чтоб задобрить компилятор. Теперь извлечём знак:
Денормализованные числа считаем равными нулю. Это не совсем правильно, но мне не хотелось с ними возиться. Резервируем один символ в буфере для округления:
Вычислить требуемый множитель за приемлемое время не получится, можно конечно, сдвигать влево-вправо, умножать-делить на 10 константу в цикле в зависимости от экспоненты, но это очень медленно. Поэтому множитель нужно брать из таблицы. Но 256 32-х битный чисел — этож целый килобайт, слишком жирно. Соседние элементы этой таблицы будут отличаться на степень 2, следовательно результат умножения для соседних элементов можно получить сдвигая его в нужную сторону на количество бит, равное расстоянию между этими элементами. Поскольку мантисса 24-х разрядная, то сдвигать её можно на 8 разрядов без потери точности. Это значит, что в таблице можно оставить только каждый 8 элемент, 32 элемента — уже приемлемо. Чтоб сохранить точность, сначала будем сдвигать мантиссу на 8 бит влево, потом умножать, потом сдвигать влево на требуемое количество бит. Так удаётся сохранить все 24 бита точности не выходя за 32-х битную арифметику (кроме умножения, конечно). Десятичную экспоненту несложно вычислить из двоичной, они связаны линейной зависимостью. Сложность только в том, чтоб подобрать масштабирующие и сдвигающие коэффициенты, при которых будет верный результат во всём диапазоне входных значений, и не вылезти при этом за пределы 16-ти битной арифметики.
Теперь удалим ведущие нули:
Удаление завершающих нулей.
В итоге в буфере знак и значащие цифры, десятичная экспонента в возвращаемом из функции значении. Теперь дело за форматированием:
Получилось достаточно многословно, но результат того стоил. На AVR в среднем 759 тактов, на STM32 — 425. Для сравнения в avr-libc есть функция dtostre, которая использует для преобразования тот-же движок, что и sprintf, но выводит всегда в экспоненциальном формате и не удаляет незначащие нули. Её среднее время в этом тесте — 1215 тактов. И это при том, что преобразование там реализовано на ассемблере.
Результат тестов с фиксированной точкой для AVR: Результат тестов с фиксированной точкой для STM32: Результат тестов с плавающей точкой для AVR: Результат тестов с плавающей точкой для STM32:
Что нужно знать про арифметику с плавающей запятой
В далекие времена, для IT-индустрии это 70-е годы прошлого века, ученые-математики (так раньше назывались программисты) сражались как Дон-Кихоты в неравном бою с компьютерами, которые тогда были размером с маленькие ветряные мельницы. Задачи ставились серьезные: поиск вражеских подлодок в океане по снимкам с орбиты, расчет баллистики ракет дальнего действия, и прочее. Для их решения компьютер должен оперировать действительными числами, которых, как известно, континуум, тогда как память конечна. Поэтому приходится отображать этот континуум на конечное множество нулей и единиц. В поисках компромисса между скоростью, размером и точностью представления ученые предложили числа с плавающей запятой (или плавающей точкой, если по-буржуйски).
Арифметика с плавающей запятой почему-то считается экзотической областью компьютерных наук, учитывая, что соответствующие типы данных присутствуют в каждом языке программирования. Я сам, если честно, никогда не придавал особого значения компьютерной арифметике, пока решая одну и ту же задачу на CPU и GPU получил разный результат. Оказалось, что в потайных углах этой области скрываются очень любопытные и странные явления: некоммутативность и неассоциативность арифметических операций, ноль со знаком, разность неравных чисел дает ноль, и прочее. Корни этого айсберга уходят глубоко в математику, а я под катом постараюсь обрисовать лишь то, что лежит на поверхности.
1. Основы
Множество целых чисел бесконечно, но мы всегда можем подобрать такое число бит, чтобы представить любое целое число, возникающее при решении конкретной задачи. Множество действительных чисел не только бесконечно, но еще и непрерывно, поэтому, сколько бы мы не взяли бит, мы неизбежно столкнемся с числами, которые не имеют точного представления. Числа с плавающей запятой — один из возможных способов предсталения действительных чисел, который является компромиссом между точностью и диапазоном принимаемых значений.
Число с плавающей запятой состоит из набора отдельных разрядов, условно разделенных на знак, экспоненту порядок и мантиссу. Порядок и мантисса — целые числа, которые вместе со знаком дают представление числа с плавающей запятой в следующем виде:
Математически это записывается так:
Основание определяет систему счисления разрядов. Математически доказано, что числа с плавающей запятой с базой B=2 (двоичное представление) наиболее устойчивы к ошибкам округления, поэтому на практике встречаются только базы 2 и, реже, 10. Для дальнейшего изложения будем всегда полагать B=2, и формула числа с плавающей запятой будет иметь вид:
Что такое мантисса и порядок? Мантисса – это целое число фиксированной длины, которое представляет старшие разряды действительного числа. Допустим наша мантисса состоит из трех бит (|M|=3). Возьмем, например, число «5», которое в двоичной системе будет равно 1012. Старший бит соответствует 2 2 =4, средний (который у нас равен нулю) 2 1 =2, а младший 2 0 =1. Порядок – это степень базы (двойки) старшего разряда. В нашем случае E=2. Такие числа удобно записывать в так называемом «научном» стандартном виде, например «1.01e+2». Сразу видно, что мантисса состоит из трех знаков, а порядок равен двум.
Допустим мы хотим получить дробное число, используя те же 3 бита мантиссы. Мы можем это сделать, если возьмем, скажем, E=1. Тогда наше число будет равно
Обратите внимание, что одно и то же число имеет несколько представлений. Это не удобно для оборудования, т.к. нужно учитывать множественность представлния при сравнении чисел и при выполнении над ними арифметических операций. Кроме того, это не экономично, поскольку число представлений — конечное, а повторения уменьшают множество чисел, которые вообще могут быть представлены. Поэтому уже в самых первых машинах начали использовать трюк, делая первый бит мантиссы всегда положительным. Такое предаставление назвали нормализованным.
Это экономит один бит, так как неявную единицу не нужно хранить в памяти, и обеспечивает уникальность представления числа. В нашем примере «2» имеет единственное нормализованное представление («1.000e+1»), а мантисса хранится в памяти как «000», т.к. старшая единица подразумевается неявно. Но в нормализованном представлении чисел возникает новая проблема — в такой форме невозможно представить ноль.
Строго говоря, нормализованное число имеет следующий вид:
Качество решения задач во многом зависит от выбора представления чисел с плавающей запятой. Мы плавно подошли к проблеме стандартизации такого представления.
2. Немного истории
В 60-е и 70-е годы не было единого стандарта представления чисел с плавающей запятой, способов округления, арифметических операций. В результате программы были крайне не портабельны. Но еще большей проблемой было то, что у разных компьютеров были свои «странности» и их нужно было знать и учитывать в программе. Например, разница двух не равных чисел возвращала ноль. В результате выражения «X=Y» и «X-Y=0» вступали в противоречие. Умельцы обходили эту проблему очень хитрыми трюками, например, делали присваивание «X=(X-X)+X» перед операциями умножения и деления, чтобы избежать проблем.
Инициатива создать единый стандарт для представления чисел с плавающей запятой подозрительно совпала с попытками в 1976 году компанией Intel разработать «лучшую» арифметику для новых сопроцессоров к 8086 и i432. За разработку взялись ученые киты в этой области, проф. Джон Палмер и Уильям Кэхэн. Последний в своем интервью высказал мнение, что серьезность, с которой Intel разрабатывала свою арифметику, заставила другие компании объединиться и начать процесс стандартизации.
Все были настроены серьезно, ведь очень выгодно продвинуть свою архитектуру и сделать ее стандартной. Свои предложения представили компании DEC, National Superconductor, Zilog, Motorola. Производители мейнфреймов Cray и IBM наблюдали со стороны. Компания Intel, разумеется, тоже представила свою новую арифметику. Авторами предложенной спецификации стали Уильям Кэхэн, Джероми Кунен и Гарольд Стоун и их предложение сразу прозвали «K-C-S».
Практически сразу же были отброшены все предложения, кроме двух: VAX от DEC и «K-C-S» от Intel. Спецификация VAX была значительно проще, уже была реализована в компьютерах PDP-11, и было понятно, как на ней получить максимальную производительность. С другой стороны в «K-C-S» содержалось много полезной функциональности, такой как «специальные» и «денормализованные» числа (подробности ниже).
В «K-C-S» все арифметические алгоритмы заданы строго и требуется, чтобы в реализации результат с ними совпадал. Это позволяет выводить строгие выкладки в рамках этой спецификации. Если раньше математик решал задачу численными методами и доказывал свойства решения, не было никакой гарантии, что эти свойства сохранятся в программе. Строгость арифметики «K-C-S» сделала возможным доказательство теорем, опираясь на арифметику с плавающей запятой.
Компания DEC сделала все, чтобы ее спецификацию сделали стандартом. Она даже заручилась поддержкой некоторых авторитетных ученых в том, что арифметика «K-C-S» в принципе не может достигнуть такой же производительности, как у DEC. Ирония в том, что Intel знала, как сделать свою спецификацию такой же производительной, но эти хитрости были коммерческой тайной. Если бы Intel не уступила и не открыла часть секретов, она бы не смогла сдержать натиск DEC.
Подробнее о баталиях при стандартизации смотрите в интервью профессора Кэхэна, а мы рассмотрим, как выглядит представление чисел с плавающей запятой сейчас.
3. Представление чисел с плавающей запятой сегодня
Разработчики «K-C-S» победили и теперь их детище воплотилось в стандарт IEEE754. Числа с плавающей запятой в нем представлены в виде знака (s), мантиссы (M) и порядка (E) следующим образом:
Замечание. В новом стандарте IEE754-2008 кроме чисел с основанием 2 присутствуют числа с основанием 10, так называемые десятичные (decimal) числа с плавающей запятой.
Чтобы не загромождать читателя чрезмерной информацией, которую можно найти в Википедии, рассмотрим только один тип данных, с одинарной точностью (float). Числа с половинной, двойной и расширенной точностью обладают теми же особенностями, но имеют другой диапазон порядка и мантиссы. В числах одинарной точности (float/single) порядок состоит из 8 бит, а мантисса – из 23. Эффективный порядок определяется как E-127. Например, число 0,15625 будет записано в памяти как
Рисунок взят из Википедии
3.1 Специальные числа: ноль, бесконечность и неопределенность
Неопределенность или NaN (от not a number) – это представление, придуманное для того, чтобы арифметическая операция могла всегда вернуть какое-то не бессмысленное значение. В IEEE754 NaN представлен как число, в котором E=Emax+1, а мантисса не нулевая. Любая операция с NaN возвращает NaN. При желании в мантиссу можно записывать информацию, которую программа сможет интерпретировать. Стандартом это не оговорено и мантисса чаще всего игнорируется.
Вернемся к примеру. Наш Emin=-1. Введем новое значение порядка, E=-2, при котором числа являются денормализованными. В результате получаем новое представление чисел:
Интервал от 0 до 0,5 заполняют денормализованные числа, что дает возможность не проваливаться в 0 рассмотренных выше примерах (0,5-0,25 и 1,5-1,25). Это сделало представление более устойчиво к ошибкам округления для чисел, близких к нулю.
Но роскошь использования денормализованного представления чисел в процессоре не дается бесплатно. Из-за того, что такие числа нужно обрабатывать по-другому во всех арифметических операциях, трудно сделать работу в такой арифметике эффективной. Это накладывает дополнительные сложности при реализации АЛУ в процессоре. И хоть денормализованные числа очень полезны, они не являются панацеей и за округлением до нуля все равно нужно следить. Поэтому эта функциональность стала камнем преткновения при разработке стандарта и встретила самое сильное сопротивление.
3.4 Очередность чисел в IEEE754
Одна из удивительных особенностей представления чисел в формате IEEE754 состоит в том, что порядок и мантисса расположены друг за другом таким образом, что вместе образуют последовательность целых чисел для которых выполняется:
4.2 Неассоциативность арифметических операций
В арифметике с плавающей запятой правило (a*b)*c = a*(b*c) не выполняется для любых арифметических операций. Например,
Допустим у нас есть программа суммирования чисел.
Некоторые компиляторы по умолчанию могут переписать код для использования нескольких АЛУ одновременно (будем считать, что n делится на 2):
Так как операции суммирования не ассоциативны, эти две программы могут выдать различный результат.
4.3 Числовые константы
Помните, что не все десятичные числа имеют двоичное представление с плавающей запятой. Например, число «0,2» будет представлено как «0,200000003» в одинарной точности. Соответственно, «0,2 + 0,2 ≈ 0,4». Абсолютная погрешность в отдельном случае может и не высока, но если использовать такую константу в цикле, можем получить накопленную погрешность.
4.4 Выбор минимума из двух значений
4.5 Сравнение чисел
Очень распространенная ошибка при работе с float-ами возникает при проверке на равенство. Например,
Ошибка здесь, во-первых, в том, что 0,2 не имеет точного двоичного представления, а во-вторых 0,2 – это константа двойной точности, а переменная fValue – одинарной, и никакой гарантии о поведении этого сравнения нет.
Лучший, но все равно ошибочный способ, это сравнивать разницу с допустимой абсолютной погрешностью:
Недостаток такого подхода в том, что погрешность представления числа увеличивается с ростом самого этого числа. Так, если программа ожидает «10000», то приведенное равенство не будет выполняться для ближайшего соседнего числа (10000,000977). Это особенно актуально, если в программе имеется преобразование из одинарной точности в двойную.
Выбрать правильную процедуру сравнения сложно и заинтересованных читателей я отсылаю к статье Брюса Доусона. В ней предлагается сравнивать числа с плавающей запятой преобразованием к целочисленной переменной. Это — лучший, хотя и не портабельный способ:
5. Проверка полноты поддержки IEE754
Думаете, что если процессоры полностью соответствуют стандарту IEEE754, то любая программа, использующая стандартные типы данных (такие как float/double в Си), будет выдавать один и тот же результат на разных компьютерах? Ошибаетесь. На портабельность и соответствие стандарту влияет компилятор и опции оптимизации. Уильям Кэхэн написал программу на Си (есть версия и для Фортрана), которая позволяет проверить удовлетворяет ли связка «архитектура+компилятор+опции» IEEE754. Называется она «Floating point paranoia» и ее исходные тексты доступны для скачивания. Аналогичная программа доступна для GPU. Так, например, компилятор Intel (icc) по умолчанию использует «расслабленную» модель IEEE754, и в результате не все тесты выполняются. Опция «-fp-model precise» позволяет компилировать программу с точным соответствием стандарту. В компиляторе GCC есть опция «-ffast-math», использование которой приводит к несоответствию IEEE754.
Заключение
Напоследок поучительная история. Когда я работал над тестовым проектом на GPU, у меня была последовательная и параллельная версия одной программы. Сравнив время выполнения, я был очень обрадован, так как получил ускорение в 300 раз. Но позже оказалось, что вычисления на GPU «разваливались» и обращались в NaN, а работа с ними в GPU была быстрее, чем с обычными числами. Интересно было другое — одна и та же программа на эмуляторе GPU (на CPU) выдавала корректный результат, а на самом GPU – нет. Позже оказалось, что проблема была в том, что этот GPU не поддерживал полностью стандарт IEEE754 и прямой подход не сработал.
Сейчас арифметика с плавающей запятой почти совершенна. Практически всегда наивный подход сработает, и программа, не учитывающая все ее особенности, выдаст правильный результат, а описанные подводные камни касаются только экзотических случаев. Но нужно всегда оставаться бдительным: в таком вопросе как компьютерная математика легко наступить на грабли.
