Find all ways to create "2016" using the numbers 1,2,3,4,5 and basic operations (parentheses, multiplication, addition, subtraction, exponentiation).

2016.py

from __future__ import division, print_function
from itertools import permutations, combinations_with_replacement, product

def guarded_exponentiation(a,b):
    """
    Raise a to the power of b, but only if the result won't be silly large, as
    silly largeness slows us way down.
    """
    if (a > 1 or a < -1) and b > 1000:
        raise OverflowError
    return a**b

op_table = {
  #'%': lambda a,b: a % b,
  '*': lambda a,b: a * b,
  '/': lambda a,b: a / b,
  '+': lambda a,b: a + b,
  '-': lambda a,b: a - b,
  '^': guarded_exponentiation
}
seeds = [1,2,3,4,5]
goal = 2016

def generate_replacement_positions(size):
    """
    Generate all permutations of an array, where the first position only has
    numbers (0 to size), the second has (0 to size-1), the third has (0 to
    size-2), ... (0,).
    e.g., for size 3:
        [0, 0, 0],
        [1, 0, 0],
        [2, 0, 0],
        [0, 1, 0],
        [1, 1, 0],
        [2, 1, 0]
    """
    return product(*[range(i) for i in range(size-1, 0, -1)])

def find_paths(operations, seeds, goal):
    i = 0
    wins = set()
    for ops in product(*([operations]*(len(seeds) - 1))):
        for nums in permutations(seeds):
            for replacement_positions in generate_replacement_positions(len(seeds)):
                if calculate(nums, ops, replacement_positions) == goal:
                    wins.add(report(nums, ops, replacement_positions))
    for win in wins:
        print(win[1:-1])
    print("TOTAL:", len(wins))

def report(nums, ops, replacement_positions):
    pool = list(nums[:])
    for op, replacement_pos in zip(ops, replacement_positions):
        a,b = pool[0:2]
        pool = pool[2:]
        result = "({}{}{})".format(a, op, b)
        pool.insert(replacement_pos, result)
    return pool[0]

def calculate(nums, ops, replacement_positions, printout=False):
    pool = list(nums[:])
    if printout:
        print()
        print(nums, ops, replacement_positions)
    for op, replacement_pos in zip(ops, replacement_positions):
        a,b = pool[0:2]
        pool = pool[2:]
        try:
            result = op_table[op](a, b)
            if result > 1000000000:
                break
            pool.insert(replacement_pos, result)
            if printout:
                print(a, op, b, "=", result, "->", replacement_pos, pool)
        except (ZeroDivisionError, ValueError, OverflowError):
            break
    else:
        return pool[0]
    return None

def find_paths_naive(operations, seeds, goal):
    for order in permutations(seeds):
        for ops in product(*([operations]*(len(seeds) - 1))):
            operand = order[0]
            for op,next_num in zip(ops, order[1:]):
                operand = op_table[op](operand, next_num)
            if operand == goal:
                print(operand, [[order[0]] + zip(ops, order[1:])])

if __name__ == "__main__":
    find_paths(op_table.keys(), seeds, goal)