gcd.hpp

#pragma once
 
#include <cassert>
#include <tuple>
#include <type_traits>
 
#include "cplib/port/bit.hpp"
 
namespace cplib {
 
template <typename T, std::enable_if_t<std::is_unsigned_v<T>>* = nullptr>
constexpr T gcd(T x, T y) {
  if (x == 0) {
    return y;
  } else if (y == 0) {
    return x;
  }
  int kx = port::countr_zero(x);
  int ky = port::countr_zero(y);
  x >>= kx;
  while (y != 0) {
    y >>= port::countr_zero(y);
    if (x > y) {
      std::swap(x, y);
    }
    y -= x;
  }
  return x << std::min(kx, ky);
}
 
template <typename T, std::enable_if_t<std::is_unsigned_v<T>>* = nullptr>
constexpr std::tuple<std::make_signed_t<T>, std::make_signed_t<T>, T> bezout(T x, T y) {
  bool swap = x < y;
  if (swap) {
    std::swap(x, y);
  }
  if (y == 0) {
    if (x == 0) {
      return {0, 0, 0};
    } else if (swap) {
      return {0, 1, x};
    } else {
      return {1, 0, x};
    }
  }
  using S = std::make_signed_t<T>;
  S s0 = 1, s1 = 0, t0 = 0, t1 = 1;
  while (true) {
    T q = x / y, r = x % y;
    if (r == 0) {
      if (swap) {
        return {t1, s1, y};
      } else {
        return {s1, t1, y};
      }
    }
    S s2 = s0 - S(q) * s1, t2 = t0 - S(q) * t1;
    x = y;
    y = r;
    s0 = s1;
    s1 = s2;
    t0 = t1;
    t1 = t2;
  }
}
 
template <typename T, std::enable_if_t<std::is_unsigned_v<T>>* = nullptr>
constexpr T mod_inverse(T x, T m) {
  auto [s, t, g] = bezout(x, m);
  assert(g == 1);
  return s < 0 ? T(s) + m : s;
}
 
}  // namespace cplib