Блог

у clang и GCC есть int __builtin_ctz(unsigned) функция. Это подсчитывает конечные нули в целом числе. В статье Википедии об этом семействе функций упоминается, что двоичный алгоритм GCD можно ускорить с помощью __builtin_ctz , но я не понимаю, как.

Пример реализации двоичного GCD выглядит следующим образом:

 unsigned int gcd(unsigned int u, unsigned int v)
{
    // simple cases (termination)
    if (u == v)
        return u;

    if (u == 0)
        return v;

    if (v == 0)
        return u;

    // look for factors of 2
    if (~u amp; 1) // u is even
        if (v amp; 1) // v is odd
            return gcd(u >> 1, v);
        else // both u and v are even
            return gcd(u >> 1, v >> 1) << 1;

    if (~v amp; 1) // u is odd, v is even
        return gcd(u, v >> 1);

    // reduce larger argument
    if (u > v)
        return gcd(u - v, v);

    return gcd(v - u, u);
}
  

Мое подозрение заключается в том, что я мог бы использовать __builtin_ctz следующим образом:

 constexpr unsigned int gcd(unsigned int u, unsigned int v)
{
    // simplified first three ifs
    if (u == v || u == 0 || v == 0)
        return u | v;

    unsigned ushift = __builtin_ctz(u);
    u >>= ushift;

    unsigned vshift = __builtin_ctz(v);
    v >>= vshift;

    // Note sure if max is the right approach here.
    // In the if-else block you can see both arguments being rshifted
    // and the result being leftshifted only once.
    // I expected to recreate this behavior using max.
    unsigned maxshift = std::max(ushift, vshift);

    // The only case which was not handled in the if-else block before was
    // the odd/odd case.
    // We can detect this case using the maximum shift.
    if (maxshift != 0) {
        return gcd(u, v) << maxshift;
    }

    return (u > v) ? gcd(u - v, v) : gcd(v - u, u);
}

int main() {
    constexpr unsigned result = gcd(5, 3);
    return resu<
}
  

К сожалению, это пока не работает. Программа выдает результат 4, когда он должен быть равен 1. Так что же я делаю не так? Как я могу правильно использовать __builtin_ctz здесь? Смотрите пока мой код на GodBolt.

Вот моя итеративная реализация из комментариев:

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

 unsigned ugcd (unsigned u, unsigned v)
{
    unsigned t = u | v;

    if (u == 0 || v == 0)
        return t; /* return (v) or (u), resp. */

    int g = __builtin_ctz(t);

    while (u != 0)
    {
        u >>= __builtin_ctz(u);
        v >>= __builtin_ctz(v);

        if (u >= v)
            u = (u - v) / 2;
        else
            v = (v - u) / 2;
    }

    return (v << g); /* scale by common factor. */
}
  

Как уже упоминалось, |u - v| / 2 шаг обычно реализуется как очень эффективный, безусловный сдвиг вправо, например, shr r32 для деления на (2) — поскольку оба (u) , (v) нечетны и, следовательно, |u - v| должны быть четными.

Это не строго необходимо, так как шаг ‘oddifying’: u >>= __builtin_clz(u); эффективно выполнит эту операцию на следующей итерации.

Предположим, что (u) или (v) имеют «случайное» распределение битов, вероятность (n) конечных нулей через tzcnt равна ~ (1/(2^n)) . Эта инструкция является улучшением по сравнению с bsf реализацией для __builtin_clz , предшествующей Haswell, IIRC.

Благодаря полезным комментаторам, я нашел критическую ошибку: Я должен был использовать min вместо max

Это окончательное решение:

 #include <algorithm>

constexpr unsigned gcd(unsigned u, unsigned v)
{
    if (u == v || u == 0 || v == 0)
        return u | v;

    // effectively compute min(ctz(u), ctz(v))
    unsigned shift = __builtin_ctz(u | v);
    u >>= __builtin_ctz(u);
    v >>= __builtin_ctz(v);

    const auto amp;[min, max] = std::minmax(u, v);

    return gcd(max - min, min) << shift;
}

int main() {
    constexpr unsigned g = gcd(25, 15); // g = 5
    return g;
}
  

Это решение также имеет очень хороший, почти безветвленный результат компиляции.

Вот некоторые результаты тестирования всех ответов на данный момент (мы на самом деле превзошли std::gcd ):