Project Euler　-Problem 51~60- - プログラミングの備忘録

こんにちは。
今回は「Project Euler」の第51問から第60問まで解いてみます。

50問目まではこちらに記事があります。

taq.hatenadiary.jp

51問目からはのんびりと、月1で増えたらいいねくらいの感覚でやっていきます。

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

Problem 51 -Prime digit replacements-
Problem 52 -Permuted multiples-
Problem 53 -Combinatoric selections-
Problem 54 -Poker hands-
Problem 55 -Lychrel numbers-
Problem 56 -Powerful digit sum-
Problem 57 -Square root convergents-
Problem 58 -Spiral primes-
Problem 59 -XOR decryption-
Problem 60 -Prime pair sets-

Problem 51 -Prime digit replacements-

2桁の数 *3 の1桁目を置き換えると、9種ある数の内6種 13, 23, 43, 53, 73, 83 が素数となる。
56**3 の3, 4桁目を同じ数字で置き換えると、この5桁の数は10種ある数の内7種 56003, 56113, 56333, 56443, 56663, 56773, 56993 が素数となる最初の例である。そして、56003 はその中で最も小さい素数となる。
数字の一部を置き換える (隣接した桁でなくても良い) ことによって8種の素数を生み出すことができるもので、生成される素数の内で最も小さいものを答えよ。

いろいろ考えてもわからなかったので、総当たり的に求めることにしました。

結果、コードがごちゃごちゃになり、実行時間も長くなってしまったので、他の方のものを参考にされた方が良いかもしれません…

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

def permutations(it):
    n = len(it)
    a = [i for i in it]
    
    if n <= 1:
        return [a]
    
    perm = []
    for i in range(n):
        a1 = [a[i]]
        a2 = a[:i] + a[i+1:]

        for j in permutations(a2):
            if a1 + j not in perm:
                perm.append(a1 + j)
    return perm

def combsWithRep(it, r):
    n = len(it)
    a = [i for i in it]
    
    if n == 0 or r == 0:
        yield []
    elif n == 1:
        yield a
    else:
        for i in range(n):
            a1 = [a[i]]
            a2 = a[i:]

            for j in combsWithRep(a2, r-1):
                pre = a1 + j
                if len(pre) < r:
                    for k in range(r-len(pre)):
                        pre.append(a[n-1])
                yield pre


target = 6

primes_list = []
for length in range(2, 10):
# length = 6;
    print("\n" + f"{length} digits")

    for s in combsWithRep('01', length):
        if s.count('0') == length:
            continue
        # print(s)
        
        for m in permutations(s):
            # print("\n" + str(m))

            residues = []
            for i in range(10**(length-1), 10**(length)):
                primes = []

                predigits = [str(k) for k in str(i)]

                residue = []
                count1 = 0
                for n in m:
                    if n == '1':
                        residue += '*'
                    else:
                        residue += predigits[count1]
                    count1 += 1
                if residue not in residues:
                    residues += [residue]
                else:
                    continue
                # print(residues)
                print("\r" + f"current: {residue}", end="")

                for j in range(0, 10):
                    # print(predigits)

                    digits = []
                    count2 = 0
                    for n in m:
                        if n == '1':
                            digits += str(j)
                        else:
                            digits += predigits[count2]
                        count2 += 1 
                    if digits[0] == '0':
                        continue

                    # print(digits)
                    # print("\r" + str(digits), end="")

                    num = 0
                    for k in range(0, len(digits)):
                        num += int(digits[-(k+1)])*(10**k)
                    # print(num)
                    # print("\r" + str(num), end="")

                    if sieve(num) == True:
                        primes += [num]

                if primes not in primes_list:
                    primes_list += primes

                    if len(primes) == target:
                        print("\n" + f">> {primes} << ")
                    # else:
                    #     print(primes)
    print("\n" + "end")

一応考え方を書いておくと、

0 と 1 から、combsWithRep() で置き換える場所の数を決める（0: 置き換えない、1: 置き換える）
得られた配列について permutations() で並び替えて、実際に置き換える場所を決める
全ての length 桁の自然数について 1 の桁を 0 ~ 9 に置き換える
素数か判定し、素数の数が target と同じだった場合に出力する

ここで、combsWithRep() は重複を許した組み合わせを、permutations() は順列を返す関数で、ライブラリを参考にしてつくっています。

taq.hatenadiary.jp

注意点として、ひとつ答えが出てもその桁の数全てについて処理し終わるまで答えが確定しないという点があります。

また、今何桁の数をこういう風に置き換えているよ、ということを出力して進捗がわかるようにすることで、動いてるんだかわからない状態でただ待たされることがないようにはしてみましたが、長いです。

実行すると、答えとして [121313, 222323, 323333, 424343, 525353, 626363, 828383, 929393] のみが出力され、この内で最小のもの 121313 を提出すると正解でした。
答えが確定するまで5時間半かかりましたが、これでも速くなった方なんです…

おそらく、全ての自然数でやっているから時間がかかるのだと思いますが、改善策が思いつかないのでこれで良しとします。

Problem 52 -Permuted multiples-

125874 とその2倍の 251748 は同じ数字を異なる順番で含む。
x, 2x, 3x, 4x, 5x, 6x が同じ数字を含むような、最小の x を求めよ。

総当たり的に求めることにしました。

m = 6 # 何倍まで考えるか

 n = 1
