最近、素数を生成するために Atkin のふるい ( http://en.wikipedia.org/wiki/Sieve_of_atkin ) を使用する C++ 素数ジェネレーターに取り組んでいます。私の目標は、任意の 32 ビット数を生成できるようにすることです。主にプロジェクトのオイラー問題に使用します。ほとんどは夏のプロジェクトです。
プログラムはビットボードを使用して素数を格納します。つまり、一連の 1 と 0 で、たとえば 11 番目のビットが 1、12 番目が 0、13 番目が 1 などになります。メモリを効率的に使用するために、これは実際にはおよび char の配列で、各 char には 8 ビットが含まれます。フラグとビットごとの演算子を使用して、ビットを設定および取得します。アルゴリズムの要旨は単純です。数値が「素数」と見なされるかどうかを定義するために、私が理解しているふりをしていないいくつかの方程式を使用して最初のパスを実行します。これでほとんどの場合正しい答えが得られますが、いくつかの非素数が素数としてマークされます。したがって、リストを反復処理するときは、見つけたばかりの素数のすべての倍数を「素数ではない」に設定します。これには、素数が大きくなるほど必要なプロセッサ時間が少なくて済むという便利な利点があります。
90% 完成しましたが、問題が 1 つあります。いくつかの素数が欠落しています。
ビットボードを調べたところ、最初のパスでこれらの素数が省略されていることがわかりました。これは基本的に、いくつかの方程式の解ごとに数値を切り替えます (ウィキペディアのエントリを参照)。私はこのコードのチャンクを何度も何度も調べてきました。ウィキペディアの記事に示されている範囲を広げてみましたが、これは効率的ではありませんが、何らかの形で省略したいくつかの数字にヒットする可能性があると考えました. 何も機能していません。これらの数値は、単純に素数ではないと評価されます。私のテストのほとんどは、128 未満のすべての素数で行われました。この範囲のうち、省略された素数は次のとおりです。
23と59。
より高い範囲では、より多くのものが欠けていることに疑いの余地はありません(それらすべてを数えたくないだけです). これらが欠落している理由はわかりませんが、欠落しています。これらの 2 つの素数について特別なことはありますか? 私は二重、三重にチェックし、間違いを見つけて修正しましたが、それでも私が見逃しているのはおそらくばかげたものです。
とにかく、ここに私のコードがあります:
#include <iostream>
#include <limits.h>
#include <math.h>
using namespace std;
const unsigned short DWORD_BITS = 8;
unsigned char flag(const unsigned char);
void printBinary(unsigned char);
class PrimeGen
{
public:
unsigned char* sieve;
unsigned sievelen;
unsigned limit;
unsigned bookmark;
PrimeGen(const unsigned);
void firstPass();
unsigned next();
bool getBit(const unsigned);
void onBit(const unsigned);
void offBit(const unsigned);
void switchBit(const unsigned);
void printBoard();
};
PrimeGen::PrimeGen(const unsigned max_num)
{
limit = max_num;
sievelen = limit / DWORD_BITS + 1;
bookmark = 0;
sieve = (unsigned char*) malloc(sievelen);
for (unsigned i = 0; i < sievelen; i++) {sieve[i] = 0;}
firstPass();
}
inline bool PrimeGen::getBit(const unsigned index)
{
return sieve[index/DWORD_BITS] & flag(index%DWORD_BITS);
}
inline void PrimeGen::onBit(const unsigned index)
{
sieve[index/DWORD_BITS] |= flag(index%DWORD_BITS);
}
inline void PrimeGen::offBit(const unsigned index)
{
sieve[index/DWORD_BITS] &= ~flag(index%DWORD_BITS);
}
inline void PrimeGen::switchBit(const unsigned index)
{
sieve[index/DWORD_BITS] ^= flag(index%DWORD_BITS);
}
void PrimeGen::firstPass()
{
unsigned nmod,n,x,y,xroof, yroof;
//n = 4x^2 + y^2
xroof = (unsigned) sqrt(((double)(limit - 1)) / 4);
for(x = 1; x <= xroof; x++){
yroof = (unsigned) sqrt((double)(limit - 4 * x * x));
for(y = 1; y <= yroof; y++){
n = (4 * x * x) + (y * y);
nmod = n % 12;
if (nmod == 1 || nmod == 5){
switchBit(n);
}
}
}
xroof = (unsigned) sqrt(((double)(limit - 1)) / 3);
for(x = 1; x <= xroof; x++){
yroof = (unsigned) sqrt((double)(limit - 3 * x * x));
for(y = 1; y <= yroof; y++){
n = (3 * x * x) + (y * y);
nmod = n % 12;
if (nmod == 7){
switchBit(n);
}
}
}
xroof = (unsigned) sqrt(((double)(limit + 1)) / 3);
for(x = 1; x <= xroof; x++){
yroof = (unsigned) sqrt((double)(3 * x * x - 1));
for(y = 1; y <= yroof; y++){
n = (3 * x * x) - (y * y);
nmod = n % 12;
if (nmod == 11){
switchBit(n);
}
}
}
}
unsigned PrimeGen::next()
{
while (bookmark <= limit)
{
bookmark++;
if (getBit(bookmark))
{
unsigned out = bookmark;
for(unsigned num = bookmark * 2; num <= limit; num += bookmark)
{
offBit(num);
}
return out;
}
}
return 0;
}
inline void PrimeGen::printBoard()
{
for(unsigned i = 0; i < sievelen; i++)
{
if (i % 4 == 0)
cout << endl;
printBinary(sieve[i]);
cout << " ";
}
}
inline unsigned char flag(const unsigned char bit_index)
{
return ((unsigned char) 128) >> bit_index;
}
inline void printBinary(unsigned char byte)
{
unsigned int i = 1 << (sizeof(byte) * 8 - 1);
while (i > 0) {
if (byte & i)
cout << "1";
else
cout << "0";
i >>= 1;
}
}
私はそれをきれいにして読みやすくするために最善を尽くしました。私はプロのプログラマーではありませんので、ご容赦ください。
pgen という名前の PrimeGen オブジェクトを初期化し、その初期ビットボードを pgen.printBoard() で出力し (next() 反復の前に 23 と 59 が欠落していることに注意してください)、次に next() を反復して得られる出力を次に示します。返された素数をすべて出力します。
00000101 00010100 01010000 01000101
00000100 01010001 00000100 00000100
00010001 01000001 00010000 01000000
01000101 00010100 01000000 00000001
5
7
11
13
17
19
29
31
37
41
43
47
53
61
67
71
73
79
83
89
97
101
103
107
109
113
127
DONE
Process returned 0 (0x0) execution time : 0.064 s
Press any key to continue.