Project Euler與Python：Problem 31~40

Project Euler針對數學和程式愛好者，設計了一系列的問題，至今已有五百多道，供大家挑戰，每一題都要求應在1分鐘內計算出答案。

我將以Python語言來解題，原始碼放在GitHub（yehnan/project_euler_python），基本上都是先採用暴力法，若運算量太大，無法在1分鐘內算出答案，才著手研究想出更快的演算法。當題目的數字很小時，隨便寫就能在1分鐘內算出答案，但是，也就只有幾題而已，大部分的題目都需要大量計算，可沒那麼好混。

我使用硬體Raspberry Pi 2與軟體Raspbian，作為執行速度的依據。

注意：不提供詳細的說明與分析。

Problem 31: Coin sums，湊硬幣問題，英國硬幣面額有1、2、5、10、20、50、100、200，請問湊出總值200的湊法共有幾種。

若使用遞迴形式，也可以，但較慢；此處使用動態規劃法（dynamic programming）。

底下函式傳入200與coins_england，即可算出答案。

coins_england = (1, 2, 5, 10, 20, 50, 100, 200)

def coin_sum(total, coins):
    ways = [1] + ([0] * total)
    for coin in coins:
        for i in range(coin, total+1):
            ways[i] += ways[i - coin]
    return ways[total]

Problem 32: Pandigital products，a × b = c，而a、b、c的位數恰恰含有數字1到9，譬如39 × 186 = 7254。請找出所有c，求總和。

使用Python內建程式庫的排列函式。

from itertools import permutations as perm

底下函式判斷某數字（字串形式）是否含有數字1到9，而且各只有一個。

ds = '123456789'

def is_pandigital(num):
    ins = sum(1 for d in ds if d in num)
    return ins == len(ds)

分析後，a、b、c的位數只有兩種可能，分別是1、4、4，以及2、3、4，採用暴力法處理。底下函式會算出答案。
           
def pd():
    result = set()
    answers = set()
    for i in ds: # 1 4 4
        for jli in perm(ds[:int(i)] + ds[int(i)+1:], 4):
            j = ''.join(jli)
            k = str(int(i) * int(j))
            if len(k) == 4 and is_pandigital(i + j + k):
                result.add((i, j, k))
                answers.add(int(k))
    for ili in perm(ds, 2): # 2 3 4
        for jli in perm(''.join(set(ds)-set(ili)), 3):
            i = ''.join(ili)
            j = ''.join(jli)
            k = str(int(i) * int(j))
            if len(k) == 4 and is_pandigital(i + j + k):
                result.add((i, j, k))
                answers.add(int(k))
    return sum(answers)

Problem 33: Digit cancelling fractions，找出誤算也正確的約分。

沒什麼好說的，看清楚題目，努力寫出程式解題。

使用Python內建程式庫的gcd，求最大公因數。

from fractions import gcd

判斷分子n與分母d，是否符合題目的條件。

