Modular Arithmetic Utilities
(weilycoder/number_theory/mod_utility.hpp)

• View this file on GitHub
• Last update: 2025-11-08 07:03:41+08:00
• Include: #include "weilycoder/number_theory/mod_utility.hpp"

Required by

• Functions for factorizing numbers using Pollard's Rho algorithm (weilycoder/number_theory/factorize.hpp)
• Modular Integer Class (weilycoder/number_theory/modint.hpp)
• Prime Number Utilities (weilycoder/number_theory/prime.hpp)
• Functions to find primitive roots modulo a prime (weilycoder/number_theory/primitive_root.hpp)
• Polynomial Derivative functions. (weilycoder/poly/elementary_func/derivative.hpp)
• Polynomial Exponential Functions (weilycoder/poly/elementary_func/exponential.hpp)
• Polynomial Exponential Functions (weilycoder/poly/elementary_func/exponential.hpp)
• Polynomial Exponential Functions (weilycoder/poly/elementary_func/exponential.hpp)
• Polynomial Integration functions. (weilycoder/poly/elementary_func/integrate.hpp)
• Polynomial Inversion functions using Newton's method and NTT. (weilycoder/poly/elementary_func/inverse.hpp)
• Polynomial Inversion functions using Newton's method and NTT. (weilycoder/poly/elementary_func/inverse.hpp)
• Polynomial Inversion functions using Newton's method and NTT. (weilycoder/poly/elementary_func/inverse.hpp)
• Polynomial Logarithm Functions (weilycoder/poly/elementary_func/logarithm.hpp)
• Polynomial Logarithm Functions (weilycoder/poly/elementary_func/logarithm.hpp)
• Polynomial Logarithm Functions (weilycoder/poly/elementary_func/logarithm.hpp)
• Polynomial Multiplication Functions (weilycoder/poly/elementary_func/multiply.hpp)
• Polynomial Multiplication Functions (weilycoder/poly/elementary_func/multiply.hpp)
• Polynomial Multiplication Functions (weilycoder/poly/elementary_func/multiply.hpp)
• Multiply two polynomials under modulo using MTT with 3 NTT-friendly moduli. (weilycoder/poly/mtt_convolve.hpp)
• Multiply two polynomials under modulo using MTT with 3 NTT-friendly moduli. (weilycoder/poly/mtt_convolve.hpp)
• Number Theoretic Transform (NTT) implementation (weilycoder/poly/ntt.hpp)
• Number Theoretic Transform (NTT) implementation (weilycoder/poly/ntt.hpp)
• Multiplying polynomials using Number Theoretic Transform (NTT) (weilycoder/poly/ntt_convolve.hpp)
• Multiplying polynomials using Number Theoretic Transform (NTT) (weilycoder/poly/ntt_convolve.hpp)

Verified with

• test/convolution_mod.test.cpp
• test/convolution_mod.test.cpp
• test/convolution_mod_1000000007.test.cpp
• test/convolution_mod_1000000007.test.cpp
• test/exp_of_formal_power_series.test.cpp
• test/exp_of_formal_power_series.test.cpp
• test/exp_of_formal_power_series.test.cpp
• test/factorize.test.cpp
• test/inv_of_formal_power_series.test.cpp
• test/inv_of_formal_power_series.test.cpp
• test/inv_of_formal_power_series.test.cpp
• test/inv_of_formal_power_series.test.cpp
• test/log_of_formal_power_series.test.cpp
• test/log_of_formal_power_series.test.cpp
• test/log_of_formal_power_series.test.cpp
• test/matrix_product.test.cpp
• test/primality_test.test.cpp
• test/primitive_root.test.cpp

Code

#ifndef WEILYCODER_MOD_UTILITY_HPP
#define WEILYCODER_MOD_UTILITY_HPP

/**
 * @file mod_utility.hpp
 * @brief Modular Arithmetic Utilities
 */

#include <cstdint>

