Edited at

C#:エラトステネスの篩を使って素数を求める

More than 1 year has passed since last update.

エラトステネスの篩 (エラトステネスのふるい)を使った素数を求めるC#のプログラムです。

エラトステネスのふるいについては、Wikipediaの こちらをみてください。


ノーマルバージョン

Primesが素数を求めるメソッドです。2から引数で与えた整数(maxnum)の間の素数を列挙しています。

この引数の値で、内部でふるい用の配列(sieve)の要素数を決定します。戻り値は、IEnumerable<T>です。


SieveOfEratosthenes.cs

using System;

using System.Collections.Generic;
using System.Linq;

namespace Gushwell.Etude {

class Program {
static void Main(string[] args) {
foreach(var n in Primes(200))
Console.Write("{0,3} ",n);
Console.WriteLine();
Console.ReadLine();
}

static IEnumerable<int> Primes(int maxnum) {
int[] sieve = Enumerable.Range(0, maxnum + 1).ToArray();
sieve[1] = 0; // 0 : 素数ではない
int squareroot = (int)Math.Sqrt(maxnum);
for (int i = 2; i <= squareroot; i++) {
if (sieve[i] <= 0)
continue;
for (int n = i * 2; n <= maxnum; n += i)
sieve[n] = 0;
}
return sieve.Where(n => n > 0);
}
}
}



省メモリバージョン

もう少し、メモリ効率を良くしたバージョンも載せておきます。intの配列ではなくて、boolの配列にしています。


PrimesSavingMemory.cs

 static IEnumerable<int> Primes(int maxnum) {

// var sieve = new BitArray(maxnum + 1, true);
bool[] sieve = Enumerable.Repeat(true, maxnum + 1).ToArray();
int squareroot = (int)Math.Sqrt(maxnum);
for (int i = 2; i <= squareroot; i++) {
if (sieve[i] == false)
continue;
for (int n = i * 2; n <= maxnum; n += i)
sieve[n] = false;
}
for (int i = 2; i <= maxnum; i++)
if (sieve[i] == true)
yield return i;
}

コメントのようにBitArray使えば、もっとメモリ効率は良くなります。

GitHubでこのソースコードを公開しています。


列挙開始を早くする

ところで、ここまで書いて思ったんですが、このメソッドは、列挙が始まるのが遅すぎるのが難点ですね。全部の素数が求まらないと列挙が始まりません。

最初のfor文の中で、yield return を書ければ、列挙を開始するのがもっと早くなるはずです。

いくつまでの素数が必要かわからない場合は、列挙開始が速くなれば、トータルの処理時間も短くなる可能性もあります。

で、書いたのが次のコードです。


PrimesStartDash.cs

 static IEnumerable<int> Primes(int maxnum) {

var sieve = new BitArray(maxnum + 1, true);
int squareroot = (int)Math.Sqrt(maxnum);
for (int i = 2; i <= squareroot; i++) {
if (sieve[i] == false)
continue;
for (int n = i * 2; n <= maxnum; n += i)
sieve[n] = false;
if (sieve[i])
yield return i; // ここで列挙してしまう
}
for (int i = squareroot + 1; i <= maxnum; i++)
if (sieve[i] == true)
yield return i;
}

うーーん、でも、Math.Sqrt(maxnum) の値までしか列挙開始を早められないから、ほとんどメリット無いかな。

一瞬、良いアイデアと思ったんですが...

もっと突き詰めて考えれば、もう少し実用的なコードになるような気もするのですが、とりあえずこれで止めておきます。


追記 ("列挙開始を早くする"の改良バージョン)

もう少し、突き詰めて考えてみました。

良く考えてみると、


  • 整数2で篩を掛けた後の配列の状態を見ると、2,3,5,7 までの素数が求まっています。

  • その後3で篩をかけると、11,13の素数が求まります。

  • さらに4で篩をかけると、17, 19, 23の素数が求まります。(実際には、4で篩を掛けても配列の状態は変わらないのだが、説明の都合上、こう書いておく)

  • さらに5で篩をかけると、29, 31の素数が求まっています。

つまり、整数nで篩を掛け終わった時には、(n+1) * (n+1) - 1 以下の素数が求まっているということです。実質的に、n+1で篩を掛けることができるのは、(n+1)*(n+1)からであり、これより小さい値は、n以下の値ですでに篩を掛けてあるからです。

そこで書き換えたのが以下のコードとなります。このコードでは、sieveのtrueとfalseの意味を反対にしています。


PrimesStartDash2.cs

 static IEnumerable<int> Primes(int maxnum) {

yield return 2;
yield return 3;
var sieve = new BitArray(maxnum + 1);
int squareroot = (int)Math.Sqrt(maxnum);
for (int i = 2; i <= squareroot; i++) {
if (sieve[i] == false) {
for (int n = i * 2; n <= maxnum; n += i)
sieve[n] = true;
}
for (int n = i * i + 1; n <= maxnum && n < (i + 1) * (i + 1); n++) {
if (!sieve[n])
yield return n;
}
}
}


GitHubでこのコードを公開しています。

2つ目のfor文がなくなり、yield returnは、最初のfor文の中に全て移動できた関係で、最初のバージョンから比べると、スムーズに列挙が進むようになりました。

それと、「列挙開始を早くする」で示したコードよりも、コードもそれほど複雑にならずに実現できました。

とは言え、引数に最大値を指定しないといけないという欠点は残ります。「1000までの素数を求めたい」という要求には答えられますが、「小さい順に1000個の素数を求めたい」という要求には、うまく対応できません。

yone64さんがコメントしてくれたような工夫が必要になりそうです。これについては、時間ができた時に、再度考えてみたいと思います。


この記事は、Gushwell's C# Programming Pageで公開したものを加筆・修正したものです。