指定した範囲内における自然数の試行回数を割り出し、試行回数の最大値を持つ自然数を算出するプログラムを作る 【コラッツ予想】

(使用する言語はpythonのみ　専門的な数学の知識は一切使わず、ライブラリも極力使用しません)

【コラッツ予想とは】

適当な自然数を決める。
自然数が偶数の場合は2で割り、奇数の場合は3倍して1を足す　それをひたすら繰り返していくとどんな自然数も1になるというのがコラッツ予想

例えば"5"で計算するとどうなるか
5*3 + 1 = 16
16/2 = 8
8/2 = 4
4/2 = 2
2/2 = 1

結果は5回試行を繰り返したのち1に収束した。
このように、どれだけ大きい自然数でも必ず1に収束するらしいが、その証明方法は未だに見つかっていない。

そこで今回は、コラッツ予想の証明は一旦諦め、適当な数字の試行回数を算出するプログラムを作ることにした。

次のようなシンプルなコードで実行することができる。

#Collatz予想

ntnum = int(input('自然数を入力してください >>'))
count = 0

while ntnum != 1:
    if (ntnum % 2 == 0):
        ntnum = ntnum/2
    else:
        ntnum = 3*ntnum + 1
    count += 1
    print('{}\n'.format(int(ntnum)))

print('試行回数は{}回'.format(count))

試しに37と入力してみると試行回数は21となり、正しい計算が行われていることがわかった。　

コラッツ予想の基本的なプログラムが実行できたのでここで本題に入る。
テーマは、指定した自然数の範囲内でもっとも試行回数の多い自然数を算出するというもの、とはいえ試行回数の最大値を算出するだけではすぐに終わってしまうので、異なる計算方法をいくつか試しどの方法が一番早く計算が終わるか比較する。

まずは先ほど書いたコードを直接利用したやり方で計算してみる　
計算方法①

#Collatz予想
import time

ntnum = int(input('自然数を入力してください >>'))

start = time.perf_counter() # 計算時間の測定開始

max_list = []
for n in range(1,ntnum):
    count2 = 0
    n2 = n
    while n != 1:
        if (n % 2 == 0):
            n = n/2
        else:
            n = 3*n + 1
        count2 += 1
    max_list.append(count2)

end = time.perf_counter() # 計算時間の測定終了

print('1~{}までの範囲内で最も試行回数の多い自然数は'.format(ntnum))
print('{}の{}回です。'.format((max_list.index(max(max_list))+1),max(max_list)))
print('計算にかかった時間は{}秒です'.format(end - start))

(ちなみにntnumというのはNatural numberの略　略語はかなり適当なので大目に見てください)
試しに1〜100000の範囲内で計算してみる。
結果は...

1~100000までの範囲内で最も試行回数の多い自然数は
77031の350回です。
計算にかかった時間は2.5451462310000004秒です

小数点第4位くらいでカットすると、およそ2.5451秒となった。　
さらに早く計算できそうなプログラムを作る。

奇数だけで計算してその中で最大値の試行回数を持つ奇数から、本来の自然数の範囲内の最大値を求めるやり方を考えたので、下記のようなコードで試してみる。
計算方法②

#Collatz予想
import time

max_list = []
ntnum = int(input('1~n までの自然数で、試行回数の最大値を算出します、nに代入する自然数を入力してください >>'))

start = time.perf_counter() # 計算時間の測定開始
m = 0

for n in range(1,ntnum+1,2): # 偶数を飛ばして計算する
    count2 = 0
    n2 = n
    while n != 1:
        if (n % 2 == 0):
            n = n/2
        else:
            n = 3*n + 1
        count2 += 1
    max_list.append(count2)

maxnum1 = (2*(max_list.index(max(max_list))+1)-1) # [2*(max_listの要素数) + 1]で最も試行回数の多い奇数を求めることができる
maxtimes1 = max(max_list)

if (maxnum1*2 <= n): # (24行目)最も試行回数の多い偶数が存在する場合のif文
    while n >= m:
        maxnum1 = maxnum1*2
        maxtimes1 += 1

