Double factorial and primorial

[Torbjorn Granlund]

nisse at lysator.liu.se (Niels Möller) writes:

  Torbjorn Granlund <tg at gmplib.org> writes:
  
  > then one could sort the indexes using a minheap.
  
  Clever!
  
  That would give the following algorithm for generating primes up to n:
  
  1. Generate a list of all primes p_k up to sqrt(n).
  
  2. Create a minheap which for each k contains the smallest odd multiple
     of p_k which is larger than sqrt(n), together with a link to p_k (can
     this be optimized? I think we should be able to get away without lots
     of divisions. E.g, for those k for which p_k^2 > sqrt(n), we can put
     these squares into the heap rather than the *smallest* multiple.
  
  3. Repeatedly extract and modify the smallest element of the heap. This
     gives us all (odd) composite numbers in the range sqrt(n), n, and any
     missing numbers are the primes, in order.
  
  Is this faster than doing the sieving blockwise using a bitmap? We get a
  very different memory access pattern. It might lose big if the heap (of
  size sqrt(n)) doesn't fit in the cache. I don't know how bad locality
  you typically get with a heap.
  
What I am suggesting is not exactly what you describe.

I am suggesting a blockwise sieving with blocksize S representing a
range of 2S oss numbers.

For primes < 2*S, we sieve as usually, i.e., good old Eratosthenes.  For
primes > 2*S we insert the squares in a minheap.

S should be small, almost surely a fraction of the L1D size.  That means
we we usually have the majority of the Pi(sqrt(n)) primes in a heap.

Representation is not clear.  Perhaps 19 bits for the heap keys, that
would alloow for the remaining 45 bits to be used by the prime.

-- 
Torbjörn