Showing posts with label project_euler. Show all posts
Showing posts with label project_euler. Show all posts

2009-08-24

Project Euler: Problem 20

I finished up Problem 20 a while back. This problem was really straightforward, especially given my work on Problem 16.


import math


if __name__ == '__main__':
    print sum((int(digit) for digit in str(math.factorial(100))))

Back to flipping out...

2009-08-13

Project Euler: Problem 19

I recently solved Problem 19 from Project Euler. Since I don't really enjoy date math, I decided to brute force it and use the datetime module from the standard library.


from datetime import date


def next_first_of_month_in_20th():
    """Generator to list every first of the month during the 20th century."""
    first = date(1901, 1, 1)
    yield first
    while first.year < 2001:
        if first.month == 12:
            first = first.replace(year=first.year + 1)
            first = first.replace(month=1)
        else:
            first = first.replace(month=first.month + 1)
        yield first


def main():
    """
    Solve `Problem 19`_ from Project Euler.
    
    .. _`Problem 19`: http://projecteuler.net/index.php?section=problems&id=19
    """
    return len([first for first in next_first_of_month_in_20th() if first.weekday() == 6])


if __name__ == '__main__':
    print main()

Back to flipping out...

2009-08-07

Project Euler, Problems 18 and 67

I recently finished up Problem 18 (and shortly thereafter, Problem 67). This was a fairly interesting problem, because I like graph theory even if I'm bad at it, and I immediately thought graph theory when I read this problem.

I started creating Vertex and Edge classes, but quickly decided that was too heavyweight for such a small problem. Instead, I chose to use a list of lists where the left-hand child of a node is stored in the next row on the same column and the right-hand child is stored in the next row and the next column.

This made representing and using the data fairly straightforward, but I kept running into the problem of doing an exhaustive search. I knew this was suboptimal (and completely unsuited for the related Problem 67), but I was so locked into thinking of the triangle from the top-down I made things over-complicated.

Eventually, I decided I needed to reset my thinking and Googled for hints. Once I started thinking of the problem from the bottom up (since the graph, after all, is not directed), the solution became clear. Since I didn't need the actual path, I modified the existing data structure to minimize memory requirements.


def find_max_sum(triangle):
    """
    Find the maximum sum for a path from the top to the bottom of ``triangle``.

    >>> test = [[3,], \
            [7, 5], \
            [2, 4, 6], \
            [8, 5, 9, 3]]
    >>> find_max_sum(test)
    23
    """
    while len(triangle) > 1:
        _reduce_triangle(triangle)
    return triangle[0][0]


def _reduce_triangle(to_reduce):
    """
    Reduce ``to_reduce`` in place by rolling up the maximum path info one row.

    Don't return anything to emphasize the in-place nature of this function.

    >>> test = [[3,], \
            [7, 5], \
            [2, 4, 6], \
            [8, 5, 9, 3]]
    >>> _reduce_triangle(test)
    >>> test
    [[3], [7, 5], [10, 13, 15]]
    >>> _reduce_triangle(test)
    >>> test
    [[3], [20, 20]]
    >>> _reduce_triangle(test)
    >>> test
    [[23]]
    """
    last_row = to_reduce[-1]
    for index in xrange(len(to_reduce) - 1):
        to_reduce[-2][index] += max(last_row[index:index + 2])
    del to_reduce[-1]


def main():
    """
    Solve `Problem 18`_ of Project Euler.

    .. _Problem 18: http://projecteuler.net/index.php?section=problems&id=18
    """
    actual = [[75,],
             [95, 64],
             [17, 47, 82],
             [18, 35, 87, 10],
             [20, 4, 82, 47, 65],
             [19, 1, 23, 75, 3, 34],
             [88, 2, 77, 73, 7, 63, 67],
             [99, 65, 4, 28, 6, 16, 70, 92],
             [41, 41, 26, 56, 83, 40, 80, 70, 33],
             [41, 48, 72, 33, 47, 32, 37, 16, 94, 29],
             [53, 71, 44, 65, 25, 43, 91, 52, 97, 51, 14],
             [70, 11, 33, 28, 77, 73, 17, 78, 39, 68, 17, 57],
             [91, 71, 52, 38, 17, 14, 91, 43, 58, 50, 27, 29, 48],
             [63, 66, 4, 68, 89, 53, 67, 30, 73, 16, 69, 87, 40, 31],
             [4, 62, 98, 27, 23, 9, 70, 98, 73, 93, 38, 53, 60, 4, 23]]
    print find_max_sum(actual)


if __name__ == '__main__':
    main()

Starting with this, the solution for Problem 67 is trivial: just write some code to parse the triangle out of a text file instead of coding it directly into the source code:


def find_max_sum(triangle):
    """
    Find the maximum sum for a path from the top to the bottom of ``triangle``.

    >>> test = [[3,], \
            [7, 5], \
            [2, 4, 6], \
            [8, 5, 9, 3]]
    >>> find_max_sum(test)
    23
    """
    while len(triangle) > 1:
        _reduce_triangle(triangle)
    return triangle[0][0]


def _reduce_triangle(to_reduce):
    """
    Reduce ``to_reduce`` in place by rolling up the maximum path info one row.

    Don't return anything to emphasize the in-place nature of this function.

    >>> test = [[3,], \
            [7, 5], \
            [2, 4, 6], \
            [8, 5, 9, 3]]
    >>> _reduce_triangle(test)
    >>> test
    [[3], [7, 5], [10, 13, 15]]
    >>> _reduce_triangle(test)
    >>> test
    [[3], [20, 20]]
    >>> _reduce_triangle(test)
    >>> test
    [[23]]
    """
    last_row = to_reduce[-1]
    for index in xrange(len(to_reduce) - 1):
        to_reduce[-2][index] += max(last_row[index:index + 2])
    del to_reduce[-1]


def _parse_triangle_from_file(data_file):
    """
    Parse out the triangle data from ``data_file``.

    >>> _parse_triangle_from_file('test_triangle.txt')
    [[3], [7, 5], [10, 13, 15]]
    """
    triangle = []
    with open(data_file, 'r') as triangle_file:
        for line in triangle_file:
            triangle.append([int(x) for x in line.split()])
    return triangle


def main():
    """
    Solve `Problem 67`_ of Project Euler.

    .. _Problem 67: http://projecteuler.net/index.php?section=problems&id=67
    """
    print find_max_sum(_parse_triangle_from_file('triangle.txt'))


