В прошлой статье мы разработали следующее:
Алгоритмы псевдослучайной генерации: xorshift64, lehmer64, xoshiro256pp;
Алгоритмы Фибоначчи и конвертеры из миль в километры через них;
Алгоритмы быстрого обратного квадратного корня;
Алгоритм бинарного возведения в степень.
В этой статье я буду рассматривать более глубокие и интересные трюки на языке C. Вам не обязательно читать первую часть, статьи в этой серии независимы друг от друга.
❯ Быстрое вычисление приближенного значения степени
В прошлой статье мне предложили разобрать быстрое вычисление приближённого значения степени:
inline double fastPow(double a, double b) {
union {
double d;
int x[2];
} u = { a };
u.x[1] = (int)(b * (u.x[1] - 1072632447) + 1072632447);
u.x[0] = 0;
return u.d;
}
Данный алгоритм — это аппроксимация функции aᵇ, основанная на особенностях представления чисел с плавающей точкой в стандарте IEEE 754. Алгоритм использует тот факт, что возведение в степень можно выразить через экспоненту и логарифм: aᵇ = exp(b × ln(a)).
В стандарте IEEE 754 тип double имеет 53 бита мантиссы. Структура числа двойной точности: знак — 1 бит, порядок — 11 бит, мантисса — 52+1 бит. Число представляется как
value = (-1)^sign × 2^(exponent - 1023) × (1 + mantissa/2^52)
Ошибка составляет примерно 1–10% в зависимости от значений a и b. Алгоритм работает значительно быстрее стандартной функции pow().
Сама функция сначала разбирает double на целочисленные компоненты через union. В старших 32 битах содержится порядок числа, который примерно равен логарифму от a по второй степени. Потом идет преобразование u.x[1], оно вычисляет u.x[1] - 1072632447 ≈ log₂a, затем умножение на b даёт b × log₂a, и в конце прибавление 1072632447 преобразует обратно в формат double. После преобразования идет обнуление младших 32 бит, устанавливая мантиссу в 0. Магическое число 1072632447 = 0x3FF00000 - 60701 — это смещённое представление логарифма.
Эта функция имеет свои ограничения: работает только для положительных a, зависит от порядка байтов и размера double.
Но у этого алгоритма есть более быстрая версия (из этой статьи):
float fastestPow(float a, float b) {
union
{
float d;
int x;
} u = { a };
u.x = (int)(b * (u.x - 1065353210) + 1065353210);
return u.d;
}
Этот алгоритм работает за счёт манипуляции битовым представлением чисел float.
Число с плавающей точкой хранится в виде знака, порядка и мантиссы. Порядок по сути является логарифмом числа по основанию 2. Алгоритм использует это свойство: вычитание магического числа 1065353210 (которое соответствует числу 1.0 в битовом представлении) приближённо вычисляет логарифм исходного числа a. Умножение на b даёт b × log₂a. Затем добавление магического числа обратно преобразует результат в битовое представление числа с плавающей точкой.
Алгоритм быстрее предыдущего метода, но и менее точен из-за использования float вместо double. По сути, float-версия делает только одну операцию с битовым представлением вместо двух и работает с данными вдвое меньшего размера, что даёт выигрыш в производительности за счёт потери точности по сравнению с double.
❯ Быстрое деление на 3
Эта функция заменяет медленное целочисленное деление на 3 быстрым умножением и сдвигом. Она использует магическую константу 0xAAAAAAAB, которая представляет собой приближение 2³²/3.
uint32_t div3(uint32_t x) {
return (uint32_t)(((uint64_t)x * 0xAAAAAAABULL) >> 33);
}
При умножении 32-битного числа на эту 64-битную константу старшие 32 бита результата содержат приблизительное значение x/3. Сдвиг вправо на 33 бита извлекает этот результат.
Почему 33 бита, а не 32? При сдвиге на 32 мы получаем 2/3 от x, а не x/3. Изначально я использовал сдвиг на 32 бита, и гадал, почему 12 / 3 = 8...
Такая оптимизация работает в 10–20 раз быстрее обычного деления, поскольку умножение и битовые операции выполняются процессором значительно быстрее операции деления. Этот метод широко используется компиляторами для оптимизации деления на известные константы.
❯ Алгоритм XOR Swap
Классический алгоритм, который позволяет обменивать значения двух переменных с помощью трёх операций XOR без создания временной третьей.
void xor_swap(int *a, int *b) {
if (a != b) {
*a ^= *b;
*b ^= *a;
*a ^= *b;
}
}
Классический, но всегда впечатляющий новичков трюк: обмен значений двух переменных с помощью трёх операций XOR. Это демонстрирует симметрию и свойства операции XOR. Хотя на современных процессорах это часто медленнее, чем использование временной переменной (из-за параллелизма), трюк остаётся красивым и поучительным.
Первая операция сохраняет разность a и b в a, вторая, используя исходное b и новое a, получает исходное a и сохраняет в b, и последняя, используя новое a и новое b, получает исходное b и сохраняет в a.
❯ Трюк Кернигана x & -x
Данный хак выделяет младший установленный бит числа. Работает это благодаря представлению отрицательных чисел в дополнительном коде. Отрицательное число -x вычисляется как ~x + 1 (инверсия всех битов плюс 1). При инверсии все нули становятся единицами. Когда добавляется 1, происходит цепная реакция переносов: все младшие нули становятся единицами, пока перенос не дойдет до первой единицы в исходном числе. Эта единица остается единственной общей между x и -x.
Пример для 12 (1100):
x = 1100
-x = 0100 (в дополнительном коде)
x & -x = 0100
Ниже предоставлено два варианта, которые дают один результат, но работают принципиально по-разному.
Первая функция работает с исходным числом: постепенно сдвигает исходное число x вправо, затем проверяет младший бит (x & 1) == 0, работает за O(n), где n — количество trailing zeros.
Вторая функция (от Кернигана) работает с выделенным битом, сначала вычисляет x & -x — получает число с одной единицей в позиции первого ненулевого бита, после сдвигает это число до тех пор, пока оно не станет равным 1. Работает за O(m), где m — позиция первого ненулевого бита (что равно n).
Эффективность также может отличаться, первая функция работает быстрее при малом количестве trailing zeros, а вторая может быть немного быстрее при работе с большим количеством trailing zeros.
int count_trailing_zeros(unsigned int x) {
if (x == 0) return sizeof(x) * 8;
int count = 0;
while ((x & 1) == 0) {
count++;
x >>= 1;
}
return count;
}
int count_trailing_zeros_kernighan(unsigned int x) {
if (x == 0) return sizeof(x) * 8;
unsigned int lowest_set_bit = x & -x;
int count = 0;
while (lowest_set_bit > 1) {
count++;
lowest_set_bit >>= 1;
}
return count;
}
❯ ГПСЧ PCG
PCG (Permuted Congruential Generator) — это семейство современных ГПСЧ, сочетающих простоту линейного конгруэнтного генератора (LCG) с качественным выходным преобразованием (permutation). PCG быстр, требует мало памяти, имеет огромный период и проходит строгие статистические тесты. Он считается одним из лучших «лёгких» ГПСЧ.
В основе — простой LCG для обновления состояния. Магия происходит в выходной функции: она берёт старшие биты состояния, применяет XOR и сдвиг, а затем циклический сдвиг (ротацию) на случайное количество бит, определяемое другими старшими битами. Это ломает линейные зависимости, присущие простым LCG.
typedef struct { uint64_t state; uint64_t inc; } pcg32_random_t;
uint32_t pcg32_random_r(pcg32_random_t* rng) {
uint64_t oldstate = rng->state;
rng->state = oldstate * 6364136223846793005ULL + (rng->inc | 1);
uint32_t xorshifted = ((oldstate >> 18u) ^ oldstate) >> 27u;
uint32_t rot = oldstate >> 59u;
return (xorshifted >> rot) | (xorshifted << ((-rot) & 31));
}
Структура typedef struct { uint64_t state; uint64_t inc; } pcg32_random_t; имеет два параметра — первый это state (текущее состояние генератора), и inc, константу смещения (increment). Этот параметр фиксируется при инициализации генератора и определяет его поток (stream). Разные inc дают абсолютно разные, непересекающиеся последовательности чисел. Это очень мощная фича для параллельных симуляций.
Вторым шагом идет классический шаг LCG: new_state = old_state * multiplier + increment. Сам по себе LCG — довольно слабый генератор. Его младшие биты имеют очень короткие периоды (например, младший бит просто чередуется 0, 1, 0, 1...), а в многомерных пространствах точки, сгенерированные LCG, выстраиваются в предсказуемые гиперплоскости (парадокс Марсальи). Поэтому мы не можем просто взять младшие 32 бита от 64-битного state — они будут выглядеть очень неслучайно.
Парадокс Марсальи (теорема Марсальи) в теории случайных чисел описывает недостатки псевдослучайных чисел, получаемых в результате линейного конгруэнтного генератора. Это утверждение связано с работой американского математика Джорджа Марсальи, который предложил в 1965 году генератор Макларена — Марсальи — криптографически стойкий генератор псевдослучайных чисел.
В конце и происходит вся магия, превращающая слабый LCG в мощный генератор. Первый этап это XSH (Xor Shift) — старшие биты состояния сдвигаются и объединяются операцией XOR с младшими битами. Это «ломает» линейные зависимости LCG, перемешивая биты. Второй шаг это RR (Random Rotate) — другие старшие биты состояния определяют, на сколько бит нужно циклически сдвинуть результат первого шага. Это добавляет нелинейность и усложняет предсказание.
Простой LCG обеспечивает скорость и большой период, а функция вывода маскирует его недостатки, создавая на выходе статистически безупречные числа. Это делает PCG быстрым, компактным и надёжным генератором.
❯ Трюк Боба Дженкинса для быстрого хеширования — lookup3
Это один из классических, быстрых и качественных некриптографических хеш-алгоритмов. Он использует смесь сложения, XOR и битовых сдвигов для тщательного перемешивания битов.
#define rot(x, k) (((x) << (k)) | ((x) >> (32 - (k))))
void jenkins_mix(uint32_t *a, uint32_t *b, uint32_t *c) {
*a -= *c; *a ^= rot(*c, 4); *c += *b;
*b -= *a; *b ^= rot(*a, 6); *a += *c;
*c -= *b; *c ^= rot(*b, 8); *b += *a;
*a -= *c; *a ^= rot(*c, 16); *c += *b;
*b -= *a; *b ^= rot(*a, 19); *a += *c;
*c -= *b; *c ^= rot(*b, 4); *b += *a;
}
Lookup3 предназначен для одной главной цели: быстро преобразовать произвольные входные данные (ключ) в целое число (хеш) с равномерным распределением. Он не является криптографическим, то есть его нельзя использовать для паролей, цифровых подписей или где важна стойкость к взлому. Его сфера применения — это в первую очередь хеш-мапы, и прочие структуры и алгоритмы, где нужна «случайность» на основе данных.
Вращение битов — сердце перемешивания, оно реализовано в макросе rot, это циклический сдвиг. Функция jenkins_mix отвечает за первоначальный микс чисел. Она берет три числа (a, b, c) и за 12 шагов так их перемешивает, что исходные данные становятся неразличимы.
Полный исходный код от Боба Дженкинса можно посмотреть по этой ссылке.
Быстрый модуль 2^n и быстрая проверка на степень двойки
Под конец статьи стоит рассказать о двух быстрых коротких функциях.
Начнем с быстрой проверки на степень двойки:
int8_t is_power_of_two(uint32_t x) {
return x && !(x & (x - 1));
}
Эта функция возвращает true (истина), если число x является степенью двойки (1, 2, 4, 8, 16, 32, …), и false (ложь) в противном случае.
Функция работает на наблюдении, что степень двойки в двоичной системе — это всегда число, в котором есть ровно одна единица (1), а все остальные биты — нули (0). Степени двойки в двоичном виде имеют ровно одну единицу. Выражение x & (x-1) сбрасывает младшую единицу числа. Если результат равен нулю и само число не ноль — значит, была всего одна единица, что соответствует степени двойки.
Разберем теперь более простую функцию быстрого модуля 2^n:
uint32_t fast_mod(uint32_t x, uint32_t mod) {
return x & (mod - 1);
}
Работает только когда mod — степень двойки (2^n). Остаток от деления на 2^n определяется последними n битами числа. Выражение mod-1 создает битовую маску из n единиц. Операция x & (mod-1) выделяет эти последние n битов, что эквивалентно взятию модуля.
Пример: 13 % 8 = 5, так как 13 (1101) & 7 (0111) = 5 (0101).
❯ Бенчмаркинг
С учетом функционала из прошлой части, я разработал бенчмарк:
PRNG Performance (10,000,000 iterations):
-----------------------------------------
xorshift64: 28.92 ms (345.76M numbers/s)
lehmer64: 43.29 ms (231.00M numbers/s)
xoshiro256pp: 33.13 ms (301.87M numbers/s)
pcg32: 27.99 ms (357.32M numbers/s)
-----------------------------------------
Conversion Methods Performance (each method called 10000 times per point):
----------------------------------------------------------------------
Basic : 0.58 ms ( 0.003 us/call)
Fibonacci Interpolation : 3.38 ms ( 0.017 us/call)
Fibonacci Cache : 2.90 ms ( 0.014 us/call)
Golden Ratio : 33.72 ms ( 0.169 us/call)
Golden Ratio (Binary) : 7.30 ms ( 0.036 us/call)
----------------------------------------------------------------------
Accuracy Comparison (5 sample points):
Miles | Basic | Interpol | Cache | Golden | GoldenBin
------+-----------+----------+----------+----------+-----------
5 | 8.05 | 0.58% | 0.58% | 0.58% | 0.58%
30 | 48.28 | 0.53% | 0.53% | 0.53% | 0.53%
55 | 88.51 | 0.55% | 0.55% | 0.55% | 0.55%
80 | 128.75 | 0.54% | 0.54% | 0.54% | 0.54%
100 | 160.93 | 0.54% | 0.54% | 0.54% | 0.54%
---------------------------------------------------------------
Math Algorithms Performance (1000000 iterations):
--------------------------------------------
fast_pow: 2.87 ms ( 0.003 us/call)
fastest_pow: 2.41 ms ( 0.002 us/call)
fast_mod: 2.45 ms ( 0.002 us/call)
is_power_of_two: 2.32 ms ( 0.002 us/call)
jenkins_hash: 58.44 ms ( 0.058 us/call)
jenkins_mix+final: 60.88 ms ( 0.061 us/call)
xor_swap: 4.51 ms ( 0.005 us/call)
div3: 2.79 ms ( 0.003 us/call)
--------------------------------------------
Как мы видим, pcg32 реально оказался самым быстрым генератором на данный момент.
И вот график по математическим и прочим алгоритмам:
❯ Заключение
Спасибо за прочтение статьи! Я надеюсь, вы узнали что-то новенькое, или может, какой-нибудь трюк натолкнул вас на другой интересный алгоритм. Если нашли нюанс в самой статье — пишите в комментарии. Возможно, я продолжу серию таких статей на тему трюков на разных языках программирования.
Если вам понравился изложенный материал, могу предложить подписаться на мой блог в телеграме. Если, конечно, вам статья понравилась и вы хотите видеть чуть больше.
Примеры работы кода вы можете увидеть в моем репозитории The Art Of Fun C.
❯ Источники
Новости, обзоры продуктов и конкурсы от команды Timeweb.Cloud — в нашем Telegram-канале ↩
Перед оплатой в разделе «Бонусы и промокоды» в панели управления активируйте промокод и получите кэшбэк на баланс.