Быстрая математика с фиксированной точкой для финансовых приложений на Java

в 16:53, , рубрики: java, высокая производительность, математика, Программирование, фиксированная точка

Не секрет, что финансовая информация (счета, проводки и прочая бухгалтерия) не очень дружит с числами с плавающей точкой, и множество статей рекомендует использовать фиксированную точку (fixed point arithmetic). В Java этот формат представлен, по сути, только классом BigDecimal, который не всегда можно использовать по соображениям производительности. Приходится искать альтернативы. Эта статья описывает самописную Java библиотеку для выполнения арифметических операций над числами с фиксированной точностью. Библиотека была создана для работы в высокопроизводительных финансовых приложениях и позволяет работать с точностью до 9 знаков после запятой при сохранении приемлемой производительности. Ссылка на исходники и бенчмарки приведены в конце статьи.

Арифметика с плавающей точкой

Cовременные компьютеры могут выполнять арифметические операции только с ограниченной точностью. Это дискретные устройства, которые могут работать не со всеми возможными числами, а только с некоторым их счетным подмножеством. Самым распространённым форматом работы с вещественными числами в памяти компьютера является плавающая (двоичная) точка — floating (binary) point, когда числа хранятся в виде M*2^E, где M и E — целые мантисса и порядок числа. Но некоторые числа, например 0.1, невозможно точно представить в этом формате. Поэтому в ходе сложных вычислений неизбежно накапливается некоторая ошибка. То есть результат машинного вычисления, скажем 0.1 + 0.1 + 0.1, не совпадает с математически правильным 0.3. Учитывая вышесказанное, при программировании сложной арифметики можно придерживаться нескольких стратегий:

Стратегия 1 — игнорировать. Не обращать внимания на погрешность, считать все операции идеально-математическими и надеятся, что имеющейся точности хватит для приемлемых результатов. Самый распространённый вариант.

Стратегия 2 — скрупулёзно подсчитать. Формулы для подсчета машинных погрешностей известны не одно десятилетие. Они позволяют оценить сверху относительную погрешность любой арифметической операции. Наверное, так и приходится делать для серьёзного численного моделирования. Проблема в том, что это очень трудоемко. По сути, каждый символ + — * / в коде должен сопровождаться вычислением погрешности. Нужно учесть все зависимости между вычислениями и повторять процедуру каждый раз при изменении кода.

Стратегия 3 — использовать десятичную точку (floating decimal point) вместо двоичной. То есть хранить числа в виде M*10^E. Это не решает проблем с погрешностью (мантисса по-прежнему округляется до конечного числа значащих цифр), но по крайней мере все «простые» для человека числа (вроде 1.1) теперь представлены в памяти точно. Расплатой будет производительность. Любая нормализация чисел (то есть эквивалентное уменьшение мантиссы и увеличение порядка) требует деления на степень 10, что очень не быстро, в отличие от деления на степень 2. А нормализовывать приходится много — при каждом сложении или вычитании с разными порядками.

Стратегия 4 — использовать фиксированную точку (fixed decimal point). Упрощение стратегии 3, когда мы фиксируем порядок E. В этом случае для сложения/вычитания не нужна нормализация. Кроме того, все вычисления будут иметь одинаковую абсолютную погрешность. Именно этой стратегии посвящена статья.

Арифметика с фиксированной точкой

В отличие от физики, где важна относительная погрешность, в финансах нужна как раз абсолютная. Если после проведения сложной финансовой транзакции клиенту выставить счёт в $1000000.23 в то время как он ожидает $1000000.18, то могут возникнуть некоторые трудности. Объяснения типа «да зачем вам точность в 8 значащих цифр??» могут не прокатить. И дело тут не в 5 центах убытка (ошибиться наоборот, «в пользу» клиента, не сильно лучше), а в нестыковках бухгалтерского учёта. Поэтому правила вычислений и округлений четко оговариваются между сторонами, и артефакты от использования double и float переменных порой усложняют жизнь.