if __name__ == '__main__':
    main()

2009-06-15

Project Euler: Problem 11, Redux

I've already solved Problem 11, but I didn't really do a great job. As BlueRaja pointed out in a comment, I was doing twice as much work as I needed to. Since I was revisiting this code anyway, I decided to play around a little bit with zip and itertools. Here's the new and improved solution:


"""Solve Problem 11 from Project Euler."""
import itertools
import operator


OFFSET = 3
GRID = [
    [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]
]
INVALID_COORDS_PRODUCT = 0


def _calc_down(row_index, col_index):
    """Calculate the product of the sequence going straight down from ``row_index``."""
    cols = itertools.repeat(col_index, OFFSET + 1)
    last_row = min(row_index + OFFSET + 1, len(GRID))
    coords = zip(xrange(row_index, last_row), cols)
    values = (GRID[row][col] for row, col in coords)
    return reduce(operator.mul, values)


def _calc_right(row_index, col_index):
    """Calculate the product of the sequence going straight right from ``col_index``."""
    rows = itertools.repeat(row_index, OFFSET + 1)
    last_col = min(col_index + OFFSET + 1, len(GRID[row_index]))
    coords = zip(rows, xrange(col_index, last_col))
    values = (GRID[row][col] for row, col in coords)
    return reduce(operator.mul, values)


def _calc_up_right(row_index, col_index):
    """Calculate the product of the sequence going up and to the right from (``row_index``, ``col_index``)."""
    last_row = max(row_index - OFFSET - 1, -1)
    last_col = min(col_index + OFFSET + 1, len(GRID[row_index]))
    coords = zip(xrange(row_index, last_row, -1), xrange(col_index, last_col))
    values = (GRID[row][col] for row, col in coords)
    return reduce(operator.mul, values)


def _calc_down_right(row_index, col_index):
    """Calculate the product of the sequence going down and to the right from (``row_index``, ``col_index``)."""
    last_row = min(row_index + OFFSET + 1, len(GRID))
    last_col = min(col_index + OFFSET + 1, len(GRID[row_index]))
    coords = zip(xrange(row_index, last_row), xrange(col_index, last_col))
    values = (GRID[row][col] for row, col in coords)
    return reduce(operator.mul, values)


def _get_max_product_for_coordinate(row, col):
    """Find the maximum product achievable from (``row``, ``col``)."""
    return max(_calc_up_right(row, col), _calc_right(row, col), _calc_down_right(row, col), _calc_down(row, col))


def problem_11():
    """Find the largest product of four adjacent elements."""
    products = set()
    for row in range(len(GRID)):
        for col in range(len(GRID[row])):
            products.add(_get_max_product_for_coordinate(row, col))
    return max(products)


if __name__ == '__main__':
    print problem_11()

Back to flipping out...

2009-04-13

Project Euler: Problem 17

I recently finished up Problem 17, and since the problem was fairly straightforward, I used it as an opportunity to explore some of Python's language features that I've been meaning to look into.


# problem17.py
"""
Find the solution to `Problem 17`_ at `Project Euler`_.

.. _Problem 17: http://projecteuler.net/index.php?section=problems&id=17
.. _Project Euler: http://projecteuler.net/
"""
__docformat__ = "restructuredtext en"

import re

_NUMBER_NAMES = {1: "one", 2: "two", 3: "three", 4: "four", 5: "five",
                6: "six", 7: "seven", 8: "eight", 9: "nine", 10: "ten",
                11: "eleven", 12: "twelve", 13: "thirteen", 14: "fourteen",
                15: "fifteen", 16: "sixteen", 17: "seventeen",
                18: "eighteen", 19: "nineteen", 20: "twenty",
                21: "twenty-one", 22: "twenty-two", 23: "twenty-three",
                24: "twenty-four", 25: "twenty-five", 26: "twenty-six",
                27: "twenty-seven", 28: "twenty-eight", 29: "twenty-nine",
                30: "thirty", 31: "thirty-one", 32: "thirty-two",
                33: "thirty-three", 34: "thirty-four", 35: "thirty-five",
                36: "thirty-six", 37: "thirty-seven", 38: "thirty-eight",
                39: "thirty-nine", 40: "forty", 41: "forty-one",
                42: "forty-two", 43: "forty-three", 44: "forty-four",
                45: "forty-five", 46: "forty-six", 47: "forty-seven",
                48: "forty-eight", 49: "forty-nine", 50: "fifty",
                51: "fifty-one", 52: "fifty-two", 53: "fifty-three",
                54: "fifty-four", 55: "fifty-five", 56: "fifty-six",
                57: "fifty-seven", 58: "fifty-eight", 59: "fifty-nine",
                60: "sixty", 61: "sixty-one", 62: "sixty-two",
                63: "sixty-three", 64: "sixty-four", 65: "sixty-five",
                66: "sixty-six", 67: "sixty-seven", 68: "sixty-eight",
                69: "sixty-nine", 70: "seventy", 71: "seventy-one",
                72: "seventy-two", 73: "seventy-three", 74: "seventy-four",
                75: "seventy-five", 76: "seventy-six", 77: "seventy-seven",
                78: "seventy-eight", 79: "seventy-nine", 80: "eighty",
                81: "eighty-one", 82: "eighty-two", 83: "eighty-three",
                84: "eighty-four", 85: "eighty-five", 86: "eighty-six",
                87: "eighty-seven", 88: "eighty-eight", 89: "eighty-nine",
                90: "ninety", 91: "ninety-one", 92: "ninety-two",
                93: "ninety-three", 94: "ninety-four", 95: "ninety-five",
                96: "ninety-six", 97: "ninety-seven", 98: "ninety-eight",
                99: "ninety-nine"}

_CHARACTERS_WE_CARE_ABOUT = re.compile("\w")

def _words_from_num(num):
    """
    Convert ``num`` to its (British) English phrase equivalent.

    If ``num`` is greater than 9,999 then raise an ``Exception``.
    >>> _words_from_num(115)
    'one hundred and fifteen'
    """
    if num >= 10000:
        raise Exception, 'This function only supports numbers less than 10000.'

    parts_list = []
    if num >= 1000:
        thousands = num // 1000
        parts_list.append(_NUMBER_NAMES[thousands])
        parts_list.append(" thousand")
        num -= thousands * 1000
    if num >= 100:
        hundreds = num // 100
        parts_list.append(_NUMBER_NAMES[hundreds])
        parts_list.append(" hundred")
        num -= hundreds * 100
    if num:
        if parts_list:
            parts_list.append(" and")
        parts_list.extend([" ", _NUMBER_NAMES[num]])

    return "".join(parts_list)

