エラトステネスの篩(決定的素数判定アルゴリズム)｜Excel VBAで学ぶ数学とアルゴリズム

今回は決定的素数判定アルゴリズムの「エラトステネスの篩ふるい」について解説していきます。
素数判定アルゴリズムは数多く存在しますがその中でも非常に単純で理解しやすいアルゴリズムです。

それもそのはず、このアルゴリズムが誕生したのは今よりはるか昔の紀元前まで遡ります。アルゴリズム名のとおり、古代ギリシアの科学者であるエラストテネスによって考案されたとされています。

本ページではこの「エラトステネスの篩」のアルゴリズムがどういった仕組みなのかを解説し、最終的にExcel VBAのコードとして実装していきます。

1. エラトステネスの篩
2. サンプルコード
3. まとめ

エラトステネスの篩

素数判定アルゴリズムの多くは入力された数値が素数かそうでないかを判定するものですが、エラストテネスの篩は”1から入力された数値までに存在する素数”をすべて取得することのできるアルゴリズムとなっており、他のものと比べると少しだけニュアンスが違います。

エラトステネスの篩のアルゴリズムは非常に単純で下記の処理を行います。

① 要素数が \(N\) 個のBoolean型の配列を用意する
※このとき0番目と1番目の要素は[False]、2番目以降の要素はすべて[True]とする
② 配列の先頭から順に要素を確認し[True]の場合、下記の処理を行う
　・その要素は素数であると判定、要素は[True]のままにする
　・その要素の倍数にあたる要素は合成数になるためすべて[False]にする
③ 確認した要素数が \( \sqrt{N} \) 回を超えるまで②を続ける

エラストテネスの篩はその名の通りで、すべての自然数(実際には入力された範囲内での自然数)をいったんすべて素数ということにしておき合成数とわかったものから順に除外していく、つまりは素数でない数字をどんどん”篩ふるい“にかけていくアルゴリズムです。

Animation Sieb des Eratosthenes.gif ©SKopp (Licensed under CC BY 3.0)

処理の流れは上のgif画像を見ればある程度イメージできると思います。
ここではさらに具体的な例で処理の流れを見ていきましょう。

―「1～30」の範囲内にある素数を判定 ―

Array(30)のBoolean型の配列を用意
0番目と1番目の要素は[False]、2番目以降の要素はすべて[True]とする
　
① 0番目の要素を確認 ⇒ [False]なので素数ではない
② 1番目の要素を確認 ⇒ [False]なので素数ではない
③ 2番目の要素を確認 ⇒ [True]なので素数、”2の倍数”番目の要素はすべて[False]にする
④ 3番目の要素を確認 ⇒ [True]なので素数、”3の倍数”番目の要素はすべて[False]にする
⑤ 4番目の要素を確認 ⇒ [False]なので素数ではない
⑥ 5番目の要素を確認 ⇒ [True]なので素数、”5の倍数”番目の要素はすべて[False]にする
⑦ 6番目の要素を確認 ⇒ [False]なので素数ではない　※ここで終了してOK
　　　　　　　　　　　：
㉚ 29番目の要素を確認 ⇒ [True]なので素数、”29の倍数”番目の要素はすべて[False]にする
㉛ 30番目の要素を確認 ⇒ [False]なので素数ではない

このように小さい数字から順に見ていき、その数字の倍数となる数字をどんどん除外していくことで素数のみが残るという考えのアルゴリズムです。ここではプログラミング的に[True][False]で判定していますが、実際は紙に数字を羅列して合成数となる数字にどんどん「×(バツ)」を付けていくという手法のため非常にアナログ的な手法です。

上記の場合、実際は30回すべてを確認する必要はなく、\( \sqrt{30} = 5.477 \) を超える回数、すなわち6番目の要素を確認した時点で終了したとしてもすべての要素を確認したことになります。
　

判定が√Nまでで十分な理由

このアルゴリズムは任意の範囲の自然数の中から合成数を除外していくアルゴリズムです。
約数の「N以下の合成数は必ず 2以上√N以下の約数を持つ」という性質から、√Nを超えるまでの数値の倍数を除外していけば「√N～N」までの範囲の合成数もすべて除外することが出来ます。

たとえば先ほどの「1～30」を例に考えると、\( \sqrt{30} = 5.477 \) を下回る数字である「2,3,(4),5」それぞれの倍数に色を付けると下図のように、1～30までの範囲ですべての合成数に色を付けることができます。(※2,3,5は除く)

つまり、このアルゴリズムは全ての数値を確認する必要はなく「√N」まで確認できればそこで終了したとしても、それ以降の合成数もまとめて確認したことになるという訳です。(※実際には \( \sqrt{30} = 5.477 \) のような場合は繰り上げた数値「6」まで確認する必要があります)
　

ちなみに「N以下の合成数は必ず 2以上√N以下の約数を持つ」という事実は、背理法を使って証明することができます。背理法とは簡単にいえば「ある事実が間違っていると仮定して、その仮定から矛盾を見つけ出すことにより、ある事実は正しい」と証明するような手法です。

今回の場合は前提を「Nが合成数の場合、最小の約数Aは√Nを超える」と仮定します。この前提で考えて何かしらの矛盾があれば、この前提は間違っている、つまりは「N以下の合成数は必ず 2以上√N以下の約数を持つ」という事実が正しいということを証明できます。

―「N が合成数の場合、2以上 √N 以下の約数が存在する」の証明 ―