f = False
# for n in range(1, 10000):
while(f == False):
    print("\r" + str(n), end="")
    
    digits_list = []
    for j in range(1, m+1):
        digits = []
        for i in str(j*n):
            digits += str(i)
        # print(digits)
        digits.sort()
        digits_list += [digits]

    # print(digits_list)

    count = 0
    for i in range(1, len(digits_list)):
        if digits_list[0] == digits_list[i]:
            count += 1

    if count == len(digits_list)-1:
        # print("\n" + str(n))
        f = True
    
    n += 1

考え方としては、

1 ~ m 倍の数について、digits_list にその各桁の数字を (sort() で並び替えて) 保存しておく
digits_list の各要素について digits_list[0] (= x) と比較し、全てが同じ (count == len(digits_list)-1) ならば、そのときの n を答えとする

これをwhile文で繰り返し、条件を満たす n が決定された時点で抜けるようにしています。

実行すると、142857 で処理が止まり、これを提出すると正解でした。
（さすがはDifficulty Ratingが5%なだけあって、51問目と比べると非常に簡単です…）

Problem 53 -Combinatoric selections-

5個の中から3個を選ぶ方法は10通りある。
（12345 なら、123, 124, 125, 134, 135, 145, 234, 235, 245, 345 の10通り）
組み合わせの表記として $\displaystyle{ 5 \choose 3 } = 10$ を使う。
一般に $\displaystyle{ n \choose r } = \frac{ n ! }{ r ! (n - r) ! }$ である。
ここで、$r \leq n$、$n ! = n \times (n - 1) \times \ldots \times 3 \times 2 \times 1$、$0 ! = 1$ である。
$n = 23$ で初めて値が100万を超える。$\displaystyle{ 23 \choose 10 } = 1144066$
$1 \leq n \leq 100$ について、$\displaystyle{ n \choose r }$ が100万を超えるものはいくつあるか。

全てのパターンについて計算して解きます。

def factorial(n):
    if n == 0:
        return 1
    else:
        a = 1
        for i in range(n, 0, -1):
            a *= i
        return a

count = 0
for n in range(1, 101):
    for r in range(1, n+1):
        value = factorial(n)/(factorial(r)*factorial(n-r))
        
        if value > 1000000:
            count += 1

print(count)

階乗を計算する関数 factorial() をつくり、全通りの組み合わせ数を計算し、それが100万を超えたときに count を1増やします。
全ての n、r について計算が終わったところで count を出力し、答えとします。

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

Problem 54 -Poker hands-

1000戦分のデータが入った「poker.txt」というファイルがある。
各行に（スペースで区切られた）10枚のカードが書かれ、初めの5枚は Player 1 の、残りの5枚が Player 2 のカードを表す。
全ての手札が有効（無効な文字や重複したカードが無い）であり、各試合に必ず勝者がいるとする。
何回 Player 1 が勝つか。

（長いので本題の部分だけ書きました。）

競プロ味のある問題が出てきました。

役の判定

まずは、手札の役を判定するコードを書いてみます。

def judgeHand(s):
    # 手札の情報を扱いやすい形に変更する
    hands = []
    for i in s:
        if i == ' ':
            continue
        hands += i
    # print(hands)
    
    numbers = []
    suits = []
    for i in range(0, 10, 2):
        numbers += hands[i]
        suits += hands[i+1]
    # print(numbers)    
    # print(suits)
    
    number = {'A': 0, '2': 0, '3': 0, '4': 0, '5': 0, '6': 0, '7': 0, '8': 0, '9': 0, 'T': 0, 'J': 0, 'Q': 0, 'K': 0}
    number_data = {}
    for i in number:
        number_data[i] = numbers.count(i)
    # print(number_data)
        
    suit = {'H': 0, 'D': 0, 'S': 0, 'C': 0}
    suit_data = {}
    for i in suit:
        suit_data[i] = suits.count(i)
    # print(suit_data)
    
    value_order = {'2': 13, '3': 12, '4': 11, '5': 10, '6': 9, '7': 8, '8': 7, '9': 6, 'T': 5, 'J': 4, 'Q': 3, 'K': 2, 'A': 1}
    value_data = []
    for i in numbers:
        value_data += [value_order[i]]
    value_data.sort()
    # print(value_data)
    
    # 役の判定
    hand = ""
    val = []
    
    count = 0
    for i in number_data:
        if number_data[i] == 2:
            val += [value_order[i]]
            count += 1
    val.sort()
    if count == 1:
        hand = "One Pair"
    elif count == 2:
        hand = "Two Pairs"
    
    for i in number_data:
        if number_data[i] == 3:
            hand = "Three of a Kind"
            val = [value_order[i]]
    
    count2 = 0
    for i in range(0, 4):
        if value_data[i]+1 == value_data[i+1]:
            count2 += 1
    if count2 == 4:
        hand = "Straight"
        val = []
        
    for i in suit_data:
        if suit_data[i] == 5:
            hand = "Flush"
            val = []
    
    for i in number_data:
        if number_data[i] == 3:
            for j in number_data:
                if number_data[j] == 2:
                    hand = "Full House"
                    val = [value_order[i]]
    
    for i in number_data:
        if number_data[i] == 4:
            hand = "Four of a Kind"
            val = [value_order[i]]
    
    for i in suit_data:
        if suit_data[i] == 4:
            count3 = 0
            for i in range(0, 4):
                if value_data[i]+1 == value_data[i+1]:
                    count3 += 1
            if count3 == 4:
                hand = "Straight Flush"
                val = []
    
    for i in suit_data:
        if suit_data[i] == 5:
            count4 = 0
            for j in ['T', 'J', 'Q', 'K', 'A']:
                if number_data[j] == 1:
                    count4 += 1
            if count4 == 5:
                hand = "Royal Flush"
                val = []
             
    if hand == "":
        hand = "High Card"
        val = []
        
    return [hand, val, value_data]
    