def _count_characters_we_care_about(string_to_count):
    """
    Count the characters in ``string_to_count``, excluding things like hyphens and spaces.

    >>> _count_characters_we_care_about("one hundred and twenty-three")
    24
    """
    return len(_CHARACTERS_WE_CARE_ABOUT.findall(string_to_count))

def problem_17(upper_bound = 1000):
    """
    Find the solution to `Problem 17`_ at `Project Euler`_.

    .. _Problem 17: http://projecteuler.net/index.php?section=problems&id=17
    .. _Project Euler: http://projecteuler.net/

    >>> problem_17(2)
    6
    """
    converted_nums = (_words_from_num(num) for num in xrange(1, upper_bound + 1))
    lengths = (_count_characters_we_care_about(phrase) for phrase in converted_nums)
    return sum(lengths)

if __name__ == '__main__':
    print problem_17()

Generator Expressions

The quick version: They're like list comprehensions, only lazy. The longer version resides in PEP 289. The verdict: They rock. It didn't make a big difference in this problem, but in general, it's nice that you don't have to allocate enough memory to contain an entire list. One problem with them is that once you've consumed an element of the generator expression, you can't get it back, so they aren't well-suited for problems where you need to iterate across the data more than once.

Regexes

Regexes are hardly unique to Python, but this is the first time I've ever used them in Python. One feature I was excited to try—and that Python pioneered—is named capturing groups, but the regex I used in this example didn't need them.

Docstrings

I've written docstrings before, but this is the first time I've tried to generate stand-alone documentation from them. To do the generation, I used epydoc. The original PEP 257 defines docstrings in terms of plaintext, but PEP 287 establishes reStructuredText as a richer alternative. Since I was generating documentation from the docstrings, I decided to go with reStructuredText. To avoid typing --docformat restructuredtext every time I invoked epydoc again:

__docformat__ = "restructuredtext en"

