欧拉工程第十一至二十题解题思路分析

十一、网格中的最大乘积(largest product in a grid)

在下面这个20×20的网格中,对角线上相邻的四个元素已经用红色标记出来了,这四个数的乘积为26 × 63 × 78 × 14 = 1788696:

IMG

求在这个网格中,同一方向(上、下、左、右以及对角线)相邻四个元素最大乘积。

分析:此题似乎没有什么讨巧的简便算法,先将数据存入到一个嵌套列表当中,然后依次遍历每个元素,求它每个方向上相邻四个元素(如果有的话,需要条件判断)的乘积,存入到结果列表中,然后求最大值。

with open('euler/ep11.txt','r') as f:
    data = [map(int, s.split()) for s in f.readlines()]

def grid_largest_prod(data):
    import numpy as np
    length = len(data)
    res = np.zeros((length,length))
    right = down = down_right = left_down = 0
    for i in range(length):
        for j in range(length):
            right_exist = all([x<length for x in [i,i+1,i+2,i+3]])
            down_exist = all([x<length for x in [j,j+1,j+2,j+3]])
            left_exist = all([x>=0 for x in [i,i-1,i-2,i-3]])
            if right_exist:
                right = data[i][j]*data[i+1][j]*data[i+2][j]*data[i+3][j]
            if down_exist:
                down = data[i][j]*data[i][j+1]*data[i][j+2]*data[i][j+3]
            if right_exist and down_exist:
                down_right = data[i][j]*data[i+1][j+2]*data[i+2][j+2]*data[i+3][j+3]
            if left_exist and down_exist:
                left_down = data[i][j]*data[i-1][j+1]*data[i-2][j+2]*data[i-3][j+3]
            res[i][j] = max(right,down,down_right,left_down)
    max_res = np.max(res)
    return max_res

print(grid_largest_prod(data))

十二、高度可除的三角数(highly divisible triangular number)

三角数即由依次排列的自然数的和构成,所以第7个三角数是\(1 + 2 + 3 + 4 + 5 + 6 + 7 = 28\),前十个三角数是:\(1,3, 6, 10, 15, 21, 28, 36, 45, 55,\cdots\),让我们列出前七个三角数的因子:

$$ \begin{equation} \begin{aligned} 1&:1\\ 3&:1,3\\ 6&: 1,2,3,6\\ 10&: 1,2,5,10\\ 15&: 1,3,5,15\\ 21&: 1,3,7,21\\ 28&: 1,2,4,7,14,28 \end{aligned} \end{equation} $$
可以看出28是第一个因子超过5的三角数,求第一个因子超过五百万的三角数。

分析:首先根据初等数论的知识,我们知道一个数的因子个数等于其各质因数的幂分别加一并相乘,如\(20=2^2\times5^1\),则其因子数量为\((2+1)(1+1)=6\)。其次,我们知道三角数的通项是\(T_n=n(n+1)/2\),因为相邻的两个数\(n\)\(n+1\)必然是互质的,即两者没有共同的质因数,则有:

$$ n=p_1^{e_1}p_2^{e_2}\cdots p_i^{e_i},\ n+1=q_1^{f_1}q_2^{f_2}\cdots q_j^{f_j} $$

因此\(n\times(n+1)\)的质因数分解形式即为:

$$ n(n+1)=p_1^{e_1}p_2^{e_2}\cdots p_i^{e_i}\times q_1^{f_1}q_2^{f_2}\cdots q_j^{f_j} $$

因此其因子数为:

$$ (e_1+1)(e_2+1)\cdots(e_i+1)\times(f_1+1)(f_2+1)\cdots(f_j+1) $$

考虑到三角数还需要在\(n(n+1)\)的基础上再除以二,而\(n\)\(n+1\)中必然有一个为偶数,则二的幂应该减去一,则三角数的因子数应该为:

$$ e_1(e_2+1)\cdots(e_i+1)\times(f_1+1)(f_2+1)\cdots(f_j+1) $$

因此为了求三角数的因子数,我们只需要求\(n\)\(n+1\)的因子数再相乘即可,因为这两个数明显要小于它们的乘积形成的三角数,所以可以更快的进行质因数分解。我们从\(n=7\)开始迭代,每次加一,计算相应的三角数的因子数。因为这次迭代中的\(n+1\)会成为下次迭代中的\(n\),因此可以把中间值保存下来避免重复计算。当迭代计算出来的因子数量超过了五百时,迭代终止,此时返回相对应的三角数即为所求。

# time cost = 127 ms ± 1.09 ms

from sympy import factorint

def number_of_divisor(n):
    d = factorint(n)
    res = 1
    for i in d:
        res *= (d[i]+1) if i%2==1 else d[i]
    return res

def main(): 
    n = 7
    nd,nnd = number_of_divisor(n),number_of_divisor(n+1)
    while nd*nnd <= 500:
        n += 1
        nd,nnd = nnd,number_of_divisor(n+1)
    return n*(n+1)//2

十三、大数之和(large sum)

计算如下一百个50位数的和(见数据文件ep13.txt),求这个和的前十位数。

分析:这道题的解法是比较直白,将这些数存入到TXT文件中,用python依次读取并加总,然后将该和转化为字符串并取前十位,即为所求。

# time cost = 159 µs ± 858 ns

def main():
    large_sum = sum(map(int,open('data/ep13.txt')))
    ans = str(large_sum)[:10]
    return ans

十四、最长的考拉兹序列(longest Collatz sequence)

对于所有正整数,考虑如下迭代序列:

$$ n\rightarrow \frac{n}{2}\ if\ n\ is\ even,n\rightarrow3n+1\ if\ n\ is\ odd $$
从13开始并使用以上迭代规则,我们可以得到以下序列:
$$ 13\rightarrow40\rightarrow20\rightarrow10\rightarrow5\rightarrow16\rightarrow8\rightarrow4\rightarrow2\rightarrow1 $$
可以看出这个序列从13开始并到1结束总共包含10个数。考拉兹猜想指出使用以上迭代规则,所有正整数都会最终回到一,虽然这个猜想仍未得到证明。求在一百万以下,那个起始数可以产生最长的考拉兹序列?(注意:序列中包含的数的个数可以超过一百万。)

分析:此题的思路相对直接,只需要把一至一百万之间所有数的考拉兹序列的长度都算出来再求最大值对应的数字即可。可以优化的地方有两点:一是当\(n\)是奇数时,\(3n+1\)必然是偶数,则接下来必定需要除以二,则可以一步走完,也就是说当\(n\)是奇数时可以多算一步,直接计算\((3n+1)/2\),这样可以更快的到达序列的终点,但注意在计算累计步数时,应该加二而不是加一。另一个可以优化的点是在计算各个数考拉兹序列时,中间过程可能出现之前已经计算过的数字,为了避免重复计算,可以把之前数字对应的序列长度都缓存下来,之后碰到同样的数字,直接引用其值并加上之前已经走过的步数即可。代码如下:

# time cost = 1.26 s ± 105 ms

def main(N=10**6):
    d = {}
    for x in range(2,N):
        i,n = x,0
        while x != 1:
            if x < i:
                n = n + d[x]
                break
            elif x % 2 == 0:
                x = x // 2
                n += 1
            else:
                x = (3*x+1) // 2
                n += 2
        d[i] = n
    return max(d,key=d.get)

十五、格子路径(lattice paths)

在一个2×2的栅格中,从左上角出来,只能向右或向下移动,总共有6条路径可以到达栅格的右下角:

IMG

求在一个20×20的栅格中,有多少条移动路径?

分析:此题至少有两种解法,一种是应用排列组合知识的方法,另一种是动态规划的方法。下面先看排列组合的方法。如题目所述,要从(0,0)到达(2,2)我们需要向下走两步,向右走两步,总共需要走四步,这四步中有两步是向右走的,两步是向下走的。因此问题变成从四步任选出两步向右走,剩下的两步则向下走。比如你可以第一步和第二步向右走,剩下两步向下走,也就是图一;或者也可以第一步和第三步向右走,第二和第四步向下走,也就是图二。因此,这只是一个简单的组合问题,即从四步中任选出两步向右走共有多少种选法,显然有\(C(4,2)=6\)种选法。一般地,要从(0,0)走到(n,n)需要向下走n步,向右走n步,总共走2n步,因此总的可选路径共有\(C(2n,n)\)条,则题目要求答案即为\(C(40,20)\)。这个可以用组合数公式轻易计算出来。

第二种是动态规划的方法。首先考虑子问题,如在题目中一步可以到达(2,2)点的点只有两个,分别是(2,1)和(1,2),则到达(2,2)点的路径数也就是到达(2,1)点的路径数加上到达(1,2)点的路径数,而到达(2,1)点的路径数等于到达(1,1)点的路径数加上到达(2,0)点的路径数。因此一般地,到达点(x,y)的路径数等于到达(x-1,y)的路径数加上到达(x,y-1)的路径数,显然可以用自下而下的动态规划轻易解决。此题的边界条件是,一旦x或者y的坐标为零,则到达该点的路径只有一条。因此,我们将该动态规划问题的状态转移方程表述如下:

$$ p(x,y)= \begin{cases} p(x-1,y)+p(x,y-1) & x,y>0\\ 1 & x=0 \or y=0 \end{cases} $$

利用该状态转移方程我们可以轻易写出自下而上动态规划的实现,我们可以利用python中的lru_cache装饰器实现自下而下动态规划的记忆化,提高代码执行的效率。两种方法的代码如下:

# approach 1, time cost = 140 ns ± 11.8 ns

from functools import lru_cache

@lru_cache(maxsize=128)
def path(x,y):
    if x == 0 or y == 0:
        return 1
    else:
        return path(x-1,y) + path(x,y-1)

# approach 2, time cost = 1.76 µs ± 83.9 ns

from math import factorial as fac

def comb_num(n,k):
    num = fac(n)//(fac(n-k)*fac(k))
    return num

def main(n,k):
    return comb_num(n+k,k)

十六、指数各位数之和(power digit sum)

\(2^{15}=32768\)且各位数之和为\(3+2+7+6+8=26\),求\(2^{1000}\)各个位数之和。

分析:求出\(2^{1000}\),将其转化字符串,并将字符串转化为整数再求和。

# time cost = 51.6 µs ± 629 ns

def main():
    ans = sum(map(int,str(2**1000)))
    return ans

十七、数词字母数量(number letter counts)

如果把数字一到五用单词写出来,即one, two, three, four, five,那么总共使用的字母数量:

$$ 3+3+5+4+4=19 $$
如果把从一到一千的数字用对应的单词写出来,总共需要多少个字母?(注意:不要计算空格与连字符,如342(three hundred and forty-two)包含23个字母,而115(one hundred and fifteen)包含20个字母,and的使用遵守英式英语的规则。)

分析:将所有需要用到的数词都存入到一个CSV文件中并用pandas导入,并将其转化为一个词典。分别考\([1:20],[21:100],[101:1000]\)不同的数词对应规则,编写一个将数字转化为相应的词语的函数,但忽略空格与连字符。列出所有1000以下的数字,转化为对应的词语并统计总的字母数量。

# time cost = 4.31 ms ± 58.6 µs

import pandas as pd

def num_to_words(n):
    if n in words_dict.keys():
        return words_dict[n]
    else:
        if 21<=n<=100:
            ten_digit = n//10*10
            digit = n%10
            return words_dict[ten_digit] + words_dict[digit]
        elif 101<=n<=1000:
            hundred_dict = n//100
            remainder = n%100
            if remainder == 0:
                return words_dict[hundred_dict] + 'hundred'
            else:
                return words_dict[hundred_dict] + 'hundredand' + num_to_words(remainder)

def main():
    df = pd.read_csv('data/ep17.csv',header=None)
    words_dict = df.set_index(0).to_dict()[1]
    ans = sum([len(num_to_words(x)) for x in range(1,1001)])
    return ans

十八、最大和的路径(maximum path sum I)

从以下这个三角形的顶部开始,向相邻的下一行的数字移动,经过之数所能得到的最大的和为23,即:\(3+7+4+9=23\)

$$ 3\\ 7\quad4\\ 2\quad4\quad6\\ 8\quad5\quad9\quad3 $$
对于以下三角形(见数据文件ep18.txt),求期从上到下所有路径中最大的和。【注:由于总共只有16384条路径,因此可以尝试每条路径并求最大值。然而,问题67与此问题类似但有三角形共有100行,从而无法用暴力枚举方法解决,这需要一个聪明的方法。】

分析:这是一道经典的动态规划问题,为求三角形从上到下的最大和,先从最下一行开始倒推,即:

$$ max(8+2,5+2)=10,max(5+4,9+4)=13,max(9+6,3+6)=15 $$

这样可以将最下二行合为一行,即:

$$ 3\\ 7\quad4\\ 10\quad13\quad15 $$

依次类推,可以继续倒推出:

$$ max(10+7,13+7)=20,max(13+4,15+4)=19 $$

即:

$$ 3\\ 20\quad19 $$

容易得到23是最大值,而路径也是题目中所给出的路径。对于题目中给出的三角形,与之类似可以得出答案。代码中同时实现了从上而下和从下而上的两种动态规划方法:

with open('data/ep18.txt') as f:
    table = [list(map(int,s.split())) for s in f.readlines()]

# approach 1: bottom-up method, time cost = 62.1 µs ± 5.71 µs

def main():
    for row in range(len(table)-1, 0, -1):
        for col in range(0, row):
            table[row-1][col] += max(table[row][col], table[row][col+1])
    return table[0][0]

# approach 2: top-down method, time cost = 131 ns ± 8.97 ns

from functools import lru_cache

@lru_cache(maxsize=128)
def mps(r=0,c=0):
    if r == len(table)-2:
        return table[r][c] + max(table[r+1][c],table[r+1][c+1])
    else:
        return table[r][c] + max(mps(r+1,c),mps(r+1,c+1))

十九、计算星期天的数量(counting Sundays)

以下信息是题目直接给予你的,当然你也可以自己做一些研究:1900年1月1日是星期一,一月有30天的月份分别为九月、四月、六月与十一月,除开二月以外其它月份都有31天。二月在闰年有29天,其它年份有28天。闰年是年份可以整除四的年份,除非这个年份跨世纪的年份,但在跨世纪的年份中如果该年份可以整除400则也为闰年。在20世纪(1900年1月1日至2000年12月31日)每个月的首日是星期天的次数是多少?

分析:本题可以使用python中datetime模块解决,从1990年1月1日开始循环至2000年12月31日结束,如果某天同时为一个月的第一天且是星期天,则计数加一,这样便可以统计出整个二十世纪符合条件的天数。事实上存在一个计算任意日期是星期几的算法,叫作Zeller’s congruence,感兴趣的读者可以研究一下。使用这个算法,也可以解决这里的问题。

# time cost = 6.46 ms ± 10.9 µs

import datetime as dt

def count_sundays():
    delta = dt.timedelta(days=1)
    start_date = dt.datetime(1901, 1, 1)
    end_date = dt.datetime(2000, 12, 31)
    count = 0
    while start_date <= end_date:
        if start_date.day==1 and start_date.isoweekday()==7: count+=1
        start_date += delta
    return count

二十、阶乘各位数的和(factorial digit sum)

\(n!\)等于\(n\times(n-1)\times(n-2)\cdots\times3\times2\times1\),例如\(10!=10\times9\times8\cdots\times2\times1=3628800\),并且其各位数的和为\(3+6+0+8+8+0+0=27\),求\(100!\)各位数之和。

分析:这道题没有什么可深入分析的地方,只需要直接直接求出\(100!\),转化为字符串得到各位数,并计算其各位数之和即可。

# time cost = 27.2 µs ± 149 ns

from math import factorial as fac

def main():
    ans = sum(map(int,str(fac(100))))
    return ans