Thursday, March 27, 2014

Problem 78 - R

So, after solving problem 45 in Rexx, I had R open and free for use. Problem 78 was a bit annoying for me. Problem 76, I decided to solve using dynamic programming as it was pretty easy to do so, and that solution generalized really well to problem 77. There was another way to solve problem 76, and that was by directly calculating the partition function called p(n) in this problem. I realized that this was what I had to do, but mathworld and wikipedia are not the best place to learn this things (and my copy of TAOPC is up at college, far from my hands), and they were not very clear what conventions were being used or what bounds were reasonable for the sum that must be calculated. My first solution took 67 minutes in Python, and then after the change to the loop bounds, in took 9 seconds. After I ported my solution to R, it took about 4 minutes to calculate the answer, but I am okay with that, according to the first estimate I looked at, about 20x slowdowns are expected moving to R from Python, and Python is already about 20x slower than C, so 4 minutes in a language like R feels like it is fair to consider it within Project Euler's one minute rule. Anyway, here is my solution, which took a while to debug off-by-one errors resulting from arrays starting with index 1.
N = 100000
p <- array(0, dim=c(N))

p[1] <- 1
p[2] <- 1
for (n in 2:N) {
    n1 <- 1
    n2 <- 1
    k  <- 1
    while (n1 > 0 || n2 > 0) {
        n1 <- n - (k*(3*k-1)) %/% 2
        n2 <- n - (k*(3*k+1)) %/% 2
        c  <- 0
        if (n1 >= 0) {
            c <- p[n1+1]
        }
        if (n2 >= 0) {
            c <- c + p[n2+1]
        }
        if (k %% 2 == 0) {
            c <- c * -1
        }
        p[n+1] <- (p[n+1] + c) %% 1000000
        k <- k + 1
    }
    if (p[n+1] == 0) {
        print(n)
        break
    }
}

No comments:

Post a Comment