namespace weilycoder {
using u128 = unsigned __int128;

/**
 * @brief Perform modular addition for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param a The first addend.
 * @param b The second addend.
 * @param modulus The modulus.
 * @return (a + b) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_add(uint64_t a, uint64_t b, uint64_t modulus) {
  if constexpr (bit32) {
    uint64_t res = a + b;
    if (res >= modulus)
      res -= modulus;
    return res;
  } else {
    u128 res = static_cast<u128>(a) + b;
    if (res >= modulus)
      res -= modulus;
    return res;
  }
}

/**
 * @brief Perform modular addition for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param a The first addend.
 * @param b The second addend.
 * @return (a + b) % Modulus
 */
template <uint64_t Modulus> constexpr uint64_t mod_add(uint64_t a, uint64_t b) {
  if constexpr (Modulus <= UINT32_MAX) {
    uint64_t res = a + b;
    if (res >= Modulus)
      res -= Modulus;
    return res;
  } else {
    u128 res = static_cast<u128>(a) + b;
    if (res >= Modulus)
      res -= Modulus;
    return res;
  }
}

/**
 * @brief Perform modular subtraction for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param a The minuend.
 * @param b The subtrahend.
 * @param modulus The modulus.
 * @return (a - b) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_sub(uint64_t a, uint64_t b, uint64_t modulus) {
  if constexpr (bit32) {
    uint64_t res = (a >= b) ? (a - b) : (modulus + a - b);
    return res;
  } else {
    u128 res = (a >= b) ? (a - b) : (static_cast<u128>(modulus) + a - b);
    return res;
  }
}

/**
 * @brief Perform modular subtraction for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param a The minuend.
 * @param b The subtrahend.
 * @return (a - b) % Modulus
 */
template <uint64_t Modulus> constexpr uint64_t mod_sub(uint64_t a, uint64_t b) {
  if constexpr (Modulus <= UINT32_MAX) {
    uint64_t res = (a >= b) ? (a - b) : (Modulus + a - b);
    return res;
  } else {
    u128 res = (a >= b) ? (a - b) : (static_cast<u128>(Modulus) + a - b);
    return res;
  }
}

/**
 * @brief Perform modular multiplication for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param a The first multiplicand.
 * @param b The second multiplicand.
 * @param modulus The modulus.
 * @return (a * b) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_mul(uint64_t a, uint64_t b, uint64_t modulus) {
  if constexpr (bit32)
    return a * b % modulus;
  else
    return static_cast<u128>(a) * b % modulus;
}

/**
 * @brief Perform modular multiplication for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param a The first multiplicand.
 * @param b The second multiplicand.
 * @return (a * b) % Modulus
 */
template <uint64_t Modulus> constexpr uint64_t mod_mul(uint64_t a, uint64_t b) {
  if constexpr (Modulus <= UINT32_MAX)
    return a * b % Modulus;
  else
    return static_cast<u128>(a) * b % Modulus;
}

/**
 * @brief Perform modular exponentiation for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param base The base number.
 * @param exponent The exponent.
 * @param modulus The modulus.
 * @return (base^exponent) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_pow(uint64_t base, uint64_t exponent, uint64_t modulus) {
  uint64_t result = 1 % modulus;
  base %= modulus;
  while (exponent > 0) {
    if (exponent & 1)
      result = mod_mul<bit32>(result, base, modulus);
    base = mod_mul<bit32>(base, base, modulus);
    exponent >>= 1;
  }
  return result;
}

/**
 * @brief Perform modular exponentiation for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param base The base number.
 * @param exponent The exponent.
 * @return (base^exponent) % Modulus
 */
template <uint64_t Modulus>
constexpr uint64_t mod_pow(uint64_t base, uint64_t exponent) {
  uint64_t result = 1 % Modulus;
  base %= Modulus;
  while (exponent > 0) {
    if (exponent & 1)
      result = mod_mul<Modulus>(result, base);
    base = mod_mul<Modulus>(base, base);
    exponent >>= 1;
  }
  return result;
}

/**
 * @brief Compute the modular inverse of a 64-bit integer using Fermat's Little
 *        Theorem.
 * @tparam Modulus The modulus (must be prime).
 * @param a The number to find the modular inverse of.
 * @return The modular inverse of a modulo Modulus.
 */
template <uint64_t Modulus> constexpr uint64_t mod_inv(uint64_t a) {
  return mod_pow<Modulus>(a, Modulus - 2);
}

/**
 * @brief Compute the modular inverse of a compile-time 64-bit integer using
 *        Fermat's Little Theorem.
 * @tparam Modulus The modulus (must be prime).
 * @tparam a The number to find the modular inverse of.
 * @return The modular inverse of a modulo Modulus.
 */
template <uint64_t Modulus, uint64_t a> constexpr uint64_t mod_inv() {
  return mod_pow<Modulus>(a, Modulus - 2);
}
} // namespace weilycoder

#endif

#line 1 "weilycoder/number_theory/mod_utility.hpp"



/**
 * @file mod_utility.hpp
 * @brief Modular Arithmetic Utilities
 */

#include <cstdint>

