Project Euler Problem 3

This post tackles Problem 3 in the Project Euler series.

Do you remember making Factor Trees in maths lessons? They are a neat way of showing how a composite number (a number that is not prime) is composed of prime factors. Suppose we wanted to factor the number 60. We might start by noticing that it is divisible by 2: 2 x 30 = 60. So we write

        60
       /  \
      2   30

2 is prime, so we're finished with that branch of the tree. Look at the other branch: we've now got to solve the same problem, but with a smaller number. 30 is also divisible by 2, so we add another layer to the tree:

        60
       /  \
      2   30
         /  \
        2    15

Applying the same procedure to 15 gives the final tree:

        60
       /  \
      2   30
         /  \
        2    15
            /  \
           3    5

The procedure stops when every branch of the tree has a prime number hanging from it. I'm not one for fortune telling, but reading these leaves will tell you something useful: the prime factorisation of the root number. In this case 60 = 2 x 2 x 3 x 5.

Problem 3 wants us to find the largest prime factor of a composite number. The Factor Trees methods suggests a nice recursive algorithm, and recursion has a very Functional feel to it. The algorithm is quite straight-forward.

We start with the number that we want to factorise, and feeling friendly, we call it n. We then try to find the first number bigger than 1 which divides n: call this number p. p must be prime: if it wasn't we'd have found another number first that divides n. We next cut n down to size by dividing it by p as may times as we can. Now for the recursive part: we start the procedure over again, this time with the smaller number that we found when dividing n by p. Any prime factor of this smaller number will be larger than p, and is also a factor of n so this trick will give us the answer we need.

There are just two special cases: first, there might not be a number p that divides n. In that case, n must be prime, so n is the answer we are looking for. Second, it might be that when we have divided n by p, we end up with 1, which won't factor any further: but in this case p is the largest prime that divides n, so p is the answer we need.

There's one other thing to note: when we are hunting for p we don't need to try all the numbers as far as n. We only need to go as far as √n. This is because n = √n x √n, so if it has a factor greater than √n, it will have a corresponding factor smaller than √n which we will find by searching only as far as √n.

We could be even cleverer and say that we only need to try dividing n by odd numbers (because if it is divisible by an even number, it will be divisible by 2, and we're chopping out 2 as a factor right at the beginning of the search). I want to keep the code simplish though, so for now I won't bother.

Now for some code:

public void Solve()
{
    long factor = FindGreatestPrimeFactor(1, 600851475143);
    Console.WriteLine(factor);
    Console.Read();
}

public long FindGreatestPrimeFactor(long factorGreaterThan, long number)
{
    long upperBound = (long)Math.Ceiling(Math.Sqrt(number));

    // find next factor of number
    long nextFactor = Range(factorGreaterThan + 1, upperBound)
        .SkipWhile(x => number % x > 0).FirstOrDefault();

    // if no other factor was found, then the number must be prime
    if (nextFactor == 0)
    {
        return number;
    }
    else
    {
        // find the multiplicity of the factor
        long multiplicity = Enumerable.Range(1, Int32.MaxValue)
            .TakeWhile(x => number % (long)Math.Pow(nextFactor, x) == 0)
            .Last();

        long quotient = number / (long)Math.Pow(nextFactor, multiplicity);

        if (quotient == 1)
        {
            return nextFactor;
        }
        else
        {
            return FindGreatestPrimeFactor(nextFactor, quotient);
        }
    }
}

private IEnumerable<long> Range(long first, long last)
{
    for (long i = first; i <= last; i++)
    {
        yield return i;
    }
}

A few words of explanation might be called for.

First, you'll recognise the Range method as being an iterator. I created this one, rather than using Enumerable.Range primarily because I needed to be able to specify first and last as long rather than int. I also took the opportunity to make it more convenient to use by making it take an last parameter, rather than a count as Enumerable.Range does.

In order to find the next factor, we generate a sequence of all the numbers greater than any factor we found previously, up to √number. SkipWhile is another Extension method that works on a sequence (it kind of does the opposite of TakeWhile which we used last time). It ignores the first part of a sequence, and only includes the part after the condition stops being met. In this case, it skips all numbers that don't divide number, and then returns the rest of the input sequence. FirstOrDefault(), as you might expect from its name takes a sequence, and returns the first element of it; or if there are no elements in the sequence, it returns a default value - the default for the Type of element in the sequence, in this case 0.

Finding the multiplicity of the prime factor is similarly straightforward (multiplicity in this case means, how many times the factor divides into number). We generate the sequence of numbers, m, where nextFactor^m divides number, and then take the last one in that sequence.

And there we have it. Not as compact a solution as for previous problems, but interesting non-the-less.