This one isn’t even funny…

Starting in the top left corner of a 2×2 grid, there are 6 routes (without backtracking) to the bottom right corner.

How many routes are there through a 20×20 grid?

Your first thought would be to generate the routes, but for a 20×20 grid, those amount to BILLIONS and you’d try to do it recursively too! Forget it.

But if you have some CompSci-level math background, though, you’ll remember this one—reading The Art of Computer Programming, Vol. 4 also helps. It’s a matter of combinatorics and if we take w for the width and h for the height, all we need to do is calculate,

[latex size=“3”]frac{(w + h)!}{w!h!}[/latex]

and that’s it. I didn’t even write a program for this, instead using Python’s interactive shell, but just for the same of completeness, here’s a “program” to calculate the possibles paths in a 20×20 grid.

import math
w = h = 20
print math.factorial(w + h) / (math.factorial(w) * math.factorial(h))

That’s it.

Project Euler’s problem #11 statement goes:

In the 20×20 grid below, four numbers along a diagonal line have been marked in red.

08 02 22 97 38 15 00 40 00 75 04 05 07 78 52 12 50 77 91 08
49 49 99 40 17 81 18 57 60 87 17 40 98 43 69 48 04 56 62 00
81 49 31 73 55 79 14 29 93 71 40 67 53 88 30 03 49 13 36 65
52 70 95 23 04 60 11 42 69 24 68 56 01 32 56 71 37 02 36 91
22 31 16 71 51 67 63 89 41 92 36 54 22 40 40 28 66 33 13 80
24 47 32 60 99 03 45 02 44 75 33 53 78 36 84 20 35 17 12 50
32 98 81 28 64 23 67 10 26 38 40 67 59 54 70 66 18 38 64 70
67 26 20 68 02 62 12 20 95 63 94 39 63 08 40 91 66 49 94 21
24 55 58 05 66 73 99 26 97 17 78 78 96 83 14 88 34 89 63 72
21 36 23 09 75 00 76 44 20 45 35 14 00 61 33 97 34 31 33 95
78 17 53 28 22 75 31 67 15 94 03 80 04 62 16 14 09 53 56 92
16 39 05 42 96 35 31 47 55 58 88 24 00 17 54 24 36 29 85 57
86 56 00 48 35 71 89 07 05 44 44 37 44 60 21 58 51 54 17 58
19 80 81 68 05 94 47 69 28 73 92 13 86 52 17 77 04 89 55 40
04 52 08 83 97 35 99 16 07 97 57 32 16 26 26 79 33 27 98 66
88 36 68 87 57 62 20 72 03 46 33 67 46 55 12 32 63 93 53 69
04 42 16 73 38 25 39 11 24 94 72 18 08 46 29 32 40 62 76 36
20 69 36 41 72 30 23 88 34 62 99 69 82 67 59 85 74 04 36 16
20 73 35 29 78 31 90 01 74 31 49 71 48 86 81 16 23 57 05 54
01 70 54 71 83 51 54 69 16 92 33 48 61 43 52 01 89 19 67 48

The product of these numbers is [pmath]26 * 63 * 78 * 14 = 1788696[/pmath].

What is the greatest product of four adjacent numbers in any direction (up, down, left, right, or diagonally) in the 20×20 grid?

This one is remarkably easy but also was quite fun. I think it’s because it reminds me of the kind of work we’d do during our Algorithms classes during my first year in college. And so this one goes to my Algorithms teacher, Ricardo Vargas Dornelles—best teacher I’ve ever had too.

