随笔 - C++ 基于标签分发的线性筛

发表于 更新于 分类于 ， 本文字数： 29 阅读时长 ≈ 1 分钟

小练习

需 C++20 及以上

代码

euler_seive.hppview raw
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
#ifndef EULER_SEIVE_HPP
#define EULER_SEIVE_HPP 1

#include <vector>
#include <type_traits>
#include <cmath>

namespace euler_seive {
namespace type_traits__ {
struct tag_base__ {};
struct null_tag: tag_base__ {};
struct min_pfactor_tag: tag_base__ {};
struct euler_phi_tag: tag_base__ {};
struct mobius_tag: tag_base__ {};
}  // namespace type_traits__
using namespace type_traits__;

namespace detail__ {
template <class Tp, bool IsNull, class Fp, class Fij, class Mulf>
requires std::integral<Tp>
constexpr auto seive_impl_(Tp n, Fp f_p, Fij f_ij, Mulf mulf)
  -> std::pair<std::vector<Tp>, std::vector<Tp>> {
  std::vector<Tp> ret, prime;
  std::vector<bool> vis((size_t)n);

  if constexpr (!IsNull) {
    ret.resize((size_t)n);
    ret[1] = 1;
  }
  prime.reserve(n <= 55 ? 16 : n / (log((long double)n) - 4));

  for (Tp i = 2; i < n; ++i) {
    if (!vis[i]) {
      vis[i] = 1;
      if constexpr (!IsNull) ret[i] = f_p(i);
      prime.push_back(i);
    }
    for (auto &&j : prime) {
      Tp ij = i * j;
      if (ij >= n) break;
      vis[ij] = 1;
      if (i % j == 0) {
        if constexpr (!IsNull) ret[ij] = f_ij(ret, i, j);
        break;
      } else if constexpr (!IsNull) ret[ij] = mulf(ret, i, j);
    }
  }
  return {prime, ret};
}

#define PARAMS_(v, i, j) (std::vector<Tp> const &v, Tp i, Tp j)

template <class Tp>
constexpr auto seive_helper(Tp n, null_tag) {
  return seive_impl_<Tp, true>(n, nullptr, nullptr, nullptr);
}
template <class Tp>
auto seive_helper(Tp n, min_pfactor_tag) {
  return seive_impl_<Tp, false>(
    n,
    [](Tp x) { return x; },
    [] PARAMS_(v, i, j) { return j; },
    [] PARAMS_(v, i, j) { return j; });
}
template <class Tp>
constexpr auto seive_helper(Tp n, euler_phi_tag) {
  return seive_impl_<Tp, false>(
    n,
    [](Tp x) { return x - 1; },
    [] PARAMS_(v, i, j) { return v[i] * j; },
    [] PARAMS_(v, i, j) { return v[i] * v[j]; });
}
template <class Tp>
constexpr auto seive_helper(Tp n, mobius_tag) {
  return seive_impl_<Tp, false>(
    n,
    [](Tp x) { return -1; },
    [] PARAMS_(v, i, j) { return 0; },
    [] PARAMS_(v, i, j) { return v[i] * v[j]; });
}

#undef PARAMS_
}  // namespace detail__

template <class Tp, class Tag>
requires std::derived_from<Tag, tag_base__>
constexpr auto seive(Tp n) {
  return detail__::seive_helper<Tp>(n, Tag{});
}
}  // namespace euler_seive

#endif

测试

euler_seive_test.cppview raw
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
#include "euler_seive.hpp"
#include <iostream>
#include <iterator>
#include <tuple>

template <class Ch, class Tr, class Ct>
requires std::random_access_iterator<typename Ct::iterator>
std::basic_ostream<Ch, Tr> &operator<<(std::basic_ostream<Ch, Tr> &os,
                                       Ct const &x) {
  if (x.begin() == x.end()) return os;
  for (auto it = x.begin(); it != x.end() - 1; ++it) os << *it << ' ';
  os << *(x.back() - 1);
  return os;
}

int main() {
  std::vector<int> x, y;

#define SINGLE_TEST__(tag, n)                                    \
  std::tie(x, y) = euler_seive::seive<int, euler_seive::tag>(n); \
  std::cout << "tag: " << #tag << std::endl                      \
            << "n: " << n << std::endl                           \
            << "x: " << x << std::endl                           \
            << "y: " << y << std::endl                           \
            << std::endl

  SINGLE_TEST__(null_tag, 100);
  SINGLE_TEST__(min_pfactor_tag, 100);
  SINGLE_TEST__(euler_phi_tag, 100);
  SINGLE_TEST__(mobius_tag, 100);

#undef SINGLE_TEST__

  return 0;
}

// clang-format off
/* Output:
tag: null_tag
n: 100
x: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
y: 

tag: min_pfactor_tag
n: 100
x: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
y: 0 1 2 3 2 5 2 7 2 3 2 11 2 13 2 3 2 17 2 19 2 3 2 23 2 5 2 3 2 29 2 31 2 3 2 5 2 37 2 3 2 41 2 43 2 3 2 47 2 7 2 3 2 53 2 5 2 3 2 59 2 61 2 3 2 5 2 67 2 3 2 71 2 73 2 3 2 7 2 79 2 3 2 83 2 5 2 3 2 89 2 7 2 3 2 5 2 97 2 3

tag: euler_phi_tag
n: 100
x: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
y: 0 1 1 2 2 4 2 6 4 6 4 10 4 12 6 8 8 16 6 18 8 12 10 22 8 20 12 18 12 28 8 30 16 20 16 24 12 36 18 24 16 40 12 42 20 24 22 46 16 42 20 32 24 52 18 40 24 36 28 58 16 60 30 36 32 48 20 66 32 44 24 70 24 72 36 40 36 60 24 78 32 54 40 82 24 64 42 56 40 88 24 72 44 60 46 72 32 96 42 60

tag: mobius_tag
n: 100
x: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
y: 0 1 -1 -1 0 -1 1 -1 0 0 1 -1 0 -1 1 1 0 -1 0 -1 0 1 1 -1 0 0 1 0 0 -1 -1 -1 0 1 1 1 0 -1 1 1 0 -1 -1 -1 0 0 1 -1 0 0 0 1 0 -1 0 1 0 1 1 -1 0 -1 1 0 0 1 -1 -1 0 1 -1 -1 0 -1 1 0 0 1 -1 -1 0 0 1 -1 0 1 1 1 0 -1 0 1 0 1 1 1 0 -1 0 0

*/