min_of_mod_linear.hpp

SECTIONMath INCLUDEnoya/min_of_mod_linear.hpp

求有限区间内 \((ax+b)\bmod m\) 的最小值及取得它的下标。

\[ \displaystyle \min_{0\le x<n}(ax+b)\bmod m \]

Complexity: Time: O(log modulus). Space: O(log modulus).

AC 记录：min_of_mod_of_linear。

跳到代码 · GitHub ↗

Implementation

ReferenceSourceExpanded

当前头文件，省略 include guard；依赖见 #include。

/// @complexity Time: O(log modulus).  Space: O(log modulus).

#include <algorithm>
#include <cassert>
#include <cstdint>
#include <numeric>
#include <utility>
#include <vector>

namespace noya {

namespace min_of_mod_linear_detail {

inline std::pair<std::vector<std::int64_t>, std::vector<std::int64_t>>
record_segments(std::int64_t mul, std::int64_t add, std::int64_t mod) {
  std::vector<std::int64_t> bou{0};
  std::vector<std::int64_t> ste;
  std::int64_t div = std::gcd(mul, mod);
  mul /= div;
  add /= div;
  mod /= div;
  // Neighbors ln/ld, rn/rd; dl/dr are their determinants against mul/mod.
  std::int64_t ln = 0;
  std::int64_t ld = 1;
  std::int64_t rn = 1;
  std::int64_t rd = 1;
  std::int64_t dl = mod - mul;
  std::int64_t dr = mul;
  std::int64_t idx = 0;
  std::int64_t val = add;

  while (val != 0) {
    std::int64_t quo = dr / dl;
    dr %= dl;
    if (dr == 0) {
      --quo;
      dr = dl;
    }
    rn += quo * ln;
    rd += quo * ld;
    while (true) {
      quo = std::max<std::int64_t>(0, (dl - val + dr - 1) / dr);
      if (dl - quo * dr <= 0) {
        break;
      }
      dl -= quo * dr;
      ln += quo * rn;
      ld += quo * rd;
      quo = val / dl;
      val -= quo * dl;
      idx += ld * quo;
      bou.push_back(idx);
      ste.push_back(ld);
    }
    quo = dl / dr;
    dl -= quo * dr;
    ln += quo * rn;
    ld += quo * rd;
  }
  return {bou, ste};
}

} // namespace min_of_mod_linear_detail

/// @brief Return an index attaining the minimum of
/// `(mul * x + add) mod mod` for `0 <= x < cnt`, together
/// with that minimum.  Continued-fraction neighbors describe every index at
/// which the prefix minimum decreases; those indices form arithmetic runs, so
/// only one run per Euclidean-algorithm step is inspected.
inline std::pair<std::uint64_t, std::uint64_t>
min_of_mod_linear_argument(std::uint64_t cnt, std::uint64_t mod,
                           std::uint64_t mul, std::uint64_t add) {
  assert(cnt >= 1 && mod >= 1);
  assert(mul < mod && add < mod);
  assert(cnt <= 1'000'000'000ULL && mod <= 1'000'000'000ULL);
  auto [bou, ste] = min_of_mod_linear_detail::record_segments(
      std::int64_t(mul), std::int64_t(add), std::int64_t(mod));
  std::int64_t idx = 0;
  for (std::size_t run = 0; run + 1 < bou.size(); run++) {
    std::int64_t l = bou[run];
    std::int64_t r = bou[run + 1];
    if (std::uint64_t(r) < cnt) {
      idx = r;
      continue;
    }
    idx = l +
          std::int64_t((cnt - 1 - std::uint64_t(l)) / std::uint64_t(ste[run])) *
              ste[run];
    break;
  }
  std::uint64_t val = std::uint64_t(
      (static_cast<unsigned __int128>(mul) * std::uint64_t(idx) + add) % mod);
  return {std::uint64_t(idx), val};
}

/// @brief Return the minimum residue of a linear function on
/// `0 <= x < cnt` using the continued-fraction prefix-minimum decomposition.
inline std::uint64_t min_of_mod_linear(std::uint64_t cnt, std::uint64_t mod,
                                       std::uint64_t mul, std::uint64_t add) {
  return min_of_mod_linear_argument(cnt, mod, mul, add).second;
}

} // namespace noya

