Skip to content

Commit 197559a

Browse files
committed
Update
1 parent d6a1ddf commit 197559a

File tree

1 file changed

+29
-37
lines changed

1 file changed

+29
-37
lines changed

cp-algo/math/sieve.hpp

Lines changed: 29 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -42,6 +42,12 @@ namespace cp_algo::math {
4242
return primes;
4343
}();
4444

45+
constexpr size_t sqrt_threshold = 1 << 16;
46+
const auto sqrt_primes = std::span(
47+
begin(base_primes),
48+
std::ranges::upper_bound(base_primes, sqrt_threshold)
49+
);
50+
4551
constexpr size_t max_wheel_size = 1 << 20;
4652
struct wheel_t {
4753
dynamic_bit_array mask;
@@ -62,16 +68,28 @@ namespace cp_algo::math {
6268
return wheel;
6369
}
6470

65-
constexpr size_t sqrt_threshold = 1 << 16;
66-
const auto sqrt_primes = std::span(
67-
begin(base_primes),
68-
std::ranges::upper_bound(base_primes, sqrt_threshold)
69-
);
70-
71-
std::span<uint32_t const> medium_primes = sqrt_primes;
72-
std::vector<wheel_t> wheels;
73-
74-
void gen_wheels();
71+
auto medium_primes = sqrt_primes;
72+
const auto wheels = []() {
73+
uint32_t product = 64;
74+
std::vector<uint32_t> current;
75+
std::vector<wheel_t> wheels;
76+
for(size_t i = 1; i < size(sqrt_primes); i++) {
77+
uint32_t p = sqrt_primes[i];
78+
if (product * p > max_wheel_size) {
79+
if (size(current) == 1) {
80+
medium_primes = sqrt_primes.subspan(i - 1);
81+
return wheels;
82+
}
83+
wheels.push_back(make_wheel(current, product));
84+
current = {p};
85+
product = 64 * p;
86+
} else {
87+
current.push_back(p);
88+
product *= p;
89+
}
90+
}
91+
assert(false);
92+
}();
7593

7694
void sieve_dense(auto &prime, uint32_t l, uint32_t r, wheel_t const& wheel) {
7795
if (l >= r) return;
@@ -127,35 +145,11 @@ namespace cp_algo::math {
127145
if (l >= r) return;
128146
while (l < r) {
129147
prime.reset(l);
130-
l += gap210[state] * p;
148+
l += p * gap210[state];
131149
state = state == 47 ? 0 : state + 1;
132150
}
133151
}
134152

135-
void gen_wheels() {
136-
static bool initialized = false;
137-
if (initialized) return;
138-
uint32_t product = 64;
139-
std::vector<uint32_t> current;
140-
for(size_t i = 1; i < size(sqrt_primes); i++) {
141-
uint32_t p = sqrt_primes[i];
142-
if (product * p > max_wheel_size) {
143-
if (size(current) == 1) {
144-
medium_primes = sqrt_primes.subspan(i - 1);
145-
initialized = true;
146-
return;
147-
}
148-
wheels.push_back(make_wheel(current, product));
149-
current = {p};
150-
product = 64 * p;
151-
} else {
152-
current.push_back(p);
153-
product *= p;
154-
}
155-
}
156-
assert(false);
157-
}
158-
159153
// Primes smaller than N
160154
dynamic_bit_array sieve(uint32_t N) {
161155
dynamic_bit_array prime(N);
@@ -164,8 +158,6 @@ namespace cp_algo::math {
164158
prime.word(i) = base_prime_bits.word(i);
165159
}
166160
cp_algo::checkpoint("init");
167-
gen_wheels();
168-
cp_algo::checkpoint("gen wheels");
169161
static constexpr uint32_t dense_block = 1 << 23;
170162
for(uint32_t start = base_threshold; start < N; start += dense_block) {
171163
uint32_t r = std::min(start + dense_block, N);

0 commit comments

Comments
 (0)