В Java есть стандартный класс для fixed point арифметики — BigDecimal. Проблемы с ним две: он медленный (из-за своей универсальности) и он немутабельный. Немутабельность означает что любая операция выделяет объект в куче. Выделение и освобождение в пересчете на объект занимает немного времени, но интенсивные вычисления в «горячем» коде создают приличную нагрузку на GC, неприемлемую в некоторых случаях. Можно понадеяться на escape-analysis и скаляризацию, но они очень нестабильны в том смысле, что даже незначительное изменение в коде или в JIT (типа ленивой загрузки новой реализации интерфейса) может перевернуть вверх ногами всю структуру инлайна, и метод, минуту назад нормально работавший, вдруг начнёт бешено выделять память.

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

Сразу покажу пример использования, прежде чем переходить к деталям реализации:

public class Sample {
    private final Decimal margin;
    private final Quantity cumQuantity = new Quantity();
    private final Quantity contraQuantity = new Quantity();
    private final Quantity cumContraQuantity = new Quantity();
    private final Price priceWithMargin = new Price();
    private final Price avgPrice = new Price();

    public Sample(int marginBp) {
        // 1 + margin / 10000
        this.margin = Decimal.create(marginBp).divRD(10000L).add(1);
    }

    public Price calculateAvgPrice(Quantity[] quantities, Price[] prices) {
        cumQuantity.set(0);
        contraQuantity.set(0);

        // avg = sum(q * p * margin) / sum(q)
        for (int i = 0; i < quantities.length; i++) {
            cumQuantity.add(quantities[i]);
            priceWithMargin.set(prices[i]).mulRD(margin);
            contraQuantity.set(quantities[i]).mulRD(priceWithMargin);
            cumContraQuantity.add(contraQuantity);
        }

        return avgPrice.quotientRD(cumContraQuantity, cumQuantity);
    }

    public static void main(String[] args) throws ParseException {
        Price p1 = Price.create("1.5");
        Price p2 = Price.create(1.6);

        Quantity q1 = Quantity.create("100");
        Quantity q2 = Quantity.create(200);

        // apply 0.05% margin to the prices
        Sample sample = new Sample(5); 
        System.out.println(sample.calculateAvgPrice(new Quantity[]{q1, q2}, new Price[]{p1, p2}));
    }
}

Идея реализации

Итак, нам нужна мутабельная обертка целочисленного примитива, если точнее — long’а, который даст нам почти 19 значащих цифр (хватит и на целую и на дробную часть). В long'е мы подразумеваем N десятичных знаков после запятой. Например, при N=2, число 2.56 хранится как 256 (двоичное 100000000). Отрицательные числа хранятся стандартно, в дополнительном коде:

-2.56
-256

(Здесь и далее курсивом обозначены «математические» числа и вычисления, а жирным – их внутреннее представление)

Также мне показалось полезным ввести NaN отдельным значением, которое возвращается в случае арифметических ошибок (вместо исключения или мусора). NaN представлен внутри как Long.MIN_VALUE, «распространяется» (propagated) через все операции и позволяет определить инвертирование знака для всех оставшихся чисел.

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

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

1.20 + 2.30 = 3.50
120 + 230 = 350

Умножение и деление требуют дополнительной нормализации, то есть умножения/деления на 10^N (на 100 в нашем примере)

1.20 * 2.00 = 2.40
120 * 200 / 100 = 240

1.20 / 2.00 = 0.60
100 * 120 / 200 = 60

Дополнительное деление — не самая быстрая операция. Но в данном случае это деление на константу, ведь мы заранее зафиксировали N=2 и 10^N=100. Деление на константу, особенно на «красивую» (типа 10), интенсивно оптимизируется в CPU и сильно быстрее деления на случайное число. Мы делаем кучу делений на 10 каждый раз, когда преобразовываем любое число в строку (например в логах), и производители CPU об этом знают (подробнее про оптимизации см "Division by a constant").

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

1.00 / 2.00 = 0.50
100 * 100 / 200 = 50

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

Округление

Попробуем обратить другое число:

