Что такое плавающая точка в процессоре

Обновлено: 07.07.2024

Важной частью архитектуры микропроцессоров Intel является наличие устройства для обработки числовых данных в формате с плавающей точкой, называемого математическим сопроцессором . Архитектура компьютеров на базе микропроцессоров вначале опиралась исключительно на целочисленную арифметику. С ростом мощи стали появляться устройства для обработки чисел с плавающей точкой. В архитектуре семейства микропроцессоров Intel 8086 устройство для обработки чисел с плавающей точкой появилось в составе компьютера на базе микропроцессора i8086/88 и получило название математический сопроцессор или просто сопроцессор. Выбор такого названия был обусловлен тем, что,

  • во-первых, это устройство было предназначено для расширения вычислительных возможностей основного процессора;
  • во-вторых, оно было реализовано в виде отдельной микросхемы, то есть его присутствие было необязательным. Микросхема сопроцессора для микропроцессора i8086/88 имела название i8087.

С появлением новых моделей микропроцессоров Intel совершенствовались и сопроцессоры, хотя их программная модель осталась практически неизменной. Как отдельные (а, соответственно, необязательные в конкретной комплектации компьютера) устройства, сопроцессоры сохранялись вплоть до модели микропроцессора i386 и имели название i287 и i387 соответственно. Начиная с модели i486, сопроцессор исполняется в одном корпусе с основным микропроцессором и, таким образом, является неотъемлемой частью компьютера.

Основные возможности математического сопроцессора:

Форма представления чисел с плавающей точкой описана здесь .

Общая форма представления вещественных чисел предполагает возможность размещения в разрядной сетке следующих типов.

Числа простой и двойной точности ( float ( DD ) и double ( DQ ) соответственно) могут быть представлены только в нормированной форме. При этом бит целой части числа является скрытым и подразумевает логическую 1. Остальные 23 (52) разряда хранят двоичную мантиссу числа.

Числа двойной расширенной точности ( long double ( DT )) могут быть представлены как в нормированной, так и в ненормированной форме, поскольку бит целой части числа не является скрытым и может принимать значения как 0, так и 1.

Основным типом данных, которыми оперирует математический сопроцессор, являются 10-байтные данные ( DT ).

Программная модель сопроцессора

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

Математический сопроцессор

В программной модели сопроцессора можно выделить три группы регистров:

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

Регистр состояния сопроцессора SWR

Регистр состояния swr – отражает текущее состояние сопроцессора после выполнения последней команды. В регистре swr содержатся поля, позволяющие определить: какой регистр является текущей вершиной стека сопроцессора, какие исключения возникли после выполнения последней команды, каковы особенности выполнения последней команды (некий аналог регистра флагов основного процессора).

Структурно регистр swr состоит из:

    6 флагов исключительных ситуаций PE, OE, UE, ZE, DE, IE.
    Исключения — это разновидность прерываний, с помощью которых процессор информирует программу о некоторых особенностях ее реального исполнения. Сопроцессор также обладает способностью возбуждения подобных прерываний при возникновении определенных ситуаций (не обязательно ошибочных). Все возможные исключения сведены к 6и типам, каждому из которых соответствует 1 бит в регистре swr . Программисту не обязательно писать обработчик для реакции на ситуацию, приведшую к некоторому исключению. Сопроцессор умеет самостоятельно реагировать на многие из них. Это так называемая обработка исключений по умолчанию. Для того чтобы вызвать обработку определенного типа исключения по умолчанию, необходимо это исключение оставить не маскированным. Такое действие выполняется с помощью установки в 1 соответствующего бита в управляющем регистре сопроцессора cwr . Типы исключений, фиксируемые с помощью регистра swr:

  • IE (Invalide operation Error) — недействительный код операция;
  • DE (Denormalized operand Error) — ненормированный операнд;
  • ZE (divide by Zero Error) — ошибка деления на нуль;
  • ОЕ (Overflow Error) — ошибка переполнения. Возникает в случае выхода порядка числа за максимально допустимый диапазон;
  • UE (Underflow Error) — ошибка антипереполнения. Возникает, когда результат слишком мал (близок к нулю);
  • РЕ (Precision Error) — ошибка точности. Устанавливается, когда сопроцессору приходится округлять результат из-за того, что его точное представление невозможно. Так, сопроцессору никогда не удастся точно разделить 10 на 3.

Регистр управления сопроцессора CWR

Регистр управления работой сопроцессора cwr – определяет особенности обработки числовых данных. С помощью полей в регистре cwr можно регулировать точность выполнения численных вычислений, управлять округлением, маскировать исключения.

Он состоит из:

  • шести масок исключений PM, UM, OM, ZM, DM, IM ;
  • поля управления точностью PC (Precision Control);
  • поля управления округлением RC (Rounding Control).

2-битовое поле управления точностью PC предназначено для выбора длины мантиссы. Возможные значения в этом поле означают:

  • PC =00 — длина мантиссы 24 бита;
  • PC =10 — длина мантиссы 53 бита;
  • PC =11 — длина мантиссы 64 бита.