#ifndef NOYA_MIN_OF_MOD_LINEAR_HPP
#define NOYA_MIN_OF_MOD_LINEAR_HPP 1

/// @complexity Time: O(log modulus).  Space: O(log modulus).

#include <algorithm>
#include <cassert>
#include <cstdint>
#include <numeric>
#include <utility>
#include <vector>

namespace noya {

namespace min_of_mod_linear_detail {

inline std::pair<std::vector<std::int64_t>, std::vector<std::int64_t>>
record_segments(std::int64_t mul, std::int64_t add, std::int64_t mod) {
  std::vector<std::int64_t> bou{0};
  std::vector<std::int64_t> ste;
  std::int64_t div = std::gcd(mul, mod);
  mul /= div;
  add /= div;
  mod /= div;
  // Neighbors ln/ld, rn/rd; dl/dr are their determinants against mul/mod.
  std::int64_t ln = 0;
  std::int64_t ld = 1;
  std::int64_t rn = 1;
  std::int64_t rd = 1;
  std::int64_t dl = mod - mul;
  std::int64_t dr = mul;
  std::int64_t idx = 0;
  std::int64_t val = add;

  while (val != 0) {
    std::int64_t quo = dr / dl;
    dr %= dl;
    if (dr == 0) {
      --quo;
      dr = dl;
    }
    rn += quo * ln;
    rd += quo * ld;
    while (true) {
      quo = std::max<std::int64_t>(0, (dl - val + dr - 1) / dr);
      if (dl - quo * dr <= 0) {
        break;
      }
      dl -= quo * dr;
      ln += quo * rn;
      ld += quo * rd;
      quo = val / dl;
      val -= quo * dl;
      idx += ld * quo;
      bou.push_back(idx);
      ste.push_back(ld);
    }
    quo = dl / dr;
    dl -= quo * dr;
    ln += quo * rn;
    ld += quo * rd;
  }
  return {bou, ste};
}

} // namespace min_of_mod_linear_detail

/// @brief Return an index attaining the minimum of
/// `(mul * x + add) mod mod` for `0 <= x < cnt`, together
/// with that minimum.  Continued-fraction neighbors describe every index at
/// which the prefix minimum decreases; those indices form arithmetic runs, so
/// only one run per Euclidean-algorithm step is inspected.
inline std::pair<std::uint64_t, std::uint64_t>
min_of_mod_linear_argument(std::uint64_t cnt, std::uint64_t mod,
                           std::uint64_t mul, std::uint64_t add) {
  assert(cnt >= 1 && mod >= 1);
  assert(mul < mod && add < mod);
  assert(cnt <= 1'000'000'000ULL && mod <= 1'000'000'000ULL);
  auto [bou, ste] = min_of_mod_linear_detail::record_segments(
      std::int64_t(mul), std::int64_t(add), std::int64_t(mod));
  std::int64_t idx = 0;
  for (std::size_t run = 0; run + 1 < bou.size(); run++) {
    std::int64_t l = bou[run];
    std::int64_t r = bou[run + 1];
    if (std::uint64_t(r) < cnt) {
      idx = r;
      continue;
    }
    idx = l +
          std::int64_t((cnt - 1 - std::uint64_t(l)) / std::uint64_t(ste[run])) *
              ste[run];
    break;
  }
  std::uint64_t val = std::uint64_t(
      (static_cast<unsigned __int128>(mul) * std::uint64_t(idx) + add) % mod);
  return {std::uint64_t(idx), val};
}

/// @brief Return the minimum residue of a linear function on
/// `0 <= x < cnt` using the continued-fraction prefix-minimum decomposition.
inline std::uint64_t min_of_mod_linear(std::uint64_t cnt, std::uint64_t mod,
                                       std::uint64_t mul, std::uint64_t add) {
  return min_of_mod_linear_argument(cnt, mod, mul, add).second;
}

} // namespace noya

