Problem 12 has us calculating divisors, but this time not of the prime variety.

The sequence of triangle numbers is generated by adding the natural numbers. So the 7th triangle number would be 1 + 2 + 3 + 4 + 5 + 6 + 7 = 28. The first ten terms would be:

1, 3, 6, 10, 15, 21, 28, 36, 45, 55, …

Let us list the divisors of the first seven triangle numbers:

 1: 1
 3: 1,3
 6: 1,2,3,6
10: 1,2,5,10
15: 1,3,5,15
21: 1,3,7,21
28: 1,2,4,7,14,28

We can see that 28 is the first triangle number to have over five divisors.

What is the value of the first triangle number to have over five hundred divisors?

In fact, we do need to examine the prime factors. Rather than going to horrible brute force way of checking every number from 1 to sqrt(n) for divisibility and counting factors that way, we can actually work it out.

The prime factors of a number have both a value and an exponent. Multiply the numbers raised to the power of their exponents and you get the original number. For example, the number 20 has six divisors; 1, 2, 4, 5, 10, 20. The prime factors of 20 are 2, 2 and 5. That is to say, 2^2 * 5^1 = 20.

If we take product of the exponents incremented by 1 we get the number of divisors of the original number. In this case 3 * 2 = 6.

If we calculate an infinite series of triangle numbers and have this quick way of finding a divisor count the rest should be easy.

    from collections import Counter
    from itertools import count
    from euler.prime import prime_factors
    from euler.math import product


    def triangle_numbers():
        x = 0
        for y in count(1):
            x = x + y
            yield x


    def num_divisors(n):
        # The number of divisors is calculated from the prime factors
        # We know 20 has divisors of 1, 2, 4, 5, 10 and prime factors
        # of 2^2 and 5^1. By taking each exponent, adding 1, and
        # multiplying them we get the total number of divisors. So:
        #     (2+1) * (1+1) = 3 * 2 = 6
        fs = dict(Counter(prime_factors(n)))
        return product(x + 1 for x in fs.values())


    def first_with_n_divisors(seq, n):
        return next(x for x in seq if num_divisors(x) >= n)
    

That would be spot on if it didn’t take 1.5 seconds to run. Another tactic, then.

An interesting thing about triangle numbers is that you can calculate the number of divisors without calculating the actual number itself. This is because of the way they are constructed; the nth triangle number is the sum of every positive integer up to and including n.

To calculate the number of divisors for the nth triangle number we need to be able to calculate prime factors. We can already do that, so the formulae are:

• for even values of n: num_divisors(n // 2) * num_divisors(n + 1)
• for odd values of n: num_divisors(n) * num_divisors(n // 2 + 1)

Once we know the index of the triangle number that has 500 divisors, we can just sum every positive integer up to that number. Note that num_divisors seemed so useful that I moved it into euler.math.

    from itertools import count
    from euler.math import num_divisors, divides_by
    from euler.iter import first, get_at


    def triangle_numbers():
        x = 0
        for y in count(1):
            x = x + y
            yield x


    def num_divisors_for_nth_triangle_number(n):
        if divides_by(n, 2):
            return num_divisors(n // 2) * num_divisors(n + 1)
        else:
            return num_divisors(n) * num_divisors(n // 2 + 1)


    def first_with_n_divisors(seq, t):
        i = first(n for n in count(1) if num_divisors_for_nth_triangle_number(n) >= t)
        return get_at(seq, i)


    def test_0012_highly_divisible_triangle_number():
        assert first_with_n_divisors(triangle_numbers(), 500) == 76576500
    

This gets the answer in 430ms, which is a significant speed increase. Not good enough, yet, but a lot better. The get_at function is just:

    def get_at(n, seq):
        return last(take(seq, n))
    

An interesting thing to note about the more optimised solution is that we are actually calculating almost every set of divisors at least twice. Consider which sets of divisors are calculated for each value of n:

    values of n   divisors calculated for
    2             1, 3
    3             3, 2
    4             2, 5
    5             5, 3
    

Even in this short range it’s obvious that each iteration uses one of the calculations from the previous iteration. We can halve the number of calculations by keeping track of the previous result and using it again.

    from itertools import count
    from euler.math import divides_by, num_divisors


    def first_with_n_divisors(t):
        x = 1
        y = 1
        for n in count(2):
            if divides_by(n, 2):
                y = num_divisors(n + 1)
            else:
                x = num_divisors(n // 2 + 1)
            if x * y >= t:
                return sum(range(1, n+1))
    

This version results the correct result of 76576500 in 254ms. I’ve had a sneaky look at the Project Euler thread for this one and the only ones that claim to be faster in a high level language like Python are either nonsense (e.g. claiming Delphi to be faster than C and ASM) or cheating (starting with a value of n very close to the actual answer). If I start my counter with n = 12000 then it finishes in 10ms, but one would have to know the answer in order to know where to start the counter. It is therefore cheating.