namespace weilycoder {
using u128 = unsigned __int128;

/**
 * @brief Perform modular addition for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param a The first addend.
 * @param b The second addend.
 * @param modulus The modulus.
 * @return (a + b) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_add(uint64_t a, uint64_t b, uint64_t modulus) {
  if constexpr (bit32) {
    uint64_t res = a + b;
    if (res >= modulus)
      res -= modulus;
    return res;
  } else {
    u128 res = static_cast<u128>(a) + b;
    if (res >= modulus)
      res -= modulus;
    return res;
  }
}

/**
 * @brief Perform modular addition for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param a The first addend.
 * @param b The second addend.
 * @return (a + b) % Modulus
 */
template <uint64_t Modulus> constexpr uint64_t mod_add(uint64_t a, uint64_t b) {
  if constexpr (Modulus <= UINT32_MAX) {
    uint64_t res = a + b;
    if (res >= Modulus)
      res -= Modulus;
    return res;
  } else {
    u128 res = static_cast<u128>(a) + b;
    if (res >= Modulus)
      res -= Modulus;
    return res;
  }
}

/**
 * @brief Perform modular subtraction for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param a The minuend.
 * @param b The subtrahend.
 * @param modulus The modulus.
 * @return (a - b) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_sub(uint64_t a, uint64_t b, uint64_t modulus) {
  if constexpr (bit32) {
    uint64_t res = (a >= b) ? (a - b) : (modulus + a - b);
    return res;
  } else {
    u128 res = (a >= b) ? (a - b) : (static_cast<u128>(modulus) + a - b);
    return res;
  }
}

/**
 * @brief Perform modular subtraction for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param a The minuend.
 * @param b The subtrahend.
 * @return (a - b) % Modulus
 */
template <uint64_t Modulus> constexpr uint64_t mod_sub(uint64_t a, uint64_t b) {
  if constexpr (Modulus <= UINT32_MAX) {
    uint64_t res = (a >= b) ? (a - b) : (Modulus + a - b);
    return res;
  } else {
    u128 res = (a >= b) ? (a - b) : (static_cast<u128>(Modulus) + a - b);
    return res;
  }
}

/**
 * @brief Perform modular multiplication for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param a The first multiplicand.
 * @param b The second multiplicand.
 * @param modulus The modulus.
 * @return (a * b) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_mul(uint64_t a, uint64_t b, uint64_t modulus) {
  if constexpr (bit32)
    return a * b % modulus;
  else
    return static_cast<u128>(a) * b % modulus;
}

/**
 * @brief Perform modular multiplication for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param a The first multiplicand.
 * @param b The second multiplicand.
 * @return (a * b) % Modulus
 */
template <uint64_t Modulus> constexpr uint64_t mod_mul(uint64_t a, uint64_t b) {
  if constexpr (Modulus <= UINT32_MAX)
    return a * b % Modulus;
  else
    return static_cast<u128>(a) * b % Modulus;
}

/**
 * @brief Perform modular exponentiation for 64-bit integers.
 * @tparam bit32 If true, won't use 128-bit arithmetic. You should ensure that
 *         all inputs are small enough to avoid overflow (i.e. bit-32).
 * @param base The base number.
 * @param exponent The exponent.
 * @param modulus The modulus.
 * @return (base^exponent) % modulus
 */
template <bool bit32 = false>
constexpr uint64_t mod_pow(uint64_t base, uint64_t exponent, uint64_t modulus) {
  uint64_t result = 1 % modulus;
  base %= modulus;
  while (exponent > 0) {
    if (exponent & 1)
      result = mod_mul<bit32>(result, base, modulus);
    base = mod_mul<bit32>(base, base, modulus);
    exponent >>= 1;
  }
  return result;
}

/**
 * @brief Perform modular exponentiation for 64-bit integers with a compile-time
 *        modulus.
 * @tparam Modulus The modulus.
 * @param base The base number.
 * @param exponent The exponent.
 * @return (base^exponent) % Modulus
 */
template <uint64_t Modulus>
constexpr uint64_t mod_pow(uint64_t base, uint64_t exponent) {
  uint64_t result = 1 % Modulus;
  base %= Modulus;
  while (exponent > 0) {
    if (exponent & 1)
      result = mod_mul<Modulus>(result, base);
    base = mod_mul<Modulus>(base, base);
    exponent >>= 1;
  }
  return result;
}

/**
 * @brief Compute the modular inverse of a 64-bit integer using Fermat's Little
 *        Theorem.
 * @tparam Modulus The modulus (must be prime).
 * @param a The number to find the modular inverse of.
 * @return The modular inverse of a modulo Modulus.
 */
template <uint64_t Modulus> constexpr uint64_t mod_inv(uint64_t a) {
  return mod_pow<Modulus>(a, Modulus - 2);
}

/**
 * @brief Compute the modular inverse of a compile-time 64-bit integer using
 *        Fermat's Little Theorem.
 * @tparam Modulus The modulus (must be prime).
 * @tparam a The number to find the modular inverse of.
 * @return The modular inverse of a modulo Modulus.
 */
template <uint64_t Modulus, uint64_t a> constexpr uint64_t mod_inv() {
  return mod_pow<Modulus>(a, Modulus - 2);
}
} // namespace weilycoder

Back to top page