Project Euler 30 -- Digit fifth powers

Question

Here's my solution for Project Euler Problem 30. Any optimizations in time, space, or style would be greatly appreciated.

from timeit import default_timer as timer
def digit_powers(exponent):
 def power(k):
 return int(k) ** exponent
 if exponent <= 1:
 return "The exponent must be at least 2."
 else:
 total_sum = 0
 upper_bound = (exponent + 1) * (9 ** exponent)
 for number in range(10, upper_bound + 1):
 digits = [x for x in str(number)]
 if number == sum(map(power, digits)):
 total_sum += number
 return total_sum
start = timer()
ans = digit_powers(5)
elapsed_time = (timer() - start) * 1000 # s --> ms
print "Found %d in %r ms." % (ans, elapsed_time)

Answer 1

In looking at your code, it seems well written and solves the problem in a straight forward manner. The one thing I would do is to remove the else after your argument check and remove a level of indent. That just helps to keep things flat... Flat is better than nested.

Now onto optimization and potentially making it more complex than it needs to be in order to make it faster. In these problems it seems like a reasonable thing to do even though....Simple is better than complex.

Here are a couple things I looked at as potential issues:

1. When doing the sum(map(power, digits)), you are calculating the power of each digit many times. If this was a much more expensive calculation then only calculating 1 time is even more important.
2. In this problem there is a conversion from int->str->int. If you can breakdown the digits without converting them it saves some time. It seems like a good thing to know how to do just in case your ever trying a different language where skipping those transitions may be even faster.

To address these, I did two things:

1. Pre-calculate 1-9 powers and store the results in a dictionary so they don't get recalculated (Memoization). It lets me look them up very quickly without calculating each time.
2. Separate the digits using math.log10 as to avoid the int->str->int conversion.

The first item is fairly straight forward, calculating for each digit:

powers = {}
for a in range(10):
 powers[a] = a ** exponent

The second part is a little trickier, what it does is uses math.log10 which gives a count of the digits, then it does a mod 10 to get the right most digit, looks up that power to add to the total and then finally does integer division to divide by 10 effectively moving the decimal place 1 position to the left (dropping the right most digit we just dealt with). This is probably a little different between Python 2/3 due to changes in int/float division. It is basically doing what you are using [x for x in str(number)] to do but without converting to string and adding the power of 5 dict lookup as it goes along:

savei = i
for _ in range(int(math.log10(i)) + 1):
 digit = i % 10
 total += powers[digit]
 i //= 10

Applying those changes to the code make it about twice as fast. Here is what I end up with for the function:

import math
def digit_powers(exponent):
 if exponent <= 1:
 return "The exponent must be at least 2."
 powers = {}
 answer = 0
 # Get the powers
 for a in range(10):
 powers[a] = a ** exponent
 limit = (exponent + 1) * (9 ** exponent)
 # Search for them
 for i in range(10, limit):
 savei = i
 total = 0
 for _ in range(int(math.log10(i)) + 1):
 digit = i % 10
 total += powers[digit]
 i //= 10
 if (total == savei):
 answer += total
 return answer

Answer 2

You can make this much faster by precalculating the powers of each digit. Also note that the order of processing for digits doesn't matter so this code starts with the ones place. Here's what I came up with, and on my machine it ran 4x faster than the original code:

def digit_powers(exponent): 
 if exponent <= 1:
 return "The exponent must be at least 2."
 powdigits = [i**exponent for i in range(10)] 
 total_sum = 0
 upper_bound = (exponent + 1) * powdigits[9]
 for number in range(10, upper_bound + 1):
 partialsum = tempnum = number
 while tempnum:
 partialsum -= powdigits[tempnum%10]
 tempnum /= 10
 if not partialsum:
 total_sum += number
 return total_sum

• clutton: I like how you avoided the math import... for some reason didn't think of doing that. Also using a comprehension for the powdigits is nice.

Answer 3

I doubt you will see my answer. Hopefully it can help someone else in the future. First I want to compare the answers above speedwise. As a tease I will include my own attempt

Comparison

-------------------------------
| Name | Time |
| Joshua (OP) | 2000 ms |
| clutton | 500 ms |
| Edward | 300 ms |
| Nebuchadnezzar | 20 ms |
-------------------------------

So the answers above give a rough 4x-5x speed improvement. So how can we increase this to a 100x speed improvement? Well naively you check every number up to the limit \$ 5*9^5 \$. This means your code does 443839 checks, and for every check it has to compute the digit sum. In order to get a significant speed improvement we need a new algorithm.

New algorithm

One key observation is to note that once we have checked \0145ドル\$ we do not need to check \1450ドル\$ or any other combination of \$(0, 1, 4, 5)\$ since they all the the same power sum

$$ P_5[(0, 1, 4, 5)] = 0^5 + 1^5 + 4^5 + 5^5 = 4150 $$

To check whether if some permutation of \$(0, 1, 4, 5)\$ is equal to it's power sum all we need to do is to sort the digits in the power sum and compare (assuming the tuple is already sorted).

Implementation

One way to find the unique tuples is with combinations_with_replacement() from the itertools libary.

The code uses two checks to further decrease the number of combinations to check. Naively with combinations_with_replacements() there are a total of 5005 numbers to check. The check

if power_sum % 10 in perm:

is pretty intuitively. It checks whether the last number in the power_sum is in the permutations. A necessity if some permutation of perm is equal to the digit power sum. This reduces the number to check to 730. The final check is

if len_power > len_perm: return False

This cuts the number needed to check down to 591. Not this is not a huge improvement, however it is a very cheap test to perform and it is more useful for higher powers. I have attached a small table below

 P S Sum # 0 # 1 # 2 time
--------------------------------------------------------------
 3 4 1301 220 60 50 0.8685 ms
 4 3 19316 715 224 165 5.0320 ms
 5 5 248860 2002 730 591 10.7200 ms 
 6 1 548834 5005 1974 1603 32.0400 ms
 7 4 25679675 11440 5040 4220 77.1300 ms
 8 3 137949578 24310 12504 10812 136.8000 ms
 9 4 2066327172 48620 24380 21224 428.4000 ms
10 1 4679307774 92378 48706 42362 858.0000 ms
11 8 418030478906 167960 92512 82456 2502.0000 ms
12 0 0 293930 175452 159471 7021.0000 ms

P is the power, S is the number of solutions. # 0 is the total number of permutations to iterate over. # 1 is the remaining numbers after the first check, and # 2 is the remaining numbers after the second check.

Code

from itertools import combinations_with_replacement as CnR
def digit_power(power):
 total = -1
 for perm in CnR(xrange(10), power):
 power_sum = sum(i**power for i in perm)
 if power_sum % 10 in perm:
 power_sum_list = map(int, str(power_sum))
 if is_fift_power(power_sum_list, perm, power):
 total += power_sum
 return total
def is_fift_power(power_lst, perm, len_perm):
 len_power = len(power_lst)
 if len_power > len_perm: return False
 power_sort = [0]*(len_perm - len_power) + sorted(power_lst)
 return power_sort == list(perm)
if __name__ == '__main__':
 import timeit
 print digit_power(5)
 times = 100
 result = timeit.timeit(
 "digit_power(5)", number=times, setup="from __main__ import digit_power")
 print 1000*result/float(times), 'ms'