While it took a while to get used to it (I've gotten pretty set in my Markdown ways), I really like it so far. In fact, I liked it enough to write this entire blog post in it and simply post the generated HTML into Blogger.

Back to flipping out…

2009-03-30

Project Euler: Problem 16

I've just finished solving Problem 16. It was a fairly straightforward problem, and at first I couldn't decide how I wanted to solve it, so I decided to do it both ways and compare them. Here's the script w/ both versions of the solution:


import math

def pure_numeric_sum_of_digits(num):
    """
    Finds the sum of the digits in num (presumed to be base-10).

    This version of the method operates purely on numbers.

    >>> pure_numeric_sum_of_digits(2 ** 15)
    26
    """
    running_total = 0
    for n in xrange(int(math.log10(num)) + 1):
        running_total += num % 10
        num //= 10
    return running_total

def string_based_sum_of_digits(num):
    """
    Finds the sum of the digits in num (presumed to be base-10).

    This version of the method works by using a string representation
    of num and converting individual characters to numbers so we can
    operate on it as a sequence.

    >>> string_based_sum_of_digits(2 ** 15)
    26
    """
    sum([int(digit) for digit in str(num)])

def problem_16():
    """Solves Problem 16 at Project Euler."""
    return string_based_sum_of_digits(2 ** 1000)

if __name__ == '__main__':
    print problem_16()

Clearly, the string-based solution is more compact, but it has to be slower, right? Strings are slower than ints and all that. Thanks to the optimization efforts of the Python implementers, the string-based solution is actually ~167% faster on my machine (I used cProfile on a 10000-iteration loop to make these measurements).

Back to flipping out...

2009-03-04

Project Euler: Problem 15

This one took me a while, largely because I didn't understand what they meant by backtracking. I was thinking of backtracking to a particular point, so I kept trying to use acyclic graphs, sometimes called forests. That's not really what they meant, though. In the context of Problem 15, no backtracking just means you can't go up and you can't go left. After working a 3x3 version of the problem, I thought the values looked a little familiar; I was right. It's a counting problem we covered in combinatorics oh so many years ago: Permutations with Repeated Elements. There are 40 steps to any given path (we have to go 20 steps down and 20 steps right). Since the definition of backtracking makes any given move down indistinguishable from any other move down, and any given move right indistinguishable from any other move right, the answer is 40!/(20!20!).

Back to flipping out...

2009-02-22

Project Euler: Problem 14

Problem 14 was a fairly straightforward problem to solve, but I liked it a lot. One reason I liked it is that I enjoy playing with mathematical sequences, but also because it marked my first foray into TDD in Python. I'd heard about doctest, and I really like the concept, so I decided to use it when solving this problem. Here's the code I used to find the solution:


"""Solves Problem 14 from Project Euler."""

def _gen_sequence(seed, cur_sequence = None):
    """Generates a sequence following the rules from Problem 14.
    
    >>> _gen_sequence(13)
    [13, 40, 20, 10, 5, 16, 8, 4, 2, 1]"""
    if not cur_sequence:
        cur_sequence = []
    cur_sequence.append(seed)
    if seed == 1:
        return cur_sequence
    elif seed % 2:
        return _gen_sequence(3 * seed + 1, cur_sequence)
    else:
        return _gen_sequence(seed / 2, cur_sequence)

def problem_14(upper_bound):
    """Finds the seed (< upper_bound and > 500000) for the longest sequence.

    >>> problem_14(502873)
    502137
    """
    cur_answer = 0, 0  # Track the length too, just out of curiosity.
    for seed in reversed(xrange(500001, upper_bound, 2)):
        cur_length = len(_gen_sequence(seed))
        if cur_length > cur_answer[1]:
            cur_answer = seed, cur_length
    return cur_answer[0]

if __name__ == '__main__':
    print problem_14(1000000)

You can execute these tests (assuming Python 2.6) from the command-line with python -m doctest -v problem14.py. The output will resemble:


Trying:
    _gen_sequence(13)
Expecting:
    [13, 40, 20, 10, 5, 16, 8, 4, 2, 1]
ok
Trying:
    problem_14(502873)
Expecting:
    502137
ok
1 items had no tests:
    problem14
2 items passed all tests:
   1 tests in problem14._gen_sequence
   1 tests in problem14.problem_14
2 tests in 3 items.
2 passed and 0 failed.
Test passed.

Not only does doctest make it easy to do TDD, it reminds me somewhat of literate programming. It's obviously not the same, but I think it encourages the sort of documentation that Knuth would approve. I know it definitely made me think more about documenting these simple functions.

Back to flipping out...

2009-02-18

Project Euler: Problem 13

Problem 13 was very straightforward to brute-force.


"""Solves Problem 13 of Project Euler."""

NUMBERS = [
                37107287533902102798797998220837590246510135740250,
                46376937677490009712648124896970078050417018260538,
                74324986199524741059474233309513058123726617309629,
                91942213363574161572522430563301811072406154908250,
                23067588207539346171171980310421047513778063246676,
                89261670696623633820136378418383684178734361726757,
                28112879812849979408065481931592621691275889832738,
                44274228917432520321923589422876796487670272189318,
                47451445736001306439091167216856844588711603153276,
                70386486105843025439939619828917593665686757934951,
                62176457141856560629502157223196586755079324193331,
                64906352462741904929101432445813822663347944758178,
                92575867718337217661963751590579239728245598838407,
                58203565325359399008402633568948830189458628227828,
                80181199384826282014278194139940567587151170094390,
                35398664372827112653829987240784473053190104293586,
                86515506006295864861532075273371959191420517255829,
                71693888707715466499115593487603532921714970056938,
                54370070576826684624621495650076471787294438377604,
                53282654108756828443191190634694037855217779295145,
                36123272525000296071075082563815656710885258350721,
                45876576172410976447339110607218265236877223636045,
                17423706905851860660448207621209813287860733969412,
                81142660418086830619328460811191061556940512689692,
                51934325451728388641918047049293215058642563049483,
                62467221648435076201727918039944693004732956340691,
                15732444386908125794514089057706229429197107928209,
                55037687525678773091862540744969844508330393682126,
                18336384825330154686196124348767681297534375946515,
                80386287592878490201521685554828717201219257766954,
                78182833757993103614740356856449095527097864797581,
                16726320100436897842553539920931837441497806860984,
                48403098129077791799088218795327364475675590848030,
                87086987551392711854517078544161852424320693150332,
                59959406895756536782107074926966537676326235447210,
                69793950679652694742597709739166693763042633987085,
                41052684708299085211399427365734116182760315001271,
                65378607361501080857009149939512557028198746004375,
                35829035317434717326932123578154982629742552737307,
                94953759765105305946966067683156574377167401875275,
                88902802571733229619176668713819931811048770190271,
                25267680276078003013678680992525463401061632866526,
                36270218540497705585629946580636237993140746255962,
                24074486908231174977792365466257246923322810917141,
                91430288197103288597806669760892938638285025333403,
                34413065578016127815921815005561868836468420090470,
                23053081172816430487623791969842487255036638784583,
                11487696932154902810424020138335124462181441773470,
                63783299490636259666498587618221225225512486764533,
                67720186971698544312419572409913959008952310058822,
                95548255300263520781532296796249481641953868218774,
                76085327132285723110424803456124867697064507995236,
                37774242535411291684276865538926205024910326572967,
                23701913275725675285653248258265463092207058596522,
                29798860272258331913126375147341994889534765745501,
                18495701454879288984856827726077713721403798879715,
                38298203783031473527721580348144513491373226651381,
                34829543829199918180278916522431027392251122869539,
                40957953066405232632538044100059654939159879593635,
                29746152185502371307642255121183693803580388584903,
                41698116222072977186158236678424689157993532961922,
                62467957194401269043877107275048102390895523597457,
                23189706772547915061505504953922979530901129967519,
                86188088225875314529584099251203829009407770775672,
                11306739708304724483816533873502340845647058077308,
                82959174767140363198008187129011875491310547126581,
                97623331044818386269515456334926366572897563400500,
                42846280183517070527831839425882145521227251250327,
                55121603546981200581762165212827652751691296897789,
                32238195734329339946437501907836945765883352399886,
                75506164965184775180738168837861091527357929701337,
                62177842752192623401942399639168044983993173312731,
                32924185707147349566916674687634660915035914677504,
                99518671430235219628894890102423325116913619626622,
                73267460800591547471830798392868535206946944540724,
                76841822524674417161514036427982273348055556214818,
                97142617910342598647204516893989422179826088076852,
                87783646182799346313767754307809363333018982642090,
                10848802521674670883215120185883543223812876952786,
                71329612474782464538636993009049310363619763878039,
                62184073572399794223406235393808339651327408011116,
                66627891981488087797941876876144230030984490851411,
                60661826293682836764744779239180335110989069790714,
                85786944089552990653640447425576083659976645795096,
                66024396409905389607120198219976047599490197230297,
                64913982680032973156037120041377903785566085089252,
                16730939319872750275468906903707539413042652315011,
                94809377245048795150954100921645863754710598436791,
                78639167021187492431995700641917969777599028300699,
                15368713711936614952811305876380278410754449733078,
                40789923115535562561142322423255033685442488917353,
                44889911501440648020369068063960672322193204149535,
                41503128880339536053299340368006977710650566631954,
                81234880673210146739058568557934581403627822703280,
                82616570773948327592232845941706525094512325230608,
                22918802058777319719839450180888072429661980811197,
                77158542502016545090413245809786882778948721859617,
                72107838435069186155435662884062257473692284509516,
                20849603980134001723930671666823555245252804609722,
                53503534226472524250874054075591789781264330331690
]

def problem_13(nums_to_sum):
    """Find the first 10 digits of the sum of nums_to_sum."""
    return int(str(sum(nums_to_sum))[:10])

if __name__ == '__main__':
    print problem_13(NUMBERS)

Back to flipping out...

2009-02-17

Project Euler: Problem 12

Problem 12 brings back an interesting mathematical concept: triangular numbers. Triangular numbers (and the formula for finding one) are the slightly more sophisticated approach I mentioned in my writeup for Problem 6. Writing a generator for triangular numbers is easy enough, as is writing some logic to factorize (not prime factorize) a number. Given those two, the solution is easy enough.


"""Solves Problem 12 of Project Euler."""

import math

def factors(to_factor):
    """Find the factors of to_factor."""
    factors = []
    divisor = 1
    while (divisor <= int(math.sqrt(to_factor))):
        if not to_factor % divisor:
            quotient = to_factor / divisor
            factors.append(divisor)
            factors.append(quotient)
        divisor += 1

    return factors

def triangular_numbers():
    """Generate the triangular numbers."""
    current = 0
    position = 1
    while True:
        current += position
        position += 1
        yield current

def problem_12(min_divisors):
    """Finds the first triangular number to have more than 500 divisors."""
    for triangular in triangular_numbers():
        cur_factors = factors(triangular)
        if len(cur_factors) > min_divisors:
            return triangular

if __name__ == '__main__':
    print problem_12(500)

Back to flipping out...

2009-02-16

Project Euler: Problem 11

Problem 11 is fairly straightforward to brute-force, although my solution uses more code than is strictly necessary. Since one of the things I enjoy most about Python is its support for functional programming concepts, I went with a reduce-based solution.


"""Solves Problem 11 from Project Euler."""

import operator

OFFSET = 3
GRID = [
    [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]
]
INVALID_COORDS_TUPLE = 0,

def _build_up(row_index, col_index):
    """Build the sequence going straight up from row_index."""
    if row_index < OFFSET:
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index - 1][col_index],
            GRID[row_index - 2][col_index],
            GRID[row_index - 3][col_index])