s = ["5H 5C 6S 7S KD", 
     "2C 3S 8S 8D TD", 
     "5D 8C 9S JS AC", 
     "2C 5C 7D 8S QH", 
     "2D 9C AS AH AC", 
     "3D 6D 7D TD QD", 
     "4D 6S 9H QH QC", 
     "3D 6D 7H QD QS", 
     "2H 2D 4C 4D 4S", 
     "3C 3D 3S 9S 9D"]

for i in s:
    print(judgeHand(i))

役の条件をそのまま書いたのでなかなか煩雑になってしまいましたが、引数に手札のデータを与えると役と手札の数字の強さを返す関数 judgeHand() をつくりました。
例えば、手札が「5H 5C 6S 7S KD」なら、['One Pair', [10], [2, 8, 9, 10, 10]] と返されます。
これは、役が10番目に強い数字 (= 5) のワンペアで、手札の強さは [2, 8, 9, 10, 10] (= K, 7, 6, 5, 5) であることを表しています。

数字の強さは、A が一番強く (= 1)、2 が13番目に強い (= 13)、という形で扱っており、カードの数字そのものではないので注意してください。

勝敗の判定

役の判定ができれば、後は各試合について適用し、勝敗を決めるだけです。

初めに、1試合分だけ勝敗判定のコードを書きました。

# def judgeHand(s):
    # 省略

# 試合
s = "5H 5C 6S 7S KD 2C 3S 8S 8D TD"
# s = "5D 8C 9S JS AC 2C 5C 7D 8S QH"
# s = "2D 9C AS AH AC 3D 6D 7D TD QD"
# s = "4D 6S 9H QH QC 3D 6D 7H QD QS"
# s = "2H 2D 4C 4D 4S 3C 3D 3S 9S 9D"

# データを扱いやすいように変形
battle = []
for i in s:
    if i == ' ':
        continue
    battle += i
# print(battle)

hands1 = []
hands2 = []
for i in range(0, 10):
    if battle[i] == ' ':
        continue
    hands1 += battle[i]
    hands2 += battle[i+10]

# 役の判定
hand1 = judgeHand(hands1)
print(f"player 1: {hand1}")
hand2 = judgeHand(hands2)
print(f"player 2: {hand2}")

# 役の強さ
hand_rank = {"High Card": 10, 
             "One Pair": 9, 
             "Two Pairs": 8, 
             "Three of a Kind": 7,
             "Straight": 6, 
             "Flush": 5, 
             "Full House": 4, 
             "Four of a Kind": 3, 
             "Straight Flush": 2, 
             "Royal Flush": 1}
            
# 勝敗の判定
win = False

# 役の強さを比較
if hand_rank[hand1[0]] < hand_rank[hand2[0]]:
    winner = 1
    win = True
    win_count += 1
elif hand_rank[hand1[0]] > hand_rank[hand2[0]]:
    winner = 2
    win = True

# 役を構成する数字の強さを比較
if win == False:
    for i in range(0, len(hand1[1])):
        if hand1[1][i] < hand2[1][i]:
            winner = 1
            win = True
            win_count += 1
            break
        elif hand1[1][i] > hand2[1][i]:
            winner = 2
            win = True
            break

# 手札の数字の強さを比較      
if win == False:
    for i in range(0, 5):
        if hand1[2][i] < hand2[2][i]:
            winner = 1
            win = True
            win_count += 1
            break
        elif hand1[2][i] > hand2[2][i]:
            winner = 2
            win = True
            break

# 勝者の出力
if winner == 1:
    print("winner: player 1")
else:
    print("winner: player 2")

例えば、s = "5H 5C 6S 7S KD 2C 3S 8S 8D TD" という試合データが与えられれば、

player 1: ['One Pair', [10], [2, 8, 9, 10, 10]]
player 2: ['One Pair', [7], [5, 7, 7, 12, 13]]
winner: player 2

が出力され、player 2 が勝つと判定されます。

ファイルの読み込み

では続いて、ファイルを読み込みます。

f = open("p054_poker.txt") # ファイルを読み込む
l = f.readlines() # 一行ごとに読み取って配列に
# print(l)
battles = []
for i in range(0, len(l)):
    battles += [l[i].rstrip("\n")] # 末尾の"\n" (改行記号) を取って試合データに

print(battles)

答えの出力

長くなりましたが、いよいよ答えを求めます。

これまでのコードを組み合わせれば良いですが、player 1 の勝つ回数を問われているので win_count の計算と出力のために一部書き加えています。

