Skip to content

Instantly share code, notes, and snippets.

@dc1394
Last active July 7, 2024 02:37
Show Gist options
  • Save dc1394/32f7604cd9ee456cba0d649e3c7517aa to your computer and use it in GitHub Desktop.
Save dc1394/32f7604cd9ee456cba0d649e3c7517aa to your computer and use it in GitHub Desktop.
C++20で素数を求めるコードその2
// コンパイルオプション:g++ -std=c++20 -Wall -O3 -march=native -mtune=native my_primes_2.cpp
#include <chrono> // for std::chrono
#include <cmath> // for std::log
#include <cstdint> // for std::int64_t
#include <format> // for std::format
#include <iostream> // for std::iostream
#include <vector> // for std::vector
namespace {
static auto constexpr MAX = 10000000000LL;
void seek_prime();
} // namespace
int main()
{
std::cout << std::format("max = {}\n", MAX);
seek_prime();
return 0;
}
namespace {
void seek_prime()
{
using namespace std::chrono;
auto const start = system_clock::now();
auto const max = MAX + 1;
auto const pmax = static_cast<std::int64_t>(static_cast<double>(max) / std::log(static_cast<double>(max)) * 1.2);
std::vector<char> sieve(max / 2 + 16);
std::vector<std::int64_t> prime(pmax);
auto * const p = reinterpret_cast<std::int64_t *>(sieve.data());
for (auto i = 0; i <= max / 16; i += 3) {
p[i + 0] = 0x0001010001010001; //0x01010001;
p[i + 1] = 0x0101000101000101; //0x00010100;
p[i + 2] = 0x0100010100010100; //0x01000101;
}
sieve[1] = 1;
for (auto i = 5LL; i * i < max; i += 2LL) {
if (sieve[i / 2]) {
auto j = i * i / 2;
sieve[j] = 0;
j += i * (i % 3 - 1);
for (; j < max / 2 - i * 2; j += i) {
sieve[j] = 0;
j += i * 2;
sieve[j] = 0;
}
if (j < max / 2) {
sieve[j] = 0;
}
}
}
prime[0] = 2;
auto pnum = 1LL;
for (auto i = 1LL; i < max / 2; i++) {
prime[pnum] = i * 2 + 1;
pnum += sieve[i];
}
auto const end = system_clock::now();
auto const elapsed = std::chrono::duration<double>(end - start).count();
std::cout << std::format("{:.3f} (sec)\n素数が{}個見つかりました\n", elapsed, pnum);
}
} // namespace
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment