Add primes function from projecteuler

Author: doctaphred

maths.py

@@ -1,4 +1,35 @@
-from itertools import product
+from itertools import count, product
+
+
+def primes():
+    """Yield prime numbers, starting with 2.
+
+    Algorithm by David Eppstein, Alex Martelli, and Tim Hochberg:
+    http://code.activestate.com/recipes/117119/#c2
+
+    >>> p = primes()
+    >>> [next(p) for _ in range(10)]
+    [2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
+    """
+    yield 2
+    # Map composites to primes witnessing their compositeness
+    composites = {}
+    # Skip even numbers
+    for n in count(3, step=2):
+        # We won't see n again, so we can delete its witness, if any
+        prime_divisor = composites.pop(n, None)
+        if prime_divisor is None:
+            # n is prime
+            yield n
+            # Record n as a divisor of its square
+            composites[n ** 2] = 2 * n
+        else:
+            # n is composite
+            # Move the witness to a new multiple
+            x = n + prime_divisor
+            while x in composites:
+                x += prime_divisor
+            composites[x] = prime_divisor
 
 
 def reordered_digit_map(exponents, base=2):