HiveBrain v1.2.0
Get Started
← Back to all entries
patterncppModerate

Compile-time sieve of Eratosthenes

Submitted by: @import:stackexchange-codereview··
0
Viewed 0 times
sievetimecompileeratosthenes

Problem

There are many instances of prime number sieve implementation both here and other places on the web, but I wanted something a little different. In particular, I wanted to create a static array of the first 1024 prime numbers at compile time as a simple reusable structure. I also wanted to allow for relatively simple creation of larger or smaller arrays.

I combined the code from this question with a little bit of macro magic from this other question to come up with code that meets my goal.

Questions:

I'd like to be able to use something more elegant than the stack of macros (templates perhaps), and would prefer to have the result be something more like an STL container than a raw array. Comments on how to do those things, or comments generally to clean up the code would be welcome.

primeconsttest.cpp

#include "primeconst.h"    
#include 
#include 
#include 

int main() {
    std::vector prime{std::begin(primes), std::end(primes)};
    for (unsigned i=0; i < prime.size(); ++i)
        std::cout << i << '\t' << prime[i] << '\n';
}


primeconst.h

extern const int primes[1024];


primeconst.cpp

```
constexpr bool isPrimeLoop(int i, int k) {
return (k*k > i)?true:(i%k == 0)?false:isPrimeLoop(i, k + 1);
}

constexpr bool isPrime(int i) {
return isPrimeLoop(i, 2);
}

constexpr int nextPrime(int k) {
return isPrime(k)?k:nextPrime(k + 1);
}

constexpr int getPrimeLoop(int i, int k) {
return (i == 0)?k:
(i % 2)?getPrimeLoop(i-1, nextPrime(k + 1)):
getPrimeLoop(i/2, getPrimeLoop(i/2, k));
}

constexpr int getPrime(int i) {
return getPrimeLoop(i, 2);
}

static_assert(getPrime(511) == 3671, "computed incorrectly");

#define K(x) J(x) J(x + 512)
#define J(x) I(x) I(x + 256)
#define I(x) H(x) H(x + 128)
#define H(x) G(x) G(x + 64)
#define G(x) F(x) F(x + 32)
#define F(x) E(x) E(x + 16)
#define E(x) D(x) D(x + 8)
#define D(x) C(x) C(x + 4)
#define C(x) B(x) B(x + 2)
#define B(x) A(x) A(x + 1)
#define A(x) getPr

Solution

I believe I can help with these things:

  • Obliterate macros



  • Obliterate raw arrays



Macros can be replaced with std::integer_sequence (C++14), and raw arrays can be replaced with std::array (C++11).

#include 
#include 

template 
constexpr auto gen_primes_helper(std::integer_sequence) {
    return std::array{{getPrime(Is)...}};
}

// T: integer type
// N: number of elements
// return: std::array
template 
constexpr auto gen_primes() {
    return gen_primes_helper(std::make_integer_sequence());
}

constexpr auto primes = gen_primes();


However, we have a problem: Compilers currently implement std::make_integer_sequence in O(N) complexity.

GCC limits N to 900, unless you add the -ftemplate-depth= compiler option. See this bug.

Clang limits N to a pathetic 256, but also provides the same -ftemplate-depth= option.

To solve the problem, std::make_integer_sequence can be implemented in O(log(N)). I have done this as follows:

#include 
#include 

namespace logseq {

// alias for std::integer_sequence, for brevity
template 
using intseq_t = std::integer_sequence;

template 
struct concat;

// A: intseq_t
// B: intseq_t
// return: typename intseq_t
//
// Example:
//   concat_t, intseq_t >
//     => intseq_t
template 
using concat_t = typename concat::type;

template 
struct concat, intseq_t> {
    using type = intseq_t;
};

template 
struct logseq;

// T: integer type
// N: number of elements
// return: typename std::make_integer_sequence
// https://gcc.gnu.org/bugzilla/show_bug.cgi?id=66059
template 
using logseq_t = typename logseq::type;

template 
struct logseq {
    using type = concat_t,logseq_t>;
};

template 
struct logseq::type> {
    using type = intseq_t;
};

template 
struct logseq::type> {
    using type = intseq_t;
};

} // namespace logseq

template 
constexpr auto gen_primes_helper(std::integer_sequence) {
    return std::array{{getPrime(Is)...}};
}

// T: integer type
// N: number of elements
// return: std::array
template 
constexpr auto gen_primes() {
    using logseq::logseq_t;
    return gen_primes_helper(logseq_t());
}

constexpr auto primes = gen_primes();


Namespace logseq provides a metafunction logseq_t that returns a typename std::integer_sequence. It is a drop-in replacement for std::make_integer_sequence, which returns the same type.

logseq_t allows for infinite (read: really freaking big) values of N.

Code Snippets

#include <array>
#include <utility>

template <typename T, T... Is>
constexpr auto gen_primes_helper(std::integer_sequence<T, Is...>) {
    return std::array<T, sizeof...(Is)>{{getPrime(Is)...}};
}

// T: integer type
// N: number of elements
// return: std::array<T,N>
template <typename T, T N>
constexpr auto gen_primes() {
    return gen_primes_helper(std::make_integer_sequence<T,N>());
}

constexpr auto primes = gen_primes<int,1024>();
#include <array>
#include <utility>

namespace logseq {

// alias for std::integer_sequence, for brevity
template <typename T, T... Is>
using intseq_t = std::integer_sequence<T, Is...>;

template <typename A, typename B>
struct concat;

// A: intseq_t
// B: intseq_t
// return: typename intseq_t
//
// Example:
//   concat_t< intseq_t<int,0,1,2>, intseq_t<int,0,1,2,3> >
//     => intseq_t<int,0,1,2,3,4,5,6>
template <typename A, typename B>
using concat_t = typename concat<A,B>::type;

template <typename T, T... As, T... Bs>
struct concat<intseq_t<T, As...>, intseq_t<T, Bs...>> {
    using type = intseq_t<T, As..., (sizeof...(As) + Bs)...>;
};

template <typename T, T N, typename = void>
struct logseq;

// T: integer type
// N: number of elements
// return: typename std::make_integer_sequence<T,N>
// https://gcc.gnu.org/bugzilla/show_bug.cgi?id=66059
template <typename T, T N>
using logseq_t = typename logseq<T,N>::type;

template <typename T, T N, typename>
struct logseq {
    using type = concat_t<logseq_t<T,N/2>,logseq_t<T,N-N/2>>;
};

template <typename T, T N>
struct logseq<T,N,typename std::enable_if<N==0>::type> {
    using type = intseq_t<T>;
};

template <typename T, T N>
struct logseq<T,N,typename std::enable_if<N==1>::type> {
    using type = intseq_t<int,0>;
};

} // namespace logseq


template <typename T, T... Is>
constexpr auto gen_primes_helper(std::integer_sequence<T, Is...>) {
    return std::array<T, sizeof...(Is)>{{getPrime(Is)...}};
}

// T: integer type
// N: number of elements
// return: std::array<T,N>
template <typename T, T N>
constexpr auto gen_primes() {
    using logseq::logseq_t;
    return gen_primes_helper(logseq_t<T,N>());
}

constexpr auto primes = gen_primes<int,1024>();

Context

StackExchange Code Review Q#93775, answer score: 16

Revisions (0)

No revisions yet.