end = time.perf_counter() # 計算時間の測定終了

print('1~{}までの範囲内で最も試行回数の多い自然数は'.format(ntnum))
print('{}の{}回です。'.format((maxnum1),maxtimes1))
print('計算にかかった時間は{}秒です'.format(end - start))

内容について軽く説明すると、10行目では、2の冪乗以外の自然数は計算途中で必ず奇数に当たるので、最初から奇数だけ計算することによって計算過程をショートカットするというやり方である。
最も試行回数の多い奇数を求めた時、その奇数の倍数(2倍(偶数),4倍(偶数),8倍(偶数)...)も指定した範囲内の自然数になった場合はそれが最も試行回数の多い自然数となるので24行目はそのためのif文である。

それではプログラムを実行してみる　結果は...　

1~100000までの範囲内で最も試行回数の多い自然数は
77031の350回です。
計算にかかった時間は1.4652075470000003秒です

計算結果は①のプログラムと一致、タイムは1.4652秒となり、①よりも1秒と少し速くなった。

もう少し計算過程を減らせないか考えてみる。
②と同様、奇数だけを計算するというやり方は残した状態で

計算方法③

#Collatz予想
import time

max_list = [0] # 3から順に奇数の試行回数を入れていく
ntnum = int(input('1~n までの自然数で、試行回数の最大値を算出します、nに代入する自然数を入力してください >>'))

start = time.perf_counter() # 計算時間の測定開始
m = 0

if (ntnum >= 3):
    for n in range(3,ntnum+1,2):
        n2 = n
        n = (3*n + 1)/2
        count2 = 3  # 上記の(3*n + 1)で1回、1/2で2回、下記のwhile文で3回
        while n != 2:
            if (n % 2 == 0):
                n = n/2
            else:
                n = 3*n + 1
            count2 += 1
        max_list.append(count2)

    maxnum1 = (2*(max_list.index(max(max_list))+1)-1) # [2*(max_listの要素数) + 1]で最も試行回数の多い奇数を求めることができる
    maxtimes1 = max(max_list)

    if (maxnum1*2 <= n): # 最も試行回数の多い偶数が存在する場合のif文
        while n >= m:
            maxnum1 = maxnum1*2
            maxtimes1 += 1

    print('1~{}までの範囲内で最も試行回数の多い自然数は'.format(ntnum)) # ntnumが3以上の場合
    print('{}の{}回です。'.format((maxnum1),maxtimes1))
else:
    print('1~{}までの範囲内で最も試行回数の多い自然数は'.format(ntnum)) # ntnumが3未満の場合
    print('{}の{}回です。'.format((ntnum),ntnum - 1))

end = time.perf_counter() # 計算時間の測定終了

print('計算にかかった時間は{}秒です'.format(end - start))

重要な点だけ解説　
まず10行目のif文は、1と2を飛ばし3から計算が始まるようになっている。すべての自然数は1に収束する前にかならず2を通るのでその過程を省略しただけである。
次に14行目のcount2 = 3の部分だが、最初から試行回数が3回分カウントされている状態から始まっている。これは13行目の計算(2回分)、そして15行目のwhile文が、1ではなく2に収束した状態で計算が終了するため2から1に収束するのに行われる計算(1回分)を足し合わせた合計3回分のことである。　
他のコードは基本的に②と同様である。

それでは①、②と同じ範囲内に自然数を設定し、プログラムを実行する。　

1~100000までの範囲内で最も試行回数の多い自然数は
77031の350回です。
計算にかかった時間は1.3793225369999993秒です

計算結果は一致、そしてタイムは1.3793秒と、②よりもさらに速く計算を行うことができた。　

