Number_Theory/OsakDivisors.hpp

View this file on GitHub
Last update: 2023-02-11 22:28:05+08:00
Include: #include "Number_Theory/OsakDivisors.hpp"

Depends on

Code

#pragma once

#include "Prime_Sieve.hpp"

template <class T>
std::vector<T> OsakDivisors(T x) {
    auto primes = prime_sieve(sqrt(x) + 10);

    std::vector<std::pair<T, int>> lt;
    for (auto p : primes) {
        if (1LL * p * p > x) break;
        int cnt = 0;
        while (x % p == 0) {
            x /= p;
            cnt++;
        }
        if (cnt >= 1) lt.emplace_back(p, cnt);
    }
    if (x > 1) lt.emplace_back(x, 1);
    std::vector<T> div;
    auto dfs = [&](auto&& rec, int i, T c) {
        if (i == int(lt.size())) {
            div.push_back(c);
            return;
        }
        for (int j = 0; j <= lt[i].second; j++) {
            rec(rec, i + 1, c);
            c *= lt[i].first;
        }
    };
    dfs(dfs, 0, 1);
    return div;
}

#line 2 "Number_Theory/OsakDivisors.hpp"

#line 2 "Number_Theory/Prime_Sieve.hpp"

#line 2 "Template/Template.hpp"

#include <bits/stdc++.h>

using i64 = std::int64_t;
#line 4 "Number_Theory/Prime_Sieve.hpp"

std::vector<int> prime_sieve(int N) {
    std::vector<bool> sieve(N / 3 + 1, 1);
    for (int p = 5, d = 4, i = 1, sqn = sqrt(N); p <= sqn;
         p += d = 6 - d, i++) {
        if (!sieve[i]) continue;
        for (int q = p * p / 3, r = d * p / 3 + (d * p % 3 == 2), s = 2 * p,
                 qe = sieve.size();
             q < qe; q += r = s - r)
            sieve[q] = 0;
    }
    std::vector<int> ret{2, 3};
    for (int p = 5, d = 4, i = 1; p <= N; p += d = 6 - d, i++)
        if (sieve[i]) ret.push_back(p);
    while (!ret.empty() && ret.back() > N) ret.pop_back();
    return ret;
}
#line 4 "Number_Theory/OsakDivisors.hpp"

template <class T>
std::vector<T> OsakDivisors(T x) {
    auto primes = prime_sieve(sqrt(x) + 10);

    std::vector<std::pair<T, int>> lt;
    for (auto p : primes) {
        if (1LL * p * p > x) break;
        int cnt = 0;
        while (x % p == 0) {
            x /= p;
            cnt++;
        }
        if (cnt >= 1) lt.emplace_back(p, cnt);
    }
    if (x > 1) lt.emplace_back(x, 1);
    std::vector<T> div;
    auto dfs = [&](auto&& rec, int i, T c) {
        if (i == int(lt.size())) {
            div.push_back(c);
            return;
        }
        for (int j = 0; j <= lt[i].second; j++) {
            rec(rec, i + 1, c);
            c *= lt[i].first;
        }
    };
    dfs(dfs, 0, 1);
    return div;
}