1.00 / 3.00 = 0.33
100 * 100 / 300 = 33

Честный математический результат лежит между 0.33 и 0.34, но мы не можем его точно представить. В какую сторону округлять? Обычно округляют к 0, и это самый быстрый способ (поддерживается аппаратно). Но, возвращаясь к реальным финансовым задачам, это не всегда так. Обычно при обработке транзакций с клиентом округление идёт «в пользу клиента». То есть цена округляется вверх, если клиент продаёт, и вниз, если клиент покупает. Но могут потребоваться и другие варианты, например арифметическое округление к ближайшему числу с подтипами (half-up, half-down, half-even) для минимизации бухгалтерских нестыковок. Или округление к ±бесконечности для отрицательных цен (у некоторых финансовых инструментов). Java BigDecimal уже содержит список стандартных режимов округления, и описываемая библиотека их все поддерживает. Режим UNNECESSARY возвращает NaN, если операция неожиданно потребует округления.

В режиме округления вверх наше вычисление должно давать:

1.00 / 3.00 = 0.34
100 * 100 / 300 + 1 = 34

Как узнать, что нужно добавить единицу? Нужен остаток от деления 10000 % 300 = 100. Который такой же медленный, как и само деление. К счастью, если написать подряд в коде "a/b; a%b", то JIT сообразит что 2 деления не нужно, достаточно одной ассебмлерной команды div возвращающей 2 числа (частное и остаток).

Другие варианты округления чуть сложнее, но тоже могут быть вычислены на основе остатка и делителя.

В API я намеренно сделал упоминание округления везде, где оно происходит, либо в виде параметра, либо в виде суффикса RoundDown в методах, где оно по умолчанию происходит к нулю.

Переполнение

Мы подходим к самой сложной части. Вспомним ещё раз наше умножение:

1.20 * 2.00 = 2.40
120 * 200 / 100 = 240

Теперь представим, что мы в 1980-х и процессоры у нас 16-битные. То есть нам доступен только short с максимальным значением 65535. Первое умножение переполнится и будет равно 240000 & 0xFFFF = 44392 (это если без знака, со знаком оно будет ещё и отрицательным), что поломает нам результат.

Так не пойдёт. У нас 2 нормальных (влезающих в наш диапазон значений) аргумента, и такой же нормальный ожидаемый результат, но мы переполняемся на полдороге. Точно такая же ситуация возможна и с 64-битным long’ом, просто числа нужны побольше.

В 1980-х нам потребовалось бы умножение, дающее 32-битный результат. Сегодня нам требуется умножение с 128-битным результатом. Самое обидное то, что оба умножения доступны в ассемблерах 8086 и x86-64 соответственно, но мы не можем использовать их из Java! JNI, даже в случае хака с быстрым JavaCritical, даёт оверхед в десятки наносекунд, привносит сложности с деплоем и совместимостью, замораживает GC на время вызова. К тому же нам каким-то образом пришлось бы возвращать 128-битный результат из native метода, а запись по ссылке в массив (в память) — это дополнительная задержка.

В общем пришлось мне писать ручное умножение и деление. Столбиком. Мне требовались 2 вспомогательные операции:

  1. A(64) * B(64) = T(128); T(128) / N(32)= Q(64),R(32) — как часть fixed point умножения A*B
  2. N(32) * A(64) = T(96); T(96) / B(64) = Q(64),R(64) — как часть fixed point деления A/B
    (в скобках указана размерность данных в битах, T — временная переменная, которая не должна переполняться)

Обе операции возвращают частное и остаток (одно – как результат метода, второе — в поле объекта). Они тоже могут переполняться, но только на последнем шаге, когда это неизбежно. Вот пример (из 1980-х):

500.00 / 0.50 = 1000.00
100 * 50000 / 50 = 100000 — переполнение!

Деление столбиком а-ля Кнут — не самый простой алгоритм. Плюс это все должно быть ещё и относительно быстрым. Поэтому код обоих операций — сотни строк достаточно суровой битовой магии, у меня самого уйдет много времени, чтобы снова вспомнить что конкретно там происходит. Я вытащил их в отдельный класс и откомментировал подробно как мог.

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