さらに計算過程を減らしてみる。
今度は、過去に計算した自然数の、計算途中で通った自然数全てとその試行回数を記録しておき、過去に出た自然数が再び出た場合にその自然数の計算を飛ばして次の計算に移行させるというやり方である。
例えば、7で計算すると、7→22→11→34→17→52→26→13→40→20→10→5→16→8→4→2→1となり、次に8、その次に9を計算することになるので9→28→14→7→...といった流れで計算が行われる。見ての通り9は計算途中で7を通ることになるので、ここで前回計算した7の計算過程と試行回数を全て記録しておけば、7以降の計算の過程を全て省略して、9の計算を終わらせることができる。
7と9以外の数字も同様である　7を計算したときに34が既に出てきているので、次に34が出てきた場合もそれ以降の計算を省略することができる。　
それではコードを書いてみる。

④

#Collatz予想
import time

max_dict = {'2':1}         #1〜nのすべての自然数の試行回数を順番に格納するもの
yobinum = []

ntnum2 = int(input('3~n までの自然数で試行を行います、nに代入する自然数を入力してください >>'))

start = time.perf_counter() # 計算時間の測定開始

for n2 in range(3,ntnum2+1):
    count2 = 0          #試行回数のカウント数           
    n = n2              
    if str(n2) in max_dict.keys():         #nがmax_dictリストにすでにあったら終了
            continue
    while n != 1:       #nが1にならない限り試行を続ける

        yobinum.append(n)
        count2 += 1
        if (n % 2 == 0):
            n = n/2
            if str(int(n)) in max_dict.keys():
                for a1 in range(0,len(yobinum)):
                    max_dict[str(int(yobinum[a1]))] = (count2  - a1 + max_dict[str(int(n))])         #maxディクショナリの数字の試行回数+今までの試行回数+1、試行回数をリストに格納  
                yobinum.clear()  
                break

        else:
            n = 3*n + 1
            if str(int(n)) in max_dict.keys():
                for a2 in range(0,len(yobinum)):
                    max_dict[str(int(yobinum[a2]))] = (count2  - a2 + max_dict[str(int(n))])        #maxディクショナリの数字の試行回数+今までの試行回数+1、試行回数をリストに格納
                yobinum.clear() 
                break

maxnum1 = [k for k, v in max_dict.items() if v == max(max_dict.values())]

end = time.perf_counter() # 計算時間の測定終了

print('1~{}までの範囲内で試行回数が最も多かったのは'.format(ntnum2))
print('{}の{}回です。'.format((maxnum1),max(max_dict.values())))
print('計算にかかった時間は{}秒です'.format(end - start))

重要な部分だけ解説　
②、③と違い、特定の試行回数が格納されているリストの番号を、その試行回数の自然数であると判定することができないので、計算が終わった自然数とその試行回数はディクショナリに入れて2つ1組で管理する。
yobinumというのは、計算途中で通った全ての自然数を一時的に入れておくリストであり、35行目と43行目のif文によって計算を省略することができるという仕組みになっている。

それではプログラムを実行してみる　結果は...

1~100000までの範囲内で試行回数が最も多かったのは
['77031']の350回です。
計算にかかった時間は774.5985557089999秒です

計算結果は一致したが、タイムは、約774秒と①〜③とは桁違いに時間がかかってしまった。

1〜10000に範囲を狭めてもう一度実行してみる。

1~10000までの範囲内で試行回数が最も多かったのは
['6171']の261回です。
計算にかかった時間は7.588531691秒です

やはり計算に時間がかかることがわかった。

【結論・考察】
今回試した中では、1〜10万の範囲内で計算を行った場合、③の計算方法が1.3793秒と最も早く試行回数の最大値を算出することができた。　
④の方法による計算は失敗したが、コードの書き方の改善と、より性能の良いコンピュータの使用によって結果が変わる可能性があるので、今後の課題にしていきたい(が、多分2度とやらない)。

【反省点】　
②と③において、指定した範囲内において最も試行回数の多い自然数が偶数だった場合、
((最も試行回数の多い奇数)×(2の冪乗))の試行回数　＜　((2番目以降に試行回数の多い奇数)×(2の冪乗))の試行回数
のパターンを想定していなかったため、答えがこのようになった場合でも　((最も試行回数の多い奇数)×(2の冪乗))の試行回数　が答えとして算出されてしまうというミスがあった。

【参考文献】　
Pythonのtimeモジュールまとめ！処理時間を計測する方法