①「Nが合成数の場合、最小の約数Aは√Nを超える」と仮定する
② 約数の性質から「A × B = N」となる自然数Bが存在する (BもNの約数)
③ この数式を並び替えると「B = N/A」となる
④ このとき仮定が「A > √N」のため、③の式にあてはめると「B < √N」が成り立つ
⑤「B < √N」となる約数Bが存在することは仮定の「最小の約数Aは√Nを超える」に矛盾する
⑥ よってこの仮定は不成立のため「N以下の合成数は必ず 2以上√N以下の約数を持つ」は正しい

サンプルコード

エラトステネスの篩をVBAで書くと下記の通りです。
「1～N」までの範囲で素数となるすべての数字が入った配列が返されます。
処理時間としては「N = 1000000」くらいまでであれば一瞬で取得できますが、
「N = 100000000」あたりからは分単位の時間がかかるようになります。

Sub main()

    Dim N As Long
    Dim iPrimes() As Long
    N = 1000

    '2～Nまでの範囲内の素数を返す
    iPrimes = SieveOfEratosthenes(N)
    
End Sub

'------------------------------------------------------------------
'   エラトステネスの篩
'       N       ：素数を判定の最大値
'       戻り値  ：2～Nまで素数の入った配列
'------------------------------------------------------------------
Function SieveOfEratosthenes(N As Long) As Long()

    Dim iNums() As Long
    Dim iPrimes() As Long
    Dim flgPrime() As Boolean
    Dim i As Long
    
    ReDim iNums(N)
    ReDim flgPrime(N)
    ReDim iPrimes(0)
    
    
    '要素がすべて[True]の素数判定配列を用意(※[0][1]の要素は使わない)
    For i = 2 To N
        flgPrime(i) = True
    Next i
    
    '素数判定(各素数の倍数を小さいものから順に[False]にしていく)
    i = 2
    Do While i * i <= N
    
        '[True]のうち最も小さい数は常に素数
        If flgPrime(i) = True Then
        
            'iの倍数をすべて[False]にする
            x = 2 * i
            Do While x <= N
                flgPrime(x) = False
                x = x + i
            Loop
            
        End If
        i = i + 1
    Loop
    
    '素数のみを配列に格納
    For i = 2 To N
        If flgPrime(i) = True Then
            
            '配列を動的に拡張しながら素数を詰め込んでいく
            iPrimes(UBound(iPrimes)) = i
            ReDim Preserve iPrimes(UBound(iPrimes) + 1)
            
        End If
    Next i
    
    '1番後ろの要素は不要なので削除
    ReDim Preserve iPrimes(UBound(iPrimes) - 1)
    

    '戻り値
    SieveOfEratosthenes = iPrimes
    

End Function

Sub main()

Dim N As Long

Dim iPrimes() As Long

N = 1000

'2～Nまでの範囲内の素数を返す

iPrimes = SieveOfEratosthenes(N)

End Sub

'------------------------------------------------------------------

' エラトステネスの篩

' N ：素数を判定の最大値

' 戻り値：2～Nまで素数の入った配列

'------------------------------------------------------------------

Function SieveOfEratosthenes(N As Long) As Long()

Dim iNums() As Long

Dim iPrimes() As Long

Dim flgPrime() As Boolean

Dim i As Long

ReDim iNums(N)

ReDim flgPrime(N)

ReDim iPrimes(0)

'要素がすべて[True]の素数判定配列を用意(※[0][1]の要素は使わない)

For i = 2 To N

flgPrime(i) = True

Next i

'素数判定(各素数の倍数を小さいものから順に[False]にしていく)

i = 2

Do While i * i <= N

'[True]のうち最も小さい数は常に素数

If flgPrime(i) = True Then

'iの倍数をすべて[False]にする

x = 2 * i

Do While x <= N

flgPrime(x) = False

x = x + i

Loop

End If

i = i + 1

Loop

'素数のみを配列に格納

For i = 2 To N

If flgPrime(i) = True Then

'配列を動的に拡張しながら素数を詰め込んでいく

iPrimes(UBound(iPrimes)) = i

ReDim Preserve iPrimes(UBound(iPrimes) + 1)

End If

Next i

'1番後ろの要素は不要なので削除

ReDim Preserve iPrimes(UBound(iPrimes) - 1)

'戻り値

SieveOfEratosthenes = iPrimes

End Function

まとめ

今回は素数判定アルゴリズムの「エラトステネスの篩」の解説でした。
エラトステネスの篩のアルゴリズムは下記の通りです。

エラトステネスの篩は指定した範囲内にあるすべての素数の入った配列を返すという、素数判定アルゴリズムの中では少し変わった性質を持つアルゴリズムです。ただ、素数の配列さえ取得できてしまえば、入力された数値がその配列内にあるかを確認するだけで他の素数判定アルゴリズムと同じようなことができます。

特にエラトステネスの篩は決定的素数判定アルゴリズムのため、そこまで大きな値の判定をしないのであれば、単純な処理で実装も簡単なのでオススメのアルゴリズムです。判定する数値が巨大であれば「フェルマーテスト」のような確率的アルゴリズムを使ったほうが速いですが、1000000くらいの桁数であればこのアルゴリズムでも十分許容できる速さで判定することができます。

素数判定アルゴリズム数多くありますがどれが優れている劣っているというのはあまりなく、判定したい数値の大きさや処理速度、実装のしやすさや可読性など、どれを優先したいかを考えて使い分けるというのが基本です。