По умолчанию устанавливается значение поля PC =11.

Поле управления округлением RC позволяет управлять процессом округления чисел в процессе работы сопроцессора. Необходимость операции округления может появиться в ситуации, когда после выполнения очередной команды сопроцессора получается не представимый результат, например, периодическая дробь. Установив одно из значений в поле RC , можно выполнить округление в необходимую сторону.
Значения поля RC с соответствующим алгоритмом округления:

  • 00 — значение округляется к ближайшему числу, которое можно представить в разрядной сетке регистра сопроцессора;
  • 01 — значение округляется в меньшую сторону;
  • 10 — значение округляется в большую сторону;
  • 11 — производится отбрасывание дробной части числа. Используется для приведения значения к форме, которая может использоваться в операциях целочисленной арифметики.

Бит 12 в регистре cwr физически отсутствует и считывается равным 0.

Регистр тегов twr – представляет собой совокупность двухбитовых полей. Каждое поле соответствует определенному физическому регистру стека и характеризует его текущее состояние. Команды сопроцессора используют этот регистр, например, для того, чтобы определить возможность записи значений в эти регистры. Изменение состояния любого регистра стека отражается на содержимом соответствующего этому регистру 2-битового поля регистра тега. Возможны следующие значения в полях регистра тега:

  • 00 — регистр стека сопроцессора занят допустимым ненулевым значением;
  • 01 — регистр стека сопроцессора содержит нулевое значение;
  • 10 — регистр стека сопроцессора содержит одно из специальных численных значений, за исключением нуля;
  • 11 — регистр пуст и в него можно производить запись. Это значение в двухбитовом поле регистра тегов не означает, что все биты соответствующего регистра стека должны быть обязательно нулевыми.
Принцип работы сопроцессора

Принцип работы сопроцессора совместно с центральным процессором
Процессор и сопроцессор имеют свои раздельные системы команд и форматы обрабатываемых данных. Несмотря на то, что сопроцессор архитектурно представляет собой отдельное вычислительное устройство, он не может существовать отдельно от основного процессора. Процессор и сопроцессор, являясь двумя самостоятельными вычислительными устройствами, могут работать параллельно. Но это распараллеливание распространяется только на выполнение команд. Оба процессора подключены к общей системной шине и имеют доступ к одной и той же информации. Инициирует процесс выборки очередной команды всегда основной процессор. После выборки команда попадает одновременно в оба процессора. Любая команда сопроцессора имеет код операции, первые пять бит, которого имеют значение 11011. Когда код операции начинается этими битами, то основной процессор по дальнейшему содержимому кода операции выясняет, требует ли данная команда обращения к памяти. Если это так, то основной процессор формирует физический адрес операнда и обращается к памяти, после чего содержимое ячейки памяти выставляется на шину данных. Если обращение к памяти не требуется, то основной процессор заканчивает работу над данной командой (не делая попытки ее исполнения) и приступает к декодированию следующей команды из текущего входного командного потока. Выбранная команда попадает в сопроцессор одновременно с основным процессором. Сопроцессор, определив по первым пяти битам, что очередная команда принадлежит его системе команд, начинает ее исполнение. Если команда требует операнды из памяти, то сопроцессор обращается к шине данных за чтением содержимого ячейки памяти, которое к этому моменту предоставлено основным процессором.

В определенных случаях необходимо согласовывать работу обоих устройств. К примеру, если во входном потоке сразу за командой сопроцессора следует команда основного процессора, использующая результаты работы предыдущей команды, то сопроцессор не успеет выполнить свою команду за то время, пока основной процессор, пропустив сопроцессорную команду, выполнит свою. При этом что логика работы программы будет нарушена. Возможна и другая ситуация. Если входной поток команд содержит последовательность из нескольких команд сопроцессора, то процессор пропустит их очень быстро, но он должен обеспечить внешний интерфейс для сопроцессора. Эти и другие, более сложные ситуации, приводят к необходимости синхронизации между собой работы двух процессоров. В первых моделях микропроцессоров это делалось путем вставки перед или после каждой команды сопроцессора специальной команды wait или fwait . Работа данной команды заключалась в приостановке работы основного процессора до тех пор, пока сопроцессор не закончит работу над последней командой. В моделях микропроцессора (начиная с i486) подобная синхронизация выполняется автоматически. Но для некоторых команд из группы команд управления сопроцессором оставлена возможность выбора между командами с синхронизацией (ожиданием) и без нее.

Цель статьи - дать читателю представление о работе с числами с плавающей точкой в программировании. Предполагается дать базовый объём информации и предоставить уровень понимания, достаточный для самостоятельного проведения оценок в задачах, возникающих на практике.

Основная причина создания статьи - наблюдаемый недостаток понимания чисел с плавающей точкой среди программистов - как опасностей так и гарантий. В отсутствие понимания гарантий борьба с опасностями ведётся методами, которые часто можно охарактеризовать как "шаманские".

