Простые числа. Решето Эратосфена

Это ста­тья по­свя­ще­на про­стым чис­лам и эф­фек­тив­ным спо­со­бам их вы­чис­ле­ния. Сра­зу ска­жу, что те ал­го­рит­мы, ко­то­рые тут при­ве­де­ны, яв­ля­ют­ся весь­ма и весь­ма нетри­ви­аль­ны­ми и са­мо­му мне не да­ва­лись до­воль­но дол­го, но в ито­ге я их всё-таки при­ду­мал, на­пи­сал и спе­шу по­де­лить­ся со все­ми вами. Ис­ход­ные тек­сты в ста­тье бу­дут при­ве­де­ны на язы­ках C# и Haskell.

Про­стое чис­ло – это на­ту­раль­ное чис­ло боль­ше еди­ни­цы, ко­то­рое име­ет ров­но два де­ли­те­ля: еди­ни­цу и само это чис­ло.

Ре­ше­то Эра­то­сфе­на – древ­ний, но при этом весь­ма эф­фек­тив­ный и до сих пор ши­ро­ко ис­поль­зу­е­мый ал­го­ритм по­ис­ка всех про­стых чи­сел не пре­вос­хо­дя­щих неко­то­ро­го N.

За­пи­шем под­ряд все чис­ла от 2 до N. Даль­ше вы­черк­нем из это­го спис­ка все чис­ла крат­ные 2, ис­клю­чая саму двой­ку, по­том вы­черк­нем все чис­ла крат­ные 3, ис­клю­чая само чис­ло 3, чис­ло 4 уже вы­черк­ну­то, вы­чер­ки­ва­ем чис­ла крат­ные 5 и т.д. Про­дол­жа­ем этот про­цесс, пока квад­рат оче­ред­но­го чис­ла не пре­вы­сит N.

Са­мая про­стая про­грамм­ная ре­а­ли­за­ция это­го ал­го­рит­ма вы­гля­дит сле­ду­ю­щим об­ра­зом:

bool[] table = new bool[100];
int i,j;
// Отмечаем все числа как простые
for ( i = 0; i < table.Length; i++ )
    table[i] = true;
// Вычеркиваем лишнее
for ( i = 2; i*i < table.Length; i++ )
    if ( table[i] )
        for ( j = 2*i; j < table.Length; j += i )
            table[j] = false;
// Выводим найденное
for ( i = 2; i < table.Length; i++ )
    if ( table[i] )
        Console.WriteLine(i);

Итак, про­грамм­ка очень про­стая, но что в ней мож­но улуч­шить? При­ве­дён­ная про­грам­ма вы­ве­дет спи­сок всех про­стых чи­сел из пер­вой сот­ни. Если необ­хо­ди­мо вы­чис­лить все чис­ла до 1000 до­ста­точ­но про­сто из­ме­нить раз­мер мас­си­ва ta­ble. Од­на­ко в C# объ­ем па­мя­ти, за­ни­ма­е­мый пе­ре­мен­ной типа bool ра­вен од­но­му бай­ту. Сле­до­ва­тель­но, для на­хож­де­ния всех чи­сел до од­но­го мил­ли­ар­да нам по­тре­бу­ет­ся па­мя­ти чуть мень­ше чем ги­га­байт. А нам ведь хо­чет­ся най­ти очень мно­го про­стых чи­сел, не толь­ко до жал­ко­го мил­ли­ар­да, прав­да? :)