def judgeHand(s):
    # 手札の情報を扱いやすい形に変更する
    hands = []
    for i in s:
        if i == ' ':
            continue
        hands += i
    # print(hands)
    
    numbers = []
    suits = []
    for i in range(0, 10, 2):
        numbers += hands[i]
        suits += hands[i+1]
    # print(numbers)    
    # print(suits)
    
    number = {'A': 0, '2': 0, '3': 0, '4': 0, '5': 0, '6': 0, '7': 0, '8': 0, '9': 0, 'T': 0, 'J': 0, 'Q': 0, 'K': 0}
    number_data = {}
    for i in number:
        number_data[i] = numbers.count(i)
    # print(number_data)
        
    suit = {'H': 0, 'D': 0, 'S': 0, 'C': 0}
    suit_data = {}
    for i in suit:
        suit_data[i] = suits.count(i)
    # print(suit_data)
    
    value_order = {'2': 13, '3': 12, '4': 11, '5': 10, '6': 9, '7': 8, '8': 7, '9': 6, 'T': 5, 'J': 4, 'Q': 3, 'K': 2, 'A': 1}
    value_data = []
    for i in numbers:
        value_data += [value_order[i]]
    value_data.sort()
    # print(values)
    
    # 役の判定
    hand = ""
    val = []
    
    count = 0
    for i in number_data:
        if number_data[i] == 2:
            val += [value_order[i]]
            count += 1
    val.sort()
    if count == 1:
        hand = "One Pair"
    elif count == 2:
        hand = "Two Pairs"
    
    for i in number_data:
        if number_data[i] == 3:
            hand = "Three of a Kind"
            val = [value_order[i]]
    
    count2 = 0
    for i in range(0, 4):
        if value_data[i]+1 == value_data[i+1]:
            count2 += 1
    if count2 == 4:
        hand = "Straight"
        val = []
        
    for i in suit_data:
        if suit_data[i] == 5:
            hand = "Flush"
            val = []
    
    for i in number_data:
        if number_data[i] == 3:
            for j in number_data:
                if number_data[j] == 2:
                    hand = "Full House"
                    val = [value_order[i]]
    
    for i in number_data:
        if number_data[i] == 4:
            hand = "Four of a Kind"
            val = [value_order[i]]
    
    for i in suit_data:
        if suit_data[i] == 4:
            count3 = 0
            for i in range(0, 4):
                if value_data[i]+1 == value_data[i+1]:
                    count3 += 1
            if count3 == 4:
                hand = "Straight Flush"
                val = []
    
    for i in suit_data:
        if suit_data[i] == 5:
            count4 = 0
            for j in ['T', 'J', 'Q', 'K', 'A']:
                if number_data[j] == 1:
                    count4 += 1
            if count4 == 5:
                hand = "Royal Flush"
                val = []
             
    if hand == "":
        hand = "High Card"
        val = []
        
    return [hand, val, value_data]


f = open("p054_poker.txt") # ファイルを読み込む
l = f.readlines() # 一行ごとに読み取って配列に
# print(l)
battles = []
for i in range(0, len(l)):
    battles += [l[i].rstrip("\n")] # 末尾の"\n" (改行記号) を取って試合データに
# print(battles)

win_count = 0
for s in battles:
    # データを扱いやすいように変形
    battle = []
    for i in s:
        if i == ' ':
            continue
        battle += i
    # print(battle)

    hands1 = []
    hands2 = []
    for i in range(0, 10):
        if battle[i] == ' ':
            continue
        hands1 += battle[i]
        hands2 += battle[i+10]

    # 役の判定
    hand1 = judgeHand(hands1)
    # print(f"player 1: {hand1}")
    hand2 = judgeHand(hands2)
    # print(f"player 2: {hand2}")

    # 役の強さ
    hand_rank = {"High Card": 10, 
                 "One Pair": 9, 
                 "Two Pairs": 8, 
                 "Three of a Kind": 7,
                 "Straight": 6, 
                 "Flush": 5, 
                 "Full House": 4, 
                 "Four of a Kind": 3, 
                 "Straight Flush": 2, 
                 "Royal Flush": 1}

# 勝敗の判定
win = False

# 役の強さを比較
if hand_rank[hand1[0]] < hand_rank[hand2[0]]:
    winner = 1
    win = True
    win_count += 1
elif hand_rank[hand1[0]] > hand_rank[hand2[0]]:
    winner = 2
    win = True

# 役を構成する数字の強さを比較
if win == False:
    for i in range(0, len(hand1[1])):
        if hand1[1][i] < hand2[1][i]:
            winner = 1
            win = True
            win_count += 1
            break
        elif hand1[1][i] > hand2[1][i]:
            winner = 2
            win = True
            break

# 手札の数字の強さを比較      
if win == False:
    for i in range(0, 5):
        if hand1[2][i] < hand2[2][i]:
            winner = 1
            win = True
            win_count += 1
            break
        elif hand1[2][i] > hand2[2][i]:
            winner = 2
            win = True
            break

if winner == 1:
    print("winner: player 1")
else:
    print("winner: player 2")
    
print(win_count)

実行すると player 1 の勝利数として 376 が出力され、これを提出すると正解でした。

Problem 55 -Lychrel numbers-

47 は、反転させたものとの和をとる (47 + 74 = 121) と回文数になる。
全ての数字が早く回文数になるとは限らない。例えば、
349 + 943 = 1292
1292 + 2921 = 4213
4213 + 3124 = 7337
これより、349 は回文数になるまで3回繰り返す必要がある。

未だ誰も証明していないが、196 のような一部の数字は回文数にはならない。
このような数を「リクレル数」と呼ぶ。
その性質やこの問題の目的から、そうでないことが示されるまではその数はリクレル数であると仮定する。
加えて 10000 以下の全ての数については、50回以下の繰り返しで回文数になるかどうかで判断するものとする。
実際、10677 が回文数になるまで50回以上の繰り返しを要する最小の数であることが示されている。
(4668731596684224866951378664 (53 回, 28 桁))

驚くべきことに、リクレル数でありながら自身が回文数となっている数もあり、最小のものは 4994 である。

10000 以下の数の中で、リクレル数はいくつあるか。

# nを反転して和をとる関数
def calc(n):
    rev = []
    for i in str(n):
        rev.insert(0, i)

    n_rev = 0
    for i in range(0, len(rev)):
        n_rev += int(rev[i])*(10**(len(rev)-i-1))
    
    # print(n, n_rev)
    return n + n_rev