Второй важной причиной стало отсутствие (насколько мне известно) на русском языке статьи по данному вопросу, сравнимой по полноте изложения с классической статьёй Голдберга. Ближайшая известная альтернатива - всё же намного менее подробна.

Статья может быть использована как учебник или как справочник ("или" - не исключающее).

Статья получилась достаточно объёмной, причём материал структурирован в логической последовательности (для облегчения использования статьи в качестве справочника), а не "от простого к сложному". В связи с этим приношу свои извинения новичкам, желающим краткого введения в работу с числами с плавающей точкой.

В качестве языков примеров используются C/C++.

Все отрывки кода в статье, если явно не оговорено обратное, находятся в общем доступе (public domain).

В статье везде используется термин "числа с плавающей точкой", вместо распространённого в русскоязычной математической литературе "числа с плавающей запятой" (а также - соответствующая типография).

74-12-25 | Детали плавающей точки

Рис. 1. Hydrodynamic number-crunching, автор: Guy L. Steele, Jr.

Вычисления с вещественными числами были одним из основных применений компьютеров с момента их появления. Уже механическая машина Z1, построенная в 1938 г. Конрадом Цузе проводила операции с числами с плавающей точкой в формате, удивительно напоминающем современный.

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

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

Процесс стандартизации занял 9 лет (1976-1985), окончательный текст был основан на "K-C-S proposal" (авторы: William Kahan, Jerome Coonen, Harold Stone).

Стандарт предоставил единый для всех систем набор форматов и удивительно высокие гарантии (с т. з. разработчиков реализаций - требования) точности.

На данный момент большинство процессоров, вообще поддерживающих вычисления с плавающей точкой, поддерживают именно IEEE 754.

Таким образом, эту часть изначальной задачи стандарта во многом можно считать выполненной.

В 1989 г. Уильям Кэхэн получил премию Тьюринга с формулировкой "for his fundamental contributions to numerical analysis" (рус. "за его фундаментальный вклад в вычислительную математику").

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

Куда уплывает точка

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

Ман­тисса — это, по сути, чис­ло, записан­ное без точ­ки. Экспо­нен­та — это сте­пень, в которую нуж­но воз­вести некое чис­ло N (как пра­вило, N = 2), что­бы при перем­ножении на ман­тиссу получить иско­мое чис­ло (с точ­ностью до раз­ряднос­ти ман­тиссы). Выг­лядит это при­мер­но так:

где m и e — целые чис­ла, записан­ные в бинар­ном виде в выделен­ных под них битах. Что­бы избе­жать неод­нознач­ности, счи­тает­ся, что 1 <= |m| < N, то есть чис­ло записа­но в том виде, как если бы оно было с одним зна­ком перед запятой, но запятую злос­тно стер­ли, и чис­ло прев­ратилось в целое.

Ман­тисса — это, по сути, чис­ло, записан­ное без точ­ки. Экспо­нен­та — это сте­пень, в которую нуж­но воз­вести некое чис­ло N (как пра­вило, N = 2), что­бы при перем­ножении на ман­тиссу получить иско­мое чис­ло.

Сов­ремен­ные вещес­твен­ные чис­ла, под­держан­ные аппа­рат­но на уров­не про­цес­сора, так­же раз­биты на ман­тиссу и экспо­нен­ту. Разуме­ется, все опе­рации, при­выч­ные для ариф­метики целых чисел, так­же под­держа­ны коман­дами про­цес­сора для вещес­твен­ных чисел и выпол­няют­ся мак­сималь­но быс­тро.

Все так неп­росто потому, что такой фор­мат записи, во‑пер­вых, поз­воля­ет про­изво­дить опе­рации умно­жения и деления с такими чис­лами дос­таточ­но эффектив­но, кро­ме того, получить исходное вещес­твен­ное чис­ло, пред­став­ленное таким фор­матом, так­же нес­ложно. Дан­ное пред­став­ление чисел называ­ется чис­лом с пла­вающей точ­кой.

Стандарт точечного плаванья

Ве­щес­твен­ные чис­ла с пла­вающей точ­кой, под­держан­ные на уров­не про­цес­сора, опи­саны спе­циаль­ным меж­дународ­ным стан­дартом IEEE 754. Основны­ми дву­мя типами для любых вычис­лений явля­ются single-precision (оди­нар­ной точ­ности) и double-precision (двой­ной точ­ности) floating-point (чис­ла с пла­вающей точ­кой). Наз­вания эти нап­рямую отра­жают раз­рядность бинар­ного пред­став­ления чисел оди­нар­ной и двой­ной точ­ности: под пред­став­ление с оди­нар­ной точ­ностью выделе­но 32 бита, а под двой­ную, как ни стран­но, 64 бита — ров­но вдвое боль­ше.

