https://taq.hatenadiary.jp/entry/2022/05/01/090000

こんにちは。
今回は「Project Euler」の第2問から第10問まで解いてみようと思います。

前回1問目を解いてみて、結構楽しかったのでとりあえずキリ良く10問くらいはやってみようということで、第2回です。

taq.hatenadiary.jp

以下、ネタバレになる（かもしれない）ので、嫌だという方は先に解いてみてください。

Problem 2 -Even Fibonacci numbers-

フィボナッチ数列の各項は、その前2つの数字を足すことで計算される。
1 2から始めると、初めの10項は1 2 3 5 8 13 21 34 55 89 となる。
400万を超えない範囲で考えたとき、偶数のフィボナッチ数を全て足すといくらになるか。

まずはフィボナッチ数列を計算するプログラムを書きます。

def fib(n):
    if n == 1:
        return 1
    elif n == 2:
        return 2
    else:
        return fib(n-2) + fib(n-1)

n = 5
print(fib(n))

定義をそのまま入れた感じです。
n = 5の場合は8が出力されます。

自分の環境では35くらいでちょっと時間がかかりました。
ちなみにそのときの答えは14930352なので、400万を超えており十分です。
ですが、本番ではこれが30回ほど繰り返されることになるため、このままでは計算に時間がかかりそうです。

そこで「メモ化」という手法を使って計算量を減らす努力をします。
メモ化の考え方は動的計画法と似ていますが、こちらの方が単純です。
（動的計画法はAtCoder Beginner Contest 248 - プログラミングの備忘録の問題Cのところで触れています。）

memo = {} #メモ
def fib(n):
    if n in memo: #結果がメモにあれば
        return memo[n] #メモから返す
    if n == 1:
        return 1
    elif n == 2:
        return 2
    else:
        f = fib(n-2) + fib(n-1)
        memo[n] = f #メモに書く
        return f

n = 50
print(fib(n))

{}で囲った場合は「辞書」と呼ばれ、[]で囲った「配列」と比べると、何故か分かりませんが辞書の方が扱いやすい。

強気にn = 50にしていますが速いです。
簡単な割に計算が爆速になるので便利です。

では、400万以下で偶数となるものだけ足していって答えを出してみます。

memo = {}
def fib(n):
    if n in memo:
        return memo[n]
    if n == 1:
        return 1
    elif n == 2:
        return 2
    else:
        f = fib(n-2) + fib(n-1)
        memo[n] = f
        return f

n = 0
ans = 0
while fib(n) < 4000000:
    if fib(n) % 2 == 0:
        ans += fib(n)
    n += 1
    
print(ans)

条件をそのまま適用した形になっています。

これを実行すると、エラーが出ます。
「RecursionError: maximum recursion depth exceeded in comparison 」とのことで、再帰回数の上限を超えてしまうようです。

memo = {0: 1, 1: 2}
def fib(n):
    if n in memo:
        return memo[n]
    else:
        f = fib(n-2) + fib(n-1)
        memo[n] = f
        return f

n = 0
ans = 0
while fib(n) < 4000000:
    if fib(n) % 2 == 0:
        ans += fib(n)
    n += 1
    
print(ans)

初期条件も辞書の中に入れただけですが、これだとうまくいくようです。
何故かはわかりません…
（手前にif文があるせいで再帰回数が割り増しになったりするのでしょうか。）

これはおまけ程度ですが、フィボナッチ数の計算部分をさらに短くしてみました。

memo = {0: 1, 1: 2}
def fib(n):
    if n not in memo:
        memo[n] = fib(n-2) + fib(n-1)
    return memo[n]

nがmemoの中にあろうがなかろうが結局返すのはmemo[n]なのでif n in memoの部分を除いた、というだけです。

以上で答えが4613732と出て、これを提出すると正解でした。

Problem 3 -Largest prime factor-

13195の素因数は5 7 13 29である。
600851475143の素因数の中で最大のものは何か。

RSA暗号チックな問題ですね。

最終的な答えは自分で入力するので出力のレギュレーションはありません。
とりあえずは素因数を片っ端から探していき、列挙してみます。

n = 6

for i in range(2, n):
    if n % i == 0:
        n /= i
        print(i)

iを2からn-1まで変えていき、nをiで割った余りが0（割り切れた）ならそのiが素因数なので、print()で出力します。
同時に、次の計算のためにnを先程求めた素因数で割っておきます。
n = 6の場合は2と3が出力されます。

とりあえず書いてみて遅かったりエラーが出たりしたら変えようと思っていたのですが、あっさり終わってしまい、答えらしき値を提出すると正解となりました。

というわけでn = 600851475143のときは71、839、1471、6857が出力され、その中で最大のもの6857を提出すると正解でした。

あまりにもあっさりしているので、割る数iを素数に限定したりして計算量を減らす努力をしても良いかもしれません。

Problem 4 -Largest palindrome product-

