Pull to refresh

Алгоритм нахождения N первых простых чисел

Algorithms *
Sandbox
В процессе обучения программированию я столкнулся с задачей поиска 1M первых простых чисел.

На Хабре уже не раз писали об этом. Но, чаще всего, речь шла о поиске всех простых чисел до числа N. Я же расскажу о том, как решал задачу нахождения N первых простых чисел.

На первый взгляд может показаться, что разницы в этих задачах практически нет. И это действительно так, пока мы не поймем, что перебор делителей справляется с задачей крайне плохо.

В результате проведенной работы на поиск 1 миллиона простых чисел у меня уходит всего 0,262 секунды, что, конечно же, не предел, но уже впечатляет. Алгоритм реализован на C++.

Как и автор предыдущей статьи о простых числах, я начал решать задачу лоб в лоб. Итак, простое число p — то, что имеет два делителя — единицу и p. Первое, что приходит в голову — проверять каждое число делением на все предыдущие числа.

Идея, как оказалось, плохая. Понял я это после первых 30 минут работы программы в поисках миллиона простых. Результаты получились такие:

10 000 — 2,349 с.
100 000 — 293,552 с.
1 000 000 — не дождался

Как и в уже упомянутой выше статье я принялся за оптимизацию алгоритма. Во-первых, проверять четные числа не имеет смысла, ведь 2 — единственное простое четное число. То же самое с делителями — на четные можно не проверять. В-третьих, в качестве делителей достаточно брать числа не превышающие половину проверяемого числа. Вот, что из этого вышло:

// Listing 1<br/>
 <br/>
#include <iostream><br/>
#include <vector><br/>
#include <ctime><br/>
 <br/>
using namespace std;<br/>
 <br/>
int main()<br/>
{<br/>
    clock_t start = clock();<br/>
 <br/>
    const int AIM = 10000;<br/>
 <br/>
    vector<int> prime_nums;<br/>
    prime_nums.push_back(2);<br/>
 <br/>
    for (int num = 3; prime_nums.size() < AIM; num += 2)<br/>
    {<br/>
        bool isprime = true;<br/>
        for (int i = 3; i <= num / 2; i += 2)<br/>
        {<br/>
            if (num % i == 0)<br/>
            {<br/>
                isprime = false;<br/>
                break;<br/>
            }<br/>
        }<br/>
 <br/>
        if (isprime)<br/>
            prime_nums.push_back(num);<br/>
    }<br/>
 <br/>
    cout << prime_nums[AIM - 1];<br/>
 <br/>
    clock_t end = clock();<br/>
    cout << "\n----------------" << endl;<br/>
    cout << "Time: " << double(end - start) / CLOCKS_PER_SEC << " sec" << endl;<br/>
}


Стало лучше, но миллиона я так и не смог дождаться
10 000 — 0,589 с.
100 000 — 72,662 с.
1 000 000 — после 20 минут ждать перестал

Понятно, что продолжать оптимизировать этот алгоритм можно еще долго, чем я в общем-то и занимался. Можно не проверять числа заканчивающиеся на 5 и 0. Можно составить список самых распространенных делителей и проверять каждое число сначала на них. В какой-то момент мне показалось, что я смогу решить проблему, проверяя новые числа только на уже найденные простые. Значительных результатов это не принесло. Проблема крылась в том, что с каждым новым числом, количество проверок возрастало, да и операция деления — одна из самых тяжелых.

Отправился в поиск. Узнал про решето Эратосфена. Но использовать алгоритм в том виде, как он есть для моей задачи нельзя. Заранее ведь неизвестно*, сколько надо просеять чисел, чтобы найти N простых. Поэтому я решил модифицировать алгоритм под свои нужды следующим образом.

Просеивать заранее большой набор натуральных чисел не выгодно. Значит надо просеивать некоторые наборы натуральных чисел до тех пор, пока не достигнем нужного результата. Причем каждый новый набор надо сначала просеять на уже найденные простые числа. Вот код, полученный уже после ряда оптимизаций:

// Listing 2<br/>
 <br/>
#include <iostream><br/>
#include <ctime><br/>
#include <vector><br/>
 <br/>
using namespace std;<br/>
 <br/>
int main()<br/>
{<br/>
    clock_t start = clock();<br/>
 <br/>
    const int AIM = 1000000;<br/>
 <br/>
    int startSize = AIM; // стартовый размер массива натуральных чисел<br/>
    int addSize = AIM; // размер дополнительного массива натуральных чисел<br/>
 <br/>
    vector<bool> nums(startSize);<br/>
    vector<int> primeNums(AIM);<br/>
 <br/>
    int foundPrimes = 0;<br/>
 <br/>
    for (int i = 2; i < startSize; i++)<br/>
        nums[i] = true;<br/>
 <br/>
    bool addition = false;<br/>
    int adder = 0;<br/>
    while (true)<br/>
    {<br/>
        if (addition)<br/>
        {<br/>
            nums.resize(nums.size() + addSize, true);<br/>
 <br/>
            // исключим составные числа простыми из предыдущих итераций<br/>
            for (int i = 0; i < foundPrimes; i++)<br/>
            {<br/>
                int cur_num = primeNums[i];<br/>
                if ((addSize + ((nums.size() - addSize) % cur_num)) < cur_num)<br/>
                    continue;<br/>
                for (int j = ((nums.size() - addSize) / cur_num) * cur_num; j < nums.size(); j += cur_num)<br/>
                    nums[j] = false;<br/>
            }<br/>
        }<br/>
        else<br/>
            addition = true;<br/>
 <br/>
 <br/>
        int iter;<br/>
        if (foundPrimes == 0)<br/>
            iter = 2;<br/>
        else<br/>
            iter = primeNums[foundPrimes - 1] + 2;<br/>
 <br/>
        for ( ; iter < nums.size(); iter++)<br/>
        {<br/>
            // выбираем очередное простое число<br/>
            if (nums[iter])<br/>
            {<br/>
                primeNums[foundPrimes] = iter;<br/>
                foundPrimes++;<br/>
                if (foundPrimes == AIM)<br/>
                    break;<br/>
                // просеиваем<br/>
                for (int j = iter + iter; j < nums.size(); j += iter)<br/>
                    nums[j] = false;<br/>
            }<br/>
            else<br/>
                continue;<br/>
 <br/>
        }<br/>
        if (foundPrimes == AIM)<br/>
            break;<br/>
    }<br/>
 <br/>
    cout << "Last prime: " << primeNums[AIM -1];<br/>
 <br/>
    clock_t end = clock();<br/>
    cout << "\n----------------" << endl;<br/>
    cout << "Time: " << double(end - start) / CLOCKS_PER_SEC << " sec" << endl;<br/>
    return 0;<br/>
}


Результаты даже слегка удивили:
10 000 — 0,001 с.
100 000 — 0,021 с.
1 000 000 — 0,262 с.

10 000 000 — 2,999 с.

На этом с задачей поиска 1M первых простых чисел я закончил. Результат меня вполне удовлетворил. Но что можно сделать ещё? Можно реализовать какой-нибудь из известных алгоритмов на основе решета Эратосфена, например, решето Сундарама. Но вряд ли результаты изменяться значительно. Возможно, решето Аткина сможет показать себя лучше — именно этот алгоритм чаще всего используется в современной криптографии.

* — неизвестно мне. halyavin рассказал в первом комментарии, что оказывается давно известно.

UPD
Относительно потребления памяти, описанные алгоритмы ведут себя следующим образом:

1 млн. простых чисел, максимальное ОЗУ
— перебор делителей: ~4Мб
— модифицированное решето Эратосфена: ~6Мб

10 млн. простых чисел, максимальное ОЗУ
— перебор делителей: ~39Мб
— модифицированное решето Эратосфена: ~61Мб
Tags:
Hubs:
Total votes 27: ↑16 and ↓11 +5
Views 26K
Comments Comments 29