Кро­ме оди­нар­ной и двой­ной точ­ности, в новой редак­ции стан­дарта IEEE 754—2008 пре­дус­мотре­ны так­же типы рас­ширен­ной точ­ности, чет­верной и даже половин­ной точ­ности. Одна­ко в C/C++, кро­ме float и double, есть раз­ве что еще тип long double, упор­но не под­держи­ваемый ком­пани­ей Microsoft, которая в Visual C++ под­став­ляет вмес­то него обыч­ный double. Ядром про­цес­сора в нас­тоящий момент так­же, как пра­вило, пока не под­держи­вают­ся типы половин­ной и чет­верной точ­ности. Поэто­му, выбирая пред­став­ления с пла­вающей точ­кой, при­ходит­ся выбирать лишь из float и double.

В качес­тве осно­вания для экспо­нен­ты чис­ла по стан­дарту берет­ся 2, соот­ветс­твен­но, при­веден­ная выше фор­мула сво­дит­ся к сле­дующей:

Рас­клад в битах в чис­лах оди­нар­ной точ­ности выг­лядит так:

1 бит под знак 8 бит экспо­нен­ты 23 бита ман­тиссы

Для двой­ной точ­ности мы можем исполь­зовать боль­ше битов:

1 бит под знак 11 бит экспо­нен­ты 52 бита ман­тиссы

В обо­их слу­чаях если бит зна­ка равен 0, то чис­ло положи­тель­ное и 1 уста­нав­лива­ется для отри­цатель­ных чисел. Это пра­вило ана­логич­но целым чис­лам с той лишь раз­ницей, что в отли­чие от целых, что­бы получить чис­ло, обратное по сло­жению, дос­таточ­но инверти­ровать один бит зна­ка.

Пос­коль­ку ман­тисса записы­вает­ся в дво­ичном виде, под­разуме­вает­ся целая часть, уже рав­ная 1, поэто­му в записи ман­тиссы всег­да под­разуме­вает­ся один бит, который не хра­нит­ся в дво­ичной записи. В битах ман­тиссы хра­нит­ся имен­но дроб­ная часть нор­мализо­ван­ного чис­ла в дво­ичной записи.

Ман­тисса записы­вает­ся в дво­ичном виде, и отбра­сыва­ется целая часть, заведо­мо рав­ная 1, поэто­му никог­да не забыва­ем, что ман­тисса на один бит длин­нее, чем в она хра­нит­ся в дво­ичном виде.

Не нуж­но иметь док­тор­скую сте­пень, что­бы вычис­лить точ­ность в десятич­ных зна­ках чисел, которые мож­но пред­ста­вить этим стан­дартом: 2 23 + 1 = 16 777 216; это явно ука­зыва­ет нам на тот факт, что точ­ность пред­став­ления вещес­твен­ных чисел с оди­нар­ной точ­ностью дос­тига­ет чуть более семи десятич­ных зна­ков. Это зна­чит, что мы не смо­жем сох­ранить в дан­ном фор­мате, нап­ример, чис­ло 123 456,78 — неболь­шое, в общем‑то, чис­ло, но уже начиная с сотой доли мы получим не то чис­ло, что хотели. Ситу­ация усложня­ется тем, что для боль­ших чисел вида 1 234 567 890, которое прек­расно помеща­ется даже в 32-раз­рядное целое, мы получим пог­решность уже в сот­нях еди­ниц! Поэто­му, нап­ример, в C++ для вещес­твен­ных чисел по умол­чанию исполь­зует­ся тип double. Ман­тисса чис­ла с двой­ной точ­ностью уже пре­выша­ет 15 зна­ков: 2 52 = 4 503 599 627 370 496 и спо­кой­но вме­щает в себя все 32-раз­рядные целые, давая сбой толь­ко на дей­стви­тель­но боль­ших 64-раз­рядных целых (19 десятич­ных зна­ков), где пог­решность в сот­нях еди­ниц уже, как пра­вило, несущес­твен­на. Если же нуж­на боль­шая точ­ность, то мы в дан­ной статье обя­затель­но в этом поможем.

Те­перь что каса­ется экспо­нен­ты. Это обыч­ное бинар­ное пред­став­ление целого чис­ла, в которое нуж­но воз­вести 10, что­бы при перем­ножении на ман­тиссу в нор­мализо­ван­ном виде получить исходное чис­ло. Вот толь­ко в стан­дарте вдо­бавок вве­ли сме­щение, которое нуж­но вычитать из бинар­ного пред­став­ления, что­бы получить иско­мую сте­пень десят­ки (так называ­емая biased exponent — сме­щен­ная экспо­нен­та). Экспо­нен­та сме­щает­ся для упро­щения опе­рации срав­нения, то есть для оди­нар­ной точ­ности берет­ся зна­чение 127, а для двой­ной 1023. Все это зву­чит край­не слож­но, поэто­му мно­гие про­пус­кают гла­ву о типе с пла­вающей точ­кой. А зря!

Примерное плаванье

