Category: Project Euler

Problem 12 on the Project Euler site asks:

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

A triangular number T(n) is defined as \(T(n) = \frac{n(n+1)}{2}\). The R code below consists of a solution, which involves the fact that the number of proper divisors of an integer n can be calculated by first computing a prime factorisation of the number n, e.g. if <\(n = p^aq^b\), where p,q are prime, then the number of proper divisors of n can be calculated as \(d(n) = (a+1)(b+1)\). This solution is extremely slow (mainly due to the naive prime sieving algorithm), and could be speeded up dramatically with a little effort. # Sieve of Eratosthenes
prime.sieve <- function(n) {
 a <- seq.int(1,n)
 p <- 1
 M <- as.integer(sqrt(n))
 while ((p <- p + 1) <= M) {
  if (a[p] != 0)
    a[seq.int(p*p, n, p)] <- 0
 }
 a[a>1 & a>0]
}
 

# Trial Division
# Returns the exponents of the prime
# factors of n
# e.g. if n = p^a*q^b
# tdiv(n) will return (a,b)
tdiv <- function(n) {
 primes <- prime.sieve(n)
 factors <- c()
 i <- 1
 curr <- 0
 
 for (p in primes) {
  while (n %% p == 0) {
   curr <- curr + 1
   n <- n %/% p
  }
  factors[i] <- curr
  i <- i + 1
  curr <- 0
 }
 
 factors[factors > 0]
}

# Compute nth triangular number
T <- function(n) {
        (n*(n+1))/2
}
 
## Problem 12
# This is a slooow solution
problem12 <- function(N) {
 n <- 0 # current triangular number Tn
 i <- 5 # \sum_{i=1}^n{i}
 
 while (TRUE) { 
  n <- T(i)
  factors <- tdiv(n) 
  if (prod(factors+1) >= N) {
   return(n)
  }
  i <- i + 1
 }
}

Problem 19 on the Project Euler website asks the user, given some initial information:

How many Sundays fell on the first of the month during the twentieth century (1 Jan 1901 to 31 Dec 2000)?

The obvious (but longer) way is to calculate the sum of the days between 1901 and 2000, given the number of days in each month, and a helper function to determine whether a year is a leap year or not:

is.leap <- function(year) {
        return (year %% 4 == 0 || (year %% 100 == 0 && year && 4 == 0))
}

# Problem 19
problem19 <- function() {
        daycount <- 1
        daylist <- list()
        i <- 1
        for (year in 1900:2000) {
                months <- c(0,31,28,31,30,31,30,31,31,30,31,30,31)
                if (is.leap(year)) {
                        months[3] <- 29
                        }
                        days <- daycount + cumsum(months)
                        daycount <- days[length(days)]
                        daylist[[i]] <- (days[–(length(days))])
                        i <- i + 1
                }
                sum(unlist(lapply(daylist[–1], function(x) {sum(x %% 7==1)} )))
}

However, with the aid of R’s chron library, there is a much easier way:

# Problem 19, method 2
library(chron)
sum(weekdays(seq.dates(”01/01/1901″, â€œ12/31/2000″, by=”months”))==”Sun”)

Problem 11 on Project Euler involves calculating the maximum product of adjacent numbers in any direction in a 20×20 matrix.

The solution below takes advantage of the symmetry of calculations to cut down on unnecessary loop operations:

problem11 < - function() {
    numbers <- scan("problem11.dat")
        m <- matrix(as.numeric(numbers), 20, byrow=TRUE)
        maxprd <- 0
        N <- 20; n <- 4
        prd1 <- 0; prd2 <- 0; prd3 <- 0
        dims <- dim(m)
        a <- (n-1)
        x <- c(0:a)
        for (i in 1:(dims[1])) {
            for (j in 1:(dims[2])) {
                prd1 <- ifelse(j <= N-a, prod(m[i,j+x]), 0) # row prod
                    prd2 <- ifelse(i <= N-a, prod(m[i+x,j]), 0) # column prod
# lower right diagonal
                    prd3 <- ifelse(i <= N-a && j <= N-a, prod(diag(m[i:(i+a),j:(j+a)])),0)
# lower left diagonal
                    prd4 <- ifelse(i <= N-a && j > a, prod(diag(m[i:(i+a),j:(j-a)])), 0)
                    maxprd < - max(prd1,prd2,prd3,prd4,maxprd)
            }
        }
    maxprd
}