# nが回文数か判定する関数
def isPalindrome(n):
    l = []
    l_rev = []
    for i in str(n):
        l += i
        l_rev.insert(0, i)

    if l == l_rev:
        return True
    else:
        return False

# nがリクレル数か判定する関数
def isLychrel(n, count=0):
    v = calc(n)
    count += 1
    # print(v, count)
    
    if isPalindrome(v) == True:
        return [False, v, count]
    else:
        if count == 50:
            return [True]
        return isLychrel(v, count)


# リクレル数か判定する数の最大値
n_max = 10000

count = 0
for i in range(1, n_max+1): # n_maxまでの全ての数について
    print("\r" + str(i), end='') # (判定中の数を出力)
    if isLychrel(i)[0] == True: # リクレル数なら
        count += 1 # カウント

print("\n" + str(count))

素直に判定をしています。
コード中のコメントが全てです。

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

Problem 56 -Powerful digit sum-

グーゴル ($10^{100}$) は巨大な数であり、1 の後に 0 が100個続く。
$100^{100}$ は考えられないほど大きく、1 の後に 0 が200個続く。
しかし、数としての大きさにもかかわらず、各桁の数字の和は 1 しかない。

$a^b$ ($a, b < 100$) について、各桁の数字の和の最大値は何か。

シンプルな計算問題がきました。

ans = 0
for a in range(0, 100):
    for b in range(0, 100):
        v = a**b
        
        s = 0
        for i in str(v):
            s += int(i)
        
        ans = max(ans, s)

print(ans)

特に解説することもありません。

実行すると 972 が出力され、これを提出すると正解でした。
ちなみにこのとき、a = 99、b = 95 でした。

Problem 57 -Square root convergents-

2 の平方根は無限に続く分数の形で表すことができる。
$$ \sqrt{ 2 } = 1 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 + \ldots } } } $$

初めの4項について計算すると、
$$ \begin{align} &1 + \frac{ 1 }{ 2 } = \frac{ 3 }{ 2 } = 1.5 \\ &1 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 } } = \frac{ 7 }{ 5 } = 1.4 \\ &1 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 } } } = \frac{ 17 }{ 12 } = 1.41666 \ldots \\ &1 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 } } } } = \frac{ 41 }{ 29 } = 1.41379 \ldots \end{align} $$

これに続く3項は $\frac{ 99 }{ 70 }$、$\frac{ 239 }{ 169 }$、$\frac{ 577 }{ 408 }$ となるが、8項目は $\frac{ 1393 }{ 985 }$ となり、初めて分子の数の桁数が分母の数の桁数を超える。

1000項目までについて、分母の数の桁数よりも分子の数の桁数が大きくなるものはいくつか。

法則性を見つけて、それを元に計算していきます。

さて、1項目について $\displaystyle 1 + \frac{ 1 }{ 2 } = \frac{ 3 }{ 2 }$ となります。
2項目は $\displaystyle 1 + \frac{ 1 }{ 2 + \frac{ 1 }{ 2 } } = 1 + \frac{ 1 }{ 1 + \frac{ 3 }{ 2 } } = 1 + \frac{ 2 }{ 2 + 3 } = \frac{ 7 }{ 5 }$ となるので、分母は「前の項の分母・分子の和」、分子は「前の項の分母・分子の和にさらに分母を足した値」で計算できることがわかります。

これを数式で書くと、

$$ \begin{align} n_{k+1} &= n_{k} + 2d_{k} \\ d_{k+1} &= n_{k} + d_{k} \end{align} $$

ここで、$n$ は分子 (numerator)、$d$ は分母 (denominator) を表します。

というわけで、これを元にコードを書いて答えを求めます。

# 分母・分子の初期値 (1項目)
n = 3
d = 2

# 分母・分子の計算
n_list = [n]
d_list = [d]
for i in range(1, 1000):
    nextN = n + 2*d
    nextD = n + d
    
    n = nextN
    d = nextD
    
    n_list += [n]
    d_list += [d]
    
# print(n_list, d_list)

# 条件を満たすものを数える
count = 0
for i in range(0, 1000):
    if len(str(n_list[i])) > len(str(d_list[i])):
        count += 1

print(count)

実行すると 153 が出力され、これを提出すると正解でした。
ちなみに1000項目は、

$$ \frac{ 72016336943533875056131468444247239328723197628440751797201898063588088312700201943482948477109536203740206649612702729920170001354454107173480483962605519493117789821758457767858986227019805650639002566946496865364666543562826303377877700877266135276209125278037204443042487153089226849544681245300260167141025277156482737568934079466850318276966893735585103845471745828701580706481 }{ 50923240208988086528630631809520139740233813238121746983184087440294876496910992032199150140268994975929379287364741267843829764622711478409393315655847391426474051205489599791876548045932140034723376284514492359420266988954982032735984900118578499266240851000814909302953919098957411314364525062171896557231165855421007859932974745573266780329830972204644845348897749854049994681209 } $$ $$ = 1.414213 562373 095145 $$

でした。
（383桁の数の分数で、mathjaxで表記するとめちゃくちゃはみ出て面白かったのでそのままにしました。）

$\sqrt{ 2 } = 1.414213 562373 095048 \ldots$ らしいので、それなりに近い値となっています。

Problem 58 -Spiral primes-

1 から始めて反時計回りに並べる。

興味深いことに右下方向に奇数の平方数が並ぶが、さらに興味深いことに斜め方向の数字13個中の8個のが素数となっている。この割合は $8/13 \approx 62\%$ である。