Что­бы ста­ло чуточ­ку понят­нее, рас­смот­рим при­мер. Закоди­руем чис­ло 640 (= 512 + 128) в бинар­ном виде как вещес­твен­ное чис­ло оди­нар­ной точ­ности:

  • чис­ло положи­тель­ное — бит зна­ка будет равен 0;
  • что­бы получить нор­мализо­ван­ную ман­тиссу, нам нуж­но поделить чис­ло на 512 — мак­сималь­ную сте­пень двой­ки, мень­шую чис­ла, получим 640 / 512 = 512 / 512 + 128 / 512 или 1 + 1/4, что дает в дво­ичной записи 1,01, соот­ветс­твен­но, в битах ман­тиссы будет 0100000 00000000 00000000;
  • что­бы получить из 1 + 1/4 сно­ва 640, нам нуж­но ука­зать экспо­нен­ту, рав­ную 9, как раз 2 9 = 512, чис­ло, на которое мы подели­ли чис­ло при нор­мализа­ции ман­тиссы, но в бинар­ном виде дол­жно быть пред­став­ление в сме­щен­ном виде, и для вещес­твен­ных чисел оди­нар­ной точ­ности нуж­но при­бавить 127, получим 127 + 9 = 128 + 8, что в бинар­ном виде будет записа­но так: 10001000.

Для двой­ной точ­ности будет поч­ти все то же самое, но ман­тисса будет содер­жать еще боль­ше нулей спра­ва в дроб­ной час­ти, а экспо­нен­та будет 1023 + 9 = 1024 + 8, то есть чуть боль­ше нулей меж­ду стар­шим битом и чис­лом экспо­нен­ты: 100 00001000.

В общем, все не так страш­но, если акку­рат­но разоб­рать­ся.

За­дание на дом: разоб­рать­ся в дво­ичной записи сле­дующих кон­стант: плюс и минус бес­конеч­ность (INF — бес­конеч­ность), ноль, минус ноль и чис­ло‑не‑чис­ло (NaN — not-a-number).

За буйки не заплывай!

Есть одно важ­ное пра­вило: у каж­дого фор­мата пред­став­ления чис­ла есть свои пре­делы, за которые зап­лывать нель­зя. При­чем обес­печивать невыход за эти пре­делы при­ходит­ся самому прог­раммис­ту, ведь поведе­ние прог­раммы на С/С++ — это сде­лать невоз­мутимое лицо при выдаче в качес­тве сло­жения двух боль­ших положи­тель­ных целых чисел малень­кое отри­цатель­ное. Но если для целых чисел нуж­но учи­тывать толь­ко мак­сималь­ное и минималь­ное зна­чение, то для вещес­твен­ных чисел в пред­став­лении с пла­вающей точ­кой сле­дует боль­ше вни­мания обра­щать не столь­ко на мак­сималь­ные зна­чения, сколь­ко на раз­рядность чис­ла. Бла­года­ря экспо­нен­те мак­сималь­ное чис­ло для пред­став­ления с пла­вающей точ­кой при двой­ной точ­ности пре­выша­ет 10 308 , даже экспо­нен­та оди­нар­ной точ­ности дает воз­можность кодиро­вать чис­ла свы­ше 10 38 . Если срав­нить с «жал­кими» 10 19 , мак­симумом для 64-бит­ных целых чисел, мож­но сде­лать вывод, что мак­сималь­ные и минималь­ные зна­чения вряд ли ког­да‑либо при­дет­ся учи­тывать, хотя и забывать про них не сто­ит.

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

Со­ответс­твен­но, любые мно­гомил­лиар­дные сум­мы будут давать зна­читель­ную пог­решность в дроб­ной час­ти. При боль­шой интенсив­ности обра­бот­ки таких чисел могут про­падать мил­лиар­ды евро, прос­то потому, что они «не помес­тились», а пог­решность дроб­ной час­ти сум­мирова­лась и накопи­ла огромный оста­ток неуч­тенных дан­ных.

Что же делать нам, прог­раммис­там на C++, если перед нами сто­ит задача обра­ботать чис­ла очень боль­шой раз­ряднос­ти, при этом не исполь­зуя высоко­уров­невые язы­ки прог­рамми­рова­ния? Да то же, что и обыч­но: запол­нить про­бел, соз­дав один неболь­шой тип дан­ных для работы с десятич­ными дро­бями высокой точ­ности, ана­логич­ный типам Decimal высоко­уров­невых биб­лиотек.

Добавим плавающей точке цемента

По­ра зафик­сировать пла­вающую точ­ку. Пос­коль­ку мы решили изба­вить­ся от типа с пла­вающей точ­кой из‑за проб­лем с точ­ностью вычис­лений, нам оста­ются целочис­ленные типы, а пос­коль­ку нам нуж­на мак­сималь­ная раз­рядность, то и целые нам нуж­ны мак­сималь­ной раз­ряднос­ти в 64 бита.