Обычно (за исключением особых случаев), обе операции содержат 4 умножения и 2 деления. Операция 1 существенно быстрее 2, так как в ней эти деления — на константу.

Кстати, если кто заметил, N(32) — это наша 10^N для нормализации. Она 32-битная, из чего следует, что N может быть максимум 9. В реальных виденных мной приложениях использовалось 2, 4 или 8 знаков после запятой. Больше 9 я не встречал, так что должно хватить. Если делать 10^N 64-битным, код усложняется (и замедляется) ещё сильнее.

Несколько разных точностей

Иногда необходимо выполнить операцию над аргументами с разным количеством знаков после запятой. Как минимум – ввести операции с участием обычного long.

К примеру:

2.0000(N=4) + 3.00(N=2) = 5.0000(N=4)
20000 + 300 * 100 = 50000

3.00 (N=2) + 2.0000(N=4) = 5.00(N=2)
300 + 20000 / 100 = 500

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

Количество знаков после запятой НЕ хранится в объекте. Вместо этого, предполагается наличие отдельного подкласса для каждой точности. Имена классов могут быть бизнес-ориентированными, например Price (N=8), Quantity (N=2). А могут быть обобщенными: Decimal1, Decimal2, Decimal3,… Чем больше точность, тем меньше диапазон хранимых значений, минимальный диапазон имеет Decimal9: ±9223372036. Предполагается, что одного-двух классов будет достаточно, чтобы покрыть необходимую функциональность, и в этом случае абстрактный метод getScale скорее всего будет девиртуализирован и заинлайнен. Подклассы (вместо дополнительного поля) позволяют строго типизировать точность аргументов и результата, а также сигнализировать о возможном округлении на этапе компиляции.

Библиотека позволяет производить операции в которых участвуют максимум 2 (но не 3) разных точности. То есть должны совпадать либо точности двух аргументов, либо точность одного из аргументов и результата. Опять-таки, поддержка 3-х разных точностей сильно замедлила бы код и усложнила бы API. В качестве аргументов можно передавать обычный long, для которого предполагается точность N=0.

2.0000 / 3.0 = 0.6667 — ok (2 разных точности)
2 / 3= 0.6667 — ok (long аргументы, decimal результат)
2 / 3.0 = 0.6667 — невозможно! (3 разных точности)

Достоинства и недостатки

Очевидно, вычисления повышенной разрядности, проводимые библиотекой, медленнее аппаратно поддерживаемых. Тем не менее, оверхед не настолько велик (см бенчмарки ниже).

Кроме того, из-за отсутствия перегрузки операторов в Java, использование методов вместо арифметических операторов усложняет восприятие кода.

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

Сложные математические алгоритмы (моделирование, статистика, прогнозирование) обычно проще проводить стандартно в double, так как их результат в любом случае не абсолютно точен.

Код и бенчмарки

Код

Benchmark Mode Cnt Score Error Units
DecimalBenchmark.control avgt 200 10.072 ± 0.074 ns/op
DecimalBenchmark.multiplyNative avgt 200 10.625 ± 0.142 ns/op
DecimalBenchmark.multiplyMyDecimal avgt 200 35.840 ± 0.121 ns/op
DecimalBenchmark.multiplyBigDecimal avgt 200 126.098 ± 0.408 ns/op
DecimalBenchmark.quotientNative avgt 200 70.728 ± 0.230 ns/op
DecimalBenchmark.quotientMyDecimal avgt 200 138.581 ± 7.102 ns/op
DecimalBenchmark.quotientBigDecimal avgt 200 179.650 ± 0.849 ns/op

В целом, умножение получается в 4 раза быстрее BigDecimal, деление — в 1.5. Скорость деления сильно зависит от аргументов, отсюда такой разброс значений.

Автор: tmaxx

Источник

* - обязательные к заполнению поля


https://ajax.googleapis.com/ajax/libs/jquery/3.4.1/jquery.min.js