#endif // NOYA_MIN_OF_MOD_LINEAR_HPP

#include <algorithm>
#include <cassert>
#include <cstdint>
#include <numeric>
#include <utility>
#include <vector>

/// @complexity Time: O(log modulus).  Space: O(log modulus).

namespace noya {

namespace min_of_mod_linear_detail {

inline std::pair<std::vector<std::int64_t>, std::vector<std::int64_t>>
record_segments(std::int64_t mul, std::int64_t add, std::int64_t mod) {
  std::vector<std::int64_t> bou{0};
  std::vector<std::int64_t> ste;
  std::int64_t div = std::gcd(mul, mod);
  mul /= div;
  add /= div;
  mod /= div;
  // Neighbors ln/ld, rn/rd; dl/dr are their determinants against mul/mod.
  std::int64_t ln = 0;
  std::int64_t ld = 1;
  std::int64_t rn = 1;
  std::int64_t rd = 1;
  std::int64_t dl = mod - mul;
  std::int64_t dr = mul;
  std::int64_t idx = 0;
  std::int64_t val = add;

  while (val != 0) {
    std::int64_t quo = dr / dl;
    dr %= dl;
    if (dr == 0) {
      --quo;
      dr = dl;
    }
    rn += quo * ln;
    rd += quo * ld;
    while (true) {
      quo = std::max<std::int64_t>(0, (dl - val + dr - 1) / dr);
      if (dl - quo * dr <= 0) {
        break;
      }
      dl -= quo * dr;
      ln += quo * rn;
      ld += quo * rd;
      quo = val / dl;
      val -= quo * dl;
      idx += ld * quo;
      bou.push_back(idx);
      ste.push_back(ld);
    }
    quo = dl / dr;
    dl -= quo * dr;
    ln += quo * rn;
    ld += quo * rd;
  }
  return {bou, ste};
}

} // namespace min_of_mod_linear_detail

/// @brief Return an index attaining the minimum of
/// `(mul * x + add) mod mod` for `0 <= x < cnt`, together
/// with that minimum.  Continued-fraction neighbors describe every index at
/// which the prefix minimum decreases; those indices form arithmetic runs, so
/// only one run per Euclidean-algorithm step is inspected.
inline std::pair<std::uint64_t, std::uint64_t>
min_of_mod_linear_argument(std::uint64_t cnt, std::uint64_t mod,
                           std::uint64_t mul, std::uint64_t add) {
  assert(cnt >= 1 && mod >= 1);
  assert(mul < mod && add < mod);
  assert(cnt <= 1'000'000'000ULL && mod <= 1'000'000'000ULL);
  auto [bou, ste] = min_of_mod_linear_detail::record_segments(
      std::int64_t(mul), std::int64_t(add), std::int64_t(mod));
  std::int64_t idx = 0;
  for (std::size_t run = 0; run + 1 < bou.size(); run++) {
    std::int64_t l = bou[run];
    std::int64_t r = bou[run + 1];
    if (std::uint64_t(r) < cnt) {
      idx = r;
      continue;
    }
    idx = l +
          std::int64_t((cnt - 1 - std::uint64_t(l)) / std::uint64_t(ste[run])) *
              ste[run];
    break;
  }
  std::uint64_t val = std::uint64_t(
      (static_cast<unsigned __int128>(mul) * std::uint64_t(idx) + add) % mod);
  return {std::uint64_t(idx), val};
}

/// @brief Return the minimum residue of a linear function on
/// `0 <= x < cnt` using the continued-fraction prefix-minimum decomposition.
inline std::uint64_t min_of_mod_linear(std::uint64_t cnt, std::uint64_t mod,
                                       std::uint64_t mul, std::uint64_t add) {
  return min_of_mod_linear_argument(cnt, mod, mul, add).second;
}

} // namespace noya