Categories
Coding Project Euler R

Project Euler Problem #21

This is a solution for problem 21 on the Project Euler website. It consists of finding the sum of all the amicable numbers under 10000. This was pretty easy to solve, but the solution could probably be improved quite a bit.

Solution #1 in R is as follows (it calculates the proper divisors of each number using prop.divs, and then adds up the sequence of amicable numbers in the main function).

[sourcecode language=”matlab”]
prop.divs <- function(x) {
if (x == 1) return (1)
divs <- integer(30)
j <- 1
divs[j] <- 1
j <- j + 1
for (i in 2:(floor(x/2))) {
if ((x %% i) == 0) {
divs[j] <- i
j <- j + 1
}
}
sum(divs[1:(j-1)])
}

problem.21 <- function(N) {
s <- 0
for (i in 2:N) {
da <- prop.divs(i)
if (da == i) next
db <- prop.divs(da)
if ( db==i ) {
s <- s + da + db
}
}
s/2
}
[/sourcecode]

The s/2 is needed as each factor is added twice during the calculation.

This gives the correct answer, but the implementation is a bit naive. I remembered coming across an article about prime factors and proper divisors on PlanetMath a while ago, and this seemed like potentially a more efficient way to calculate the factors involved. Specifically, the sum of proper divisors of a number n can be given by:

\[ \prod_{i=1}^k\frac{p_i^{m_i+1}-1}{p_i – 1}-\prod_{i=1}^kp_i^{m_i} \]

The second attempt at this problem looked like the following:

[sourcecode language=”r”]
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]
}

sum.proper.divisors <- function(x) {
primes <- prime.sieve( x )
primes <- primes[ x %% primes == 0]

geo.sum <- numeric(length(primes))
i <- 1

for (prime in primes) {
n <- x
curr <- 0
while (n %% prime == 0) {
curr <- curr + 1
n <- n %/% prime
}
geo.sum[i] <- ( (prime^(curr+1) – 1)/(prime – 1) )
i <- i + 1
}
prod(geo.sum)-x
}

problem.21_2 <- function(N) {
s <- 0
for (i in 2:N) {
da <- sum.proper.divisors(i)
if (da == i) next
db <- sum.proper.divisors(da)
if (db==i) s <- s + da +db
}
s/2
}
[/sourcecode]

This also gives the correct answer, but with much reduced runtime overhead:


> system.time(problem.21(10000))
user system elapsed
103.943 0.511 106.978
> system.time(problem.21_2(10000))
user system elapsed
24.834 0.160 26.565

Categories
Coding Project Euler R

Project Euler Problem #28

Problem 28 on the Project Euler website asks what is the sum of both diagonals in a 1001×1001 clockwise spiral. This was an interesting one: the relationship between the numbers on the diagonals is easy to deduce, but expressing it succinctly in R took a little bit of tweaking. I’m sure it could be compressed even further.

[sourcecode lang=”r”]
spiral.size <- function(n) {
stopifnot(n %% 2 ==1)

if (n==1) {
return(1)
}
sum(cumsum(rep(2*seq(1,floor(n/2)), rep(4,floor(n/2))))+1)+1
}

spiral.size(1001)
[/sourcecode]

Categories
Coding Project Euler R

Project Euler Problem #22

Problem 22 on Project Euler proves a text file containing a large number of comma-delimited names and asks us to calculate the numeric sum of the alphabetical score for each name multiplied by the name’s position in the original list. This is made slightly easier by the presence of the predefined LETTERS variable in R.

problem22 <- function() {
        namelist <- scan(file="c:/temp/names.txt", sep=",", what="", na.strings="")
        sum(unlist(
                lapply(namelist, 
                        function(Z) which(namelist==Z) * sum(match(unlist(strsplit(Z,"")), LETTERS)))))
}