さらに外側に並べていくことを考えると、一辺に何個並んだときに斜め方向にある素数の割合が $10\%$ を切るか。

いわゆる「ウラムの螺旋」というやつでしょうか。

（過去にprocessingでつくったことがあるので宣伝しておきます。）

taq.hatenadiary.jp

さて、宣伝だけして終わりにするつもりでしたが方法を参考にできそうだったので、全く同じではないですが似た方法を使いました。
詳しくは先程挙げた記事を参照してください。
（わかりやすいとは限りませんが…）

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


# 一辺の長さ
l = 3

ratio = 1
n = 1
count = 0
d = 2
diagonals = [n]
primes_diagonal = []
while ratio > 0.1:
    # 斜め方向の数を配列に
    newdiagonals = []
    while n < l**2:
        n = n + d
        newdiagonals += [n]
        if count % 4 == 3:
            d += 2
        count += 1 
    diagonals.extend(newdiagonals)
    # print(diagonals)
    # print(len(diagonals))

    # 斜め方向の数について素数を数える
    for i in newdiagonals:
        if sieve(i) == True:
            primes_diagonal += [i]
    # print(primes_diagonal)
    # print(len(primes_diagonal))

    # 割合の計算
    ratio = len(primes_diagonal)/len(diagonals)
    print(f"\r{str(l)} {str(ratio)}", end="")
    
    l += 2

考え方としては、斜め方向の番号を配列 diagonals に用意しておき、この配列中の素数 primes_diagonal と全体 diagonals の数を数えて割合 ratio を出す、という感じです。
これを1周増やす (l += 2) ごとに行い、割合 ratio が 0.1 より小さくなるまで繰り返します。

ただ、毎回 l**2 まで diagonals をつくったり、diagonals の全てについて素数を数えていたのでは時間がかかってしまいます。
せっかくその前のステップで計算した分があるのにそれを利用しない手はありません。

ということで、前の結果を利用するような工夫も一応入れ、それなりの速度で答えが出るようにはできました。

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

Problem 59 -XOR decryption-

コンピューター上の文字は特定のコードが割り当てられており、その規格をASCII (American Standard Code for Information Interchange) という。例えば、A = 65、* = 42、k = 107。

現代の暗号化方法は、テキストの各バイトをASCIIに変換し、秘密鍵から得られた値でXORを計算する。XOR関数の利点は、暗号文に対して暗号鍵を使えば元の文が復元されることである。例えば、65 XOR 42 = 107、107 XOR 42 = 65。

破られない暗号のために、鍵は元の文と同じ長さでランダムなバイト列からなるべきである。ユーザーは暗号文と暗号鍵を別の場所に保存しておけば、それらどちらかのみではメッセージを復号することはできない。

不運にも、この方法は多くのユーザーに対して非実用的なものであるため、改良された方法では鍵としてパスワードを使う。パスワードがメッセージよりも短ければ、鍵はメッセージに対して繰り返し適用される。この方法で、パスワードは安全性のためには十分に長いほうが良いが、覚えるためには短いほうが良い。

今回の課題は簡単なもので、暗号鍵は3個の小文字からなる。暗号化されたASCIIコードが書かれた「p059_cipher.txt」からメッセージを解読し (元の文は一般的な英単語が書かれている)、元の文のASCIIコードの和を計算せよ。

慣れない分野なので訳がぎこちない感じがしますし、表現が間違っている部分もあるかもしれません。

XORとは

さて、ここで「XOR」が何かというと、「排他的論理和」とも言う、

どちらか片方が真でもう片方が偽の時には結果が真となり、両方とも真あるいは両方とも偽の時は偽となる演算

らしいです。

真を負の数、偽を正の数に置き換えれば、これらの積を考えることでXORの演算ができますね。

問題文に例として挙げられている「65 XOR 42 = 107」や「107 XOR 42 = 65」については、各数を二進数に変換して考えれば答えが出せます。
例えば「65 XOR 42 = 107」なら、$65_{ (10) } = 1000001_{ (2) }$、$42_{ (10) } = 101010_{ (2) }$ となるので、1の位、2の位、2²の位…と演算をしていくと $1101011_{ (2) }$ が得られ、これは $107_{ (10) }$ となります。
（pythonなら十進数に対して「^」を行うだけで演算できるのであまり関係ありませんが…）

そして、XORの利点として「暗号文に対して暗号鍵を使えば元の文が復元される」ということでしたが、これは

$$ \begin{align} \mathrm{ 元の文 ~ XOR ~ 暗号鍵 } &= \mathrm{ 暗号文 } \\ \mathrm{ 暗号文 ~ XOR ~ 暗号鍵 } &= \mathrm{ 元の文 } \end{align} $$

が成り立つ、ということです。

$\mathrm{ 元の文 ~ XOR ~ 暗号鍵 } = \mathrm{ 暗号文 }$ は作問者さんがやってくれているので、我々は $\mathrm{ 暗号文 ~ XOR ~ 暗号鍵 } = \mathrm{ 元の文 }$ をやれば良いということになります。

つまり今回の課題は「暗号鍵を探してXOR暗号を解読せよ」でした。

暗号鍵の探索

まずは暗号鍵を探しますが、暗号鍵は「3個の小文字」であると与えられているので、これを全通り試して良い感じの文字列になるものを選べば良いです。

暗号鍵 key 全てを配列 keys に入れるコードは以下です。