Се­год­ня в учеб­ных целях мы рас­смот­рим, как соз­дать пред­став­ление вещес­твен­ных чисел с гаран­тирован­ной точ­ностью до 18 зна­ков пос­ле точ­ки. Это дос­тига­ется прос­тым ком­биниро­вани­ем двух 64-раз­рядных целых для целой и дроб­ной час­ти соот­ветс­твен­но. В прин­ципе, ник­то не меша­ет вмес­то одно­го чис­ла для каж­дой из ком­понент взять мас­сив зна­чений и получить пол­ноцен­ную «длин­ную» ариф­метику. Но будет более чем дос­таточ­но сей­час решить проб­лему точ­ности, дав воз­можность работать с точ­ностью по 18 зна­ков до и пос­ле запятой, зафик­сировав точ­ку меж­ду дву­мя эти­ми зна­чени­ями и залив ее цемен­том.

Отсыпь и мне децимала!

Сна­чала нем­ного теории. Обоз­начим наше две ком­понен­ты, целую и дроб­ную часть чис­ла, как n и f, а само чис­ло будет пред­ста­вимо в виде

x = n + f * 10 -18 , где n, f — целые, 0 <= f < 10 18 .

Для целой час­ти луч­ше все­го подой­дет зна­ковый тип 64-бит­ного целого, а для дроб­ной — без­зна­ковый, это упростит мно­гие опе­рации в даль­нейшем.

Такой метод является компромиссом между точностью и диапазоном представляемых значений. Представление чисел с плавающей точкой рассмотрим на примере чисел двойной точности (double precision). Такие числа занимают в памяти два машинных слова (8 байт на 32-битных системах). Наиболее распространенное представление описано в стандарте IEEE 754.

Кроме чисел двойной точности также используются следующие форматы чисел:

  • половинной точности (half precision) (16 бит),
  • одинарной точности (single precision) (32 бита),
  • четверной точности (quadruple precision) (128 бит),
  • расширенной точности (extended precision) (80 бит).

При выборе формата программисты идут на разумный компромисс между точностью вычислений и размером числа.

Недостатком такой записи является тот факт, что числа нельзя записать однозначно: [math] 0.01 = 0.001 \times 10^1 [/math] .

Число с плавающей точкой хранится в нормализованной форме и состоит из трех частей (в скобках указано количество бит, отводимых на каждую секцию в формате double):

  1. знак
  2. экспонента (показатель степени) (в виде целого числа в коде со сдвигом)
  3. мантисса (в нормализованной форме)

В качестве базы (основания степени) используется число [math] 2 [/math] . Экспонента хранится со сдвигом [math] -1023 [/math] .

Итоговое значение числа вычисляется по формуле:
[math] x = (-1)^ \times (1.mant) \times 2^ [/math] .
  1. В нормализованном виде любое отличное от нуля число представимо в единственном виде. Недостатком такой записи является тот факт, что невозможно представить число 0.
  2. Так как старший бит двоичного числа, записанного в нормализованной форме, всегда равен 1, его можно опустить. Это используется в стандарте IEEE 754.
  3. В отличие от целочисленных стандартов (например, integer), имеющих равномерное распределение на всем множестве значений, числа с плавающей точкой (double, например) имеют квазиравномерное распределение.
  4. Вследствие свойства 3, числа с плавающей точкой имеют постоянную относительную погрешность (в отличие от целочисленных, которые имеют постоянную абсолютную погрешность).
  5. Очевидно, не все действительные числа возможно представить в виде числа с плавающей точкой.
  6. Точно в таком формате представимы только числа, являющиеся суммой некоторых обратных степеней двойки (не ниже -53). Остальные числа попадают в некоторый диапазон и округляются до ближайшей его границы. Таким образом, абсолютная погрешность составляет половину величины младшего бита.
  7. В формате double представимы числа в диапазоне [math] [2.3 \times 10^, 1.7 \times 10^] [/math] .

В нормализованной форме невозможно представить ноль. Для его представления в стандарте зарезервированы специальные значения мантиссы и экспоненты.

Согласно стандарту выполняются следующие свойства:

  • [math] +0 = -0 [/math]
  • [math] \frac< \left| x \right| >= -0\,\![/math] (если [math]x\ne0[/math] )
  • [math] (-0) \cdot (-0) = +0\,\![/math]
  • [math] \left| x \right| \cdot (-0) = -0\,\![/math]
  • [math] x + (\pm 0) = x\,\![/math]
  • [math] (-0) + (-0) = -0\,\![/math]
  • [math] (+0) + (+0) = +0\,\![/math]
  • [math] \frac<-\infty>= +0\,\![/math]
  • [math] \frac<\left|x\right|>= -\infty\,\![/math] (если [math]x\ne0[/math] )

Для приближения ответа к правильному при переполнении, в double можно записать бесконечное значение. Так же, как и в случае с нолем, для этого используются специальные значение мантиссы и экспоненты.

Бесконечное значение можно получить при переполнении или при делении ненулевого числа на ноль.

В математике встречается понятие неопределенности. В стандарте double предусмотрено псевдочисло, которое арифметическая операция может вернуть даже в случае ошибки.

