Постановка задачи

Рассматривается целое (беззнаковое) число шириной бит. Требуется вычислить частное и остаток при делении

имея в распоряжении -битную арифметику.

Данная величина - назовем ее extended reciprocal (ER) - может пригодиться, например, при организации алгоритма деления "широкого" числа на "узкое"

Здесь под широким мы понимаем число разрядности , а под узким - число разрядности бит.

Решение задачи

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

Расширим число максимально возможно, умножая его на два, что соответствует битовому сдвигу числа в сторону старшего разряда "до упора", когда старшая единица займет крайний разряд -битного числа:

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

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

Данная "дельта" - это то, что скорректирует частное и даст искомый остаток:

Окончательно имеем

Заметим, что данный алгоритм актуален при работе с беззнаковыми числами фиксированной разрядности бит, поэтому знак минус работает так, как если принять аксиому , при этом рассматривать знак отдельно (в уме) нельзя - он становится неотъемлемой частью своего числа.

Расширение беззнакового числа может быть сделано, например, функцией библиотеки STL

auto i = std::countl_zero(D);

В итоге, получается такая функция

template <std::unsigned_integral T>
inline std::pair<T, T> reciprocal_and_extend(T x)
{
    assert(x != 0);
    const auto x_old = x;
    const auto i = std::countl_zero(x);
    x <<= i;
    const T &r = (-x) % x_old;
    const T &q = (-x) / x_old + (T{1} << i);
    return {q, r};
}

Алгоритм деления широкого числа на узкое

Пусть стоит задача нахождения частного и остатка от деления широкого числа на узкое:

Здесь - битовая ширина узкого числа (или полуширина широкого числа). Данное деление может пригодиться как самостоятельно, так и, например, как часть алгоритма деления широких чисел

В статье реализован итеративный алгоритм деления широкого числа на узкое, используя дробь

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

В рассматриваемом алгоритме будем вычислять дробь

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

Здесь

то есть задача сводится к новому делению (следующей итерации), и при количество итераций возрастет, так как в этом случае деление выгоднее было бы вычислять, используя дополнение:

что мы, в итоге, и сделаем, но при условии, если дополненный остаток меньше оригинального:

В этом случае удобно увеличить на единицу, потому что

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

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

Это так, потому что если при увеличении остатка на некую величину происходит переполнение, то имеем:

Здесь

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

В итоге, получается следующий алгоритм увеличения остатка.

template <std::unsigned_integral T>
inline T smart_remainder_adder(T &remainder, T delta, const T &modulo)
{
    assert(modulo != 0);
    assert(remainder < modulo);
    const auto rem_old = remainder;
    delta %= modulo;    
    const T &summ = rem_old + delta;
    remainder = summ;
    T carry = summ >= modulo ? T{1} : T{0};
    if (auto overflow = summ < std::min(rem_old, delta); overflow)
    {
        const auto& [_, r_rec] = reciprocal_and_extend(modulo);
        remainder += r_rec;
        carry = T{1};
    }   
    remainder %= modulo;
    return carry;
}

Пара для следующей итерации

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

template <class U, std::unsigned_integral T>
inline std::pair<U, T> div(U X, const T &Y)
{
    assert(Y != 0);
    if (Y == T{1})
        return {X, 0};    
    U Q{0};
    T R = 0;
    auto rcp = reciprocal_and_extend(Y);
    const auto &rcp_compl = Y - rcp.second;
    const bool make_inverse = rcp_compl < rcp.second; // Для ускорения сходимости.
    rcp.first += make_inverse ? T{1} : T{0};
    const auto X_old = X;
    for (;;)
    {
        const bool x_has_high = X.high() != 0;
        Q += x_has_high ? mult64(X.high(), rcp.first) : 0ull;
        Q += U{X.low() / Y};
        const auto& carry = smart_remainder_adder(R, X.low(), Y, rcp.second);
        Q += carry;
        X = X.high() != 0ull ? mult64(X.high(), make_inverse ? rcp_compl : rcp.second) : 0ull;
        if (X == U{0})
        {
            if (Q > X_old) // Коррекция знака.
            {
                Q = -Q;
                R = Y - R; // mod Y
                if (R == Y)
                {
                    Q.inc();
                    R = 0;
                }
            }
            break;
        }
        if (make_inverse)
        {
            Q = -Q;
            R = Y - R; // mod Y
        }
    }
    return {Q, R};
}

Использование трюка с чередованием знака (на примере -битных чисел) снижает максимальное значение количества итераций с до примерно (вероятно, до ~, ), при этом среднее значение (мат. ожидание) снижается примерно с до . Мода распределения практически не изменяется.