def _build_down(row_index, col_index):
    """Build the sequence going straight down from row_index."""
    if row_index + OFFSET >= len(GRID):
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index + 1][col_index],
            GRID[row_index + 2][col_index],
            GRID[row_index + 3][col_index])

def _build_left(row_index, col_index):
    """Build the sequence going straight left from col_index."""
    if col_index < OFFSET:
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index][col_index - 1],
            GRID[row_index][col_index - 2],
            GRID[row_index][col_index - 3])

def _build_right(row_index, col_index):
    """Build the sequence going straight right from row_index."""
    if col_index + OFFSET >= len(GRID[row_index]):
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index][col_index + 1],
            GRID[row_index][col_index + 2],
            GRID[row_index][col_index + 3])

def _build_up_right(row_index, col_index):
    """Build the sequence going up and to the right from \
    (row_index, col_index)."""
    if row_index < OFFSET:
        return INVALID_COORDS_TUPLE
    if col_index + OFFSET >= len(GRID[row_index]):
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index - 1][col_index + 1],
            GRID[row_index - 2][col_index + 2],
            GRID[row_index - 3][col_index + 3])

def _build_down_right(row_index, col_index):
    """Build the sequence going down and to the right from \
    (row_index, col_index)."""
    if row_index + OFFSET >= len(GRID):
        return INVALID_COORDS_TUPLE
    if col_index + OFFSET >= len(GRID[row_index]):
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index + 1][col_index + 1],
            GRID[row_index + 2][col_index + 2],
            GRID[row_index + 3][col_index + 3])

def _build_down_left(row_index, col_index):
    """Build the sequence going down and to the left from \
    (row_index, col_index)."""
    if row_index + OFFSET >= len(GRID):
        return INVALID_COORDS_TUPLE
    if col_index < OFFSET:
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index + 1][col_index - 1],
            GRID[row_index + 2][col_index - 2],
            GRID[row_index + 3][col_index - 3])

def _build_up_left(row_index, col_index):
    """Build the sequence going down and to the left from \
    (row_index, col_index)."""
    if row_index < OFFSET:
        return INVALID_COORDS_TUPLE
    if col_index < OFFSET:
        return INVALID_COORDS_TUPLE

    return (GRID[row_index][col_index],
            GRID[row_index - 1][col_index - 1],
            GRID[row_index - 2][col_index - 2],
            GRID[row_index - 3][col_index - 3])

def _build_all_possible_sequences():
    """Build a set of all possible 4-element sequences."""
    sequences = set()
    for row in range(len(GRID)):
        for col in range(len(GRID[row])):
            sequences.add(_build_up(row, col))
            sequences.add(_build_down(row, col))
            sequences.add(_build_left(row, col))
            sequences.add(_build_right(row, col))
            sequences.add(_build_up_right(row, col))
            sequences.add(_build_down_right(row, col))
            sequences.add(_build_down_left(row, col))
            sequences.add(_build_up_left(row, col))
    return sequences

def problem_11():
    """Find the largest product of four adjacent elements."""
    products = [reduce(operator.mul, sequence)
                for sequence in _build_all_possible_sequences()]
    return max(products)

if __name__ == '__main__':
    print problem_11()

Of course, this isn't the best fit for a functional solution. For one thing, it wastes an awful lot of memory. The more imperative approach below runs in less than half the time on my box, and it's (a little) less code, to boot. The lesson, as always: KISS.


"""Solves Problem 11 from Project Euler."""

OFFSET = 3
GRID = [
    [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]
]
INVALID_COORDS_PRODUCT = 0

def _calc_up(row_index, col_index):
    """Build the sequence going straight up from row_index."""
    if row_index < OFFSET:
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index - 1][col_index] * \
            GRID[row_index - 2][col_index] * \
            GRID[row_index - 3][col_index]

def _calc_down(row_index, col_index):
    """Build the sequence going straight down from row_index."""
    if row_index + OFFSET >= len(GRID):
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index + 1][col_index] * \
            GRID[row_index + 2][col_index] * \
            GRID[row_index + 3][col_index]

def _calc_left(row_index, col_index):
    """Build the sequence going straight left from col_index."""
    if col_index < OFFSET:
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index][col_index - 1] * \
            GRID[row_index][col_index - 2] * \
            GRID[row_index][col_index - 3]

def _calc_right(row_index, col_index):
    """Build the sequence going straight right from row_index."""
    if col_index + OFFSET >= len(GRID[row_index]):
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index][col_index + 1] * \
            GRID[row_index][col_index + 2] * \
            GRID[row_index][col_index + 3]

def _calc_up_right(row_index, col_index):
    """Build the sequence going up and to the right from \
    (row_index, col_index)."""
    if row_index < OFFSET:
        return INVALID_COORDS_PRODUCT
    if col_index + OFFSET >= len(GRID[row_index]):
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index - 1][col_index + 1] * \
            GRID[row_index - 2][col_index + 2] * \
            GRID[row_index - 3][col_index + 3]