import operator
numbers = [[8,2,22,97,38,15,0,40,0,75,4,5,7,78,52,12,50,77,91,8],
[49,49,99,40,17,81,18,57,60,87,17,40,98,43,69,48,4,56,62,0],
[81,49,31,73,55,79,14,29,93,71,40,67,53,88,30,3,49,13,36,65],
[52,70,95,23,4,60,11,42,69,24,68,56,1,32,56,71,37,2,36,91],
[22,31,16,71,51,67,63,89,41,92,36,54,22,40,40,28,66,33,13,80],
[24,47,32,60,99,3,45,2,44,75,33,53,78,36,84,20,35,17,12,50],
[32,98,81,28,64,23,67,10,26,38,40,67,59,54,70,66,18,38,64,70],
[67,26,20,68,2,62,12,20,95,63,94,39,63,8,40,91,66,49,94,21],
[24,55,58,5,66,73,99,26,97,17,78,78,96,83,14,88,34,89,63,72],
[21,36,23,9,75,0,76,44,20,45,35,14,0,61,33,97,34,31,33,95],
[78,17,53,28,22,75,31,67,15,94,3,80,4,62,16,14,9,53,56,92],
[16,39,5,42,96,35,31,47,55,58,88,24,0,17,54,24,36,29,85,57],
[86,56,0,48,35,71,89,7,5,44,44,37,44,60,21,58,51,54,17,58],
[19,80,81,68,5,94,47,69,28,73,92,13,86,52,17,77,4,89,55,40],
[4,52,8,83,97,35,99,16,7,97,57,32,16,26,26,79,33,27,98,66],
[88,36,68,87,57,62,20,72,3,46,33,67,46,55,12,32,63,93,53,69],
[4,42,16,73,38,25,39,11,24,94,72,18,8,46,29,32,40,62,76,36],
[20,69,36,41,72,30,23,88,34,62,99,69,82,67,59,85,74,4,36,16],
[20,73,35,29,78,31,90,1,74,31,49,71,48,86,81,16,23,57,5,54],
[1,70,54,71,83,51,54,69,16,92,33,48,61,43,52,1,89,19,67,48]]

maxp = 0
for row in range(0,16):
        for col in range(0,16):
                maxp = max(maxp, reduce(lambda x,y: x*y, 
                                        [n for n in numbers[row][col:col+4]]))
                maxp = max(maxp, numbers[row][col] * numbers[row + 1][col + 1] * 
                                 numbers[row+2][col+2] * numbers[row+3][col+3])
                maxp = max(maxp, reduce(lambda x,y: x*y, 
                                        [n[col] for n in numbers[row:row+4]]))
                if col > 3:
                        maxp = max(maxp, numbers[row][col] * 
                                         numbers[row + 1][col - 1] * 
                                         numbers[row+2][col-2] * 
                                         numbers[row+3][col-3])

print maxp

As well, this runs in 0.020s, so not bad at all.

I decided to take on Project Euler’s problem #10. Its statement goes like this:

The sum of the primes below 10 is [pmath]2 + 3 + 5 + 7 = 17[/pmath].

Find the sum of all the primes below two million.

My first attempt used brute force and testing primality using only previously found primes:

primes = []

def is_prime(n):
        if not (n < 2 or any(n % x == 0 
           for x in filter(lambda x: x < math.ceil(n ** 2), primes))):
                primes.append(n)
                return True
        return False

sum_primes = 0
n = 1   
while n < 2000000:
    if (is_prime(n)):
        sum_primes += n
    n += 1  
print sum_primes

It worked if you consider an execution time of over 7 minutes reasonable. Clearly a bad solution. So I decided to rethink it a bit and tried again. This time, I decided to use a little creative thinking. For every number I tested for primality, I took the time to mark its multiples as not primes (when I first test, say, 3, I already know that [pmath]6, 9, 12, … n * 3 < 2000000[/pmath] will not be a prime numbers, so I don’t have to test their primality later.)

max_primes = 2000000
numbers = [0] * max_primes

def is_prime(n):
        if n <= 2:
                return True
        if numbers[n] == 1:
                return False
        c = 2
        m = 0
        while True:
                m = n * c
                if m < max_primes:
                        numbers[m] = 1
                else:
                        break
                c+= 1

        #if any(n % x == 0 for x in xrange(2, int(n ** 0.5) + 1)):
        i = 2
        while i < int(n** 0.5 + 1):
                if n % i == 0:
                        return False
                i += 1

        numbers[n] = 1
        return True