Неопределенность можно получить в нескольких случаях. Приведем некоторые из них:

  • [math] f(NaN) = NaN [/math] , где [math] f [/math] - любая арифметическая операция
  • [math] \infty + (-\infty) = NaN [/math]
  • [math] 0 \times \infty = NaN [/math]
  • [math] \frac<\pm0><\pm0>= \frac<\pm \infty><\pm \infty>= NaN [/math]
  • [math] \sqrt = NaN [/math] , где [math] x \lt 0 [/math]

Денормализованные (denormalized numbers) - способ увеличить количество представимых числе в окрестности нуля. Каждое такое число по модулю меньше самого маленького нормализованного.< Согласно стандарту числа с плавающей точкой можно представить в следующем виде:

  • [math] (-1)^s \times 1.M \times 2^E [/math] , в нормализованном виде если [math] E_ \leq E \leq E_ [/math] ,
  • [math] (-1)^s \times 0.M \times 2^E_ [/math] , в денормализованном виде если [math] E = E_ - 1 [/math] ,

где [math] E_ [/math] - минимальное значение порядка, используемое для записи чисел (единичный сдвиг), [math] E_ - 1 [/math] - минимальное значение порядка, которое он может принимать - все биты нули, нулевой сдвиг.

Ввиду сложности, денормализованные числа обычно реализуют на программном уровне, а не на аппаратном. Из-за этого резко возрастает время работы с ними. Это недопустимо в областях, где требуется большая скорость вычислений (например, видеокарты). Так как денормализованные числа представляют числа мало отличные от нуля и мало влияют на результат, зачастую они игнорируются (что резко повышает скорость). При этом используются две концепции:

  1. Flush To Zero (FTZ) - в качестве результата возвращается нуль, как только становится понятно, что результат будет представляться в денормализованном виде.
  2. Denormals Are Zero (DAZ) - денормализованные числа, поступающие на вход, рассматриваются как нули.

Начиная с версии стандарта IEEE 754 2008 года денормализованные числа называются "субнормальными" (subnormal numbers), то есть числа, меньшие "нормальных".

Таким образом, компьютер не различает числа [math] x [/math] и [math] y [/math] , если [math] 1 \lt \frac \lt 1 + \varepsilon_m [/math] . Из свойств чисел двойной точности следует, что для них [math] \varepsilon_m = 2^<-54>[/math] .

Мера единичной точности используется для оценки точности вычислений.


Приведем пример кода на Python, который показывает, при каком значении числа [math] x [/math] компьютер не различает числа [math] x [/math] и [math] x + 1 [/math] .

То есть [math] x = 2^ [/math] , так как мантисса числа двойной точности содержит 53 бита (в памяти хранятся 52). В C++ для расчета расстояния между двумя числами двойной точности можно воспользоваться функцией [math] \mathrm [/math] .

  • [math] a \oplus b = (a + b) (1 + \delta), |\delta| \leq \varepsilon_m [/math] ,
  • [math] a \ominus b = (a - b) (1 + \delta), |\delta| \leq \varepsilon_m [/math] ,
  • [math] a \otimes b = ab (1 + \delta), |\delta| \leq \varepsilon_m [/math] .
[math] \forall a, b, c \in D^2, \tilde = (b_x \ominus a_x) \otimes (c_y \ominus a_y) \ominus (b_y \ominus a_y) \otimes (c_x \ominus a_x) [/math]

[math] \exists \tilde \in D: [/math]

  1. [math] \tilde \gt \tilde \Rightarrow (b - a) \times (c - a) \gt 0 [/math]
  2. [math] \tilde \lt -\tilde \Rightarrow (b - a) \times (c - a) \lt 0 [/math]

Обозначим [math] v = (b - a) \times (c - a) = (b_x - a_x) (c_y - a_y) - (b_y - a_y) (c_x - a_x)[/math] .

Теперь распишем это выражение в дабловой арифметике.

[math]\tilde = (b_x \ominus a_x) \otimes (c_y \ominus a_y) \ominus (b_y \ominus a_y) \otimes (c_x \ominus a_x) = \\ = [ (b_x - a_x) (c_y - a_y) (1 + \delta_1) (1 + \delta_2) (1 + \delta_3) - \\ - (b_y - a_y) (c_x - a_x) (1 + \delta_4) (1 + \delta_5) (1 + \delta_6) ] (1 + \delta_7),[/math]

[math] |\delta_i| \leq \varepsilon_m [/math]

Заметим, что [math] v \approx \tilde [/math]

Теперь оценим абсолютную погрешность [math] \epsilon = |v - \tilde|. [/math]