Пер­вая и оче­вид­ная оп­ти­ми­за­ция – это ис­поль­зо­вать для хра­не­ния ин­фор­ма­ции об од­ном чис­ле од­но­го бита. В этом слу­чае в один ги­га­байт вой­дут все чис­ла не пре­вос­хо­дя­щие 8`589`934`592. Для ра­бо­ты с та­ки­ми чис­ла­ми уже не под­хо­дят 32-х бит­ные пе­ре­мен­ные. Очень хо­ро­шо! Пе­ре­хо­дим на 64-х бит­ные. Сама про­грам­ма:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66

using System;
using System.Collections.Generic;
using System.Text;

namespace sieve
{
    class Program
    {
        static uint[] table;

        static void reset_bit(long bit_number)
        {
            long i = bit_number / 32;
            long p = bit_number % 32;
            table[i] = table[i] & (uint)(~(1 << (int)p));
        }

        static bool get_bit(long bit_number)
        {
            long i = bit_number / 32;
            long p = bit_number % 32;
            if ( (table[i] & (uint)(1 << (int)p)) != 0 )
                return true;
            return false;
        }

        static void Main(string[] args)
        {
            DateTime dt = DateTime.Now;
            // One Gb!
            // table = new uint[1024*1024*1024 / 4];

            table = new uint[1024*1024*10 / 4];
            long bit_count = (long)table.Length * (long)32;

            long i,j;

            for ( i = 0; i < table.Length; i++ )
                table[i] = 0xffffffff;

            for ( i = 2; i*i < bit_count; i++ )
                if ( get_bit(i) )
                    for ( j = i*2; j < bit_count; j += i )
                        reset_bit(j);

            TimeSpan calc_time = DateTime.Now - dt;
            dt = DateTime.Now;

            long prime_count = 0;
            long max = 0;

            for ( i = 2; i < bit_count; i++ )
                if ( get_bit(i) )
                {
                    //Console.WriteLine(i);
                    prime_count++;
                    max = i;
                }

            Console.WriteLine("Всего найдено: {0}", prime_count);
            Console.WriteLine("Максимальное: {0}", max);
            Console.WriteLine("Время на вычисление: {0}", calc_time);
            Console.WriteLine("Время на вывод: {0}", DateTime.Now - dt);
        }
    }
}

Кому хва­тит для жиз­ни, до­стиг­ну­тых вось­ми с по­ло­ви­ной мил­ли­ар­дов, мо­гут даль­ше не чи­тать. А со все­ми кому мало, пе­ре­хо­дим к са­мо­му ин­те­рес­но­му. Если вни­ма­тель­но по­смот­реть, мож­но уви­деть, что сра­зу по­сле пер­вой ите­ра­ции ока­зы­ва­ют­ся вы­черк­ну­ты­ми все чет­ные чис­ла. Та­ким об­ра­зом, бук­валь­но по­сле пер­вых несколь­ких ите­ра­ций бу­дет вы­черк­ну­то боль­ше по­ло­ви­ны мас­си­ва. Су­ще­ству­ет ли спо­соб не хра­нить эти чис­ла в мас­си­ве, тем са­мым су­ще­ствен­но сэко­но­мив па­мять? Я не буду при­во­дить про­грам­му, хра­ня­щую толь­ко нечет­ные чис­ла, я сра­зу обоб­щу за­да­чу сле­ду­ю­щим об­ра­зом: ре­а­ли­зо­вать ал­го­ритм ре­ше­та, та­ким об­ра­зом, что­бы мас­сив не вклю­чал чи­сел, крат­ных неко­то­ро­му ко­ли­че­ству пер­вых про­стых (на­при­мер: 2,3,5,7).

Для ис­сле­до­ва­ния раз­ных за­ко­но­мер­но­стей рас­пре­де­ле­ния чи­сел, не де­ля­щих­ся на пер­вые несколь­ко про­стых, я ис­поль­зо­вал язык Haskell. Здесь я при­ве­ду толь­ко ту це­поч­ку рас­суж­де­ний, и по­ка­жу толь­ко тот код, ко­то­рые при­ве­ли к вер­но­му ре­ше­нию. На са­мом деле я пред­при­нял мно­го по­пы­ток.

Итак, за­пус­ка­ем hugs или ghci и на­чи­на­ем кол­до­вать :) .

Для на­ча­ла нам нуж­но опре­де­лить спи­сок чи­сел, из ко­то­ро­го ис­клю­че­ны крат­ные неко­то­ро­му мно­же­ству чи­сел:

nums a b = [x | x <- [a..], not (isDiv x b)]
    where
    isDiv num lst =
        case lst of
        []     -> False
        (h:t)   -> if num `mod` h == 0 then True else isDiv num t

> take 10 (nums 3 [2])
[3,5,7,9,11,13,15,17,19,21]
> take 10 (nums 5 [2,3])
[5,7,11,13,17,19,23,25,29,31]

Для удоб­ства я опре­де­лю несколь­ко ва­ри­ан­тов, с ко­то­ры­ми и буду ра­бо­тать:

sieve_2_3 = nums 5   [2,3]
sieve_2_5 = nums 7   [2,3,5]
sieve_2_7 = nums 11  [2,3,5,7]

Нам по­на­до­бить­ся ра­бо­тать с раз­но­стя­ми меж­ду бли­жай­ши­ми чис­ла­ми в по­сле­до­ва­тель­но­стях. Функ­ция пре­об­ра­зо­ва­ния неко­то­рой по­сле­до­ва­тель­но­сти в спи­сок раз­но­стей меж­ду её со­сед­ни­ми эле­мен­та­ми:

rzn lst = [b - a | (a,b) <- zip lst (tail lst)]

> take 10 (rzn $ nums 3 [2])
[2,2,2,2,2,2,2,2,2,2]
> take 10 (rzn sieve_2_3)
[2,4,2,4,2,4,2,4,2,4]
> take 10 (rzn sieve_2_5)
[4,2,4,2,4,6,2,6,4,2]

Ага! Уже что-то есть, мож­но за­ме­тить пе­ри­о­ди­че­ски по­вто­ря­ю­щи­е­ся по­сле­до­ва­тель­но­сти. Пи­шем функ­цию по­ис­ка по­вто­ря­ю­щих­ся по­сле­до­ва­тель­но­стей:

get_period lst =
    let
        is_all_eq lst =
            let
                fs = head lst
            in
                foldl1 (&&) (map (== fs) (tail lst))
        slice_period p lst = [take p lst] ++ (slice_period p (drop p lst))
        find_period p =
            if is_all_eq $ take 5 (slice_period p lst) then p
                else find_period (p+1)
        found_period = find_period 1
    in
        take found_period lst

Я не буду на ней по­дроб­но оста­нав­ли­вать­ся. Про­сто пред­по­ло­жим, что она есть и кор­рект­но ра­бо­та­ет :) . Кому ин­те­рес­но, я ду­маю, смо­гут в ней разо­брать­ся. Про­ве­ря­ем с по­мо­щью этой функ­ции наши по­сле­до­ва­тель­но­сти:

> get_period $ rzn $ nums 3 [2]
[2]
> get_period $ rzn $ sieve_2_3
[2,4]
> get_period $ rzn $ sieve_2_5
[4,2,4,2,4,6,2,6]
> get_period $ rzn $ sieve_2_7
[2,4,2,4,6,2,6,4,2,4,6,6,2,6,4,2,6,4,6,8,4,2,4,2,4,
8,6,4,6,2,4,6,2,6,6,4,2,4,6,2,6,4,2,4,2,10,2,10]

Те­перь раз­бе­рем­ся, в ка­кую сто­ро­ну нам пред­сто­ит дви­гать­ся даль­ше. Соб­ствен­но для пол­но­го сча­стья нам необ­хо­ди­мо на­учить­ся де­лать две вза­и­мо­об­рат­ные вещи: во-пер­вых, по но­ме­ру чис­ла в мас­си­ве опре­де­лять само это чис­ло (ведь чис­ло в мас­си­ве это все­го один бит), а во-вто­рых, по чис­лу опре­де­лять его по­зи­цию в мас­си­ве.

Са­мое вре­мя пе­рей­ти к кон­крет­ным при­ме­рам. Бу­дем раз­би­рать­ся с мас­си­вом без чи­сел крат­ных двум и трем.

Пер­вые 10 чи­сел это­го спис­ка вы­гля­дят так:

5,7,11,13,17,19,23,25,29,31

Его пе­ри­од:

2,4

Сле­до­ва­тель­но, этот мас­сив со­сто­ит из под­ряд иду­щих пар, чис­ла внут­ри пары от­ли­ча­ют­ся на 2, раз­ни­ца меж­ду пер­вы­ми чис­ла­ми под­ряд иду­щих пар рав­на 6. При этом не за­бы­ва­ем, что эле­мент с ин­дек­сом 0 ра­вен 5.

Име­ем функ­цию на­хож­де­ния ин­дек­са по чис­лу:

find_idx_by_num_2_3 num =
    let
        a = num-5
        d = a `div` 6
        m = a `mod` 6
    in
        d*2 + (if m < 2 then 0 else 1)

Об­ра­ти­те вни­ма­ние, мы про­ве­ря­ем оста­ток от де­ле­ния на 6, если он боль­ше 2 (раз­ни­цы меж­ду чис­ла­ми в паре), то, зна­чит, мы име­ем дело со вто­рым чис­лом в паре и, сле­до­ва­тель­но, долж­ны к ин­дек­су при­ба­вить еди­ни­цу.

> find_idx_by_num_2_3 5
0
> find_idx_by_num_2_3 7
1
> find_idx_by_num_2_3 25
7

Об­рат­ная функ­ция. На­хож­де­ние чис­ла по ин­дек­су:

find_num_by_idx_2_3 i =
    let
        a = i `div` 2
        b = i `mod` 2
    in
        a * 6 + (if b == 1 then 2 else 0) + 5

Тут мы уже де­лим и на­хо­дим оста­ток от де­ле­ния ар­гу­мен­та на 2 – дли­ну пе­ри­о­да для на­ше­го мас­си­ва.

> find_num_by_idx_2_3 0
5
> find_num_by_idx_2_3 1
7
> find_num_by_idx_2_3 7
25

Ну что ж, в част­ном слу­чае за­да­ча ре­ше­на. Оста­лась са­мая ма­лость – обоб­щить по­лу­чен­ный ре­зуль­тат до спис­ка с про­из­воль­ным чис­лом вы­черк­ну­тых пер­вых про­стых.

Ин­декс по чис­лу:

find_idx_by_num num period first_number =
    let
        period_sum = foldl1 (+) period
        period_len = integer_len period
        a = num - first_number
        d = a `div` period_sum
        m = a `mod` period_sum
        find_inc lst sum len n =
            case lst of
            []     -> len
            (h:t)     ->
                if n < sum then len
                else find_inc t (sum+h) (len+1) n
    in
        d*period_len + (find_inc period 0 0 m) - 1

Здесь наи­боль­ший ин­те­рес пред­став­ля­ет функ­ция find_inc, она на­хо­дит по­ло­же­ние чис­ла внут­ри пе­ри­о­да: сум­ми­ру­ем по­оче­ред­но чис­ла из пе­ри­о­да, пока по­лу­чен­ное чис­ло не ста­нет боль­ше чем ис­ко­мое чис­ло (оста­ток от де­ле­ния ар­гу­мен­та num, функ­ции find_idx_by_num на сум­му всех чи­сел пе­ри­о­да).

Чис­ло по ин­дек­су:

find_num_by_idx i period first_number =
    let
        period_sum = foldl1 (+) period
        period_len = integer_len period
        a = i `div` period_len
        b = i `mod` period_len
        sm n list
            | n == 0    = 0
            | otherwise    =
                case list of
                []    -> 0
                (h:t)    -> h + sm (n-1) t
    in
        a * period_sum + (sm b period) + first_number

Здесь функ­ция sm сум­ми­ру­ет пер­вые n чи­сел спис­ка list.

Са­мые гла­за­стые за­ме­ти­ли ис­поль­зо­ва­ние некой неопи­сан­ной функ­ции in­te­ger_len, ис­прав­ля­юсь:

integer_len l =
    case l of
    []    -> 0
    (h:t)    -> 1 + integer_len t

Это обыч­ный length, от­ли­ча­ю­щий­ся от стан­дарт­но­го ти­пом воз­вра­ща­е­мо­го зна­че­ния, в дан­ном слу­чае это In­te­ger вме­сто Int. Haskell язык со стро­гой ти­пи­за­ци­ей, и воль­но­стей меж­ду Int и In­te­ger он не поз­во­ля­ет.

Про­ве­ря­ем, про­ве­ря­ем, про­ве­ря­ем:

> take 10 sieve_2_7
[11,13,17,19,23,29,31,37,41,43]
> find_idx_by_num 19 (get_period $ rzn $ sieve_2_7) 11
3
> find_num_by_idx 3 (get_period $ rzn $ sieve_2_7) 11
19

Ну что ж, тео­ре­ти­че­скую часть на этом мож­но счи­тать за­кон­чен­ной, спо­соб ра­бо­ты с нуж­ным фор­ма­том мас­си­ва най­ден, с чем я себя и вас по­здрав­ляю! Ура! :)

Пе­ре­хо­дим к прак­ти­че­ской ре­а­ли­за­ции все­го вы­ше­на­кол­до­ван­но­го.

Сра­зу при­хо­дит в го­ло­ву при­мер­но та­кая ре­а­ли­за­ция: бе­рем оче­ред­ной нену­ле­вой бит на­ше­го мас­си­ва, по его но­ме­ру на­хо­дим чис­ло, ко­то­ро­му он со­от­вет­ству­ет, даль­ше на­чи­на­ем пе­ре­би­рать все крат­ные ему чис­ла и вы­чер­ки­вать их из мас­си­ва с по­мо­щью функ­ции по­ис­ка но­ме­ра чис­ла по са­мо­му это­му чис­лу. Зву­чит вро­де не пло­хо и так уже мож­но сде­лать ра­бо­та­ю­щую ре­а­ли­за­цию, од­на­ко нуж­но учесть, что каж­дый по­иск чис­ла по ин­дек­су и по­иск ин­дек­са по чис­лу со­дер­жит внут­ри себя пе­ре­бор пе­ри­о­да для на­ше­го фор­ма­та, хо­те­лось бы от этой опе­ра­ции из­ба­вить­ся.

Сно­ва ухо­дим в глу­би­ну Haskell’я (по­след­ний раз, кля­нусь :) ). По­смот­рим, как рас­пре­де­ле­ны чис­ла крат­ные неко­то­ро­му N в на­ших мас­си­вах. Функ­ция на­хо­дит ин­дек­сы всех чи­сел крат­ных num в спис­ке list:

div_idx list num = [b | (a,b) <- zip list [0..], a `mod` num == 0]

К при­ме­ру, ин­дек­сы чи­сел крат­ных 5 в мас­си­ве без дво­ек и тро­ек:

> take 10 (div_idx sieve_2_3 5)
[0,7,10,17,20,27,30,37,40,47]

Ну что, зна­ко­мая кар­ти­на, явно наша лю­би­мая пе­ри­о­ди­че­ски по­вто­ря­ю­ща­я­ся по­сле­до­ва­тель­ность раз­но­стей меж­ду со­сед­ни­ми чис­ла­ми:

> get_period $ rzn $ div_idx sieve_2_3 5
[7,3]

Ко­неч­ный ал­го­ритм бу­дет вы­гля­деть сле­ду­ю­щим об­ра­зом: бе­рем оче­ред­ной еди­нич­ный бит, на­хо­дим со­от­вет­ству­ю­щее ему чис­ло. На­хо­дим N+1 бли­жай­ших крат­ных ему чи­сел (при этом ис­клю­ча­ем чис­ла крат­ные про­стым из спис­ка вы­черк­ну­тых) и на­хо­дим раз­но­сти меж­ду их ин­дек­са­ми. Даль­ше на­чи­на­ем по­сле­до­ва­тель­но об­ну­лять ин­дек­сы всех крат­ных най­ден­но­му чи­сел. К при­ме­ру, мы об­на­ру­жи­ли еди­нич­ный бит с но­ме­ром 0 со­от­вет­ству­ю­щий чис­лу 5 в мас­си­ве без 2 и 3, мы зна­ем, что ин­дек­сы чи­сел крат­ных 5 из­ме­ня­ют­ся с пе­ри­о­дом [7,3] (см. выше), со­от­вет­ствен­но, об­ну­ля­ем биты с но­ме­ра­ми 7,10,17,20,27 и т.д.

В этом ме­сте сле­до­ва­ло бы ска­зать, что и это ещё не всё, од­на­ко это­го (пока?) не про­изой­дет. Даль­ней­ший путь уско­ре­ния я вижу в вы­яв­ле­нии ма­те­ма­ти­че­ско­го за­ко­на, фор­ми­ру­ю­ще­го пе­ри­од раз­но­стей ин­дек­сов чи­сел крат­ных неко­то­ро­му N. Гос­по­да ма­те­ма­ти­ки, если у вас есть ка­кие-то идеи на этот счет, ми­ло­сти про­шу пи­сать ком­мен­та­рии тут или мне на e-mail.

Я не буду тут при­во­дить пол­ный ва­ри­ант фи­наль­ной про­грам­мы на C#, оста­нов­люсь раз­ве толь­ко на C# ре­а­ли­за­ции функ­ций по­лу­че­ния ин­дек­са по чис­лу и чис­ла по ин­дек­су:

public Int64 get_index_by_number(Int64 num)
{
    Int64 a = num - first_number;
    Int64 d = a / period_sum;
    Int64 m = a % period_sum;

    int sum = 0;
    int i;

    for ( i = 0; i < period.Length; i++ )
    {
        if ( m < sum ) break;
        sum += period[i];
    }

    i--;

    return d*period.Length + i;
}

public Int64 get_number_by_index(Int64 idx)
{
    Int64 a = idx / period.Length;
    Int64 b = idx % period.Length;

    Int64 sum = 0;
    for ( int i = 0; i < b; i++ )
        sum += period[i];

    return a * period_sum + sum + first_number;
}

Здесь pe­ri­od, pe­ri­od_sum и first_num­ber опре­де­ле­ны в фай­ле gen.cs, этот файл ге­не­ри­ру­еть­ся про­грам­мой на Haskell, один из ва­ри­ан­тов это­го фай­ла:

namespace PrimeNumbers
{
    partial class PNTable
    {
        protected static readonly int[] period = {4,2,4,2,4,6,2,6};
        protected static readonly int[] base_primes = {2,3,5};
        protected const long period_sum = 30;
        protected const long first_number = 7;
    }
}

Ге­не­ра­тор для та­ких фай­лов на­хо­дит­ся в пап­ке haskell внут­ри при­ло­жен­но­го к ста­тье фай­ла. Он по­лу­ча­ет из ко­манд­ной стро­ки чис­ло пер­вых про­стых, ко­то­рые долж­ны быть ис­клю­че­ны и вы­во­дит в стан­дарт­ный по­ток вы­во­да сге­не­ри­ро­ван­ный C# код.

Те­перь по­смот­рим в циф­рах, чего мы всё-таки до­би­лись:

Клас­си­че­ский ва­ри­ант, для хра­не­ния чис­ла ис­поль­зу­ет­ся один бит. Раз­мер мас­си­ва один ги­га­байт:

Всего найдено простых чисел: 393`615`806
Максимальное простое число: 8`589`934`583
Время на вычисление: 00:41:44