def sum_primes():
        soma = 0
        for i in range (2, max_primes):
                if is_prime(i):
                        soma += i
        return soma

print sum_primes()

This is a much better algorithm and it gave the correct result in 54 seconds. Project Euler has a “one-minute rule”, which state that all its problems should be solvable “on a modestly powered computer in less than one minute,” which means the solution applies. Still, it’s a brute force solution and I wanted a better one.

So today I picked my copy of Concrete Mathematics and read about primes. Honestly, I’ve been reading more about primes lately than I even thought I would.

Anyways, Knuth—yes, I know the book has three authors, but I can’t help thinking of it as a Knuth book—talks about several strategies to test primality, but it also mentions sieve algorithms that are used to generate lists of prime numbers. Exactly what I wanted!

def prime_list(limit):
        if limit < 2:
                return []
        sieve_size = limit / 2
        sieve = [True] * sieve_size

        for i in range(0, int(limit ** 0.5) / 2):
                if not sieve[i]:
                        continue
                for j in range((i * (i + 3) * 2) + 3, sieve_size, (i * 2) + 3):
                        sieve[j] = False

        primes = [2]
        primes.extend([(i * 2) + 3 for i in range(0, sieve_size) if sieve[i]])

        return primes

print reduce(lambda x,y: x+y, prime_list(2000000))

Et voilà! Correct result in 0.456s! The code is based on the description of the Sieve of Eratosthenes found in Concrete Mathematics. Also interesting to note, my idea in the second algorithm above is actually the basis of this sieve algorithm, I was just thinking “backwards.”

So I recently wrote an ugly solution to Project Euler’s problem #9. It was written in Python and even though every other Project Euler solution I wrote ran in less than one second, this one took over 18 seconds to get to the answer. Obviously it’s my fault for using a purely brute force algorithm.

Anyway, boredom is a great motivator and I just rewrote the exact same algorithm in C, just to see how much faster it would be.

#include <stdio.h>

int is_triplet(int a, int b, int c)
{
        return (((a < b) && (b < c)) && ((a * a + b * b) == (c * c)));
}

int main(void)
{
        for (int a = 0; a < 1000; a++)
                for (int b = a + 1; b < 1000; b++)
                        for (int c = b + 1; c < 1000; c++)
                                if ((a + b + c == 1000) && is_triplet(a, b, c)) {
                                        printf("%d %d %d = %dn", a, b, c,
                                               a * b * c);
                                        return 0;
                                }
}

It’s the exact same algorithm, but instead of 18 seconds, it ran in 0.072s.

Update: And here’s the version suggested by Gustavo Niemeyer:

#include <stdio.h>

int is_triplet(int a, int b, int c)
{
        return ((a * a + b * b) == (c * c));
}

int main(void)
{
        for (int a = 1; a < 1000; a++)
                for (int b = a + 1; b < (1000 - a) / 2; b++) {
                        int c = 1000 - a - b;
                        if (is_triplet(a, b, c)) {
                                printf("%d %d %d = %dn", a, b, c,
                                       a * b * c);
                                return 0;
                        }
                }
}

It now runs in 0.002s :–)

So the other night I was a bit bored and decided to do something to pass the time. I first came across Project Euler a while ago, but had never gone further than problem #1. Boredom is a great motivator and I went through problems #2 thru #9 last night and I decided to post my solutions in search of better ones. Feel free to comment with your suggestions.

Project Euler’s Problem #9 statement is —

A Pythagorean triplet is a set of three natural numbers, $latex a < b < c$, for which,

[latex]a² + b² = c²[/latex]

For example, $latex 3² + 4² = 9 + 16 = 25 = 5²$.

There exists exactly one Pythagorean triplet for which [pmath]a + b + c = 1000[/pmath].

Find the product [pmath]abc[/pmath].

Butt ugly, super slow, brute-force solution:

from sys import exit