def _calc_down_right(row_index, col_index):
    """Build the sequence going down and to the right from \
    (row_index, col_index)."""
    if row_index + OFFSET >= len(GRID):
        return INVALID_COORDS_PRODUCT
    if col_index + OFFSET >= len(GRID[row_index]):
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index + 1][col_index + 1] * \
            GRID[row_index + 2][col_index + 2] * \
            GRID[row_index + 3][col_index + 3]

def _calc_down_left(row_index, col_index):
    """Build the sequence going down and to the left from \
    (row_index, col_index)."""
    if row_index + OFFSET >= len(GRID):
        return INVALID_COORDS_PRODUCT
    if col_index < OFFSET:
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index + 1][col_index - 1] * \
            GRID[row_index + 2][col_index - 2] * \
            GRID[row_index + 3][col_index - 3]

def _calc_up_left(row_index, col_index):
    """Build the sequence going down and to the left from \
    (row_index, col_index)."""
    if row_index < OFFSET:
        return INVALID_COORDS_PRODUCT
    if col_index < OFFSET:
        return INVALID_COORDS_PRODUCT

    return GRID[row_index][col_index] * \
            GRID[row_index - 1][col_index - 1] * \
            GRID[row_index - 2][col_index - 2] * \
            GRID[row_index - 3][col_index - 3]

def _calc_all_possible_sequences():
    """Build a set of all possible 4-element sequences."""
    products = set()
    for row in range(len(GRID)):
        for col in range(len(GRID[row])):
            products.add(_calc_up(row, col))
            products.add(_calc_down(row, col))
            products.add(_calc_left(row, col))
            products.add(_calc_right(row, col))
            products.add(_calc_up_right(row, col))
            products.add(_calc_down_right(row, col))
            products.add(_calc_down_left(row, col))
            products.add(_calc_up_left(row, col))
    return products

def problem_11():
    """Find the largest product of four adjacent elements."""
    return max(_calc_all_possible_sequences())

if __name__ == '__main__':
    print problem_11()

Back to flipping out...

2009-02-11

Project Euler: Problem 10

I just finished up Problem 10; given my earlier work on Problem 7, it was trivial to adapt it and arrive at the following program.


"""Solves Problem 10 for Project Euler."""
import math

def is_prime(candidate, known_primes):
    """Determines whether candidate is prime by trial division using \
    known_primes.

    For this function to work, known_primes *must* be accurate.
    """
    last_possible = math.sqrt(candidate)
    for current_prime in known_primes:
        if current_prime > last_possible:
            break
        if not candidate % current_prime:
            return False
    return True

def primes_generator(upper_bound):
    """A generator for all the primes < upper_bound."""
    candidates = xrange(2, upper_bound)
    primes = []
    for n in candidates:
        if is_prime(n, primes):
            primes.append(n)
            yield n

def problem_10():
    """Sum the primes less than 2 million."""
    return sum(primes_generator(2000000))

if __name__ == '__main__':
    print problem_10()

Back to flipping out...

2009-02-10

Project Euler: Problem 9

Problem 9 deals with one of the more interesting things I learned in high school geometry: Pythagorean triples. In high school, I just memorized the 2 most common ones (3, 4, 5 and 5, 12, 13) and thought to myself: Wouldn't it be cool if I could generate all of these? But that was in the before time; Wikipedia didn't exist and my text book wasn't cool enough to dwell on them. At any rate, now I know Euclid's formula for generating Pythagorean triples: m ** 2 - n ** 2, 2 * m * n, m ** 2 + n ** 2

Here's the script I used to solve the actual problem:


"""Solves Problem 9 from Project Euler."""

import operator

def triples(upper_bound):
    """Generator for Pythagorean Triples (represented as tuples).

    Uses Euclid's formula to generate Pythagorean Triples
    (see http://en.wikipedia.org/wiki/Pythagorean_triple#Generating_a_triple).
    """
    for m in xrange(2, upper_bound):
        for n in xrange(1, m):
            yield m ** 2 - n ** 2, 2 * m * n, m ** 2 + n ** 2
# Only uncomment if the triples we get from the original Euclid's are insufficient.
# for k in xrange(1, upper_bound):
#     yield k * (m ** 2 - n ** 2), k * (2 * m * n), k * (m ** 2 + n ** 2)

def problem_9():
    """Finds the product of the Pythagorean Triple where a + b + c = 1000."""
    for triple in triples(1000):
        if sum(triple) == 1000:
            return reduce(operator.mul, triple)
    return 0

if __name__ == '__main__':
    print problem_9()

Back to flipping out...

2009-02-09

Project Euler: Problem 8

Just finished up with Problem 8. Brute-forcing it was pretty straightforward, so I decided to play about with some of the more functional aspects of Python. Enter reduce. Here's the original version I used to solve the problem.


"""Solves Problem 8 from Project Euler."""

def problem_8(num_in_question):
    """Finds and returns the greatest product of 5 consecutive digits \
    of num_in_question."""
    to_process = str(num_in_question)
    offset = 0
    highest_product = 0
    last_possible_start = len(to_process) - 5
    while (offset < last_possible_start):
        digits = [int(digit) for digit in to_process[offset:offset + 5]]
        product = 1
        for n in digits:
            product *= n
 
        if product > highest_product:
            highest_product = product
 
        offset += 1
 
    return highest_product

if __name__ == '__main__':
    print problem_8("73167176531330624919225119674426574742355349194934\
96983520312774506326239578318016984801869478851843\
85861560789112949495459501737958331952853208805511\
12540698747158523863050715693290963295227443043557\
66896648950445244523161731856403098711121722383113\
62229893423380308135336276614282806444486645238749\
30358907296290491560440772390713810515859307960866\
70172427121883998797908792274921901699720888093776\
65727333001053367881220235421809751254540594752243\
52584907711670556013604839586446706324415722155397\
53697817977846174064955149290862569321978468622482\
83972241375657056057490261407972968652414535100474\
82166370484403199890008895243450658541227588666881\
16427171479924442928230863465674813919123162824586\
17866458359124566529476545682848912883142607690042\
24219022671055626321111109370544217506941658960408\
07198403850962455444362981230987879927244284909188\
84580156166097919133875499200524063689912560717606\
05886116467109405077541002256983155200055935729725\
71636269561882670428252483600823257530420752963450")

And here's the same problem_8 function using reduce:


def problem_8(num_in_question):
    """Finds and returns the greatest product of 5 consecutive digits \
    of num_in_question."""
    to_process = str(num_in_question)
    offset = 0
    highest_product = 0
    last_possible_start = len(to_process) - 5
    while (offset < last_possible_start):
        digits = [int(digit) for digit in to_process[offset:offset + 5]]
        product = reduce(operator.mul, digits)
 
        if product > highest_product:
            highest_product = product
 
        offset += 1
 
    return highest_product

Pretty similar: 1 less line of code, 1 more line of imports, almost identical performance. I guess it all comes down to taste. One note: if you're doing functional programming and need to use a function supported by the operator module, that's the recommended way of doing it. Since it's part of the standard library, it's more obvious what's going on than a comparable lambda, plus they're implemented in C to give better performance. But since we're getting all functional, we might as well do it all the way. Here's another version:


def problem_8(num_in_question):
    """Finds and returns the greatest product of 5 consecutive digits \
    of num_in_question.

    This function expects num_in_question to be a string so we can
    slice it into 5-digit sequences.
    """
    SEQUENCE_LENGTH = 5
    sequences = [num_in_question[offset:offset + SEQUENCE_LENGTH] \
        for offset in range(len(num_in_question) - SEQUENCE_LENGTH)]
    nums = []
    for sequence in sequences:
        nums.append([int(num) for num in sequence])

    return max([reduce(operator.mul, num_list) for num_list in nums])

Back to flipping out...

2009-02-06

Project Euler: Problem 7, redux

Remember how I mentioned the Sieve wasn't the most performant solution? Here's a much faster solution. On my system the time dropped from 11.661 seconds to .332 seconds. Also, it makes use of one of my favorite features in Python, so far: generators.


"""Solves Problem 7 from Project Euler."""

import math
import sys

def is_prime(candidate, known_primes):
    """Determines whether candidate is prime by trial division using \
    known_primes.

    For this function to work, known_primes *must* be accurate.
    """
    last_possible = math.sqrt(candidate)
    for current_prime in known_primes:
        if current_prime > last_possible:
            break
        if not candidate % current_prime:
            return False
    return True

def primes_generator():
    """A generator for all the primes <= sys.maxint."""
    candidates = xrange(2, sys.maxint)
    primes = []
    for n in candidates:
        if is_prime(n, primes):
            primes.append(n)
            yield n

def problem_7(n):
    """Finds the nth prime number."""
    primes = primes_generator()
    i = 0
    while i < n - 1:
        i += 1
        primes.next()

    return primes.next()

if __name__ == '__main__':
    print problem_7(10001)

Back to flipping out...

Project Euler: Problem 7

My new math trick of the day? Upper and lower bounds on the nth prime:

n * ln(n) + n * ln(ln(n - 1)) is less than p[n] is less than n * ln(n) + n * ln(ln(n)) for n greater than or equal to 6

Using this little nugget, I can give myself an upper bound on the numbers I need to test for primality and unleash the Sieve of Eratosthenes on Problem 7. BTW, I realize the Sieve isn't the most performant way to do most of these tests; I just think it's a really elegant solution for finding primes and it isn't slow enough to be an issue for most numbers of this size. At any rate, here's my script:


"""Solves Problem 7 from Project Euler."""

import math

FIRST_SIX_PRIMES = {1:2, 2:3, 3:5, 4:7, 5:11, 6:13}

def e_sieve(upper_bound):
    """Uses the Sieve of Eratosthenes to get a list of the primes up to max."""
    primes = []
    candidates = range(2, upper_bound)
    while candidates:
        head = candidates[0]
        primes.append(head)
        candidates = [n for n in candidates[1:] if n % head]

    return primes

def problem_7(n):
    """Finds the nth prime number.

    Thanks to
    http://en.wikipedia.org/wiki/Prime_number_theorem#Approximations_for_the_nth_prime_number
    we know that it will be less than n * ln(n) + n * ln(ln(n))
    (for n >= 6).
    """
    if n < 6:
        return FIRST_SIX_PRIMES[n]

    upper_bound = int(n * math.log(n) + n * math.log(math.log(n)))

    primes = e_sieve(upper_bound)
    return primes[n - 1]

if __name__ == '__main__':
    print problem_7(10001)

Back to flipping out...

2009-02-05

Project Euler: Problem 6

I just finished up Problem 6. This one was really straightforward. The interesting part was learning a new math trick: Square Pyramidal Numbers. Here's the most brute-force way (I thought of) to solve it:

"Naïve" Solution


"""Solves Problem 6 from Project Euler."""

def sum_of_squares(upper_bound):
    """Sums the squares of all the natural numbers from 1 to \
    upper_bound (inclusive)."""
    return sum(n ** 2 for n in range(1, upper_bound + 1))

def square_of_sums(upper_bound):
    """Sums all the numbers from 1 to upper_bound."""
    return (sum(range(1, upper_bound + 1)) ** 2)

def problem_6(n):
    """Finds the difference between the square of the sums and the \
    sum of the squares for the first n natural numbers."""
    return square_of_sums(n) - sum_of_squares(n)

if __name__ == '__main__':
    print problem_6(100)

Even though this solution has the worst Big O, it wasn't noticeably slower than any of the others for numbers the size of the ones here.

I used a slightly more sophisticated approach because I remembered an insight for summing the first N natural numbers, discovered by none other than Gauss himself (especially entertaining to use Gauss to solve problems on Project Euler, I know): n * (n + 1)) / 2. I fell back on brute force for the sum of squares part, though.

Solution I Used


"""Solves Problem 6 from Project Euler."""

def sum_of_squares(upper_bound):
    """Sums the squares of all the natural numbers from 1 to \
    upper_bound (inclusive)."""
    return sum(n ** 2 for n in range(1, upper_bound + 1))

def square_of_sums(upper_bound):
    """Sums all the numbers from 1 to upper_bound."""
    # Use Gauss' insight about n/2 * (n + 1).
    return (upper_bound * (upper_bound + 1) / 2) ** 2

def problem_6(n):
    """Finds the difference between the square of the sums and the \
    sum of the squares for the first n natural numbers."""
    return square_of_sums(n) - sum_of_squares(n)

if __name__ == '__main__':
    print problem_6(100)

Since there was a math trick for solving the summing problem, I suspected there was one for the sum of squares problem, as well. A quick Google revealed: Square Pyramidal Numbers. I wish I had known about this when I was doing stupid brain teasers; I always knew I was wasting too much time counting squares.

Least Brutish Solution


"""Solves Problem 6 from Project Euler."""

