On Wednesday, Mar 23, 2005, at 04:00 America/Chicago, [EMAIL PROTECTED] wrote:


Now it would be fine to have an *equally fast*
infinite prime number generator.

Has anybody any suggestions?


I think when you add the condition of it being an infinite generator, you are changing the rules and can't end up with something as fast. What makes the sieve so fast is that we are generating all the primes in a given range using a very simple "strike out" method. In the infinite generator scenario, you will get all the primes up to n with the sieve and can yield those back, but at some point you will have to stop and get "the next one." To get the next one you will have to pick a range in which a next prime is likely to be found. The previous primes can be used to clear out this range.


What follows is an attempt based on the previous tutor-evolved sieve that extends the range in which to find the next prime by a factor of 2 over the last known prime. A similar algorithm is on ASPN (I bellieve), under

Space-efficient version of sieve of Eratosthenes.
D. Eppstein, May 2004

It's slower than the sieve but about twice as fast as Eppstein's up to 1,000,000. I'll leave it to someone else to see when it finally blows up :-) The output of primes was checked through 1,000,000 against Eppstein's generator without error.

/c

###
def esieve():
'''extended sieve generator that returns primes on each call. When the end of the
existing list is reached, more primes are sought in the range that extends to
twice the last known prime. The known primes are used to clear out this extended
range using the Sieve of Eratosthenes.'''


# get started with primes less than 100; you need at least [2, 3]
primes=[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37,
41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]
i=0
while 1:
yield primes[i]
i+=1
if i==len(primes): #extend range
minn=primes[-1]+2 #this is an odd number
maxx=2*minn+1 #there should be a prime in this range; +1 makes it odd
sqrt=int(maxx**.5) #don't use primes bigger than this
length = (maxx-minn)//2 #this is used for computing the crossing out None values
nums=range(minn,maxx+1,2) #here is the range in which a primes will be found
for p in primes:
if p>sqrt:break
j=minn%p #find out where the striking out should start
if j<>0:
j=(p-j) #if evens were present, this is where to start, but...
if (minn+j)%2==0:j+=p #if it lands on an even, go to the next one
j//=2 #and now account for the fact that evens aren't present
nums[j::p]=[None]*(1+(length-j)//p) #cross out multiples of p
x=filter(None,nums) #clean up the range
assert len(x)>0 #there should be a prime in here, but check anyway
primes.extend(x) #add these to the existing primes and continue yielding
###



_______________________________________________ Tutor maillist - Tutor@python.org http://mail.python.org/mailman/listinfo/tutor

Reply via email to