Project Euler Problem 12 in F#

After a brief digression with an RPN Calculator, here we are again with a new Project Euler problem:

The sequence of triangle numbers is generated by adding the natural numbers.

So the 7^(th) 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 factors 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?

Again, this problem is relatively simple to solve and it's about finding an algorithm with decent performances.

The naive (slow) solution

#light

let triangles =
    Seq.unfold (fun (acc, state) -> Some (acc, (state + acc, state + 1))) (0, 1)
    |> Seq.skip 1 // skip the initial 0

let divides y x =
    x % y = 0

let rec all_factors_slow_rec n i acc =
    match i with
        | 1 -> 1::acc
        | x when divides x n -> all_factors_slow_rec n (x-1) (x::acc) 
        | _ -> all_factors_slow_rec n (i-1) acc
        
let all_factors_slow n =
    all_factors_slow_rec n n []

let find_index m =
    triangles
    |> Seq.map all_factors_slow
    |> Seq.tryFindIndex (fun x -> List.length x >= m)

let resolve_problem_12 =
    let i = find_index 500
    match i with
        | None -> failwith "Cannot resolve problem"
        | Some(i) ->
            let list = triangles |> Seq.take (i+1) |> List.ofSeq
            list.[i]

• Considering that triangle numbers are a mathematical sequence defined as follows:
  $$t_n = t_{n-1} + n$$ where \(t_0 = 0\), then we easily define an F# sequence that yields triangle numbers using Seq.unfold.
• We then define a function all_factors_slow_rec to calculate all factors (not just primes) of a given number. The slow solution presented in the first listing starts from a number \(n\) and then by brute force tries all the numbers down to \(1\) as a candidate factor.
• Once we have these two pieces together, it is just a matter of enumerating all possible triangle numbers, finding for each their factors (using Seq.map) and finally determining the first triangle number having at least 500 factors (using Seq.tryFindIndex). If \(i\) is the index of such a number, we then go back to the original sequence of triangles to extract the \(i^{th}\) element.

An improved version

let divides y x =
    x % y = 0

let triangles =
    Seq.unfold (fun (acc, state) -> Some (acc, (state + acc, state + 1))) (0, 1)
    |> Seq.skip 1 // skip the initial 0

let rec all_factors_quick_rec n i factors =
    if divides i n then
        let y = n / i
        if (i < y) then
            all_factors_quick_rec n (i + 1) (i::y::factors)
        elif (i = y) then
            // we have reached the square root value
            i::factors
        else
            factors
    elif i > int (sqrt (float n)) then
        // we are beyond the square root value
        factors
    else
        // try with the next number
        all_factors_quick_rec n (i + 1) factors
        
let all_factors_quick n =
    all_factors_quick_rec n 1 []

let find_index m =
    triangles
    |> Seq.map all_factors_quick
    |> Seq.tryFindIndex (fun x -> List.length x >= m)

let resolve_problem_12 =
    let i = find_index 500
        match i with
        | None -> failwith "Cannot resolve problem"
        | Some(i) ->
            let list = triangles |> Seq.take (i+1) |> List.ofSeq
            list.[i]

The bottleneck of the previous solution is the naive approach to find all the factors of a number. On my machines this solution takes several minutes.

Better performances can be obtained observing that in order to determine all the factors of a number we don't have try all smaller numbers. As a matter of fact, the following is true (considering n the number to be factorized):

• if \(n / x = y\), obviously both \(x\) and \(y\) are factors, so I don't have to probe y anymore
• If I start probing from 2 going up, I can stop probing when I have reached \(\sqrt n\). In fact, any number bigger than \(\sqrt n\) must have been already found probing a smaller number.

Suppose for example that we have to factor the number 36 (one of the triangle numbers):

• 1 is always a factor -> factors = {1, 36}
• 2 * 18 = 36 -> factors = {1, 36, 2, 18}
• 3 * 12 = 36 -> factors = {1, 36, 2, 18, 3, 12 }
• 4 * 9 = 36 -> factors = {1, 36, 2, 18, 3, 12, 4, 9}
• 5 does not divide -> factors = {1, 36, 2, 18, 3, 12, 4, 9}
• 6 * 6 = 36 -> factors = {1, 36, 2, 18, 3, 12, 4, 9, 6}

We can stop here because if 7 divided 36 evenly, we would have already found a number x < 6 such that 7 * x = 36.

The improved version runs on my machine in less than 1 second.