Алгоритм деления 2W-разрядных чисел с использованием операций с числами разрядностью W
На примере 32-битных целых чисел рассматривается масштабируемый алгоритм деления, использующий числа с двукратно меньшей (16 бит) разрядностью. Для иллюстрации работоспособности алгоритма приведен код тестового приложения на языке С++.
Описание алгоритма
Представим некоторое 
-битное число 
в виде:
где 
, 
— старшая и младшая части 
-битных компонент числа, соответственно.
Тогда деление двух чисел, 
и 
, можно записать в виде дроби:
Заметим, что если 
-битное число. Ограничимся этим случаем. Для 
в примере ниже на С++ реализован алгоритм деления «широкого» числа на «узкое», основанный на представлении делимого 
в виде произведения частного 
и делителя 
плюс слагаемого-остатка 
:
Считаем далее, что 
, иначе — результат деления равен нулю. Представим число 
в виде:
Здесь 
— целая часть от деления, а 
— остаток от деления 
на 
.
Перепишем дробь:
Пренебрежем слагаемым 
, чтобы получить нижнюю оценку результата деления плюс дополнительно упростить формулу:
Выделим слагаемое 
в качестве главной компоненты:
Сделаем замену переменных 
Дело в том, что «тяжелые» случаи соответствуют параметру 
достаточно большому, поэтому введенный параметр «дельта» при этом будет мал и формула быстрее приведет к результату:
В числителе и знаменателе дроби стоят «широкие» числа (разрядностью 
). Так как допускается использовать алгоритм деления широкого числа на узкое, то добьемся «узости» числа в знаменателе, вынеся множитель 
за скобки и выполняя деление последовательно:
Таким образом, «широкий» числитель последовательно делим на «узкие» знаменатели. Первый знаменатель иногда может равняться множителю 
, что необходимо отслеживать в алгоритме: в регистрах ЦПУ число фактически обнулится и возникнет исключение. Альтернативно, можно заранее вычесть единицу и не заботиться о граничных условиях, так как для данного алгоритма 
, что в итоге даст окончательный вариант:
Численный эксперимент показал, что достаточно одной вышеописанной итерации. Физически это объясняется тем, что алгоритм разработан специально для args = {
{51774, 28457, 50018, 10280}, {28792, 5507, 37, 64804},
{65258, 18362, 87, 35198}, {65526, 63280, 198, 52129},
{56139, 10364, 39, 36881}, {65498, 60804, 204, 20825},
{58092, 52199, 1, 57003}, {65498, 60804, 204, 20825},
{64666, 34598, 1, 60805}, {30903, 7652, 143, 48035},
{30161, 40182, 3351, 26310}, {40824, 35384, 13, 49151},
{60215, 18033, 165, 58003}, {42499, 42189, 4, 58879},
{16384, 16384, 0, 1}, {16384, 16384, 1, 0}};
auto make_test = [](const Quadrupole q) -> bool {
return test_div({.mHigh = q.A, .mLow = q.B},
{.mHigh = q.C, .mLow = q.D});
};
bool is_ok = true;
for (const auto& arg : args) {
is_ok &= make_test(arg);
}
assert(is_ok);
}
void test_division_randomly(long long N) {
auto make_test = [](const Quadrupole q, const Signess s) -> bool {
return test_div(U32{.mHigh = q.A, .mLow = q.B, .mSign = s.s1},
U32{.mHigh = q.C, .mLow = q.D, .mSign = s.s2});
};
long long counter = 0;
long long ext = 0;
bool is_ok = true;
while (true) {
++counter;
const Quadrupole q{rollu16(), rollu16(), rollu16(), rollu16()};
const Signess s{rollu16() % 2, rollu16() % 2};
if (q.C == 0 && q.D == 0) {
continue;
}
is_ok &= make_test(q, s);
assert(is_ok);
if (counter % (1 ll << 27) == 0) {
ext++;
std::cout << "... iterations: " << counter << '\n';
}
if (ext >= N) {
break;
}
}
}
int main(int argc, char* argv[]) {
long long N = 10;
if (argc > 1) {
N = std::atoi(argv[1]);
std::cout << "You set " << N << " external iterations\n";
}
{
U32 x{.mHigh = 1, .mLow = 0};
std::cout << x.value() << '\n';
}
{
U32 x{.mHigh = 1, .mLow = 65535};
std::cout << x.value() << '\n';
}
std::cout << "Test division quick\n";
test_division_quick();
std::cout << "Ok\n";
std::cout << "Test division randomly...\n";
test_division_randomly(N);
std::cout << "Ok\n";
return 0;
}
Выводы
Предложен и протестирован алгоритм деления чисел, состоящих из старшей и младшей половинок, масштабирующийся на произвольную разрядность кратно 
. Данный алгоритм в некотором смысле является вариантом «умного» деления в столбик: сначала вычисляется первое приближение, равное делению старших половинок числа, а затем — второе, равное скорректированному первому. Корректор равен последовательному делению некоторого широкого числа на два узких.
Предложенный алгоритм легко масштабируется на 128-битный вариант с использованием встроенной 64-битной арифметики. Однако, вариант с масштабированием, например, на 256-бит, требует реализации в структуре U128 полноценного умножения, что можно сделать, масштабируя реализованный оператор «половинчатого» умножения: U32 на u16. Также потребуется реализация оператора побитового сдвига. В конечном итоге, при реализации всех необходимых операторов можно реализовать шаблонную структуру с произвольной разрядностью 
, рекурсивно (делением пополам) спадающую в арифметику 64-битных встроенных чисел.