# itの中から重複を許してr個選んだときの組み合わせを返す関数
def combsWithRep(it, r):
    n = len(it)
    a = [i for i in it]
    
    if n == 0 or r == 0:
        yield []
    elif n == 1:
        yield a
    else:
        for i in range(n):
            a1 = [a[i]]
            a2 = a[i:]

            for j in combsWithRep(a2, r-1):
                pre = a1 + j
                if len(pre) < r:
                    for k in range(r-len(pre)):
                        pre.append(a[n-1])
                yield pre

# itを並べ替えたものを返す関数
def permutations(it):
    n = len(it)
    a = [i for i in it]
    
    if n <= 1:
        yield it
    
    for i in range(n):
        a1 = [a[i]]
        a2 = a[:i] + a[i+1:]

        for j in permutations(a2):
            if j != []:
                yield a1 + j


keys = []
for comb in combsWithRep('abcdefghijklmnopqrstuvwxyz', 3):
    # print(i) 
    for key in permutations(comb):
        if key not in keys:
            keys += [key]
            print(f"\r{key}", end='')
# print(keys)

しかし26個の中から3個を重複を許して選ぶとき、その組み合わせは $26^{3} = 17576$ 通り (=len(keys)) になるので、人力で探すのは大変そうです。
そこで復号後の文に、ある程度単純で出現頻度が高そうな単語である「the」が含まれるものだけを出力し、その中から良い感じのものを選ぶことにしました。

ファイルを読み込みます。

file = open("p059_cipher.txt")
encryption = file.read().split(',')
print(encryption)

とりあえずは始めの50個くらいを復号してみて、['t', 'h', 'e'] が含まれるものだけを出力します。

# itの中から重複を許してr個選んだときの組み合わせを返す関数
def combsWithRep(it, r):
    n = len(it)
    a = [i for i in it]
    
    if n == 0 or r == 0:
        yield []
    elif n == 1:
        yield a
    else:
        for i in range(n):
            a1 = [a[i]]
            a2 = a[i:]

            for j in combsWithRep(a2, r-1):
                pre = a1 + j
                if len(pre) < r:
                    for k in range(r-len(pre)):
                        pre.append(a[n-1])
                yield pre

# itを並べ替えたものを返す関数
def permutations(it):
    n = len(it)
    a = [i for i in it]
    
    if n <= 1:
        yield it
    
    for i in range(n):
        a1 = [a[i]]
        a2 = a[:i] + a[i+1:]

        for j in permutations(a2):
            if j != []:
                yield a1 + j

                
# ファイルの読み込み
file = open("p059_cipher.txt")
encryption = file.read().split(',')
# print(encryption)
         
keys = []
for comb in combsWithRep('abcdefghijklmnopqrstuvwxyz', 3):
    # print(i) 
    for key in permutations(comb):
        if key not in keys:
            keys += [key]
            print(f"\r{key}", end='')
            
            # 復号
            decryption = []
            i = 0
            for e in encryption[:50]:
                decryption += [int(e) ^ ord(key[i])]

                i += 1
                if i == 3:
                    i = 0
            # print(decryption)

            # 得られたASCIIコードを文字に変換
            original = []
            for d in decryption:
                original += chr(d)
            # print(original)
            
            # 「the」を含むものだけを選ぶ
            t = ['t','h','e']
            for p in range(0, len(original)-len(t)):
                if original[p] == t[0]:
                    count = 1
                    for q in range(1, len(t)):
                        if original[p+q] == t[q]:
                            count += 1
                    if count == 3:
                        print(*original)

これを実行すると、