# n: numerator, d: denominator
def is_it(n, d):
    ns = str(n)
    ds = str(d)
    if ns[0] == ds[0]:
        x = int(ns[1])
        y = int(ds[1])
        if (n // gcd(n, d) == x // gcd(x, y) and
           d // gcd(n, d) == y // gcd(x, y)):
           return True
    elif ns[0] == ds[1]:
        x = int(ns[1])
        y = int(ds[0])
        if (n // gcd(n, d) == x // gcd(x, y) and
           d // gcd(n, d) == y // gcd(x, y)):
           return True
    elif ns[1] == ds[0]:
        x = int(ns[0])
        y = int(ds[1])
        if (n // gcd(n, d) == x // gcd(x, y) and
           d // gcd(n, d) == y // gcd(x, y)):
           return True
    return False

判斷某數字是否為能輕易判斷的。
   
def is_dd(num): # double digit, like 11, 22...
    return (num // 10) == (num % 10)
   
def is_m10(num): # mutiple of 10, like 10, 20...
    return (num % 10) == 0

題目本身，用暴力法迭代，找出答案。呼叫底下函式即可算出答案，會回傳乘積約分後的分子與分母，以及共有幾組組合符合條件，其中分母是此題答案。

def dcf():
    prod_n = 1
    prod_d = 1
    count = 0
    for d in range(10, 99+1):  # two digits
        for n in range(10, d): # less than 1
            if (not is_dd(n) and not is_dd(d) and
                not is_m10(n) and not is_m10(d) and
                is_it(n, d)):
                prod_n *= n
                prod_d *= d
                count += 1
    gg = gcd(prod_n, prod_d)
    return prod_n//gg, prod_d//gg, count 

Problem 34: Digit factorials，某些數字，把一個個位數拿出來，作階乘運算後相加，又會等於原本的數字。請找出這些數字、求總和。

使用Python內建程式庫的階乘函式，先算出0到9的階乘，加快速度。

from math import factorial
facts = [factorial(i) for i in range(10)]

找出上限，7位數已經遠遠超過其位數階乘和。

def limit(): # 9999999 is way more than its fact-sum
    n = 2
    while True:
        if 10**(n-1) > facts[9]*n:
            return n
        n += 1

檢查某數字是否符合條件。

def check(n):
    tmp = n
    total = 0
    while tmp > 0:
        total += facts[tmp % 10]
        tmp //= 10
    return total == n

問題本身，暴力法。呼叫底下函式便可算出答案。 
 
def df():
    result = 0
    for n in range(10, facts[9] * (limit()-1)):
        if check(n):
            result += n
    return result

Problem 35: Circular primes，197、971、719都是質數，稱為環形質數（circular prime），請問一百萬以下有幾個環形質數。

使用埃氏篩法找出小於n的質數。

def prime_sieve(n):
    primes = set(range(2, n+1))
    for i in range(2, (n+1) // 2):
        if i in primes:
            m = 2
            while i*m <= n:
                primes.discard(i*m)
                m += 1
    return primes

給定某數n，找出所有的環形排列。
  
def circulars(n):
    result = []
    s = str(n)
    for i in range(len(s)):
        result.append(int(s[i:] + s[0:i]))
    return result

判斷是否都是質數。 

def all_in(iterable, primes):
    for x in iterable:
        if x not in primes:
            return False
    return True

問題本身，底下函式傳入1000000，可算出答案。
   
def cp(max):
    primes = prime_sieve(max)
    count = 0
    for p in primes:
        cs = circulars(p)
        if all_in(cs, primes):
            count += 1
    return count

Problem 36: Double-base palindromes，數字585，不論是十進位表示法還是二進位表示法1001001001，都是迴文數（palindromic ）。小於一百萬的數字中，找出符合此條件的數字，求總和。

判斷十進位、二進位是不是迴文數。

def is_palindrome_10(n):
    s = str(n)
    return s == s[::-1]
   
def is_palindrome_2(n):
    s = bin(n)[2:]
    return s == s[::-1]
   
def is_palindrome(n):
    return is_palindrome_10(n) and is_palindrome_2(n)

問題本身，暴力法，底下函式傳入1000000，可算出答案。
   
def dbp(max):
    return sum(i for i in range(1, max) if is_palindrome(i))

Problem 37: Truncatable primes，3797若從左邊逐一拿掉一個位數，會變成797、97、7，通通都是質數，從右邊逐一拿掉也是。這種數字總共只有11個，不包括2、3、5、7，請找出這些數字求總和。

首先是判斷質數的函式。

def is_prime(n):
    if n < 2:
        return False
    if n == 2:
        return True
    if n % 2 == 0:
        return False
    for x in range(3, int(n ** 0.5) + 1, 2):
        if n % x == 0:
            return False
    return True

然後從左右兩邊逐一拿掉位數，並檢查是否仍為質數。

def is_truncatable_left(n):
    s = str(n)
    s_len = len(s)
    for i in range(1, s_len):
        x = s[i:]
        if not is_prime(int(x)):
            return False
    return True
       
def is_truncatable_right(n):
    s = str(n)
    s_len = len(s)
    for i in range(1, s_len):
        x = s[:-i]
        if not is_prime(int(x)):
            return False
    return True

為了簡便，把以上的檢查通通放在一個函式內。

def check_all(n):
    return is_prime(n) and is_truncatable_left(n) and is_truncatable_right(n)

最後是問題本身，因已知道共有11個，所以使用內建程式庫itertools的count，從10開始，判斷是否符合條件。底下函式傳入11，可求出答案。
 
import itertools

def tp(cnt_max):
    result = 0
    cnt = 0
    for n in itertools.count(10):
        if check_all(n):
            result += n
            cnt += 1
            if cnt == cnt_max:
                break
    return result
   

Problem 38: Pandigital multiples，以192為例，乘上1、2、3後得到192、384、576，此三數連接起來得到192384576，稱為連接乘積，恰好是9位數，而且剛好含有1到9。請問連接乘積中最大的是哪個。

首先撰寫函式判斷某數字是否為9位數、而且剛好含有1到9。

def is_pandigital(num):
    ds = '123456789'
    s = str(num)
    ins = sum(1 for d in ds if d in s)
    return ins == len(ds)

然後，底下函式給定某數字參數num，會乘上1、2、3、等等，然後連接，若位數是9個才回傳，若不是則回傳False。回傳值有兩個，前一個是連接乘積，後一個是最大乘上了多少，譬如192最大需乘上3。
           
def prod(num):
    s = str(num)
    n = 2
    while True:
        s += str(num * n)
        if len(s) == 9:
            return s, n
        elif len(s) > 9:
            return False, False
        n += 1
    return False, False

最後是問題本身，迭代1到4位數，因為5位數以上的連接乘積不會是9位數。然後以暴力法迭代每個數字，算出連接乘積，判斷是否符合條件，找出最大的。呼叫底下函式即可算出答案。

def pm():
    num_max = 0
    for i in range(1, 4+1):
        for num in range(10**(i-1), 10**i - 1):
            s, n = prod(num)
            if s and is_pandigital(s) and int(s) > num_max:
                num_max = int(s)
    return num_max

Problem 39: Integer right triangles，給定直角三角形周長p為120的話，可算出三組整數邊長的解{20,48,52}、{24,45,51}、{30,40,50}。在周長小於1000（包含）的情況下，請問哪個周長的直角三角形的整數邊長解最多。

用暴力法，會超過1分鐘，必須進行分析，縮減需要處理的計算量。

def irt(p_upperbound):
    p_max = 0
    t_max = 0
    for p in range(p_upperbound//4*2, p_upperbound+1, 2):
        t = 0
        for a in range(2, int(p / (2**0.5+1)) + 1):
            if p * (p-2*a) % (2* (p-a)) == 0:
                t += 1
        if t > t_max:
            t_max = t
            p_max = p
   
    return p_max

Problem 40: Champernowne's constant，把從1開始的數字連接起來，變成123456789101112131415161...，請問第1、10、100、1000、10000、100000、1000000個位數的乘積為何。

把從1開始的數字不斷連接起來，直到超過給定長度dn。

def d(dn):
    li = []
    i = 1
    n = 0
    while True:
        s = str(i)
        li.append(s)
        n += len(s)
        if n > dn:
            break
        i += 1
    return ''.join(li)

然後，拿出想要的位數，相乘。底下函式傳入7，可算出答案。

def cc(pow):
    s = d(10**pow)
    result = 1
    for i in range(pow):
        result *= int(s[10**i-1])
    return result