Re: [Tutor] primes (generator)

[Gregor Lingl]

Hi all of you!
Many thanks for your remarks, ideas and experiments.
When I posted my input yesterday:
[x for x in range(2,100) if not [y for y in range(2,x) if x%y==0]]
my intention was indeed to find a "shortest" a program, disregarding
efficiency. Thus range instead of xrange etc.
(About one year ago I was involved in a "coding contest" at Linux-Tag
in Karlsruhe, where the goal was to achive a certain output with a
shortest program ;-) Strange. I used Python, but I didn't win :-(  )
Clearly there are several criteria concerning the quality of programs
as e. g. conciseness, clarity, efficiency, beauty(!) etc. which may not
be achieved optimally all together ...
I definitely intend to sum up your contibutions (for me) later, but for
now (because of lack of time) I'll only present one more idea,
induced by Danny's (unfortunately buggy - because of recursion overflow) 
contribution:

Thought it would be nice to use Python's new generator-expressions
instead of itertool module. First idea:
##################################################################
# infinite prime generator
# induced by ideas of Danny Yoo,
# which themselves were induced by material of SICP
import math
def prime(n):
    return [x for x in range(2,1+int(math.sqrt(n))) if n%x==0]==[]
def infrange(n):
    while True:
        yield n
        n+=1
primes = (p for p in infrange(2) if prime(p))
###############  The generator primes works like this:
>>> primes.next()
2
>>> primes.next()
3
>>> primes.next()
5
>>> [primes.next() for i in range(10)]
[7, 11, 13, 17, 19, 23, 29, 31, 37, 41]
>>> [primes.next() for i in range(100)]
[43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 
113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 
193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 
271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 
359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 
443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 
541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617]
>>>
>>> ### You may also use it this way:
>>> for p in primes:
	if p > 1000: break
	print p,

	
619 631 641 643 647 653 659 661 673 677 683 691 701 709 719 727 733 739 
743 751 757 761 769 773 787 797 809 811 821 823 827 829 839 853 857 859 
863 877 881 883 887 907 911 919 929 937 941 947 953 967 971 977 983 991 997

With a fresh generator
[primes.next() for p in xrange(10000)]
produces (on my machine) a list of the first 10000 primes (upto 104729)
within 4.68 secs.
Now for a bit of optimization: instead of building a list of
all divisors use a loop and get out as soon as possible:
So I replace prime by:
def prime(n):
    for d in range(2,1+int(math.sqrt(n))):
        if n%d == 0: return  False
    return True
This makes the generator nearly 4 times as fast. I get my list in
1.33 secs.
Lets optimize a bit more: test only 2 and odd numbers:
def infrange(n):
    yield 2
    n=3
    while True:
        yield n
        n+=2
and do a similar thing for searching for divisors:
def prime(n):
    if n==2: return True
    for d in [2]+ range(3,1+int(math.sqrt(n)),2):
        if n%d == 0: return  False
    return True
Now the first 10000  primes get computed in 0.67 secs.
Again double speed. Together seven times as fast as
first version. But IMO at the cost of clarity
and conciseness.
Maybe you have amendments (of course a matter of taste)
to this code. I'd be interested in them ..
Maybe sombody likes to compare these algorithms to the
ones mentioned before in this thread. I would be
interested in the results.
Regards,
Gregor

Danny Yoo schrieb:
On Thu, 17 Mar 2005, Gregor Lingl wrote:

Hi!
Who knows a more concise or equally concise but more efficient
expression, which returns the same result as
...
######
from itertools import ifilter, count
def notDivisibleBy(n):
...     def f(x):
...         return x % n != 0
...     return f
...
def sieve(iterable):
...     nextPrime = iterable.next()
...     yield nextPrime
...     restOfPrimes = sieve(ifilter(notDivisibleBy(nextPrime), iterable))
...     for prime in restOfPrimes:
...         yield prime
...
primes = sieve(count(2))
primes.next()
2
primes.next()
...
primes.next()
23
######
The neat thing about this code is that it produces an infinite iterator of
prime numbers.  The original code in Scheme is itself quite concise and
quite nice.  It is described here:
http://mitpress.mit.edu/sicp/full-text/book/book-Z-H-24.html#call_footnote_Temp_455
Best of wishes to you!

--
Gregor Lingl
Reisnerstrasse 3/19
A-1030 Wien
Telefon: +43 1 713 33 98
Mobil:   +43 664 140 35 27
Autor von "Python für Kids"
Website: python4kids.net
_______________________________________________
Tutor maillist  -  Tutor@python.org
http://mail.python.org/mailman/listinfo/tutor