def sum_of_squares(upper_bound):
    """Sums the squares of all the natural numbers from 1 to \
    upper_bound (inclusive)."""
    # Use the rule for square pyramidal numbers
    # (http://en.wikipedia.org/wiki/Square_pyramidal_number).
    return ((2 * upper_bound ** 3) + (3 * upper_bound ** 2) + upper_bound) / 6

def square_of_sums(upper_bound):
    """Sums all the numbers from 1 to upper_bound."""
    # Use Gauss' insight about n/2 * (n + 1).
    return (upper_bound * (upper_bound + 1) / 2) ** 2

def problem_6(n):
    """Finds the difference between the square of the sums and the \
    sum of the squares for the first n natural numbers."""
    return square_of_sums(n) - sum_of_squares(n)

if __name__ == '__main__':
    print problem_6(100)

All in all, I'm really glad I started working these problems. I've usually learned/remembered at least one interesting math trick or Python trick during each one I solved, or at least had good programming practice (KISS) driven home.

Back to flipping out...

2009-02-04

Project Euler: Problem 5

At last! A Project Euler problem I didn't brute-force. This problem was a fairly straight-forward LCM problem. I just reused my factorize function from before and implemented the algorithm to find LCM by prime factorization.


"""Solves Problem 5 of Project Euler."""

def factorize(to_factor):
    """Use trial division to factorize to_factor and return all the resulting \
    factors."""
    factors = []
    divisor = 2
    while (divisor < to_factor):
        if not to_factor % divisor:
            to_factor /= divisor
            factors.append(divisor)
            # Note we don't bump the divisor here; if we did, we'd have
            # non-prime factors.
        elif divisor == 2:
            divisor += 1
        else:
            # Trivial optimization: skip even numbers that aren't 2.
            divisor += 2
    if not to_factor % divisor:
        # Don't forget the last factor
        factors.append(to_factor)
    return factors

def lcm(numbers):
    """Finds the Least Common Multiple of numbers."""
    highest_degree_factors = {}
    for number in numbers:
        degrees_by_factor = {}
        for factor in factorize(number):
            # Translate the raw list of factors into a dictionary of degrees
            # keyed on the factor.
            current_degree = degrees_by_factor.setdefault(factor, 0)
            degrees_by_factor[factor] = 1 + current_degree

        # Update the top-level dict so it really is tracking the highest
        # degrees.
        for k, v in degrees_by_factor.iteritems():
            highest_degree_factors.setdefault(k, v)
            if highest_degree_factors[k] < v:
               highest_degree_factors[k] = v 

    running_product = 1
    for factor, degree in highest_degree_factors.iteritems():
        running_product *= factor ** degree

    return running_product

if __name__ == '__main__':
    print lcm(range(1,20))

Back to flipping out...

Project Euler: Problem 4

I really enjoyed solving this one; palindromes have intrigued me since learned about them in grade school. Rather than repeat my mistake of getting too cute from Problem 3, I decided to brute-force this one. That approach did lead me to learn some nifty Python tricks though:

  • [::-1] for reversing a string
  • the new reversed feature
  • xrange for lazy range generation; I'd seen this one, but this is the first time I'd written it

Without further ado, here's my solution to Problem 4.


"""Solves Problem 4 from Project Euler.
"""

def is_palindrome(to_test):
    """Determine whether to_test is a palindrome."""
    # Uses extended slice syntax to reverse a sequence.
    return to_test == to_test[::-1]

def problem_4():
    """Function to solve problem 4."""
    # I arrived at the values in the range by finding the smallest and largest
    # products of 2 3-digit numbers.
    for n in reversed(xrange(10000, 998001)):
        if is_palindrome(str(n)):
            for factor_1 in reversed(xrange(1000)):
                for factor_2 in reversed(xrange(1000)):
                    if n == (factor_2 * factor_1):
                        return n

    return 255  # Return a sentinel value to indicate failure.

if __name__ == '__main__':
    print problem_4()

Back to flipping out...

2009-02-02

Project Euler: Problem 3

I just finished playing around with Problem 3 from Project Euler. My rough idea for the problem was to build a list of the primes up to the number in question and start testing them, since I already knew the basics of the Sieve of Eratosthenes and I didn't know any good algorithms for factoring. After I finished, I looked up some of those algorithms (see Wikipedia's (incomplete) list), and I probably wouldn't have bothered with anything that complex anyway. Here's my original "solution" (I use scare quotes since I know it has bugs, and I knew it at the time).


"""Solves Problem 3 from Project Euler."""

import math

def e_sieve(upper_bound):
    """Uses the Sieve of Eratosthenes to get a list of the primes up to max."""
    primes = []
    candidates = range(2, upper_bound)
    while candidates:
        head = candidates[0]
        primes.append(head)
        candidates = [n for n in candidates[1:] if n % head]

    return primes

def find_highest_prime_factor(to_factor):
    """Find the highest prime factor of to_factor."""
    highest_prime = None
    primes = e_sieve(int(math.sqrt(to_factor)))
    for prime in primes:
        if not to_factor % prime:
            quotient = to_factor / prime
            if quotient in primes:
                return quotient
            highest_prime = prime

    if highest_prime:
        return highest_prime

    return to_factor

if __name__ == "__main__":
    print find_highest_prime_factor(600851475143)

After I coded up my solution and verified I had the right answer (after what seemed like forever but was really ~15 minutes), I just knew there had to be a faster way. A little work with cProfile quickly revealed I was spending all my time in the Sieve. It turns out that even a "naïve" algorithm for factoring is better than generating the primes for a number this size. My new, "naïve" solution ran in 44 ms.


"""Solves Problem 3 from Project Euler."""

def factorize(to_factor):
    """Use trial division to factorize to_factor and return all the resulting \
    factors."""
    factors = []
    divisor = 2
    while (divisor < to_factor):
        if not to_factor % divisor:
            to_factor /= divisor
            factors.append(divisor)
            # Note we don't bump the divisor here; if we did, we'd have
            # non-prime factors.
        elif divisor == 2:
            divisor += 1
        else:
            # Trivial optimization: skip even numbers that aren't 2.
            divisor += 2
    if not to_factor % divisor:
        # Don't forget the last factor
        factors.append(to_factor)
    return factors

if __name__ == "__main__":
    print max(factorize(600851475143))

As an added bonus, this code has fewer bugs. One of these days, I'll learn to follow KISS. One of these days….

Back to flipping out...