Elementary Math
(src/math/elementary_math.hpp)

• View this file on GitHub
• Last update: 2024-04-07 17:34:21+09:00
• Include: #include "src/math/elementary_math.hpp"

概要

基本的な数学関数の集合．

関数	効果	時間計算量
divisor(n)	$n$ の約数を列挙する．	$O(\sqrt{n})$
prime_factor(n)	$n$ の素因数を列挙する．	$O(\sqrt{n})$
osa_k(n)	$n$ 以下の各数 $i$ について，その最小素因数を記録した配列を返す．	$O(n \log \log n)$
prime_factor(min_factor, n)	osa_k で構築された最小素因数の配列 min_factor をもとに $n$ の素因数を列挙する．	$O(\log n)$
modpow(x, n, mod)	$x^n \pmod{mod}$ を計算する．	$O(\log n)$
extgcd(a, b, x, y)	$x, y$ を参照渡しすることで $ax + by = d$ を満たす $x, y$ に変更し，$a$ と $b$ の最大公約数 $d$ を返す．	$O(\log \max(a,b))$
inv_mod(a, mod)	$\mathbb{Z}/mod\mathbb{Z}$ における $a$ の逆元を返す．$a$ と $mod$ が互いに素でなく逆元が存在しない場合は -1 を返す．	$O(\log mod)$
euler_phi(n)	Euler’s totient function $\varphi(n)$ を返す．ここで，$\varphi(n)$ は $n$ 以下の $n$ と互いに素な自然数の個数である．	$O(\sqrt{n})$
euler_phi_table(n)	$n$ 以下の各数 $i$ について，$\varphi(i)$ を記録した配列を返す．	$O(n \log \log n)$
order(x, m)	$x^i \equiv 1 \pmod{m}$ を満たす最小の正整数 $i$ を返す．存在しない場合は -1 を返す．	$O(\sqrt{m}\log{m})$
quotient_ranges(n)	$\left\lfloor \frac{n}{i} \right\rfloor = k \iff i \in [l,r]$ となるとき，3 数の組 $(l,r,k)$ を $l$ の昇順に列挙する．	$O(\sqrt{n})$

Verified with

• test/aoj/NTL_1_A.test.cpp
• test/aoj/NTL_1_B.test.cpp
• test/aoj/NTL_1_D.test.cpp
• test/aoj/NTL_1_E.test.cpp

Code

#pragma once
#include <algorithm>
#include <cassert>
#include <numeric>
#include <tuple>
#include <vector>

namespace elementary_math {

template <typename T> std::vector<T> divisor(T n) {
    std::vector<T> res;
    for (T i = 1; i * i <= n; i++) {
        if (n % i == 0) {
            res.emplace_back(i);
            if (i * i != n) res.emplace_back(n / i);
        }
    }
    return res;
}

template <typename T> std::vector<std::pair<T, int>> prime_factor(T n) {
    std::vector<std::pair<T, int>> res;
    for (T p = 2; p * p <= n; p++) {
        if (n % p == 0) {
            res.emplace_back(p, 0);
            while (n % p == 0) {
                res.back().second++;
                n /= p;
            }
        }
    }
    if (n > 1) res.emplace_back(n, 1);
    return res;
}

std::vector<int> osa_k(int n) {
    std::vector<int> min_factor(n + 1, 0);
    for (int i = 2; i <= n; i++) {
        if (min_factor[i]) continue;
        for (int j = i; j <= n; j += i) {
            if (!min_factor[j]) {
                min_factor[j] = i;
            }
        }
    }
    return min_factor;
}

std::vector<int> prime_factor(const std::vector<int>& min_factor, int n) {
    std::vector<int> res;
    while (n > 1) {
        res.emplace_back(min_factor[n]);
        n /= min_factor[n];
    }
    return res;
}

long long modpow(long long x, long long n, long long mod) {
    assert(0 <= n && 1 <= mod && mod < (1LL << 31));
    if (mod == 1) return 0;
    x %= mod;
    long long res = 1;
    while (n > 0) {
        if (n & 1) res = res * x % mod;
        x = x * x % mod;
        n >>= 1;
    }
    return res;
}

long long extgcd(long long a, long long b, long long& x, long long& y) {
    long long d = a;
    if (b != 0) {
        d = extgcd(b, a % b, y, x);
        y -= (a / b) * x;
    } else
        x = 1, y = 0;
    return d;
}

long long inv_mod(long long a, long long mod) {
    assert(1 <= mod);
    long long x, y;
    if (extgcd(a, mod, x, y) != 1) return -1;
    return (mod + x % mod) % mod;
}

template <typename T> T euler_phi(T n) {
    auto pf = prime_factor(n);
    T res = n;
    for (const auto& p : pf) {
        res /= p.first;
        res *= p.first - 1;
    }
    return res;
}

std::vector<int> euler_phi_table(int n) {
    std::vector<int> res(n + 1, 0);
    std::iota(res.begin(), res.end(), 0);
    for (int i = 2; i <= n; i++) {
        if (res[i] != i) continue;
        for (int j = i; j <= n; j += i) res[j] = res[j] / i * (i - 1);
    }
    return res;
}

// minimum i > 0 s.t. x^i \equiv 1 \pmod{m}
template <typename T> T order(T x, T m) {
    T n = euler_phi(m);
    auto cand = divisor(n);
    std::sort(cand.begin(), cand.end());
    for (auto& i : cand) {
        if (modpow(x, i, m) == 1) {
            return i;
        }
    }
    return -1;
}

template <typename T> std::vector<std::tuple<T, T, T>> quotient_ranges(T n) {
    std::vector<std::tuple<T, T, T>> res;
    T m = 1;
    for (; m * m <= n; m++) res.emplace_back(m, m, n / m);
    for (; m >= 1; m--) {
        T l = n / (m + 1) + 1, r = n / m;
        if (l <= r and std::get<1>(res.back()) < l) res.emplace_back(l, r, n / l);
    }
    return res;
}

}  // namespace elementary_math