Ва­ри­ант ре­ше­та с ис­клю­че­ни­ем чи­сел крат­ных 2,3,5,7,11,13. Раз­мер мас­си­ва один ги­га­байт:

Всего найдено простых чисел: 1`907`461`679
Максимальное простое число: 44`783`981`911
Время на вычисление: 00:58:58

С по­мо­щью улуч­шен­но­го ал­го­рит­ма уда­лось най­ти по­чти в 5 раз боль­ше про­стых чи­сел, при этом за­тра­ты по вре­ме­ни уве­ли­чи­лись толь­ко в пол­то­ра раза. Очень хо­ро­ший ре­зуль­тат!

Я дол­жен за­ме­тить, что даль­ше ис­клю­че­ния по­сле­до­ва­тель­но­сти 2,3,5,7,11,13 про­дви­гать­ся нет смыс­ла. Во-пер­вых за­тра­ты вре­ме­ни на вы­чис­ле­ние мас­си­ва с раз­но­стя­ми ин­дек­сов для каж­до­го оче­ред­но­го чис­ла ста­но­вят­ся очень боль­ши­ми. Во-вто­рых эф­фект от каж­до­го оче­ред­но­го чис­ла силь­но сни­жа­ет­ся. Рас­смот­рим, как ме­ня­ет­ся диа­па­зон чи­сел до­ступ­ных для ра­бо­ты. Бу­дем брать со­тый член раз­ных по­сле­до­ва­тель­но­стей:

> (nums 1 []) !! 100
101
> (nums 3 [2]) !! 100
203
> (nums 5 [2,3]) !! 100
305
> (nums 7 [2,3,5]) !! 100
379
> (nums 11 [2,3,5,7]) !! 100
443
> (nums 13 [2,3,5,7,11]) !! 100
491
> (nums 17 [2,3,5,7,11,13]) !! 100
529
> (nums 19 [2,3,5,7,11,13,17]) !! 100
569
> (nums 23 [2,3,5,7,11,13,17,19]) !! 100
593
> (nums 29 [2,3,5,7,11,13,17,19,23]) !! 100
601

Вид­но, что наи­боль­ший эф­фект име­ет вы­чер­ки­ва­ние двой­ки, чуть худ­ший – трой­ки и чем даль­ше, тем мень­ший эф­фект от вы­чер­ки­ва­ния оче­ред­но­го про­сто­го чис­ла.

Та­к­же надо ска­зать, что ис­поль­зо­вать мас­сив с вы­черк­ну­ты­ми 2,3,5,7,11,13 мож­но толь­ко на боль­ших объ­е­мах па­мя­ти, по край­ней мере в моей ре­а­ли­за­ции, ина­че не бу­дет хва­тать чи­сел для вы­чис­ле­ния мас­си­ва с раз­но­стя­ми для оче­ред­но­го чис­ла и этот мас­сив бу­дет за­пол­нен некор­рект­но.

Мне было очень ин­те­рес­но ре­шать эту за­да­чу, на­де­юсь, вам та­к­же было ин­те­рес­но это все чи­тать. Я по­ста­рал­ся по­ка­зать на при­ме­ре очень ин­те­рес­ное свой­ство функ­ци­о­наль­ных язы­ков: код жи­вет бук­валь­но на кон­чи­ках паль­цев, про­грам­му мож­но кру­тить и раз­во­ра­чи­вать как угод­но, мож­но очень быст­ро по­ве­рять свои идеи и смот­реть на воз­ни­ка­ю­щие за­ко­но­мер­но­сти. И, кро­ме все­го, при­ят­но чув­ство­вать себя оси­лив­шим дол­го не да­вав­шу­ю­ся слож­ную за­да­чу :) .

Пол­ные ис­ход­ные тек­сты всех при­ве­ден­ных в ста­тье про­грамм на­хо­дят­ся в при­ла­га­е­мом фай­ле.