🔮 Как считать экспоненту с помощью умножения и сложения 🔮
const int L = 1 << 23; // == 2^23
const float a = L / std::log(2);
const float b = L * 127;
float fastExp(float x) {
x = a * x + b;
return std::bit_cast<float>((int)x);
}
Выше упрощенная версия кода, который считает приблизительное значение экспоненты, который я недавно встретил на GitHub.
Я был супер удивлен, что тут умножение + сложение + преобразование типов, и всё — экспонента подсчитана.
В моей голове экспонента всегда была нетривиальной операцией, на которую требуется много тактов. Как же так, давайте разбираться.
Целые и дробные числа хранятся в битах по-разному
Числа в формате
int хранятся ожидаемым образом, например, для числа
123:
123: int = 0000000000000000000000000 1111011
Числа в формате
float хранятся в виде маски битов:
s eeeeeeee mmmmmmmmmmmmmmmmmmmmmmm
s — 1 бит
знака. Определяет, положительное число или отрицательное.
e — 8 битов
экспоненты. Показывают степень двойки.
m — 23 бита
мантиссы. Хранят дробную часть числа.
Само число = (-1)^{s} * 2^{e} * 1.{m}
Например для нашего числа 123:
123: float = 0 00000110 11101100000000000000000
e = 00000110b = 6
m = 0.1110110b = 0.921875
Проверяем:
2^{e} x 1.{m} = 2^6 x 1.921875 = 123. Сходится!
Переводить представления из
int в
float и обратно несложно, и процессор умеет делать это очень эффективно.
Чтобы понять суть алгоритма, вместо e^{x} будем считать 2^{x}
Считаем быстро 2^{x} если x - int
Заметим, что если мы знаем двоичное представление
x, то если мы его засунем в биты экспоненты
float, то задача уже решена (
e =
x,
m = 0). Это ключевое наблюдение для понимания алгоритма!
В
int представлении все биты числа
x у нас хранятся в правом конце числа. Чтобы сдвинуть их в область экспоненты, мы умножим наш
int x на
2^23. Тем самым просто допишем 23 нуля справа.
В этих двух идеях и есть вся суть. Как перейти от 2^{x} к e^{х} в других деталях можно разобраться самому или заглянуть в оригинальную статью.
Cчитать экспоненту в десятичной системе мы (человеки) тоже можем с помощью одного умножения (wow!)
Например, чтобы оценить exp(36) мы можем переписать
exp(36) = 10^{36 * log10(e)} = 10^{15.634} ~ 10^16
Эта неожиданная похожесть алгоритмов, получается из-за того что научная запись чисел вроде
5.34e7 очень похожа на то, как хранятся
float в компьютере.
P. S. На самом деле степень двойки в экспоненте
float сдвинута на 127, чтобы мы могли записывать как большие числа, так и очень маленькие. Мы учитываем этот сдвиг, добавляя
b = L * 127
. Я убрал его, чтобы он не отвлекал от основной идеи.
P. P. S. Для тех кто всё понял: почему умножение на 2^32 в коде происходит в
float, а не в
int и что происходит с мантиссой в этом случае?
@lovesyuk