素数判定プログラムを改良｜Pythonで数学を学ぼう！第5回

約数をすべて表示する

前回の素数判定プログラム（prime1）は「素数ではありません」「素数です」だけの判定をする7行のコードでした。

今回はこれをもとにいくつか改良してみます。

プログラム：prime2

>>>     n = int(input('素数判定したい2以上の自然数nを入れてね n='))   # 入力されたnを整数に変換
>>>     p = 0    # 約数の個数カウンター
>>>     for k in range(1,n+1):   # k=1,...,n
>>>          if n % k == 0:     # n÷kの余りが0ならば、（kはnの約数ならば）
>>>                           print(f'{n} は {k} を約数にもつ')    # 約数kを表示
>>>                           p = p + 1        # 約数の個数カウンターpを+1
>>>     if p > 2:               # for文を抜け出した後約数の個数で条件分岐 2個よりも大きい場合
>>>          print(f'{n} は約数を{p}個もつ合成数で素数ではありません')
>>>     else:                    # そうでない場合（p=2）
>>>          print(f'{n} は約数が2個だから素数！ (k={k})')

素数判定は、与えられたnを1から順にnまで割って余りを計算することで行います（前回prime1は2からn-1まで）。

余りが0すなわちの場合にはその数（約数）を表示するようにします。

>>> print(f'{n} は {k} を約数にもつ') # 約数kを表示

前回説明したf文字列(f-string)表示方法です。

ついでに約数の個数もカウント・表示するようにします。カウンター用変数をpとします。

>>> p = 0 # 約数の個数カウンター

割る数kが約数である場合にカウントアップするようにします。

>>> p = p + 1 # 約数の個数カウンターpを+1

さて、割る数kがnまで渡り終わるとfor文を抜け出し、

>>> if p > 2:

が実行されます。

約数の個数pが2より大きい場合（if p > 2 :）、すなわち3以上の場合、nは「素数ではない」の判定を行い、

>>> print(f'{n} は約数を{p}個もつ合成数で素数ではありません')

そうでない場合（else:）すなわち、nが2である場合には「素数である」判定を行います。

>>> print(f'{n} は約数が2個だから素数！ (k={k})')

prime2をGoogle Colaboratory（グーグルコラボラトリー）に書いて実行してみると、次のように表示されます。

計算時間を表示する

前回「56547511の素数判定にどれほどの時間がかかるか、試してみてください」と最後に書きました。

そこで、素数判定にかかる計算時間を計測・表示するようにコードを改良します。

プログラム：prime3

>>>     import time                               # timeモジュールのインポート
>>>     n = int(input('素数判定したい2以上の自然数nを入れてね n='))   # 入力されたnを整数に変換
>>>     start = time.time()                       # 計測開始UNIX時間
>>>     p = 0                                     # 約数の個数カウンター
>>>     for k in range(1,n+1):                    # k=1,...,n
>>>          if n % k == 0:                        # n÷kの余りが0ならば
>>>                           print(f'{n} は {k} を約数にもつ')   # 約数kを表示
>>>                           p = p + 1                         # 約数の個数カウンターpを+1
>>>     if p > 2:                                 # 約数の個数が3以上の場合
>>>          print(f'{n} は約数を{p}個もつ合成数で素数ではない')    # 素数ではない
>>>          t = time.time() - start               # 計測終了UNIX時間 -> 計算時間算出
>>>          print(f'計算にかかった時間{t}秒')        # 計算時間表示
>>>     else:                                     # そうでない場合（p=2）
>>>          print(f'{n} は素数！')                 # 素数
>>>          t = time.time() - start               # 計測終了UNIX時間 -> 計算時間算出
>>>          print(f'計算にかかった時間{t}秒')        # 計算時間表示

1行目でtimeモジュールをインポートします。

>>> import time # timeモジュールのインポート

これで時間を扱うことができるようになります。

>>> start = time.time() # 計測開始UNIX時間

このコードが実行された時点でのUNIX時間（エポック秒）を取得します。

次のコードを実行してみましょう。

>>> import time
>>> print(time.time())

1611654943.353461

これがUNIX時間（エポック秒）で、単位は秒です。

nの入力後直後のUNIX時間をstartとしてマークします。