#line 2 "src/math/elementary_math.hpp"
#include <algorithm>
#include <cassert>
#include <numeric>
#include <tuple>
#include <vector>

namespace elementary_math {

template <typename T> std::vector<T> divisor(T n) {
    std::vector<T> res;
    for (T i = 1; i * i <= n; i++) {
        if (n % i == 0) {
            res.emplace_back(i);
            if (i * i != n) res.emplace_back(n / i);
        }
    }
    return res;
}

template <typename T> std::vector<std::pair<T, int>> prime_factor(T n) {
    std::vector<std::pair<T, int>> res;
    for (T p = 2; p * p <= n; p++) {
        if (n % p == 0) {
            res.emplace_back(p, 0);
            while (n % p == 0) {
                res.back().second++;
                n /= p;
            }
        }
    }
    if (n > 1) res.emplace_back(n, 1);
    return res;
}

std::vector<int> osa_k(int n) {
    std::vector<int> min_factor(n + 1, 0);
    for (int i = 2; i <= n; i++) {
        if (min_factor[i]) continue;
        for (int j = i; j <= n; j += i) {
            if (!min_factor[j]) {
                min_factor[j] = i;
            }
        }
    }
    return min_factor;
}

std::vector<int> prime_factor(const std::vector<int>& min_factor, int n) {
    std::vector<int> res;
    while (n > 1) {
        res.emplace_back(min_factor[n]);
        n /= min_factor[n];
    }
    return res;
}

long long modpow(long long x, long long n, long long mod) {
    assert(0 <= n && 1 <= mod && mod < (1LL << 31));
    if (mod == 1) return 0;
    x %= mod;
    long long res = 1;
    while (n > 0) {
        if (n & 1) res = res * x % mod;
        x = x * x % mod;
        n >>= 1;
    }
    return res;
}

long long extgcd(long long a, long long b, long long& x, long long& y) {
    long long d = a;
    if (b != 0) {
        d = extgcd(b, a % b, y, x);
        y -= (a / b) * x;
    } else
        x = 1, y = 0;
    return d;
}

long long inv_mod(long long a, long long mod) {
    assert(1 <= mod);
    long long x, y;
    if (extgcd(a, mod, x, y) != 1) return -1;
    return (mod + x % mod) % mod;
}

template <typename T> T euler_phi(T n) {
    auto pf = prime_factor(n);
    T res = n;
    for (const auto& p : pf) {
        res /= p.first;
        res *= p.first - 1;
    }
    return res;
}

std::vector<int> euler_phi_table(int n) {
    std::vector<int> res(n + 1, 0);
    std::iota(res.begin(), res.end(), 0);
    for (int i = 2; i <= n; i++) {
        if (res[i] != i) continue;
        for (int j = i; j <= n; j += i) res[j] = res[j] / i * (i - 1);
    }
    return res;
}

// minimum i > 0 s.t. x^i \equiv 1 \pmod{m}
template <typename T> T order(T x, T m) {
    T n = euler_phi(m);
    auto cand = divisor(n);
    std::sort(cand.begin(), cand.end());
    for (auto& i : cand) {
        if (modpow(x, i, m) == 1) {
            return i;
        }
    }
    return -1;
}

template <typename T> std::vector<std::tuple<T, T, T>> quotient_ranges(T n) {
    std::vector<std::tuple<T, T, T>> res;
    T m = 1;
    for (; m * m <= n; m++) res.emplace_back(m, m, n / m);
    for (; m >= 1; m--) {
        T l = n / (m + 1) + 1, r = n / m;
        if (l <= r and std::get<1>(res.back()) < l) res.emplace_back(l, r, n / l);
    }
    return res;
}

}  // namespace elementary_math

Back to top page