С ++ Сито Аткина упускает из виду несколько простых чисел

Недавно я работал над генератором простых чисел C ++, который использует решето Аткина (http://en.wikipedia.org/wiki/Sieve_of_atkin) для генерации простых чисел. Моя цель - создать любое 32-битное число. Я буду использовать его в основном для задач проекта Эйлера. в основном это просто летний проект.

Программа использует битовую доску для хранения простоты: то есть серии единиц и нулей, где, например, 11-й бит будет 1, 12-й - 0, 13-й - 1 и т. Д. Для эффективного использования памяти это фактически и массив символов, каждый из которых содержит 8 бит. Я использую флаги и побитовые операторы для установки и извлечения битов. Суть алгоритма проста: сделайте первый проход, используя некоторые уравнения, которые я не пытаюсь понять, чтобы определить, считается ли число «простым» или нет. По большей части это даст правильные ответы, но пара непростых чисел будет помечена как простые. Следовательно, при итерации по списку вы устанавливаете все числа, кратные только что найденному простому числу, как «не простое». Это имеет удобное преимущество: чем больше простое число, тем меньше процессорного времени.

У меня все готово на 90%, с одной загвоздкой: отсутствуют некоторые простые числа.

Изучив битовую доску, я убедился, что эти простые числа опускаются во время первого прохода, который в основном переключает число для каждого решения, которое оно имеет для ряда уравнений (см. Запись в википедии). Я снова и снова просматривал этот кусок кода. Я даже попытался увеличить границы до того, что показано в статьях в Википедии, что менее эффективно, но я подумал, что может достичь нескольких цифр, которые я как-то пропустил. Ничего не получилось. Эти числа просто оцениваются как непростые. Большая часть моего теста проводилась со всеми простыми числами меньше 128. Из этого диапазона пропущены следующие простые числа:

23 и 59.

Я не сомневаюсь, что на более высоком диапазоне будет не хватать большего (просто не хочу пересчитывать их все). Я не знаю, почему они отсутствуют, но они есть. Есть ли что-нибудь особенное в этих двух простых числах? Я дважды и трижды проверял, находил и исправлял ошибки, но все еще, вероятно, что-то глупое, чего мне не хватает.

в любом случае вот мой код:

#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;
    }
}

Я изо всех сил старался очистить его и сделать его читабельным. Я не профессиональный программист, так что будьте милосердны.

Вот результат, который я получаю, когда я инициализирую объект PrimeGen с именем pgen, распечатываю его начальную битовую панель с помощью pgen.printBoard () (обратите внимание, что 23 и 59 отсутствуют перед итерацией next ()), а затем перебираю 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.