>>> start = time.time() # 計測開始UNIX時間

2つの判定完了後それぞれで直後のUNIX時間からstartを引いて計測時間

>>> t = time.time() - start # 計測終了UNIX時間 -> 計算時間算出
>>> print(f'計算にかかった時間{t}秒') # 計算時間表示

prime3をGoogle Colaboratory（グーグルコラボラトリー）に書いて実行してみると次のように表示されます。

8桁56547511の判定にかかった計算時間は6.4秒でした。

10桁1234567890を試してみた結果が次です。195秒でした。

メルセンヌ素数とGIMPS

www.mersenne.org

素数探査はこの瞬間も世界中で行われています。

GIMPS（Great Internet Mersenne Prime Search）が活躍しています。メルセンヌ素数の発見を目的に1996年に発足したプロジェクトです。

分散コンピューティングにより、インターネットで繋がれた世界中のPCを使って素数判定に必要な計算が行われます。

2018年12月7日に判定完了した素数は、なんと2486万2048桁です。

メルセンヌ素数とは2のn乗-1（2ⁿ-1）の形をした素数です。nが3の場合2^3-1=7、5の場合2^5-1=31はメルセンヌ素数です。

メルセンヌ素数の判定法が研究されてきたこととコンピューターの発達によってGIMPSが誕生しました。

メルセンヌ素数を判定する

プログラム：prime4

>>>     import time                               # timeモジュールのインポート
>>>     import math                               # mathモジュールをインポート
>>>     n = eval(input('素数判定したい2以上の自然数nを入れてね n='))   # 入力されたnを計算してnに代入
>>>     start = time.time()                       # 計測開始UNIX時間
>>>     p = 0                                     # 約数の個数カウンター
>>>     for k in range(1,n+1):                    # k=1,...,n
>>>          if n % k == 0:                        # n÷kの余りが0ならば
>>>                           print(f'{n} は {k} を約数にもつ')   # 約数kを表示
>>>                           p = p + 1                         # 約数の個数カウンターpを+1
>>>     if p > 2:                                 # 約数の個数が3以上の場合
>>>          print(f'{n} は約数を{p}個もつ合成数で素数ではない')    # 素数ではない
>>>          t = time.time() - start               # 計測終了UNIX時間 -> 計算時間算出
>>>          print(f'計算にかかった時間{t}秒')        # 計算時間表示
>>>     else:                                     # そうでない場合（p=2）
>>>          print(f'{n} は素数！')                 # 素数
>>>          t = time.time() - start               # 計測終了UNIX時間 -> 計算時間算出
>>>          print(f'計算にかかった時間{t}秒')        # 計算時間表示

そこでメルセンヌ素数が判定できるようにプログラムを改良してみます。

最初のinput文

>>> n = int(input('素数判定したい2以上の自然数nを入れてね n=')) # 入力されたnを整数に変換

は入力された文字列nを整数型に変換するコードです。ここに2**7-1（2^7-1のこと）と入力すると、Pythonに処理できないと叱られます。

2**5-1（2^5-1のこと）を処理できるようにするには、いろいろなコードを考えることができますが、ここは一番簡単な変更を考えてみます。

プログラムprime4の変更点は1箇所だけ。

>>> n = eval(input('素数判定したい2以上の自然数nを入れてね n=')) # 入力されたnを計算してnに代入

intをevalに変えるだけです。eval(文字列)は文字列を計算して数値を返してくれる魔法のコードです。これで2**7-1を入力できます。

1772年、オイラーは2^31-1= 2147483647の素数判定に成功しました。8番目のメルセンヌ素数の発見です。

同じ問題をprime4で計算してみましょう。2**31-1を入力します。6分弱かかって計算できました。

魔法のevalを味方につけたので、nには2**7-1以外の数式も入力できます。階乗（5!=1・2・3・4・5）を入力できるようにしてみます。

>>> import math # mathモジュールをインポート

を最初に書けばOKです。math.factorial()で階乗が計算できます。

>>> import math
>>> factorial(5)

120

では、7!-1を判定してみましょう。「math.factorial(7)-1」と入力します。

結果は素数でした。

いかがでしたでしょうか。今回は素数判定プログラムを改良しました。みなさんも独自の改良をしながら数学を学んでみてください。