Двоичные вычисления для десятичной арифметики

    Продолжая исследовать проблему точности десятичных вычислений средствами двоичной арифметики, начатую в предыдущих постах [1,2,3,4], мне удалось разработать алгоритмы вычисления вещественных чисел, представленных в формате десятичных чисел с плавающей точкой, которые дают такой же точный результат, как если бы вычисления велись вручную.

    В этих алгоритмах использована двоичная арифметика, регламентированная стандартом IEEE754. Для проверки работы алгоритмов была разработана тестовая программа на C++, реализующая 18-ти разрядный десятичный калькулятор.
    Поскольку объем материала превышает формат поста, я изложил основные моменты в виде тезисов. Назовем этот пост «Майскими тезисами»:(.

    Итак.

    Известно, что


    Привычная для пользователя арифметика, это десятичная арифметика.

    Существуют также b-ичные арифметики, где b- база системы счисления, принимающая любое ненулевое значение [5].

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

    Если степень числа фиксирована и мантисса числа является целым числом, то такой формат называется форматом с фиксированной точкой. Частным случаем формата с фиксированной точкой является число, в котором степень равна нулю. Такой формат является форматом целого числа.

    Если мантисса представляет собой дробное число в b-ичной системе счисления с целой частью c≠0 и c < b, то такое число называется нормализованным.

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

    Под точными вычислениями в арифметике подразумевают получение результата с заданным количеством верных значащих цифр после точки [6].

    Все вычисления в компьютере производятся в двоичном виде. Для них база b = 2.

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

    Десятичные числа, которые имеют точный двоичный эквивалент, называют представимыми.
    Десятичные числа, которые не имеют точного двоичного эквивалента, называются непредставимыми.

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

    Чем большим количеством двоичных разрядов представлен десятичный эквивалент числа в двоичном виде, тем меньше ошибка представления. Этим объясняется стремление постоянно наращивать разрядность операционного регистра процессора.

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

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

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

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

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

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

    Аналогично, любое десятичное число можно округлить до требуемого количества десятичных цифр, отбрасывая лишние цифры.

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

    Любое вещественное число, выраженное в форме десятичной дроби, может быть точно представлено в формате числа с плавающей точкой (ЧПТ), в котором мантисса является целым числом. Экспонента в ЧПТ будет указывать положение точки в этом числе.

    Если число представлено в формате ЧПТ с целочисленной мантиссой, то мантисса и экспонента этого числа могут быть точно проконвертированы в двоичный код.

    Новое


    Формат, в котором мантисса десятичного ЧПТ представлена двоичным эквивалентом десятичной целочисленной мантиссы, а экспонента является двоичным эквивалентом степени числа 10 (база b=10), будем называть смешанным десятично — двоичным форматом (СДДФ).

    Отличие СДДФ от двоичного формата ЧПТ в том, что экспонента в СДДФ равна степени базы b=10, в то время как экспонента двоичного формата ЧПТ равна двоичной степени базы b=2. Соответственно, для СДДФ число будет представлено как $F=M_{2}10^{e}$, а для ЧПТ, в стандарте IEEE754, как $F=M_{2}2^{e}$.

    Отличие СДДФ от двоично-десятичного формата (ДДФ или BCD) ЧПТ в том, что в ДДФ мантисса и экспонента представляют собой целые десятичные числа, в которых каждая цифра записана в виде байта или тетрады, в то время, как в СДДФ все десятичные числа выражаются их целыми двоичными эквивалентами.

    Таким образом, любое десятичное вещественное число можно представить в СДДФ двоичным кодом с точностью до N значащих десятичных цифр.

    Все арифметические операции над десятичными ЧПТ в СДДФ проводятся по правилам обычной арифметики, где все аргументы являются целыми числами.

    Перед каждой арифметической операцией десятичное число представляется в формате СДДФ:[S,M,z,e]. Где S-код знака числа (0 или 1). M — двоичный целочисленный эквивалент мантиссы числа с количеством десятичных цифр N. Где N — точность вычислений. z — знак экспоненты, e -значение экспоненты. Такое представление является нормализованным представлением десятичного числа.

    Например, для вычислений с точностью до N=7 значащих цифр, число 123,456 должно быть представлено как 1234560*10^-4.

    Минимальное значение десятичной мантиссы числа для N=7 будет равно M=1000000.

    Максимальное значение десятичной мантиссы числа для N=7 будет равно M=9999999.

    Двоичный эквивалент максимального 7-ми разрядного числа 9999999 будет равен 100 110 001 001 011 001 111 111. Он содержит 24 двоичных разряда. Следовательно, для представления десятичных мантисс в диапазоне от 1000000 до 9999999 требуется двоичный 24-разрядный регистр.

    Если в 32-х разрядном двоичном машинном слове, в котором 24 разряда отвести под мантиссу, один разряд отвести под знак числа S, один разряд под знак экспоненты z, а также 6 разрядов под экспоненту, то в таком СДДФ могут быть представлены десятичные вещественные числа с точностью до N=7 значащих десятичных чисел. Абсолютные значения этих чисел лежат в диапазоне от 1000000*10^-64 до 9999999*10^64.

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

    Пример.

    Найти с точностью до N=7 результат выражения (9675,423*10^2-9,675421*10^5)*10^6 — 199992
    Вычисленное вручную, или на калькуляторе Windows, это выражение будет равно числу 8,000000
    Запишем операнды в нормализованном виде:

    A=9,675423*10^5= 9675423*10^-1
    B= 9675,421*10^2 = 9675421*10^-1
    C=1000000 = 1000000*10^0
    D = 1999920*10^-1


    В СДДФ эти операнды будут представлены как:

    A=[0, 9675423,1, 1]
    B=[0,9675421,1, 1]
    C=[0, 1000000,0, 0]
    D=[0, 1999920,1, 1]


    Найдем разность S=A-B. Поскольку экспоненты операндов A и B одинаковы, найдем мантиссу их разности:

    9675423-9675421=2

    Для нормализации мантиссы S надо умножить ее на 10^6, одновременно экспоненту надо уменьшить на 6. Тогда S =2*1000000=2000000*10^-7

    Вычислим произведение P=D*C. Для этого перемножим мантиссы сомножителей и сложим экспоненты:

    M= 2000000*10^-7*1000000*10*0=2000000000000*10^-7
    После нормализации мантиссы получим P=2000000 *10^-1.
    Результат R вычисления будет равен:
    R=P-D=2000000 *10^-1-1999920*10^-1=80*10^-1
    После нормализации получим R = 8000000*10^-6.

    Для сравнения, вычисление этого выражения в Excel дает результат R = 8,0000698E+00.

    Автором разработан алгоритм калькулятора в СДДФ, осуществляющий сложение, вычитание, умножение и деление десятичных чисел с точностью до 18-ти значащих цифр. Для подтверждения правильности алгоритма была написана программа на C++. Поскольку автор не является профессиональным программистом, разработанная программа предназначена только для исследования алгоритма вычислений.

    Ниже, для примера, представлен скриншот, демонстрирующий вычисление следующего выражения:

    1,23456789098765432*10^8 * 9,87654321234567891*10^(-9) - 1,2193263123914037*10^0≈ 3.0*10^(-17)



    Для проверки быстродействия, в цикле была запущена операция умножения двух 18-ти разрядных чисел. Программа запускалась на компьютере Intel® Core(TM) i3-4330 CPU@3.50GHz 3.50 GHz. ОЗУ 8,0 ГБ. Тип системы: 64-разрядная. Скорость получилась равной ≈ 2.4*10^6 умножений в секунду.

    Сравнить с быстродействием калькуляторов Windows и Excel я пока не могу, не хватает образования:(. Что же касается точности вычислений, то она такая же, как если бы расчеты велись вручную.

    Литература:

    1. Взгляд со стороны: Стандарт IEEE754
    2. Нужна ли нормализация в числах с плавающей точкой?
    3. Фатальные ошибки двоичной арифметики при работе с числами с плавающей точкой
    4. Снова о числах с плавающей точкой
    5. Системы счисления
    6. Основные правила приближенных вычислений
    Поделиться публикацией
    AdBlock похитил этот баннер, но баннеры не зубы — отрастут

    Подробнее
    Реклама

    Комментарии 33

      0
      Может вы не знаете но есть IEEE 754—2008 в котором есть decimal32, decimal64, decimal128.
        –1
        Я вам больше скажу, в настоящее время в IEEE рабочая группа по стандарту 754 обсуждает новую редакцию стандарта, IEEE754-2018. Но это никак не противоречит тому, что я написал.
        Что касается decimal128, это нерегламентированный формат. Не думаю, что стоит здесь в это углубляться.
          –1
          Извините, в IEEE754-2018 формат decimal128 прописан. Но алгоритм упаковки для формата обмена в decimal и арифметика, использующая BDC весьма громоздкие.
            +1
            А как с громоздскостью у Вас? Опишите, как происходит сложение чисел с разной экспонентой.
        0
        Рассмотрим сумму двух чисел 1.234567*10^3+7.654321*10^-1. Вычисление проведем с точностью до 7 значащих цифр, или с точностью до 6 знаков после точки.
        Для этого представим наши числа в нормализованном виде: 1234567*10^-3; 7654321*10^-1. В СДДФ они будут выглядеть как [0, 000100101101011010000111,1,11] и
        [0, 011101001100101110110001,1, 1].
        Будем складывать в соответствие с правилами арифметики. 1234567*10^-3+ 7654321*10^-1=(1234567*10^-2+ 7654321) *10^-1= (12345,67 + 7654321)*10^-1=7666666,67*10^-1≈7666667*10^-1. Или в двоичном виде: 1234567*10^-2= 000100101101011010000111/1100100≈ 11000000111001.1≈11000000111010=12346. Сумма чисел в скобках будет 11000000111010+011101001100101110110001=11101001111101111101011=7666667. В результате получим число 7666667*10^-1, такое же, как если бы мы вычисляли вручную с точностью до 7 значащих цифр. В СДДФ это число будет выглядеть как
        [0, 11101001111101111101011,1, 1]. В нем S=0, M=11101001111101111101011=7666667, z=1, e=1.
          0
            0
            Не совсем понял, что говорит google?
              0
              1.234567*10^3 + 7.654321*10^-1 = 1235.3324321
                0
                Так, что здесь не так?
                  +1
                  1235.3324321 не равно 7666667*10^-1
                  0
                  Если в первом слагаемом вместо 10^3 взять 10^-3, то получится 7666667*10^-1.
                    0
                    А ну норм тогда. Альтернативная математика?
                      0
                      Извините, не сразу врубился. Конечно, должно быть так:
                      1.234567*10^3=1234567*10^-3
                      7.654321*10^-1=7654321*10^-7
                      1234567*10^-3+7654321*10^-7=(1234567+7654321*10^-4)*10^-3=
                      (1234567+765,4321)*10^-3=1235332,4321*10^-3≈1235332*10^-3
                0
                >(1234567*10^-2+ 7654321) *10^-1= (12345,67 + 7654321)*10^-1

                Иными словами, приходится умножать на число, не равное основанию системы.
                Это точно менее громоздко, чем двоично-десятичная арифметика?
                  0
                  В СДДФ все операции производятся над целыми числами, представленными в двоичном виде. Мы умножаем/делим ни на число, которое «не равно основанию системы», а на целое число, равное 10 в какой-то степени, которое представлено в двоичном виде.
                    0
                    Мы умножаем/делим ни на число, которое «не равно основанию системы», а на целое число, равное 10 в какой-то степени, которое представлено в двоичном виде.

                    10 по-вашему, равно основанию двоичной системы?
                    Впрочем, это лирика. Суть в том, что Вам приходится делать настоящее умножение вместо сдвига.
                    Вы проверяли, Ваш алгоритм действительно менее громоздок, чем двоично-десятичный?
                      0
                      10 по-вашему, равно основанию двоичной системы?

                      У меня смешанный десятично-двоичный формат.

                      В современных АЛУ операция умножения/деления двух целых чисел делается за 1 такт.
                      А вот, чтобы сложить два десятичных ЧПТ в BCD, а тем более с различными экспонентами, приходится попотеть. www.e-reading.club/chapter.php/99776/140/Cvetkova_-_Informatika_i_informacionnye_tehnologii__konspekt_lekciii.html
                        0
                        У меня смешанный десятично-двоичный формат.

                        Складываете и умножаете Вы двоичные числа. Как храните — вопрос десятый.

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

                        Интересно было бы сравнить СДДФ с ДДФ на процессоре без команд умножения…
                          0
                          Складываете и умножаете Вы двоичные числа. Как храните — вопрос десятый.

                          Как сказать? Если обрабатывать большой массив чисел из памяти, то вопрос распаковки весьма актуален.
                          Интересно было бы сравнить СДДФ с ДДФ на процессоре без команд умножения…

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

                            Вы точно не путаете понятия «поддерживается аппаратное умножение/деление» и «умножение/деление выполняется за 1 такт»?
                              0
                              Конечно, я имел ввиду FPU, в котором поддерживается аппаратное умножение/деление.
                        0
                        >> В современных АЛУ операция умножения/деления двух целых чисел делается за 1 такт.

                        А можно пруф?
                        Или хотя бы просто бенчмарк, который 2 случайных числа делит быстрее чем за наносекунду (такт сильно меньше)

                        По ссылке выше просто информация о медленности BCD, что логично. Поэтому никто в здравом уме в BCD мантиссы и не хранит.
                          0
                          Речь идет об FPU
                            +1
                            Поэтому никто в здравом уме в BCD мантиссы и не хранит.

                            В IEEE754-2008 мантисса может храниться в двух видах. Как целое двоичное число и как упакованное BCD число.
                              0
                              На FPU нет деления за 1 такт, вас обманули.
                              Минимум чего обещают — 11 тактов (руководство по оптимизации от Intel). И это на плавающей точке на команде divss и в хороших условиях, скорее всего применяется метод Goldschmidt'а. На фиксированной — заметно дольше (будет SRT или что-то подобное).

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

                              Ну случайное число ему там не сильно нужно — деление на несколько фиксированных степеней числа 10 может быть сделано через умножение на вшитый обратный множитель для каждого реализованного случая (а вот их уже можно делать через степени двойки в показателе). А на умножение, IMHO, лучше писать 3 такта (если меньше, это такты слишком длинные, а не умножение лёгкое).

                              > Поэтому никто в здравом уме в BCD мантиссы и не хранит.

                              Увы, IEEE754-2008 их таки ввёл в BCD тоже. Зачем, для кого они это сделали — ХЗ.
                              Я пытался найти бенчмарки binary vs decimal в железе, где есть последнее (IBM zSeries и POWER), но пока безуспешно. Есть только такое, где BCD раза в полтора медленнее чисто двоичной, но по разным операциям разница катастрофически неровна => ещё явно надо оптимизировать и отлаживать.
                                0
                                Да, я как-то замерял switch-case на 9 константных делителей vs обычное деление — switch сильно быстрее. Но тут их не 9, а 128, так что хз. Может сработает, а может загнётся кэш инструкций и будет грузить каждое конкретное деление из ОЗУ за те же десятки наносекунд.

                                IEEE действительно странный, допускает как BCD, так и binary. Может они заложились на гипотетическую платформу, которая с BCD работает не сильно медленнее. Нормализовывать десятичную точку все-таки удобнее в BCD.

                      0

                      Я правильно понимаю, что описано вот это, просто по-русски и другими словами?
                      https://en.m.wikipedia.org/wiki/Decimal_floating_point


                      Вычисления в decimal floating point действительно более «привычные» для человека, хотя остаётся проблема с округлением (10^64 + 10^-64 — 10^64 == 0). Основная проблема с ними — производительность, любое сложение/вычитание с разными экспонентами требует деления на степень десятки, что на порядок медленнее обычного сложения.

                        0
                        Необходимо различать следующие этапы арифметической обработки данных.
                        1. Распаковка. Она состоит из декодирования числа, записанного в памяти машинного слова в так называемом формате обмена. Число хранится в упакованном виде.
                        2. Арифметические операции, которые производятся над распакованными числами с использованием операционных регистров в расширенном виде.
                        3. Запаковка результата.
                        Сравнивая формат decimal IEEE и СДДФ, можно видеть, что упаковка-распаковка для decimal производится по достаточно сложному алгоритму, состоящему из проверок множества условий.
                        Числа в СДДФ записываются в память машины без какой-либо предварительной обработки.
                        Десятичная арифметика в настоящее время реализована двоично-десятичным форматом (BCD). Для округления в этом формате требуется операция деления/умножения на 10.
                        Реализовать сложение десятичных чисел в двоичном виде с разными экспонентами без деления на 10 можно только приблизительно.
                        Основная проблема в компьютерной бинарной арифметике, это получение необходимой точности вычислений. Существует много способов повысить точность вычислений. В частности: увеличение разрядности операционных регистров и специальные алгоритмы обработки (алгоритм Кэхэна), которые приводят к существенным дополнительным программно-аппаратным затратам.
                        Алгоритм, использующий СДДФ позволяет получать результаты, такие же, как сделанные вручную с небольшими накладными расходами…

                          0
                          Вы пишите что
                          Десятичная арифметика в настоящее время реализована двоично-десятичным форматом (BCD).
                          . То есть мантиссу вы сначала превращаетесь в bcd (в тетрадах байта цифры числа) и только потом выполняете арифметическую операцию?
                            0
                            Нет, мантисса выражена, как целое двоичное число. Арифметическая операция выполняется над целыми двоичными числами- целочисленная мантисса и двоичный эквивалент характеристики (10^e).
                            Десятичная арифметика в настоящее время реализована двоично-десятичным форматом (BCD).

                            Здесь речь идет о том, как это реализовано в компьютере сегодня.
                        0
                        остаётся проблема с округлением (10^64 + 10^-64 — 10^64 == 0).

                        В СДДФ такой проблемы не существует. Как и в обычной арифметике вычисление, допустим, 1,2345-1,2343 с точностью до 4 значащих цифр даст 0.
                          0

                          Почему этот пост в хабах c#, c, c++?

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

                          Самое читаемое