- 本記事はProjectEulerの「100番以下の問題の説明は記載可能」という規定に基づいて回答のヒントが書かれていますので、自分である程度考えてみてから読まれることをお勧めします。
問題 23. 過剰数の和
原文 Problem 23: Non-abundant sums
問題の要約:過剰数の和にならない数の合計を求めよ
過剰数(Abundant number)(Wikipedia)の定義はこちら。今回は2つの過剰数の和で表せない数を列挙して合計を求めよと言うもの。問題では28123より大きい数は必ず和で表せるとありますが、Wikipediaによると20161が最大のようです。
過剰数の判定は「Problem 21: 友愛数」で作った関数sumpdivを使います。重複組合せのitertools.combinations_with_replacementを使ってすべての2つの数の和を求めて、その数をリストから除いて残った数を合計します。
import itertools
N = 20161
# create the list of abundunt numbers
abn = [i for i in range(1, N+1) if sumpdiv(i)> i]
abn2sum = [i for i in range(N+1)] # list of all natural numbers
for i, j in itertools.combinations_with_replacement(abn, 2):
if i+j<=N:
abn2sum[i+j] = 0 # erase the number by adding 2 abundant numbers combination
print(f"Answer: {sum(abn2sum)}")
過剰数の和にならない数は数が大きくなると減っていくみたいなので、グラフを書いてみました。確かにn=10000を超えるとかなり少なくなっていのが分かります。描画プログラムは下にあります。
import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(1,1,1)
abn2sumlist = [abn2sum[i] for i in range(1,N+1) if abn2sum[i]>0]
ax.hist(abn2sumlist, bins=300, range=(0,N))
ax.set_title('Numbers that can not be written as the sum of two abundant numbers ')
ax.set_xlabel('n')
ax.set_ylabel('Numbers')
fig.show()
(開発環境:Google Colab)