['t', 'h', 'a']P ~ 1 t h e c q r e 0 e p { t  0 w c  | 1 d y t 0 x  d c ~ t d r d x ~ ~ 1 ~ v 1 ~ ~ t 1  w 1 U
['b', 's', 'm']F e = b s i u j ~ s + i f ` x i + { u d p '  u b + t i  o h o h d  t h e = h m = h e x ' d { ' N
['b', 'y', 'v']F o & b y r u ` e s ! r f j c i ! ` u n k ' u n b ! o i u t h e s d u o h o & h g & h o c ' n ` ' D
['c', 'h', 'm']G ~ = c h i t q ~ r 0 i g { x h 0 { t  p & d u c 0 t h d o i t h e d t i ~ = i v = i ~ x &  { & U
['c', 'i', 'm']G  = c i i t p ~ r 1 i g z x h 1 { t ~ p & e u c 1 t h e o i u h e e t i  = i w = i  x & ~ { & T
['c', 'q', 'v']G g & c q r t h e r ) r g b c h ) ` t f k & } n c ) o h } t i m s e } o i g & i o & i g c & f ` & L
['e', 'x', 'p']A n   e x t r a c t   t a k e n   f r o m   t h e   i n t r o d u c t i o n   o f   o n e   o f   E
['n', 'i', 'q']J  ! n i u y p b  1 u j z d e 1 g y ~ l + e i n 1 h e e s d u t h e h d  ! d w ! d  d + ~ g + T
['o', 'x', 'j']K n : o x n x a y ~   n k k  d   | x o w * t r o   s d t h e d o i t s e n : e f : e n  * o | * E
['l', 'v', 'p']H `   l v t { o c } . t h e e g . f { a m ) z h l . i g z r f j u j z i f `   f h   f ` e ) a f ) K
['r', 't', 'l']V b < r t h e m  c , h v g y y , z e c q 7 x t r , u y x n x h i t x u x b < x j < x b y 7 c z 7 I
['o', 'x', 'q']K n ! o x u x a b ~   u k k d d   g x o l * t i o   h d t s e d t i t h e n ! e f ! e n d * o g * E
['z', 'z', 'z']

が出力され、この中で英文らしくなっているものは暗号鍵として ['e', 'x', 'p'] を選んだときだとわかります。
（ただの数列からこのような意味のわかる文字列が出てくると結構感動します。）

答えの出力

暗号鍵がわかったので、暗号文全てについて復号してASCIIコードの和を計算します。

# ファイルの読み込み
file = open("p059_cipher.txt")
encryption = file.read().split(',')
# print(encryption)

# 復号
key = ['e', 'x', 'p']
decryption = []
i = 0
for e in encryption:
    decryption += [int(e) ^ ord(key[i])]

    i += 1
    if i == 3:
        i = 0
        
print(sum(decryption))

実行すると 129448 が出力され、これを提出すると正解でした。
ちなみに暗号化前の英文は…と思いましたが、これは自分で解けたときに見てみてください。

Problem 60 -Prime pair sets-

3, 7, 109, 673 という4つの素数は、その中から2つ選んで結合させると素数となる。例えば、7 と 109 を結合させると 7109 や 1097 となるが、これらは素数である。
この4つの素数の和は 792 であり、このような性質を持った4つの素数の集合の中で最小となる。

同様の性質を持つ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

def combs(it, r):
    n = len(it)
    a = [i for i in it]
    
    if n == 0 or r == 0 or n < r:
        yield []
    elif n == 1 or n == r:
        yield a
    else:
        for i in range(1, 2**n):
            pre = []
            for j in range(len(bin(i))):
                if i>>j & 1:
                    pre.append(a[j])
            if len(pre) == r:
                yield pre
                
def isPrimePairSet(arr):
    count = 0
    for c in combs(arr, 2):
        num1 = int(str(c[0]) + str(c[1]))
        num2 = int(str(c[1]) + str(c[0]))
        
        if sieve(num1) == True and sieve(num2) == True:
            count += 1

    if count == len(arr)*(len(arr)-1)/2:
        return True
    else:
        return False


n_max = 10000 # 考える素数の範囲
r = 5 # 求める集合中の素数の数

# 素数を配列に入れる
primes = []
for n in range(2, n_max):
    if sieve(n) == True:
        primes += [n]      
# print(primes)

array = [i for i in range(0, len(primes)-r-1)]
# print(array)

# primesの中から選ぶr個のインデックスを決める
for pos in combs(array, r):
    # print(pos)
    
    # 選んだ5個を配列に入れる
    prime_set = []
    i = 0
    for p in pos:
        prime_set += [primes[p+i]]
        i += 1
        
        if i == 5:
            i = 0
    print(f"\r{prime_set}", end='')
    
    # 題意を満たす集合かを判定
    if isPrimePairSet(prime_set) == True:
        # print(prime_set)
        print(sum(prime_set))

案の定、やたら時間がかかってしまったので考え方を変えることにしました。

結合すると素数となる集合は、その部分数合も結合すると素数となります。
つまり、上の例で r = 2 を考え、結合すると素数となるものだけについて1個素数を追加して (r = 3 と同義) 判定すれば、考えるものが少なくなります。
これを繰り返していき、集合の要素数が5個 (r = 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
           
def isPrimePairSet(arr):
    count = 0
    for c in combs(arr, 2):
        num1 = int(str(c[0]) + str(c[1]))
        num2 = int(str(c[1]) + str(c[0]))
        
        if sieve(num1) == True and sieve(num2) == True:
            count += 1

    if count == len(arr)*(len(arr)-1)/2:
        return True
    else:
        return False
    

n_max = 10000 # 考える素数の範囲

# 素数を配列に入れる
primes = []
for n in range(3, n_max):
    if sieve(n) == True:
        primes += [n]      
# print(primes)

# 集合を端から試していく
for p1 in primes:
    # print(p1)
    
    prime_set = [p1]
    for p2 in primes:
        if p2 > max(prime_set):
            prime_set += [p2]
            # print(f"\r{p2}")
            # print(prime_set)
            
            if isPrimePairSet(prime_set) == True:
                # print(prime_set)
                # print(sum(prime_set))
                
                for p3 in primes:
                    if p3 > max(prime_set):
                        prime_set += [p3]
                        # print(f"\r{p3}")
                        # print(prime_set)

                        if isPrimePairSet(prime_set) == True:
                            # print(prime_set)
                            # print(sum(prime_set))
                            
                            for p4 in primes:
                                if p4 > max(prime_set):
                                    prime_set += [p4]
                                    # print(f"\r{p4}")
                                    # print(prime_set)

                                    if isPrimePairSet(prime_set) == True:
                                        # print(prime_set)
                                        # print(sum(prime_set))
                                        
                                        for p5 in primes:
                                            if p5 > max(prime_set):
                                                prime_set += [p5]
                                                # print(f"\r{p5}")
                                                # print(f"\r{prime_set}", end='')

                                                if isPrimePairSet(prime_set) == True:
                                                    print(prime_set)
                                                    print(sum(prime_set))
                                                else:
                                                    del prime_set[4:]
                                    else:
                                        del prime_set[3:]
                        else:
                            del prime_set[2:]
            else:
                del prime_set[1:]

似たような構造が繰り返されるので、頑張ればまとめられると思いますがやめました。

実行すると [13, 5197, 5701, 6733, 8389] とその和 26033 が出力され、これを提出すると正解でした。
n_max は適当に設定しているので値を変えたら新たな候補が出てくると思います。
とりあえず 10000 で始めてみて要素5個のものが出たので、提出してみたら正解だった、という感じでした。

まあまあ時間はかかりましたが答えが出たので良しとします。

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

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