Two Factoring Games
February 18, 2011
Over at mersenneforum.org, a small community of people interested in factoring large integers talk about their craft, share tips and accomplishments, and sometimes play factoring games. Today’s exercise is about two of those games.
The home prime of a number n is computed by factoring the number into its prime factors, concatenating the digits of the prime factors (which have been sorted into ascending order), and repeating until the result is prime; this prime number is known as the home prime. For instance, 99 factors as 3 · 3 · 11; 3311 factors as 7 · 11 · 43, and 71143 is prime, so HP(99) = 71143, computed in three steps. Sloane has the list of home primes at A037274. Note that the home prime of a prime number is the number itself. Many numbers resolve to their home prime in just a few steps, as above, but others take longer, and it not known if all numbers eventually resolve to a home prime. At mersenneforum, and at the World of Numbers, Alex Gruppa, Paul Leyland, and Nicolas Daminelli have been working on the calculation of HP49 for ten years, and are currently stalled, after 105 steps, by a 210-digit composite.
Another factoring game involves the Euclid-Mullin sequence, which starts with a1 = 2 and continues 
Your task is to write two functions to compute the home prime of a given number and to compute the Euclid-Mullin sequence. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.
Programming Praxis 
[…] today’s Programming Praxis exercise, our goal is to implement two algorithms related to prime factors; one […]
My Haskell solution (see http://bonsaicode.wordpress.com/2011/02/18/programming-praxis-two-factoring-games/ for a version with comments):
import Data.List import Data.Numbers.Primes homePrime :: Integer -> Integer homePrime = head . filter isPrime . iterate (read . concatMap show . primeFactors) euclidMullin :: [Integer] euclidMullin = unfoldr (\p -> let a = head (primeFactors $ p + 1) in Just (a, p * a)) 1My solution isn’t very long, but the machinery behind the scenes is rather
extensive. In an effort to save memory, I’ve tried to make heavy use of
iterators—sort of like lazy evaluation of lists in Haskell. Below is what I
wrote for this exercise; I also used an implementation of the Sieve of
Eratosthenes for primes and a Miller-Rabin primality test.
The full code is available here.
def times(x, n): t = 0 while not n % x: t += 1 n //= x return t from itertools import chain def prime_factors(n): pfs = () for p in primes(n): if not n % p: pfs = chain(pfs, [p for i in xrange(times(p, n))]) return pfs def home_prime(n): while not is_prime(n): n = int(''.join(str(p) for p in prime_factors(n))) return n def lpf(n): return next(prime_factors(n)) def euclid_mullin(): a = 2 p = 1 while True: yield a p *= a a = lpf(1 + p)