2

I'm can't figure out why the output of this comes out as:

>  File "<pyshell#177>", line 6, in func
>     list.append(next(PrimeGen)) 
>  StopIteration

when it makes so much sense in my head!! Anyways basically i'm trying to make a Prime Generator using an ifprime function and a generator to collect the primes in a list.

This determines if prime and returns the value if it is, otherwise nothing.

def prime(x):
    if (2**(x-1))%x ==1:
        return x

This makes the generator that should return a list full of prime numbers up to x but instead gives the error above. I started list with a 2 inside it because the above function prime(x) doesn't consider 2 as a prime (so range would start at 3)

def func(x):
  count=0
  list=[2]
  PrimeGen = (prime(X) for X in range(3,x+1))
  while count<99:
      list.append(next(PrimeGen))
      count+=1
  print list

Could anybody explain why it doesn't work? Thank you in advance! V.

Valentine Bondar
  • 109
  • 2
  • 10
  • related: [fastest way to list all primes below N in Python](http://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n-in-python) – jfs Aug 01 '12 at 19:47

3 Answers3

5

Fewer than 99 values were generated. Use itertools.islice() instead of a loop.

Ignacio Vazquez-Abrams
  • 699,552
  • 132
  • 1,235
  • 1,283
3

Note that your prime test

def prime(x):
    if (2**(x-1))%x ==1:
        return x

is wrong, that declares e.g. 341 = 11*31 and 2047 = 23*89 prime.

Also, for larger x, that generates very large intermediate values, much more efficient is

pow(2,x-1,x)

that reduces the intermediate values.

A moderately efficient implementation of a stronger primality check:

# The primes below 200 used as bases for the strong Fermat test,
# prime bases have more discriminatory power than composite bases,
# therefore we use prime bases first
bases = [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, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199]

# The strong Fermat test for base b, where n is odd and large,
# n - 1 = m * 2^s with odd m
# the strong test checks whether b**m % n == 1 or
# b**(2**j) % n == n-1 for a 0 <= j < s
# If the test returns False, n is definitely composite, otherwise probably prime
def strongFermat(b,n,m,s):
    a = pow(b,m,n)
    if a == 1:
        return True
    n1 = n-1
    for i in xrange(s):
        if a == n1:
            return True
        a = (a*a) % n
    return False

# Multiple strong Fermat tests, by default use 10 bases
# The probability that a composite passes is less than 0.25**iters
def sppTest(n, iters = 10):
    # Assumes n > 1 and with no prime divisors < 200
    m = n-1
    s = 0
    while (m & 1) == 0:
        m >>= 1
        s += 1
    pbases = iters if iters < 47 else 46
    for i in xrange(pbases):
        if not strongFermat(bases[i],n,m,s):
            return False
    if pbases < iters:
        for i in xrange(iters-pbases):
            if not strongFermat(211 + 2*i,n,m,s):
                return False
    return True

# Trial division to weed out most composites fast
def trialDivisionPrime(n):
    if n < 2:
        return 0        # Not a prime
    if n < 4:
        return 2        # definitely prime
    if n % 2 == 0 or n % 3 == 0:
        return 0        # definitely composite
    for d in xrange(5, 200, 6):
        if d*d > n:
            return 2    # definitely prime
        if n % d == 0 or n % (d+2) == 0:
            return 0    # composite
    return 1            # not yet decided

# The prime test, first do trial division by numbers < 200,
# if that doesn't decide the matter, use some strong Fermat tests
# using 20 tests is the paranoid setting for largish numbers,
# for numbers in 64-bit range, ten or fewer suffice
def prime(n):
    td = trialDivisionPrime(n)
    return td > 1 or (td == 1 and sppTest(n,20))

# just check a couple of larger numbers
for c in xrange(100000):
    if prime(c + 10**25):
        print c
Daniel Fischer
  • 174,737
  • 16
  • 293
  • 422
  • Oh thats strange, I used Fermat's little theorem I just assumed it was bulletproof... or did i translate it wrong to python? Thanks for the more efficient one though – Valentine Bondar Aug 01 '12 at 19:52
  • 3
    The translation is correct, but the theorem only states a _necessary_ condition for being prime. Composite integers having the tested property are known as "Fermat pseudoprimes". There are many more fast probabilistic prime tests, but the known definitive prime tests are still significantly slower. You can get a higher correctness rate by using more bases, and I'd recommend to switching to the strong Fermat test, that's more powerful and hardly more complex. – Daniel Fischer Aug 01 '12 at 19:58
1

This is a bit nicer (it seems kind of silly to have your "prime" function return either the number if it's prime or None otherwise, though it would actually work just fine in place of the slightly modified version below due to how bool(None) == False) and doesn't result in you getting StopIteration:

def isprime(x):
    return (2**(x-1))%x==1

def func(x):
    list=[2]
    list.extend(X for X in range(3,x+1) if isprime(X))
    print list

But going by what Daniel Fischer stated, your prime function is incorrect anyway.

JAB
  • 19,150
  • 4
  • 64
  • 78