def is_triplet(a, b, c):
        if a < b and b < c:
                return a ** 2 + b ** 2 == c ** 2

for a in range(0, 1000):
        for b in range(a + 1, 1000):
                for c in range(b + 1, 1000):
                        if a + b + c == 1000 and is_triplet(a, b, c):
                                print a * b * c
                                exit(0)

I’ll have to rethink this, as it’s really inefficient. As it is, it runs in 18.034s!

Updated: take a look at another take at this problem and algorithm.

Updated 2010/08/04: Gustavo Niemeyer pointed an obvious optimization to the algorithm above (written in C)—the innermost loop is unnecessary. I rewrote it in Python to see the difference:

from sys import exit

def is_triplet(a, b, c):
        return (a ** 2 + b ** 2 == c ** 2)

for a in range(1, 1000):
        for b in range(a + 1, (1000 - a) / 2):
                c = 1000 - a - b
                if is_triplet(a, b, c):
                        print "%d * %d * %d = %d" % (a, b, c, a * b * c)

From 18s to 0.076s. As well, following Eduardo Habkost’s suggestion, I used psyco and execution time went down to 0.023s.

So the other night I was a bit bored and decided to do something to pass the time. I first came across Project Euler a while ago, but had never gone further than problem #1. Boredom is a great motivator and I went through problems #2 thru #9 last night and I decided to post my solutions in search of better ones. Feel free to comment with your suggestions.

Project Euler’s Problem #6 statement is —

The sum of the squares of the first ten natural numbers is,

[pmath]1² + 2² + … + 10² = 385[/pmath]

The square of the sum of the first ten natural numbers is,

pmath² = 55² = 3025[/pmath]

Hence the difference between the sum of the squares of the first ten natural numbers and the square of the sum is [pmath]3025 – 385 = 2640[/pmath].

Find the difference between the sum of the squares of the first one hundred natural numbers and the square of the sum.

Can’t be more straightforward than that, really.

sqr_sum = 0
num_sum = 0
for i in range(1,100 + 1):
        num_sum += i
        sqr_sum += i**2

num_sum = num_sum**2
print sqr_sum - num_sum

I’m pretty sure Niemeyer could rewrite this in one line, but whatever. This runs in 0.015s.

So the other night I was a bit bored and decided to do something to pass the time. I first came across Project Euler a while ago, but had never gone further than problem #1. Boredom is a great motivator and I went through problems #2 thru #9 last night and I decided to post my solutions in search of better ones. Feel free to comment with your suggestions.

Project Euler’s Problem #5 statement is —

2520 is the smallest number that can be divided by each of the numbers from 1 to 10 without any remainder.

What is the smallest positive number that is evenly divisible by all of the numbers from 1 to 20?

My first attempt was using brute force, but it didn’t go very well. It took a considerable amount of time, so I decided to think about the problem a bit. In the end, it boils down to the a relation between the lesser common multiplier and the greater common divisor.

def gcd(a, b):
        if b== 0: return a
        return gcd(b, a % b)

def lcm(a, b):
        return a * b / gcd(a, b)

if __name__ == "__main__":
        b = 2
        for i in range(3,20):
                b = lcm(b, i)
        print b

While the original code took several seconds before I killed it, this runs in 0.015s.

So the other night I was a bit bored and decided to do something to pass the time. I first came across Project Euler a while ago, but had never gone further than problem #1. Boredom is a great motivator and I went through problems #2 thru #9 last night and I decided to post my solutions in search of better ones. Feel free to comment with your suggestions.

Project Euler’s Problem #4 statement is —

A palindromic number reads the same both ways. The largest palindrome made from the product of two 2-digit numbers is 9009 = 91 99.

Find the largest palindrome made from the product of two 3-digit numbers.

And here’s the solution.

def int_to_list(number):
  return map(int, str(number))

def is_palindrome(number):
  local_list = int_to_list(number)
  return local_list == list(reversed(local_list))