回文数はどちらから読んでも同じである。
2桁の数の積として表せる最大の回文数は9009 = 91 × 99である。
3桁の数の積として表せる最大の回文数は何か。

とりあえずはいろんな組み合わせで積を計算していき、回文数になったものを列挙してみます。

i = j = 999
while i > 0:
    while j > 0:
        a = i * j
        l = [k for k in str(a)]
        lr = [k for k in str(a)]
        lr.reverse()
        if l == lr:
            print(a)
        j -= 1
    i -= 1
print("end")

i、jを999から0まで変えていく間、i * jを計算し、それを1文字ずつ配列lに入れるとともにそれを逆から並べた配列lrをつくります。
このとき、lとlrが同じであったらaを出力します。
この場合はたくさん出力される中で最大のものは90909でした。
しかし、不正解。

それもそのはず、while文は条件式が正になったら終わってしまいます。
つまり、jが0になった時点で繰り返しは実質終わってしまっていたわけです。

というわけで、素直にfor文を使って繰り返すことにします。
さらに、大事なiとjの桁数の条件も加えました。

d = 3
pre = 0
for i in range(1000):
    for j in range(1000):
        a = i * j
        l = [k for k in str(a)]
        li = [k for k in str(a)]
        li.reverse()
        if len(str(i)) == len(str(j)) == d and l == li:
            print(a)

このままでは出力が多すぎてどれが最大値かわからないので、max()を使います。

d = 3
b = 0
for i in range(1000):
    for j in range(1000):
        a = i * j
        l = [k for k in str(a)]
        li = [k for k in str(a)]
        li.reverse()
        if len(str(i)) == len(str(j)) == d and l == li:
            b = max(b, a)

print(b)

aとbを比較し、より大きいほうをbに入れ直して次のステップの比較に備える、という感じです。

少し時間はかかりますが906609が出力され、これを提出すると正解でした。

Problem 5 -Smallest multiple-

2520は1から10の数で割り切ることができる最小の数である。
1から20の数で割り切ることができる自然数の中で最小のものは何か。

1~20で割り切れるということは$20!$かと思ってしまいますが、20で割り切れるということは2、4、5、10でも割り切れるということなので、その分は省いて考えなければなりません。

考え方としては、1~20まで素因数分解して、各素因数を見比べたときに数が一番多いものを採用し、それを数値に直すことで答えを得る、という感じです。（いわゆる「最小公倍数」）
1~10の場合はこんな感じの計算をすることになります。

まずは素因数分解をする関数をつくります。

def factorize(n):
    l = {}
    for i in range(2, int(n**0.5)+1):
        if n % i == 0:
            c = 0
            while n % i == 0:
                c += 1
                n //= i
            l[i] = c
    if n != 1:
        l[n] = 1
    return l

n = 100
print(factorize(n))

nについて割り切れるかどうか端から試していくわけですが、processingの備忘録 -ウラムの螺旋- - プログラミングの備忘録でエラトステネスの篩をつくったときと同様、range(2, n)とすると計算量が増えてしまうのでrange(2, int(n**0.5)+1)、つまりnの平方根（一応+1）までにしています。
int(n**0.5)まで繰り返し割ったけど割り切れずに残った分については、それも素数なのでif n != 1以下でaに追加しています。

出力の形式は、この後のことを考えて素因数とその個数という形にすることにしました。
n = 100なら、{2: 2, 5: 2}と出力されます。（$100 = 2^{2} \times 5^{2}$）

では続いて、各素因数で一番数が多いものを採用する、という操作を実装します。

def factorize(n):
    l = {}
    for i in range(2, int(n**0.5)+1):
        if n % i == 0:
            c = 0
            while n % i == 0:
                c += 1
                n //= i
            l[i] = c
    if n != 1:
        l[n] = 1
    return l

def sieve(n):
    if n == 1:
        return False
    for i in range(2, int(n**0.5)+1):
        if n % i == 0:
            return False
    return True

n = 10

ans = {}
for i in range(1, n+1):
    if sieve(i) == True:
        ans[i] = 0
    
for i in range(1, n+1):
    factors = factorize(i)
    for j in factors:
        ans[j] = max(factors[j], ans[j])

a = 1
for i in ans:
    a *= i ** ans[i]
print(a)

思ったより大作になってしまいました。

1フロックずつ見てみます。

ans = {}
for i in range(1, n+1):
    if sieve(i) == True:
        ans[i] = 0

値の比較をする上で答えを入れておく辞書ansをつくるのですが、キーとして素数が必要だったのでエラトステネスの篩を導入して素数かどうかの判定をしています。
素数であったら、辞書の中にキーとして入れておきます。

n = 10ならans = {2: 0, 3: 0, 5: 0, 7: 0}となります。

for i in range(1, n+1):
    factors = factorize(i)
    for j in factors:
        ans[j] = max(factors[j], ans[j])

素因数分解factrize()の結果をfactorsに保存しておき、先程つくったansと比較して各キーにおいてより数が大きいものを採用してansを更新する、という感じです。

