Sunday, January 15, 2012

Project Euler in R: Problem 23

I was just looking through the programming language statistics on Project Euler. It shows that only 7% of the problems have been solved in R, whereas 8% have been solved on any kind of spreadsheet.  This is outrageous!

Let's look at the solution of problem 23 in R. Here is the statement of Project Euler's problem 23:
A perfect number is a number for which the sum of its proper divisors is exactly equal to the number. For example, the sum of the proper divisors of 28 would be 1 + 2 + 4 + 7 + 14 = 28, which means that 28 is a perfect number.
A number n is called deficient if the sum of its proper divisors is less than n and it is called abundant if this sum exceeds n.
As 12 is the smallest abundant number, 1 + 2 + 3 + 4 + 6 = 16, the smallest number that can be written as the sum of two abundant numbers is 24. By mathematical analysis, it can be shown that all integers greater than 28123 can be written as the sum of two abundant numbers. However, this upper limit cannot be reduced any further by analysis even though it is known that the greatest number that cannot be expressed as the sum of two abundant numbers is less than this limit.
Find the sum of all the positive integers which cannot be written as the sum of two abundant numbers.
Find proper divisors
First we write a function to find the proper divisors (that is divisors less than the number itself) of a number n.
divisor <- function (n) {
    divs <- c(1)
    i = 2
    lim = floor (sqrt(n))
    while (i <= lim) {
        if (n %% i == 0) {
          divs <- c(divs, i)
          if (n/i != i) { divs <- c(divs, n/i) }
        i = i + 1
    return (divs)

Check if number is abundant
Next we create a simple function which checks whether a number n is abundant.
abundant <- function (n) {
    if (sum (divisor(n)) > n) return (TRUE) else return (FALSE)

The solution
With the above two functions we have the tools we need to solve the main problem: finding the sum of all positive integers which can not be written as the sum of two abundant numbers.

All the numbers greater than 28,123 can be written as sum of two abundant numbers, so we only need to check numbers smaller than 28,123 for this property.  The vector abunlist in the code below holds all abundant numbers starting from 12 to 28,123.  With every new abundant number than we find, we create all possible two abundant number sums and store them in the vector sumno.  The answer we're looking for is the sum of all numbers 1 to 28,123 except the ones in the vector sumno.

abunlist <- c(12)            # The smallest abundant number = 12
sumno <- c(24)               # The smallest abundant number sum = 12 + 12 = 24
for (i in 13:28123) {
  if (abundant(i)) {
    abunlist <- c(abunlist, i)                   # Latest abundant number = i
    tmp <- abunlist[(abunlist + i) <= 28123] + i   # New abundant sums <= 28123 using i
    sumno <-unique(c(sumno, tmp))                # Add to existing list

This is certainly not the most efficient solution, but runs in approximately 25 seconds - well within the stated 1 minute goal.

Note: The code for this post was formatted in Pretty R at