pal = 0
if __name__ == "__main__":
  for i in range(100,1000):
        for j in range(100,1000):
                if is_palindrome(i * j):
                        print "%d * %d = %d"%(i, j, i*j)
                        if i * j > pal:
                                pal = i * j

  print pal

So the other night I was a bit bored and decided to do something to pass the time. I first came across Project Euler a while ago, but had never gone further than problem #1. Boredom is a great motivator and I went through problems #2 thru #9 last night and I decided to post my solutions in search of better ones. Feel free to comment with your suggestions.

Project Euler’s Problem #3 statement is —

The prime factors of 13195 are 5, 7, 13 and 29.

What is the largest prime factor of the number 600851475143 ?

My solution is really ugly and is basically dumb-, brute-force.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41

#!/usr/bin/python

import math

top_number = 600851475143

def is_divisor(divisor):
        return top_number % divisor == 0

def is_prime(divided):
        divisor = 3
        sqrt_divided = math.sqrt(divided)
        if divided == 1:
                return true
        while divisor <= sqrt_divided:
                if divided == divisor:
                        return True
                elif divided % divisor == 0:
                        return False
                else:
                        divisor += 2
        return True

def main():
        count = 3
        go_to = top_number

        first_list =[]
        while count <= go_to:
                if is_divisor(count):
                        first_list.append(count)
                        go_to = top_number / count
                        first_list.append(go_to)

                count += 2

        second_list = map(is_prime, first_list)
        print "%s" % max(zip(second_list, first_list))[-1]

if __name__ == "__main__":
        main()

So the other night I was a bit bored and decided to do something to pass the time. I first came across Project Euler a while ago, but had never gone further than problem #1. Boredom is a great motivator and I went through problems #2 thru #9 last night and I decided to post my solutions in search of better ones. Feel free to comment with your suggestions.

Project Euler’s Problem #8 statement is —

Find the greatest product of five consecutive digits in the 1000-digit number.

73167176531330624919225119674426574742355349194934
96983520312774506326239578318016984801869478851843
85861560789112949495459501737958331952853208805511
12540698747158523863050715693290963295227443043557
66896648950445244523161731856403098711121722383113
62229893423380308135336276614282806444486645238749
30358907296290491560440772390713810515859307960866
70172427121883998797908792274921901699720888093776
65727333001053367881220235421809751254540594752243
52584907711670556013604839586446706324415722155397
53697817977846174064955149290862569321978468622482
83972241375657056057490261407972968652414535100474
82166370484403199890008895243450658541227588666881
16427171479924442928230863465674813919123162824586
17866458359124566529476545682848912883142607690042
24219022671055626321111109370544217506941658960408
07198403850962455444362981230987879927244284909188
84580156166097919133875499200524063689912560717606
05886116467109405077541002256983155200055935729725
71636269561882670428252483600823257530420752963450

This is another very straightforward problem. My solution —

number = '73167176531330624919225119674426574742355349194934
96983520312774506326239578318016984801869478851843
85861560789112949495459501737958331952853208805511
12540698747158523863050715693290963295227443043557
66896648950445244523161731856403098711121722383113
62229893423380308135336276614282806444486645238749
30358907296290491560440772390713810515859307960866
70172427121883998797908792274921901699720888093776
65727333001053367881220235421809751254540594752243
52584907711670556013604839586446706324415722155397
53697817977846174064955149290862569321978468622482
83972241375657056057490261407972968652414535100474
82166370484403199890008895243450658541227588666881
16427171479924442928230863465674813919123162824586
17866458359124566529476545682848912883142607690042
24219022671055626321111109370544217506941658960408
07198403850962455444362981230987879927244284909188
84580156166097919133875499200524063689912560717606
05886116467109405077541002256983155200055935729725
71636269561882670428252483600823257530420752963450'

max_prod = 0
for i in range(1, len(number) - 4):
        prod = reduce(lambda x,y: x*y, [int(s) for s in number[i:i+5]])
        if prod > max_prod:
                max_prod = prod

print max_prod

This ran in 0m0.023s.