a = 1
for i in ans:
    a *= i ** ans[i]
print(a)

最後にansから素因数とその数を取り出し、数値に直して出力します。

n = 10では2520が出力されたので合っているとみなし、n = 20とすると232792560と出力されたのでこれを提出すると正解でした。

Problem 6 -Sum square difference-

1から10までの自然数の二乗の和は385、和の二乗は3025である。
よって、この2つの差は3025 - 385 = 2640となる。
同じことを1から100までの自然数で行うと答えはどうなるか。

これは単純に計算すれば良いだけです。

n = 100

SUMofSQ = 0
s = 0
for i in range(1, n+1):
    SUMofSQ += i**2
    s += i
    SQofSUM = s**2
    
print(SQofSUM - SUMofSQ)

実行すると25164150と出力され、これを提出すると正解でした。

Problem 7 -10001st prime-

初めの6つの素数を並べると2 3 5 7 11 13となり、6番目の素数は13だとわかる。
10001番目の素数は何か。

Problem 5でつくったエラトステネスの篩が使えそうです。

def sieve(n):
    if n == 1:
        return False
    for i in range(2, int(n**0.5)+1):
        if n % i == 0:
            return False
    return True

n = 10
prime = []
i = 2
while len(prime) < 10001:
    if sieve(i) == True:
        prime.append(i)
    i += 1
print(prime[10000])

エラトステネスの篩を使ってiが素数かを判定し、素数なら配列primeに入れていきます。
そして、要素数が10001個になったら終わりにします。

結果、prime[10000]は104743となり、これを提出すると正解でした。

Problem 8 -Largest product in a series-

以下の1000桁の数の中で隣接した4個の数の積の最大値は9×9×8×9=5832である。
隣接する13個の数の積の最大値はいくつか。

片っ端から試してみます。

num = ~~~ #1000桁の数
l = [i for i in str(num)]

n = 13
ans = 0
for i in range(1000-n):
    product = 1
    for j in range(n):
        product *= int(l[i+j])
    ans = max(ans, product)
    
print(ans)

1000桁の数numを1文字ずつ配列lに入れ、0から1000-nまで隣接したn個分の数を取り出してその積をとる（product *= int(l[i+j])）という計算を繰り返します。
答え用にansを用意し、productとの比較でより大きい方を新たなansにしています（ans = max(ans, product)）。

n = 4では5832が出力されることが確認できたので、n = 13にすると23514624000が出力され、これを提出すると正解でした。

Problem 9 -Special Pythagorean triplet-

ピタゴラス数は以下の条件を満たす3つの数の組み合わせである。
a^2 + b^2 = c^2, a < b < c
例えば、3^2 + 4^2 = 9 + 16 = 25 = 5^2
a + b + c = 1000となるピタゴラス数が1組存在するが、その積abcは何か。

まずはピタゴラス数を探すコードを書いてみます。

n = 10
for a in range(1, n):
    for b in range(1, n):
        if a < b:
            c2 = a*a + b*b
            c = c2**0.5
            if b < c and c.is_integer() == True:
                print(a, b, c)

とりあえずaとbを0からnまで繰り返して（ただしa < b）cの二乗c2を計算します。
もしb < cかつcが整数（c.is_integer() == True）であったらそのときのa、b、cを出力します。

n = 10のときは3 4 5と6 8 10の組み合わせが出力されたので合っているということにしました。

続いて、a + b + c = 1000を満たすかどうかの判定を加えます。

n = 1000
for a in range(1, n):
    for b in range(1, n):
        if a < b:
            c2 = a*a + b*b
            c = c2**0.5
            if b < c and c.is_integer() == True:
                if a+b+c == 1000:
                    print(a, b, c)

といっても、さらにif a+b+c == 1000を加えるだけです。

和が1000になるところまで探せば良いので、とりあえずn = 1000としておいて実行してみると200 375 425のみが出力されました。
最後にこれらの積を計算して31875000となり、これを提出すると正解でした。

Problem 10 -Summation of primes-

10以下の素数の和は2 + 3 + 5 + 7 = 17である。
200万以下の素数の和はいくらになるか。

Problem 5でつくったエラトステネスの篩が役に立ちます。

def sieve(n):
    if n == 1:
        return False
    for i in range(2, int(n**0.5)+1):
        if n % i == 0:
            return False
    return True

prime = []
i = 2
ans = 0
while True:
    if sieve(i) == True:
        ans += i
    i += 1
    if i > 2000000:
        break

print(ans)

while Trueによる無限ループの間iが素数か判定し、もし素数なら答えに足していきます。
iが200万を超えたところで無限ループから抜け出し、ansを出力します。

n = 10で17が出力されたので良しとし、n = 2000000にすると少し時間がかかりましたが142913828922と出力されました。
これを提出すると正解でした。

以上、Project Eulerの問題を解いてみました。

Problem 100までならこのような形で答えを載せても良いらしいので、できるところまでやってみようかと思っています。