# from: http://r...content-available-to-author-only...e.org/wiki/Sieve_of_Eratosthenes#Fast_infinite_generator_using_a_wheel
def primes():
for p in [2,3,5,7]: yield p # base wheel primes
gaps1 = [ 2,4,2,4,6,2,6,4,2,4,6,6,2,6,4,2,6,4,6,8,4,2,4,2,4,8 ]
gaps = gaps1 + [ 6,4,6,2,4,6,2,6,6,4,2,4,6,2,6,4,2,4,2,10,2,10 ] # wheel2357
def wheel_prime_pairs():
yield (11,0); bps = wheel_prime_pairs() # additional primes supply
p, pi = next(bps); q = p * p # adv to get 11 sqr'd is 121 as next square to put
sieve = {}; n = 13; ni = 1 # into sieve dict; init cndidate, wheel ndx
while True:
if n not in sieve: # is not a multiple of previously recorded primes
if n < q: yield (n, ni) # n is prime with wheel modulo index
else:
npi = pi + 1 # advance wheel index
if npi > 47: npi = 0
sieve[q + p * gaps[pi]] = (p, npi) # n == p * p: put next cull position on wheel
p, pi = next(bps); q = p * p # advance next prime and prime square to put
else:
s, si = sieve.pop(n)
nxt = n + s * gaps[si] # move current cull position up the wheel
si = si + 1 # advance wheel index
if si > 47: si = 0
while nxt in sieve: # ensure each entry is unique by wheel
nxt += s * gaps[si]
si = si + 1 # advance wheel index
if si > 47: si = 0
sieve[nxt] = (s, si) # next non-marked multiple of a prime
nni = ni + 1 # advance wheel index
if nni > 47: nni = 0
n += gaps[ni]; ni = nni # advance on the wheel
for p, pi in wheel_prime_pairs(): yield p # strip out indexes
ps = primes() # 100k: 0.50s-9.4MB
for n in range(1, 1000000): next(ps) # 200k: 1.04s-9.4MB n^1.06 !!
print(next(ps)) # 400k: 2.07s-9.4MB n^1.00 !!!
# 1m: 5.56s-9.4MB n^1.08 !!! 15485863
IyBmcm9tOiAgaHR0cDovL3IuLi5jb250ZW50LWF2YWlsYWJsZS10by1hdXRob3Itb25seS4uLmUub3JnL3dpa2kvU2lldmVfb2ZfRXJhdG9zdGhlbmVzI0Zhc3RfaW5maW5pdGVfZ2VuZXJhdG9yX3VzaW5nX2Ffd2hlZWwKCmRlZiBwcmltZXMoKToKICAgIGZvciBwIGluIFsyLDMsNSw3XTogeWllbGQgcCAgICAgICAgICAgICAgICAgIyBiYXNlIHdoZWVsIHByaW1lcwogICAgZ2FwczEgPSBbIDIsNCwyLDQsNiwyLDYsNCwyLDQsNiw2LDIsNiw0LDIsNiw0LDYsOCw0LDIsNCwyLDQsOCBdCiAgICBnYXBzID0gZ2FwczEgKyBbIDYsNCw2LDIsNCw2LDIsNiw2LDQsMiw0LDYsMiw2LDQsMiw0LDIsMTAsMiwxMCBdICMgd2hlZWwyMzU3CiAgICBkZWYgd2hlZWxfcHJpbWVfcGFpcnMoKToKICAgICAgICB5aWVsZCAoMTEsMCk7IGJwcyA9IHdoZWVsX3ByaW1lX3BhaXJzKCkgIyBhZGRpdGlvbmFsIHByaW1lcyBzdXBwbHkKICAgICAgICBwLCBwaSA9IG5leHQoYnBzKTsgcSA9IHAgKiBwICAgICAgICAgICAgIyBhZHYgdG8gZ2V0IDExIHNxcidkIGlzIDEyMSBhcyBuZXh0IHNxdWFyZSB0byBwdXQKICAgICAgICBzaWV2ZSA9IHt9OyBuID0gMTM7IG5pID0gMSAgICAgICAgICAgICAgIyAgIGludG8gc2lldmUgZGljdDsgaW5pdCBjbmRpZGF0ZSwgd2hlZWwgbmR4CiAgICAgICAgd2hpbGUgVHJ1ZToKICAgICAgICAgICAgaWYgbiBub3QgaW4gc2lldmU6ICAgICAgICAgICAgICAgICAgIyBpcyBub3QgYSBtdWx0aXBsZSBvZiBwcmV2aW91c2x5IHJlY29yZGVkIHByaW1lcwogICAgICAgICAgICAgICAgaWYgbiA8IHE6IHlpZWxkIChuLCBuaSkgICAgICAgICAjIG4gaXMgcHJpbWUgd2l0aCB3aGVlbCBtb2R1bG8gaW5kZXgKICAgICAgICAgICAgICAgIGVsc2U6CiAgICAgICAgICAgICAgICAgICAgbnBpID0gcGkgKyAxICAgICAgICAgICAgICAgICMgYWR2YW5jZSB3aGVlbCBpbmRleAogICAgICAgICAgICAgICAgICAgIGlmIG5waSA+IDQ3OiBucGkgPSAwCiAgICAgICAgICAgICAgICAgICAgc2lldmVbcSArIHAgKiBnYXBzW3BpXV0gPSAocCwgbnBpKSAjIG4gPT0gcCAqIHA6IHB1dCBuZXh0IGN1bGwgcG9zaXRpb24gb24gd2hlZWwKICAgICAgICAgICAgICAgICAgICBwLCBwaSA9IG5leHQoYnBzKTsgcSA9IHAgKiBwICAjIGFkdmFuY2UgbmV4dCBwcmltZSBhbmQgcHJpbWUgc3F1YXJlIHRvIHB1dAogICAgICAgICAgICBlbHNlOgogICAgICAgICAgICAgICAgcywgc2kgPSBzaWV2ZS5wb3AobikKICAgICAgICAgICAgICAgIG54dCA9IG4gKyBzICogZ2Fwc1tzaV0gICAgICAgICAgIyBtb3ZlIGN1cnJlbnQgY3VsbCBwb3NpdGlvbiB1cCB0aGUgd2hlZWwKICAgICAgICAgICAgICAgIHNpID0gc2kgKyAxICAgICAgICAgICAgICAgICAgICAgIyBhZHZhbmNlIHdoZWVsIGluZGV4CiAgICAgICAgICAgICAgICBpZiBzaSA+IDQ3OiBzaSA9IDAKICAgICAgICAgICAgICAgIHdoaWxlIG54dCBpbiBzaWV2ZTogICAgICAgICAgICAgIyBlbnN1cmUgZWFjaCBlbnRyeSBpcyB1bmlxdWUgYnkgd2hlZWwKICAgICAgICAgICAgICAgICAgICBueHQgKz0gcyAqIGdhcHNbc2ldCiAgICAgICAgICAgICAgICAgICAgc2kgPSBzaSArIDEgICAgICAgICAgICAgICAgICMgYWR2YW5jZSB3aGVlbCBpbmRleAogICAgICAgICAgICAgICAgICAgIGlmIHNpID4gNDc6IHNpID0gMAogICAgICAgICAgICAgICAgc2lldmVbbnh0XSA9IChzLCBzaSkgICAgICAgICAgICAjIG5leHQgbm9uLW1hcmtlZCBtdWx0aXBsZSBvZiBhIHByaW1lCiAgICAgICAgICAgIG5uaSA9IG5pICsgMSAgICAgICAgICAgICAgICAgICAgICAgICMgYWR2YW5jZSB3aGVlbCBpbmRleAogICAgICAgICAgICBpZiBubmkgPiA0Nzogbm5pID0gMAogICAgICAgICAgICBuICs9IGdhcHNbbmldOyBuaSA9IG5uaSAgICAgICAgICAgICAjIGFkdmFuY2Ugb24gdGhlIHdoZWVsCiAgICBmb3IgcCwgcGkgaW4gd2hlZWxfcHJpbWVfcGFpcnMoKTogeWllbGQgcCAgICMgc3RyaXAgb3V0IGluZGV4ZXMKCnBzID0gcHJpbWVzKCkgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIyAxMDBrOiAwLjUwcy05LjRNQgpmb3IgbiBpbiByYW5nZSgxLCAxMDAwMDAwKTogbmV4dChwcykgICAgICAgICMgMjAwazogMS4wNHMtOS40TUIgICBuXjEuMDYgISEKcHJpbnQobmV4dChwcykpICAgICAgICAgICAgICAgICAgICAgICAgICAgICAjIDQwMGs6IDIuMDdzLTkuNE1CICAgbl4xLjAwICEhIQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICMgMW06ICAgNS41NnMtOS40TUIgICBuXjEuMDggISEhICAgMTU0ODU4NjM=