[math] |v - \tilde| = |(b_x - a_x) (c_y - a_y) - (b_y - a_y) (c_x - a_x) - \\ - (b_x - a_x) (c_y - a_y) (1 + \delta_1) (1 + \delta_2) (1 + \delta_3) (1 + \delta_7) + \\ + (b_y - a_y) (c_x - a_x) (1 + \delta_4) (1 + \delta_5) (1 + \delta_6) (1 + \delta_7)| = \\ = |(b_x - a_x) (c_y - a_y) (1 - (1 + \delta_1) (1 + \delta_2) (1 + \delta_3) (1 + \delta_7)) - \\ - (b_y - a_y) (c_x - a_x) (1 - (1 + \delta_4) (1 + \delta_5) (1 + \delta_6) (1 + \delta_7))| \leq \\ \leq |(b_x - a_x) (c_y - a_y) (1 - (1 + \delta_1) (1 + \delta_2) (1 + \delta_3) (1 + \delta_7))| + \\ + |(b_y - a_y) (c_x - a_x) (1 - (1 + \delta_4) (1 + \delta_5) (1 + \delta_6) (1 + \delta_7))| = \\ = |(b_x - a_x) (c_y - a_y)| \cdot |((1 + \delta_1) (1 + \delta_2) (1 + \delta_3) (1 + \delta_7) - 1)| + \\ + |(b_y - a_y) (c_x - a_x)| \cdot |((1 + \delta_4) (1 + \delta_5) (1 + \delta_6) (1 + \delta_7) - 1)| = \\ = |(b_x - a_x) (c_y - a_y)| \cdot |\delta_1 + \delta_2 + \delta_3 + \delta_7 + \delta_1 \delta_2 \ldots| + \\ + |(b_y - a_y) (c_x - a_x)| \cdot |\delta_4 + \delta_5 + \delta_6 + \delta_7 + \delta_4 \delta_5 \ldots| \leq \\ \leq |(b_x - a_x) (c_y - a_y)| \cdot (|\delta_1| + |\delta_2| + |\delta_3| + |\delta_7| + |\delta_1 \delta_2| \ldots) + \\ + |(b_y - a_y) (c_x - a_x)| \cdot (|\delta_4| + |\delta_5| + |\delta_6| + |\delta_7| + |\delta_4 \delta_5| \ldots) \leq \\ \leq |(b_x - a_x) (c_y - a_y)| \cdot (4 \varepsilon_m + 6 \varepsilon_m^2 + 4 \varepsilon_m^3 + \varepsilon_m^4) + \\ + |(b_y - a_y) (c_x - a_x)| \cdot (4 \varepsilon_m + 6 \varepsilon_m^2 + 4 \varepsilon_m^3 + \varepsilon_m^4) = \\ = (|(b_x - a_x) (c_y - a_y)| + |(b_y - a_y) (c_x - a_x)|)(4 \varepsilon_m + 6 \varepsilon_m^2 + 4 \varepsilon_m^3 + \varepsilon_m^4)[/math]

Пусть [math] t = (|(b_x - a_x) (c_y - a_y)| + |(b_y - a_y) (c_x - a_x)|).[/math] Получаем, что

[math] \epsilon = |v - \tilde| \leq t \cdot (4 \varepsilon_m + 6 \varepsilon_m^2 + 4 \varepsilon_m^3 + \varepsilon_m^4). [/math]

[math]\tilde = (|(b_x - a_x) (c_y - a_y) (1 + \delta_1) (1 + \delta_2) (1 + \delta_3)| + \\ + |(b_y - a_y) (c_x - a_x) (1 + \delta_4) (1 + \delta_5) (1 + \delta_6)|) (1 + \delta_7) \geq \\ \geq |(b_x - a_x) (c_y - a_y) (1 - \varepsilon_m)^3)|(1 - \varepsilon_m) + \\ + |(b_y - a_y) (c_x - a_x) (1 - \varepsilon_m)^3)|(1 - \varepsilon_m) = \\ = |(b_x - a_x) (c_y - a_y)| (1 - \varepsilon_m)^4 + |(b_y - a_y) (c_x - a_x)| (1 - \varepsilon_m)^4 = \\ = (|(b_x - a_x) (c_y - a_y)| + |(b_y - a_y) (c_x - a_x)|) (1 - \varepsilon_m)^4 = t \cdot (1 - \varepsilon_m)^4[/math]

[math] t \leq \tilde \frac <(1 - \varepsilon_m)^4>= \tilde (1 + 4 \varepsilon_m + 10 \varepsilon_m^2 + 20 \varepsilon_m^3 + \cdots) [/math]

[math] \epsilon = |v - \tilde| \leq \tilde \leq \tilde (1 + 4 \varepsilon_m + 10 \varepsilon_m^2 + 20 \varepsilon_m^3 + \cdots) (4 \varepsilon_m + 6 \varepsilon_m^2 + 4 \varepsilon_m^3 + \varepsilon_m^4) [/math]

[math] \tilde \lt 8 \varepsilon_m \tilde [/math]

Заметим, что это довольно грубая оценка. Вполне можно было бы написать [math] \tilde \lt 4.25 \varepsilon_m \tilde [/math] или [math] \tilde \lt 4.5 \varepsilon_m \tilde.[/math]

Читайте также: