Project Euler 012
三角数の数列は自然数の和で表わされ, 7番目の三角数は 1 + 2 + 3 + 4 + 5 + 6 + 7 = 28 である. 三角数の最初の10項は:
1, 3, 6, 10, 15, 21, 28, 36, 45, 55, ...
となる.
最初の7項について, その約数を列挙すると, 以下のとおり.
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
これから, 7番目の三角数である28は, 5個より多く約数をもつ最初の三角数であることが分かる.
では, 500個より多く約数をもつ最初の三角数はいくつか.
->次の問題
この問題は割と力押しでやっており、処理時間が1秒を切れていません…
考え方
自然数nの約数の数は、nを構成する素因数の乗数から求めることができます。
例えばn=60のとき、
60 = 2^2 \times 3^1 \times 5^1
となります。この乗数にそれぞれ+1したものの積が約数の数になります。
(2+1) \times (1+1) \times (1+1) = 12
実際に60の約数は[1,2,3,4,5,6,10,12,15,20,30,60]の12個で合ってますね。
関数の呼び出し回数が増えると速度が遅くなってしまいますが、可読性と再利用のため約数の数を返す関数を作ります。以前作った素因数を返すprime_factor_generatorを使っています。
from math_functions import prime_factor_generator
def divisors_count(n: int) -> int:
"""
自然数nの約数の数を返す関数
素因数の乗数+1の総積が約数の数になるという方法を使用
:param n: int 自然数
:return: int 約数の数
"""
# nの素因数のlistを作成 n=60なら[2, 2, 3, 5]
prime_list = list(prime_factor_generator(n))
# countを使ってlistの中に何個を同じ素数があるか計算しそれに+1する
counts = [prime_list.count(p) + 1 for p in set(prime_list)]
result = 1
for c in counts:
result *= c # countsの総積を計算
return result
import numpy as np
from collections import deque
def prime_factor_generator(n: int) -> int:
"""
素因数分解するgenerator。素因数を順番に返す
:param n: int
:return: int
"""
prime_que = deque(get_prime_arr(int(n ** 0.5) + 2))
prime = 0
while prime ** 2 < n:
prime = prime_que.popleft()
while n % prime == 0:
n //= prime
yield prime
if n != 1:
yield n
def get_prime_arr_bool(n: int) -> np.ndarray:
"""
自然数n未満の素数をboolean配列で返す関数
2,3,5 -> [False, False, True, True, False, True]
:param n: int
:return: numpy.ndarray
"""
boolean_arr = np.ones((n,), dtype=np.bool) # n個の要素をもつのboolean配列(0~n-1に対応)を生成、全てTrueとする
boolean_arr[0:2] = False # [0,1]をFalseとする
boolean_arr[4::2] = False # 2以外の偶数をFalseとする
boolean_arr[6::3] = False # 3以外の3の倍数をFalseとする
for i in range(5, ceil(sqrt(n)), 6): # iは初項5、公差6の等差数列とする
if boolean_arr[i]: # 6k-1に対応。iに一致する要素がTrue=素数のとき
boolean_arr[i * 2::i] = False # iの倍数の要素(非素数)をFalseとする
if boolean_arr[i + 2]: # 6k+1に対応。
boolean_arr[(i + 2) * 2::(i + 2)] = False
return boolean_arr # boolean配列を返す
def get_prime_arr(n: int) -> np.ndarray:
"""
自然数n未満の素数の配列を返す関数
:param n: int
:return: numpy.ndarray
"""
return np.arange(n)[get_prime_arr_bool(n)] # 0~n-1の配列を生成、boolean配列でTrueとなっている数のみを返す
これを使って約数の個数を出していき、500を超えるまで検証します。力技です。
三角数は階差数列なので一般項から出すこともできますが、どちらにしろ小さい数から順に検証していくことになるのでシンプルに足していきます。
ちなみにsympyには約数の数を返す関数sympy.divisor_countがあります。
コード
from math_functions import prime_factor_generator
def main():
n = 500
triangle_number = 1 # 初項 1
difference = 2 # 階差数列の差
while divisors_count(triangle_number) <= n:
triangle_number += difference
difference += 1 # 差を1つずつ増やしていく
print(triangle_number)
def divisors_count(n: int) -> int:
"""
自然数nの約数の数を返す関数
素因数の乗数+1の総積が約数の数になるという方法を使用
:param n: int 自然数
:return: int 約数の数
"""
# nの素因数のlistを作成 n=60なら[2, 2, 3, 5]
prime_list = list(prime_factor_generator(n))
# countを使ってlistの中に何個を同じ素数があるか計算しそれに+1する
counts = [prime_list.count(p) + 1 for p in set(prime_list)]
result = 1
for c in counts:
result *= c # countsの総積を計算
return result
if __name__ == '__main__':
main()
結果
76576500
main処理時間: 1.590303897857666 sec.
このコードでは2秒前後、sympy.divisor_countを使っても0.5秒くらいかかってしまいます。
現在下から順番に総当りしている部分にどうにか高速化できる余地があるはずですが、思いつかないです…