Книга: О чём не пишут в книгах по Delphi
Назад: 3.1. Неочевидные особенности целых чисел
Дальше: 3.3. Тонкости работы со строками
3.2. Неочевидные особенности вещественных чисел
Если рассмотренные в предыдущих разделах особенности целых чисел могли быть неочевидными только начинающим, то вещественные числа могут преподнести сюрпризы даже достаточно опытным программистам, т.к. их поведение существенно дальше от интуитивных представлений, и эти неожиданности не ограничиваются выходом та пределы диапазона. Существующая литература по Delphi, в основном, считает этот вопрос несущественным и обходит его стороной, в результате чего программист, впервые столкнувшийся с одним из таких сюрпризов, впадает в недоумение и испытывает желание "попрыгать вокруг компьютера с бубном". Здесь мы попытаемся восполнить этот пробел и показать, что необъяснимые на первый взгляд явления на самом деле просты и предсказуемы, если известно, как реализуется вещественная арифметика компьютером.
3.2.1. Двоичные дроби
Для начала — немного математики. В школе мы проходим два вида дробей простые и десятичные. Десятичные дроби, по сути дела, представляют собой разложение числа по степеням десяти. Так, запись 13,6704 означает число, равное 1·101 + 3·100 + 6·10-1 + 7·10-2 + 0·10-3 + 4·10-4. Но внутреннее представление всех чисел в компьютере, в том числе и вещественных, не десятичное, а двоичное. Поэтому он использует двоичные дроби. Они во многом похожи на десятичные, но основанием степени у них служит двойка. Так, двоичная дробь 101.1101 — это 1·22 + 0·21 + 1·20 + 1·2-1 + 1·2-2 + 0·2-3 + 1·2-4. В десятичном представлении это число равно 5,8125, в чем нетрудно убедиться с помощью любого калькулятора.
Теперь вспомним научный формат записи десятичного числа. Первым в этой записи идет знак числа (плюс или минус). Дальше идет так называемая мантисса (число от 1 до 10). Затем идет экспонента (степень десяти, на которую надо умножить мантиссу, чтобы получить нужное число). Итак, уже упоминавшееся число 13,6704 запишется в этом формате как 1.36704·101 (или 1.36704E1 по принятым в компьютере правилам). Если записываемое число меньше единицы, экспонента будет отрицательной. Аналогичная запись существует и в двоичной системе. Так, 101.1101 запишется в виде 1.011101*1010 (везде использована двоичная форма записи, так что 1010 означает 22). Именно такое представление реализовано в компьютере. Двоичная точка в такой записи не остается на одном месте, а сдвигается на величину, указанную в экспоненте, поэтому такие числа называются числами с плавающей точкой (floating point numbers).
3.2.2. Вещественные типы Delphi
В Delphi существует четыре вещественных типа: Single, Double, Extended и Real. Их общий формат одинаков (рис. 3.1, а).
Знак — это всегда один бит. Он равен нулю для положительных чисел и единице для отрицательных. Что же касается размеров мантиссы и экспоненты, то именно в них и заключается различие между типами.
Прежде чем перейти к конкретным цифрам, рассмотрим подробнее тип Real, сделав для этого небольшой экскурс в историю. Real — это стандартный тип языка Паскаль, присутствовавший там изначально. Когда создавался Паскаль, процессоры еще не имели встроенной поддержки вещественных чисел, поэтому все операции с данными типа Real сводились к операциям с целыми числами. Соответственно, размер полей в типе Real был подобран так, чтобы оптимизировать эти операции.
а) общий вид вещественного числа
б) Двоичное представление числа типа Single
Рис. 3.1. Хранение вещественного числа в памяти
Микропроцессор Intel 8086/88 и его улучшенные варианты — 80286 и 80386 — также не имели аппаратной поддержки вещественных чисел. Но у систем на базе этих процессоров была возможность подключения так называемого сопроцессора. Эта микросхема работала с памятью через шины основного процессора и обеспечивала аппаратную поддержку вещественных чисел. В системах средней руки гнездо сопроцессора обычно было пустым, т.к. это уменьшало цену (разумеется, вставить туда сопроцессор не было проблемой). Для каждого центрального процессора выпускались свои сопроцессоры, маркировавшиеся Intel 8087, 80287 и 80387 соответственно. Были даже сопроцессоры, выпускаемые другими фирмами. Они работали быстрее, чем сопроцессоры Intel, но появлялись на рынке позже. Тип вещественных чисел, поддерживаемый сопроцессорами, не совпадает с Real. Он определяется стандартом IEEE (Institute of Electrical and Electronics Engineers).
Чтобы обеспечить в своих системах поддержку типов IEEE, Borland вводит в Turbo Pascal типы Single, Double и Extended. Extended — это основной для сопроцессора тип, a Single и Double получаются из него очень простым усечением. Система команд сопроцессора допускает работу с этими типами: при загрузке числа типа Single или Double во внутренний регистр сопроцессора последний конвертирует их в Extended. Напротив, при выгрузке чисел этих типов из регистра в память сопроцессор усекает их до нужного размера. Внутренние же операции всегда выполняются с данными типа Extended (впрочем, из этого правила есть исключение, на котором мы остановимся позже, после детального рассмотрения формата различных типов). Single и Double позволяют экономить память. Ни один из них также не совпадает с типом Real. В системах с сопроцессорами новые типы обрабатываются заметно (в 2–3 раза) быстрее, чем Real (это с учетом того, что тип Real после соответствующего преобразования также обрабатывался сопроцессором; если же сравнивать обработку типа Extended на машине с сопроцессором и Real на машине без сопроцессора, то там на отдельных операциях достигалась разница в скорости примерно в 100 раз). Чтобы программы с этими типами можно было выполнять и в системах без сопроцессора, была предусмотрена возможность подключать к ним программный эмулятор сопроцессора. Обработка этих типов эмулятором была медленнее, чем обработка Real.
Начиная с 486-й серии Intel берет курс на интеграцию процессора и сопроцессора в одной микросхеме. Процент брака в микросхемах слишком велик, поэтому Intel идет на хитрость: если у микросхемы брак только в сопроцессорной части, то на этом кристалле прожигаются перемычки, блокирующие сопроцессор, и микросхема продается как процессор 80486SX, не имеющий встроенного сопроцессора (в отличие от полноценной версии, которую назвали 80486DX). Бывали и обратные ситуации, когда сопроцессор повреждений не имел, зато процессор был неработоспособен. Такие микросхемы превращали в "сопроцессор 80487". Но это уже из области экзотики, и, по имеющейся у нас информации, до России такой сопроцессор не дошел.
Процессор Pentium во всех своих вариантах имел встроенный блок вычислений с плавающей точкой (FPU — Floating Point Unit), и отдельный сопроцессор ему не требовался. Таким образом, с приходом этого процессора тип Real остался только для обратной совместимости, а на передний план вышли типы Single, Double и Extended. Начиная с Delphi 4, тип Real становится синонимом типа Double, а старый 6-байтный тип получает название Real48.
Здесь и далее под словом Real мы будем понимать старый 6-байтный тип.
Примечание
Существует директива компилятора {$REALCOMPATIBILITY ON/OFF}, при включении которой (по умолчанию она отключена) Real становится синонимом Real48, а не Double.
Размеры полей для различных вещественных типов указаны в табл. 3.1.
Таблица 3.1. Размеры полей в вещественных типах
| Single | 4 | 23 | 8 |
| Double | 8 | 52 | 11 |
| Extended | 10 | 64 | 15 |
| Real | 6 | 40 | 7 |
Другие параметры вещественных типов, такие как диапазон и точность, можно найти в справке Delphi.
3.2.3. Внутренний формат вещественных чисел
Рассмотрим тип Single, т.к. он самый короткий и, следовательно, самый простой для понимания. Остальные типы отличаются от него только количественно. В дальнейшем числа в формате Single мы будем записывать как s eeeeeeee mmmmmmmmmmmmmmmmmmmmmmm, где s означает знаковый бит, е — бит экспоненты, m — бит мантиссы. Порядок хранения битов в типе Single показан на рис. 3.1, б (по принятым в процессорах Intel правилам байты в многобайтных значениях переставляются так. что младший байт идет первым, а старший — последним, и вещественных чисел это тоже касается В мантиссе хранится двоичное число. Чтобы получить истинное значение мантиссы, к ней надо мысленно добавить слева единицу с точкой (т.е., например, мантисса 1010000000000000000000 означает двоичную дробь 1.101). Таким образом, имея 23 двоичных разряда, мы записываем числа с точностью до 24-х двоичных разрядов.
Экспонента — по определению всегда целое число. Но способ записи экспоненты в вещественных числах не совпадает с рассмотренным ранее способом записи чисел со знаком. Ноль в этом представлении записывается как 01111111 (в обычном представлении это равно 127). Соответственно. 10000000 (128 в обычном представлении) означает единицу, а 01111110 (126) означает -1, и т. д. (т.е. из обычного беззнакового числа надо вычесть 127, и получится число, закодированное в экспоненте). Такая запись чиста называется нормализованной.
Из описанных правил есть исключения. Так, если все биты экспоненты равны нулю (т.е. там стоит число -127), то к мантиссе перед ее началом надо добавлять не "1.", а "0." (денормализованная запись). Это позволяет увеличить диапазон вещественных чисел. Если бы этого исключения не было, то минимально возможное положительное число типа Single было бы равно примерно 5,9·10-39. А так появляется возможность использовать числа до 1,4·10-45. Побочным эффектом этого является то, что числа, меньшие чем 1,17·10-38, представляются с меньшей, чем 24 двоичных разряда, точностью. Если все биты в экспоненте равны единице, а в мантиссе — нулю, то мы получаем комбинацию, известную как INF (от англ. Infinity — бесконечность). Эта комбинация используется тогда, когда результат вычислений превышает максимально допустимое форматом число. В зависимости от значения бита s бесконечность может быть положительной или отрицательной. Если же при такой экспоненте в мантиссе хоть один бит не равен нулю, такая комбинация называется NAN (Not A Number — не число). Попытки работы с комбинациями NAN или INF приводят к ошибке времени выполнения.
Для задания нуля все биты мантиссы и экспоненты должны быть равны нулю (формально это означает 0·10-127). С учетом описанных правил, если хотя бы один бит экспоненты не будет равен нулю (т.е. экспонента будет больше -127), запись будет считаться нормализованной, и нулевая мантисса будет рассматриваться как единица. Поэтому никакие другие комбинации значений мантиссы и экспоненты не могут дать ноль.
Тип Double устроен точно так же, разница только в количестве разрядов и в том, какое значение экспоненты берется за ноль. Итак, мы имеем 11 разрядов для экспоненты. За ноль берется значение 1023.
Несколько иначе устроен Extended. Кроме количественных отличий добавляется еще и одно качественное: в мантиссе явно указывается первый разряд. Это означает, что мантисса 1010… интерпретируется как 1.01, а не как 1.101, как это было в типах Single и Float. Поэтому если 23-битная мантисса типа Single обеспечивает 24-знаковую точность, а 52-битная мантисса Double — 53-битную, то 64-битная мантисса Extended обеспечивает 64-х, а не 65-битную точность. Соответственно, при денормализованной форме записи первый разряд мантиссы явно содержит 0. За ноль экспоненты принимается значение 16 383.
Тип Real, как уже упоминалось, стоит особняком. Во-первых, в нем биты следуют в другом порядке, а во-вторых, нет денормализованной формы. Мы не будем касаться внутреннего устройства типа Real, т.к. эта информация уже перестала быть актуальной.
3.2.4. "Неполноценный" Extended
Ранее мы отметили, что FPU всегда выполняет все операции в формате Extended, оговорившись при этом, что есть исключение из этого правила. Здесь мы рассмотрим это исключение.
У FPU существует специальный двухбайтный регистр, называемый управляющим словом. Установка отдельных битов этого регистра диктует то или иное поведение при выполнении операций. Прежде всего, это связано с тем, какие исключения может возбуждать FPU. Другие биты этого регистра отвечают за то, как будут округляться числа, как FPU понимает бесконечность, — всё это можно при необходимости узнать из документации Intel. Нас же будут интересовать только два бита из этого слова: восьмой и девятый. Именно они определяют, как будут обрабатываться числа внутри сопроцессора.
Если восьмой бит содержит единицу (так установлено по умолчанию), то десять байтов внутренних регистров сопроцессора будут задействованы полностью, и мы получим "полноценный" Extended. Если же этот бит равен нулю, то все определяется значением бита 9. Если он равен единице, то используется только 53 разряда мантиссы (остальные всегда равны нулю). Если же этот бит равен нулю — только 24 разряда мантиссы. Это увеличивает скорость вычислений, но уменьшает точность. Другими словами, точность работы сопроцессора может быть понижена до типа Double или даже Single. Но это касается только мантиссы, экспонента в любом случае будет содержать 15 бит, так что диапазон типа Extended сохраняется в любом случае.
Для работы с управляющим словом сопроцессора в модуле System описана переменная Default8087CW типа Word и процедура Set8087CW(CW: Word). При запуске программы в переменную Default8087CW записывается то управляющее слово, которое установила система при запуске программы. Функция Set8087CW одновременно записывает новое значение в управляющее слово и в переменную Default8087CW.
Такое поведение этой функции не всегда удобно — иногда бывает нужно сохранить старое значение переменной Default8087CW (впрочем, это несложно сделать, заведя дополнительную переменную). С другой стороны, если значение управляющею слова изменить, не используя Set8087CW (а в дальнейшем мы увидим, что такие изменения могут происходить помимо нашей воли), то с помощью функции Default8087CW просто нет возможности узнать текущее значение управляющего слова. В Delphi 6 и выше появилась функция Get8087CW, позволяющая узнать значение именно контрольного слова, а не переменной Default8087CW. В более ранних версиях существовал единственный способ получить значение этого слова — встроенный в Delphi ассемблер.
Итак, установить значение управляющего слова можно с помощью команды FLDCW, прочитать с помощью FNSTCW. Обе эти команды имеют один аргумент — переменную типа Word. Чтобы, например, установить 53-значную точность, не изменив при этом другие биты управляющего слова нужно выполнить такую последовательность команд:
asm
FNSTCW MyCW
AND MyCW, 0FC00h
OR MyCW, 200h
FLDCW MyCW
end;
Начиная с Delphi 6, в модуле Math появилась еще одна функция, позволяющая устанавливать точность FPU без манипуляции с отдельными битами управляющего слова — SetPrecisionMode. В зависимости от значения аргумента (pmSingle, pmDouble или pmExtended) она устанавливает требуемую точность. Современные сопроцессоры обрабатывают числа с такой скоростью, что при обычных вычислениях вряд ли может возникнуть необходимость в ускорении за счет точности — выигрыш будет ничтожен. Эта возможность необходима, в основном, в тех случаях, когда вычисления с плавающей точкой составляют значительную часть программы, а высокая точность не имеет принципиального значения (например, в 3D-играх). Однако забывать об этой особенности работы сопроцессора не следует, потому что она может преподнести один неприятный сюрприз, о котором чуть позже.
3.2.5. Бесконечные дроби
Из школы мы все помним, что не каждое число может быть записано конечной десятичной дробью. Бесконечные дроби бывают двух видов: периодичные и непериодичные. Примером непериодичной дроби является число π, периодичной — число ⅓ или любая другая простая дробь, не представимая в виде конечной десятичной дроби.
Примечание
Напомним, что периодичные дроби — это такие дроби которые содержат бесконечно повторяющуюся последовательность цифр. Например, 1/9=0,11111..., 1/12=0,08333333..., 1/7=0,142857142857... Такие числа записывают со скобками — в них заключают повторяющуюся часть. Те же числа должны быть записаны так: 1/9=0,1(1), 1/12=0,08(3), 1/7=0,1(428571)
Вопрос о периодичности или непериодичности числа нас сейчас не интересует, нам достаточно знать, что не все числа можно представить в виде конечной десятичной дроби. При работе с такими числами мы всегда имеем не точное, а приближенное значение, поэтому ответ получается тоже приближенным. Это нужно учитывать в своих расчетах.
До сих пор мы говорили только о десятичных бесконечных дробях. Но двоичные дроби тоже могут быть бесконечными. Даже более того, любое число, выражаемое конечной двоичной дробью, может быть также выражено и десятичной конечной дробью. Но существуют числа (например, 1/5), которые выражаются конечной десятичной дробью, но не могут быть выражены конечной двоичной дробью. Это и есть наиболее важное отличие аппаратной реализации вещественных чисел от наших интуитивных представлений. Теперь у нас достаточно теоретических знаний, чтобы перейти к рассмотрению конкретных примеров — "подводных камней", приготовленных вещественными числами.
3.2.6. "Неправильное" значение
Самый первый "подводный камень", на котором спотыкаются новички — это то, что вещественная переменная может получить не совсем то значение, которое ей присвоено. Рассмотрим это на простом примере (листинг 3.9, примеp WrongValue на компакт-диске).
Листинг 3.9. Пример присваивания "неправильного" вещественного значения
procedure TForm1.Button1Click(Sender: TObject);
var
R: Single;
begin
R := 0.1;
Label1.Caption = FloatToStr(F);
end;
Что мы увидим, когда нажмем кнопку? Разумеется, не 0.1, иначе не было бы смысла писать этот пример. Мы увидим 0.100000001490116, т.е. расхождение в девятой значащей цифре. Из справки по Delphi мы знаем, что точность типа Single — 7–8 десятичных разряда, так что нас, по крайнем мере, никто не обманывает. В чем же причина? Просто число 0,1 не представимо в виде конечной двоичной дроби, оно равно 0,0(0011). И эта бесконечная двоичная дробь обрубается на 24-х знаках; мы получаем не 0,1, а некоторое приближенное число (какое именно — см. выше). А если мы присвоим переменной R не 0.1, а 0.5? Тогда мы получим на экране 0.5, потому что 0.5 предоставляется в виде конечной двоичной дроби. Немного поэкспериментировав с различными числами, мы заметим, что точно представляются те числа, которые выражаются в виде m/2n, где m, n — некоторые целые числа (разумеется, n не должно превышать 24, а то нам не хватит точности типа Single). В качестве упражнения предлагаем доказать, что любое целое число, для записи которого хватает 24-х двоичных разряда, может быть точно передано типом Single.
Примечание
Если в этом примере изменить тип переменной R с Single на Double или на Extended, на экран будет выведено 0.1. Но это не значит, что в переменную будет записано ровно 0.1 — это просто особенности работы функции FloatToStr, которая не учитывает столь малую разницу между 0,1 и переданным ей числом.
3.2.7. Сравнение
Теперь попробуем сравнить значение переменной и константы, которую мы ей присвоили (листинг 3.10, пример Compare1 на компакт-диске).
Листинг 3.10. Пример ошибки при сравнении вещественной переменной и константы
procedure TForm1.Button1Click(Sender: TObject);
var
R: Single;
begin
R := 0.1;
if R = 0.1 then Label1.Caption := 'Равно'
else Label1.Caption := 'He равно';
end;
При нажатии кнопки мы увидим надпись Не равно. На первый взгляд это кажется абсурдом. Действительно, мы уже знаем, что переменная R получает значение 0.100000001490116 вместо 0.1. Но ведь "0.1" в правой части равенства тоже должно преобразоваться по тем же законам, т.к. работает аналогичный алгоритм. Тут самое время вспомнить, что FPU работает только с 10-байтным типом Extended, поэтому и левая, и правая часть равенства сначала преобразуется в этот тип, и лишь потом производится сравнение. То число, которое оказалось в переменной R вместо 0.1, хотя и выглядит страшно, но зато представляется в виде конечной двоичной дроби. Информация же о том, что это на самом деле должно означать "0.1", нигде не сохранилась. При преобразовании этого числа в Extended младшие, избыточные по сравнению с типом Single разряды мантиссы просто заполняются нулями, и мы снова получим то же самое число, только записанное в формате Extended. А "0.1" из правой части равенства преобразуется в Extended без промежуточного превращения в Single. Поэтому некоторые из младших разрядов мантиссы будут содержать единицы. Другими словами, мы получим хоть и не точное представление числа 0.1, но все же более близкое к истине, чем 0.100000001490116.
Из-за таких хитрых преобразований оказывается, что мы сравниваем два близких, но все же не равных числа. Отсюда — закономерный результат в виде надписи Не равно.
Тут уместна аналогия с десятичными дробями. Допустим, в одном случае мы делим 1 на три с точностью до трех знаков и получаем 0,333. Потом мы делим 1 на три с точностью до четырех знаков и получаем 0,3333. Теперь мы хотим сравнить эти два числа. Для этого приводим их к точности в четыре разряда. Получается, что мы сравниваем 0,3330 и 0,3333. Очевидно, что это разные числа.
Если попробовать заменить число 0,1 на 0,5, то мы увидим надпись Равно. Полагаем, что читатели уже догадались, почему, но все же приведем объяснение. Число 0,5 — это конечная двоичная дробь. При прямом приведении ее к типу Extended в младших разрядах оказываются нули. Точно такие же нули оказываются в этих разрядах при превращении числа 0,5 типа Single в тип Extended. Поэтому в результате мы сравниваем два равных числа. Это похоже на процесс деления 1 на 4 с точностью до трех и до четырех значащих цифр. В первом случае получили бы 0,250, во втором — 0,2500. Приведя оба значения к точности в четыре знака, получим сравнение 0,2500 и 0,2500. Очевидно, что эти числа равны.
3.2.8. Сравнение разных типов
Теперь попытаемся сравнить переменную не с константой, а с другой переменной (листинг 3.11, пример Compare2 на компакт-диске).
Листинг 3.11. Пример ошибки при сравнении переменных разных типов
procedure TForm1.Button1Click(Sender: TObject);
var
R1: Single;
R2: Double;
begin
R1 := 0.1;
R2 := 0.1;
if R1 = R2 then Label1.Caption := 'Равно'
else Label1.Caption := 'He равно';
end;
Почему этот пример также выдаст Не равно, понять проще, чем в предыдущем случае. При R1 бесконечная дробь обрывается на 24-х разрядах, а при R2 — на 53-х. Таким образом, в дополнительных по сравнению с типом Single разрядах переменной R2 будут единицы. При дополнении значений нулями до 10-байтной точности мы получим разные числа, что и определяет результат сравнения. Это напоминает ситуацию, когда мы сравниваем 0,333 и 0,3333, приводя их к точности в пять знаков: числа 0,33300 и 0,33330 не равны.
Как и в предыдущем случае, замена 0,1 на 0,5 даст результат Равно.
3.2.9. Вычитание в цикле
Рассмотрим еще один пример, иллюстрирующий ситуацию, которая часто озадачивает начинающего программиста (листинг 3.12, пример Subtraction на компакт-диске).
Листинг 3.12. Накапливание ошибки при вычитании
procedure TForm1.Button1Click(Sender: TObject);
var
R: Single;
I: Integer;
begin
R := 1;
for I := 1 to 10 do R := R - 0.1;
Label1.Caption := FloatToStr(R);
end;
В результате выполнения этого кода на экране появится -7.3015691270939E-8 вместо ожидаемого нуля. Объяснение этому достаточно очевидно, если вспомнить то, о чем мы говорили ранее. Число 0,1 не может быть передано точно ни в одном из вещественных типов, а при каждом вычислении происходит преобразование Single в Extended и обратно, причем последнее — с потерей точности. Эти потери приводят к тому, что мы получаем в результате не ноль, а "почти ноль".
3.2.10. Неожиданная потеря точности
Изменим в предыдущем примере тип переменной R с Single на Double. Значение, выводимое программой, станет 1.44327637948555E-16. Вполне логичный и предсказуемый результат, т.к. тип Double точнее, чем Single, и, следовательно, все вычисления имеют меньшую погрешность, мы просто обязаны получить более точный результат. Хотя, разумеется, абсолютная точность (т.е. ноль) для нас остается недостижимым идеалом.
А теперь — вопрос на засыпку. Изменится ли результат, если мы заменим Double на более точный Extended? Ответ не такой однозначный, каким его хотелось бы видеть. В принципе, после такой замены вы должны получить -6.7762635780344E-20. Но в некоторых случаях от замены Double на Extended результат не изменится, и вы снова получите 1.44327637948555Е-16. Это зависит от операционной системы и версии Delphi.
Все дело в использовании "неполноценного" Extended. При запуске программы любая система устанавливает такое управляющее слово FPU, чтобы Extended был полноценным. Но затем программа вызывает много разных функций Windows API. Какая-то (или какие-то) из этих многочисленных функций некорректно работают с управляющим словом, меняя его значение и не восстанавливая при выходе. Такая проблема встречается, в основном, в Windows 95 и старых версиях Windows 98. Также имеются сведения о том, что управляющее слово может "портиться" и в Windows NT, причем эффект наблюдался не сразу после установки системы, а лишь через некоторое время, после доустановки других программ. Проблема именно в некорректности поведения системных функций; значение управляющего слова, устанавливаемое системой при запуске программы, всегда одинаково. Таким образом, приходим к неутешительному выводу: к тем проблемам с вещественными числами, которые обусловлены особенностями их аппаратной реализации, добавляются еще и ошибки операционной системы. Правда, радует то, что в последнее время эти ошибки встречаются крайне редко — видимо, новые версии системы от них избавлены. Тем не менее полностью исключать такую возможность нельзя, особенно если ваша программа будет запускаться на старой технике с устаревшими системами. Чтобы приведенный пример всегда выдавал правильное значение -6.7762635780344E-20, достаточно поставить в начале нашей процедуры Set8080CW(Get8087CW or $0100), и программа в любой системе будет устанавливать сопроцессор в режим максимальной точности.
Примечание
В версиях Delphi по 5-ю включительно, где отсутствует функция Get8087CW, можно использовать такую конструкцию : Set8087CW(Default8087CW). При этом следует учитывать, что она возвращает к начальному состоянию все флаги, а не только интересующий нас. Если это неприемлемо, управляющее слово придется изменять с помощью встроенного ассемблера.
Раз уж мы заговорили об управляющем слове, давайте немного поэкспериментируем с ним. Изменим первую строчку на Set8087CW(Get8087CW and $FCFF or $0200). Тем самым мы перевезем сопроцессор в режим 53-разрядной точности представления мантиссы. Теперь в любой системе мы увидим 1.44327637948555Е-16, несмотря на использование Extended. Если же мы изменим первую строчку на Set8087CW(Get8087CW and $FCFF), то будем работать в режиме 24-разрядной точности. Соответственно, в любой системе будет результат -7.3015691270939Е-8.
Заметим, что при загрузке в 10-байтный регистр сопроцессора числа типа Extended в режиме пониженной точности "лишние" биты не обнуляются. Только результаты математических операций представляются с пониженной точностью. Кроме того, при сравнении двух чисел также учитываются все биты, независимо от точности. Поэтому код, приведенный в листинге 3.10 при выборе любой точности даст Не равно.
3.2.11. Борьба с потерей точности в VCL
В том, что описанная проблема с потерей точности встречается все реже, есть заслуга и разработчиков VCL. Зная, вызовы каких функций могут привести к изменению управляющего слова FPU, они перед этими вызовами запоминают управляющее слово, а затем восстанавливают. В более поздних версиях Delphi количество таких "оберток" больше, чем в ранних, поэтому чем новее версия Delphi, тем меньше шанс столкнуться с описанной проблемой. Здесь мы рассмотрим несколько примеров из исходного кода стандартных модулей Delphi 2007.
Для динамической загрузки DLL предназначена API-функция LoadLibrary. В модуле SysUtils для этой функции предлагается обертка, называющаяся SafeLoadLibrary (листинг 3.13).
Листинг 3.13. Функция SysUtils.SafeLoadLibrary
{ SafeLoadLibrary calls LoadLibrary, disabling normal Win32 error message popup dialogs if the requested file can't be loaded. SafeLoadLibrary also preserves the current FPU control word (precision, exception masks) across the LoadLibrary call (in case the DLL you're loading hammers the FPU control word in its initialization, as many MS DLLs do) }
function SafeLoadLibrary(const Filename: string; ErrorMode: UINT): HMODULE;
var
OldMode: UINT;
FPUControlWord: Word;
begin
OldMode := SetErrorMode(ErrorMode);
try
asm
FNSTCW FPUControlWord
end;
try
Result := LoadLibrary(PChar(Filename));
finally
asm
FNCLEX
FLDCW FPUControlWord
end;
end;
finally
SetErrorMode(OldMode);
end;
end;
Как видно из комментария, проблема в том, что многие системные библиотеки изменяют управляющее слово FPU при своей инициализации.
В функции CreateADOObject (внутренняя функция модуля ADODB) тоже сохраняется и восстанавливается управляющее слово (листинг 3.14).
Листинг 3.14. Функция CreateADOObject модуля ADODB
function CreateADOObject(const ClassID: TGUID): IUnknown;
var
Status: HResult;
FPUControlWord: Word;
begin
asm
FNSTCW FPUControlWord
end;
Status :=
CoCreateInstance(ClassID, nil, CLSTX_INPROC_SERVER or CLSCTX_LOCAL_SERVER, IUnknown, Result);
asm
FNCLEX
FLDCW FPUControlWord
end;
if (Status = REGDB_E_CLASSNOTREG) then
raise Exception.CreateRes(@SADOCreateError)
else OleCheck(Status);
end;
Здесь восстанавливать управляющее слово приходится после вызова системной функции CoCreateInstance, создающей СОМ-объект. Но, судя по тому, что больше нигде при вызове CoCreateInstance такой код не используется, проблема не в самой функции, а в тех конкретных ADO-объектах, которые создаются здесь с ее помощью.
Аналогичную защиту можно обнаружить в модуле Dialogs, в методе TCommonDialog.TaskModalDialog. Комментарий к этой защите гласит: "Avoid FPU control word change in NETRAP.dll, NETAPI32.dll, etc".
В модуле Windows особым образом импортируются функции CreateWindow и CreateWindowEx, которые, видимо, тоже были замечены в некорректном обращении с управляющим словом FPU. Вот как, например, выглядит импорт функции CreateWindowEx (листинг 3.15).
Листинг 3.15. Импорт функции CreateWindowEx модулем Windows
function _CreateWindowEx(dwExStyle: WORD; lpClassName: PChar; lpWindowName: PChar; dwStyle: DWORD; X, Y, nWidth, nHeight: Integer; hWndParent: HWND; hMenu: HMENU; hInstance: HINST; lpParam: Pointer): HWND; stdcall; external user32 name 'CreateWindowExA';
function CreateWindowEx(dwExStyle: DWORD; lpClassName: PChar; lpWindowName: PChar; dwStyle: DWORD; X, Y, nWidth, nHeight: Integer; hWndParent: HWND; hMenu: HMENU; hInstance: HINST; lpParam: Pointer): HWND;
var
FPUCW: Word;
begin
FPUCW := Get8087CW;
Result :=
_CreateWindowEx(dwExStyle, lpClassName, lpWindowName,
dwStyle, X, Y, nWidth, nHeight, hWndParent, hMenu,
hInstance, lpParam);
Set8087CW(FPUCW);
end;
Модуль Windows импортирует функцию CreateWindowExA из библиотеки user32.dll, но дает ей измененное название и не показывает ее в своем интерфейсе. Вместо этого он экспортирует другую функцию с названием CreateWindowEx (и аналогичную с названием CreateWindowExA), которая является оберткой над настоящей CreateWindowExA и обеспечивает сохранение значения управляющего слова FPU. Аналогичным способом импортируется и Unicode-вариант функции. Таким образом, стандартные библиотеки обеспечивают вызов безопасного варианта CreateWindowEx в любой программе.
Примечание
В модуле Windows можно обнаружить еще одну интересную деталь: функции CreateWindowA и CreateWindowW из библиотеки user32.dll этим модулем вообще не импортируются. Вместо этого одноименные обертки вызывают импортированные функции _CreateWindowExA и _CreateWindowExW, передавая им 0 в качестве значения параметра dwExStyle.
3.2.12. Машинное эпсилон
Когда мы имеем дело с вычислениями с ограниченной точностью, возникает такой парадокс. Пусть, например, мы считаем с точностью до трех значащих цифр. Прибавим к числу 1,00 число 1,00·10-4. Если бы все было честно, мы получили бы 1,0001. Но у нас ограничена точность, поэтому мы вынуждены округлять до трех значащих цифр. В результате получается 1,00. Другими словами, к некоторому числу мы прибавляем другое число, большее нуля, а в результате из-за ограниченной точности мы получаем то же самое число. Наименьшее положительное число, которое при добавлении его к единице дает результат, не равный единице, называется машинным эпсилон.
Понятие машинного эпсилон у новичков нередко путается с понятием наименьшего числа, которое может быть записано в выбранном формате. Это неправильно. Машинное эпсилон определяется только размером мантиссы, а минимально возможное число оказывается существенно меньше из-за сдвига плавающей двоичной точки с помощью экспоненты.
Прежде чем искать машинное эпсилон программно, попытаемся найти его из теоретических соображений. Итак, мантисса типа Extended содержит 64 разряда. Чтобы закодировать единицу, старший бит мантиссы должен быть равен 1 (денормализованная запись), остальные биты — нулю. Очевидно, что при такой записи наименьшее из чисел, для которых выполняется условие x > 1, получается, когда самый младший бит мантиссы тоже будет равен единице, т.е. х = 1,00...001 (в двоичном представлении, между точкой и младшей единицей 62 нуля). Таким образом, машинное эпсилон равно х-1, т.е. 0.00...001. В более привычной десятичной форме записи это будет 2-63, т.е. примерно 1,084·10-19.
Листинг 3.16 показывает, как можно найти это число (пример Epsilon на компакт-диске).
Листинг 3.16. Поиск машинного эпсилон
procedure TForm1.Button1Click(Sender: TObject);
var
R: Extended;
I: Integer;
begin
R := 1;
while 1 + R/2 > 1 do R := R / 2;
Label1.Caption := FloatToStr(R);
end;
Запустив этот код, мы получим на экране 1.0842021724855Е-19 в полном соответствии с нашими теоретическими выкладками.
Примечание
В тех системах, где наблюдается описанная проблема с уменьшением точности, программа выдаст 2.22044604925031Е-16. Если вы увидели у себя это число, добавьте код, который переведет FPU в режим максимальной точности.
А теперь изменим тип переменной R с Extended на Double. Результат не изменится. На Single — опять не изменится. Но такое поведение лишь на первый взгляд может показаться странным. Давайте подробнее рассмотрим выражение 1 + R / 2 > 1. Итак, все вычисления (в том числе и сравнение) сопроцессор выполняет с данными типа Extended. Последовательность действий такова: число R загружается в регистр сопроцессора, преобразуясь при этом к типу Extended. Дальше оно делится на 2, а затем к результату прибавляется 1, и все это в Extended, никакого обратного преобразования в Single или Double не происходит. Затем это число сравнивается с единицей. Очевидно, что результат сравнения не должен зависеть от исходного типа R, т.к. диапазона даже типа Single вполне хватает, чтобы разместить машинное эпсилон.
3.2.13. Методы решения проблем
Подведем итоги сказанному. Значения, которые мы получаем, могут отличаться от ожидаемых, даже если речь идет о простом присваивании. Во многих случаях (например, в научных расчетах) это несущественно, т.к. сам метод расчета дает еще большую погрешность. Проблемы начинаются там, где мы хотим вывести число на экран или сравнить его с другим. Универсальных рецептов на все случаи жизни не существует, но во многих ситуациях помогают следующие советы:
□ Если ваша задача — просто получить "красивое" представление числа на экране, то функцию FloatToStr заменяйте на ее более мощный аналог FloatToStrF или на функцию Format — они позволяют указать желаемое количество символов после точки.
□ Сравнение вещественных чисел следует выполнять с учетом погрешности, т.е. вместо if а = b … писать if Abs(а - b) < Ерs …, где Eps — некоторая величина, задающая допустимую погрешность (в модуле Math, начиная с Delphi 6, существует функция SameValue, с помощью которой это же условие можно записать как if SameValue(a, b, Eps) …).
□ Для денежных расчетов следует выбирать тип Currency, реализующий число с фиксированной, а не плавающей, десятичной точкой. Отметим также, что не следует пытаться решить проблему неточного представления числа (0,100000001490116 вместо 0,1) с помощью функции RoundTo, поскольку эта функция не может обеспечить точность бо́льшую, чем точность аппаратного представления вещественных чисел.
