Project Euler 21: Sum of Amicable Numbers

Question

Problem 21:

Let \$d(n)\$ be defined as the sum of proper divisors of \$n\$ (numbers less than \$n\$ which divide evenly into \$n\$). If \$d(a) = b\$ and \$d(b) = a\,ドル where \$a ≠ b\,ドル then \$a\$ and \$b\$ are an amicable pair and each of \$a\$ and \$b\$ are called amicable numbers.

For example, the proper divisors of 220 are 1, 2, 4, 5, 10, 11, 20, 22, 44, 55 and 110; therefore \$d(220) = 284\$. The proper divisors of 284 are 1, 2, 4, 71 and 142; so \$d(284) = 220\$.

Evaluate the sum of all the amicable numbers under 10000.

Am I abusing reduce or itertools here? I'm concerned this solution could be more readable. Also, is there a more efficient solution rather than calculating everything up front?

from itertools import chain, count
from operator import mul
def factorize(n):
 for factor in chain((2,), count(start=3,step=2)):
 if factor*factor > n:
 break
 exp = 0
 while n % factor == 0:
 exp += 1
 n //= factor
 if exp > 0:
 yield factor, exp
 if n > 1:
 yield n, 1
def sum_of_factors(n):
 """
 >>> sum_of_factors(220)
 284
 >>> sum_of_factors(284)
 220
 """
 total = reduce(mul,
 ((fac**(exp+1)-1)/(fac-1) for fac,exp in factorize(n)),
 1)
 return total - n
if __name__ == '__main__':
 cache = {k: sum_of_factors(k)
 for k in xrange(1, 10000)
 }
 print sum(k for k, v in cache.iteritems()
 if cache.get(v, None) == k
 and v != k)

Answer 1

Why U No Doctest factorize()?

You wrote a doctest for sum_of_divisors, but why'd you forget about factorize()? Just because it's a generator doesn't mean you can't doctest it:

def factorize(n):
 """ 
 >>> list(factorize(10))
 [(2, 1), (5, 1)]
 >>> list(factorize(360))
 [(2, 3), (3, 2), (5, 1)]
 >>> list(factorize(991))
 [(991, 1)]
 """
 ...

Answer 2

Your code looks nice but I have a few suggestions to make it easier to understand according to my personal preferences.

---

It might be worth defining a function to encapsulate your call to reduce and mul (if you go forward in Project Euler problems, you'll reuse it many times) :

def mult(iterable, start=1):
 """Returns the product of an iterable - like the sum builtin."""
 return functools.reduce(operator.mul, iterable, start)

That you can easily use :

def sum_of_factors(n):
 """
 >>> sum_of_factors(220)
 284
 >>> sum_of_factors(284)
 220
 """
 total = mult((fac**(exp+1)-1)/(fac-1) for fac,exp in factorize(n))
 return total - n

---

Your function factorize generates tuples of the form (prime factor, power) which is probably worth documenting (either in the doc or in the function name).

I think it'd be easier to generate prime factors without grouping them as this can be delegated to a different function.

You'd have something like :

def factorize(n):
 """Generate prime factors (potentially repeated) in order."""
 for factor in chain((2,), count(start=3,step=2)):
 if factor*factor > n:
 break
 while n % factor == 0:
 n //= factor
 yield factor
 if n > 1:
 yield n

This function does one thing (and does it well) and can be easily tested :

for i in range(1, 1000):
 fact = list(factorize(i))
 assert mult(fact) == i # see, i've reused mult already :-)
 assert sorted(fact) == fact # there was a time where I needed that property but it's not really important now

Grouping and counting can be delected to a collection.Counter.

def sum_of_factors(n):
 """
 >>> sum_of_factors(220)
 284
 >>> sum_of_factors(284)
 220
 """
 total = mult((fac**(exp+1)-1)/(fac-1) for fac,exp in collections.Counter(factorize(n)).items())
 return total - n

---

The function sum_of_factors probably deserves some additional explanation (especially because the name says sum while the code says mult). Maybe a link to some reference would be enough.

---

You can avoid a few cache lookups. Indeed, at the moment, you consider all pairs (a, b) with a != b and add a whenever you find one. You could limit yourself to pairs (a, b) with a < b and add a+b when you find it. Also, by doing so, you know for sure that the value you are looking for will be in the cache.

sum_facts = {k: sum_of_factors(k) if k else 0 for k in range(10000) }
print(sum(k + v for k, v in sum_facts.items() if v < k and sum_facts[v] == k))

• Barry: I guess since the factors will be sorted, you could just use groupby instead of Counter. Last suggestion is solid.
• SylvainD: groupby would work indeed but it is slightly more complicated to use and I couldn't be bothered ;-)