Алгебра в программе Mathematica


Пример 2



Пример 2




Чтобы в таблице сначала шли меньшие модули s (они помогают отбраковать больше кандидатов), таблицу пришлось отсортировать, предварительно вычеркнув из нее элементы NULL, которые соответствуют тем простым модулям д, для которых сравнение 3*2n +1 = 0 (mod q) неразрешимо.

Каждая пара в этой таблице описывает арифметическую профессию из показателей, для которых выполняется некоторое тождественное сравнение. Пусть, например, в таблице имеется пара (s, t). Она была получена для некоторого простого числа q. Ее наличие в таблице означает, что 3*2n+1 +1 = 0 (mod q) для любых натуральных k.

Теперь давайте предварительно оценим количество чисел, которые может помочь отбраковать такая таблица. Для этого сначала составим программу, которая подсчитывает количество чисел, отбракованных с помощью таблицы. Для этой программы нам понадобится функция, отбраковывающая заданное число по заданной таблице. Вот код этой функции.
fx[n_,t_]:= Block[{1= Length[t],i=l,c=(i<=l) }, While[c, {s,r}=t[[i]]; answer= (Mod [M, s] !=r) ; i++; c=(i<=l)S&answer]; answer]
Эта функция принимает значение True, если число я не было отсеяно с помощью таблицы t, т.е. если оно подлежит дальнейшему испытанию. Теперь можем составить программу подсчета отбракованных чисел.
prog:= Block[{nn=l00000,Yes=0,No=0,j=l}, While[j]=nn, If[fx[j,t],Yes++,No++]; J++1; Print["Yes=",Yes,"; No=",No]]
Заметьте, что программа отличается от аналогичной программы для чисел вида 5*2n+ 1. Отличие связано с тем, что на этот раз перебор не ограничивается только нечетными показателями. Конечно, как и для чисел вида 5-2n + 1, количество отбракованных чисел и время отбраковки зависят от размера таблицы t. Сначала давайте сгенерируем совсем небольшую таблицу.

t=Sort[DeleteCases[Array[f3,10,3] ,Null]]

{{3,1},{4,3},{10,7},{12,2],{18,14],{28,9),{36,28}}

Теперь запускаем программу.

prog//Timing Yes= 33650 ; N0= 66350

{4.937 Second,Null}

Как видите, менее чем за 5 секунд удалось отбраковать более 66 % чисел! Давайте попробуем увеличить таблицу.
t=Sort[DeleteCases[Array[f3,25,3],Null]] {{3,1},{4,3},{10,7),{12,2},{18,14),{28,9},{36,28),{39,29}, {48,5},{51,4 2},{52,9},{58,37},{60,24},{66,60},{82,51},{100,81}}
Опять считаем процент отбракованных чисел и замеряем время.

prog//Timing Yes= 25645 ; N0= 74355

{7.188 Second,Null}

Время отбраковки возросло чуть более чем на 2 секунды, а отбраковали мы более 74 % чисел! Значит, таблицу стоит увеличить.

t=Sort[DeleteCases[Array[f3,100,3], Null]];

Опять измеряем время.

prog//Timing Yes = 21055 ; N0= 78945

{14.765 Second,Null}

Как видим, таблицу увеличивать стоило! Время отбраковки, конечно, возросло, притом более чем вдвое, но 15 секунд по сравнению с часами счета роли не играют. Так что увеличение таблицы окупает себя! Увеличиваем таблицу еще.

t=Sort[DeleteCases[Array[f3,1000,3],Null]];

 Теперь снова измеряем время и процент отбраковки.

prog//Timing Yes= 15388 ; N0= 84612

{84.64 Second,Null}

Конечно, время возросло: теперь отбраковка занимает более одной минуты. Но проверять по-настоящему придется чуть более 15 % чисел! Потратив дополнительно на предварительную отбраковку немногим более одной минуты, мы почти в 8 раз уменьшим количество чисел, которые надо будет проверить на простоту. Так что увеличение таблицы оправдало себя и на этот раз. Поэтому еще раз увеличиваем таблицу.

t=Sort[DeleteCasestArray[f3,l0000,3],Null]];//Timing

{90.703 Second,Null}

На этот раз сама генерация таблицы занимает более полутора минут. Но и это время окупается.

prog//Timing Yes = 12123 ; N0= 87877

{630.734 Second,Null}

Количество чисел, которые подлежат проверке на простоту, мы уменьшили более чем в 8 раз! Может быть, стоит увеличить таблицу еще? Ну, на этот раз, во всяком случае, уж точно не в 10 раз. Ведь ждать генерации таблицы несколько часов — занятие весьма неприятное! Поэтому лучше сначала попытаться увеличить таблицу, скажем, в 2,5 раза.

t=Sort[DeleteCases[Array!f3,25000,3],Null]];//Timing

{629.781 Second,Null}

На генерацию таблицы понадобилось почти десять с половиной минут. Заметно, что время генерации таблицы растет нелинейно с ростом таблицы. Наверное, нужно проверить, не пора ли остановиться.

prog//Timing Yes= 11129 ; N0= 88871

{1405.2 Second,Null}

На отбраковку потрачено почти двадцать четыре минуты. Это существенно больше, чем в предыдущем случае. Несмотря на это, дополнительно отбраковать удалось чуть более одного процента чисел. Кажется, пора остановиться. (Конечно, это довольно грубый и притом субъективный критерий останова.)

Теперь можем приступить к написанию программы. Пусть ее заголовок (имя и список параметров) будет М3n[n_]. (Параметр n— верхний конец диапазона, в котором происходит перебор.) Прежде всего нужно решить следующий вопрос: как использовать таблицу так, чтобы не пропустить ни одного простого числа вида 3*2n+1 даже в том случае, если оно использовалось при построении таблицы. Конечно, можно заблаговременно составить список n для таких чисел. Но это значит, что список п для таких чисел нужно составить по другой программе, решающей фактически ту же задачу, «решением которой мы сейчас как раз и заняты. При решении этой вспомогательной задачи, правда, можно ограничиться перебором на меньших начальных отрезках натурального ряда. Конечно, этот способ решения вполне корректен. Но он ведет к программе, которая выглядит несколько искусственно. Фактически получится не программа, а подглядывание в раздел ответов в случае небольших значений п. Ничего плохого в этом нет. Просто этот способ реализовать нельзя, если мы не знаем нужную нам часть ответа для решаемой задачи. Поэтому мы поступим иначе. Мы вообще сделаем вид, что не знаем даже начальных значений п. Сначала найдем в программе все те начальные значения n, при которых числа вида 3*2n+1 не превосходят тех простых чисел q, которые использовались для построения таблицы. Хотя это и будет чистый перебор, перебирать придется недолго, поскольку числа вида 3*2n +1 с ростом и растут гораздо быстрее, чем n-е простое число ра. Например, если мы будем строить таблицу, используя первые 25 000 простых чисел, то перебирать без дополнительной отбраковки придется лишь три-четыре десятка значений п. Затем можно будет построить таблицу и запустить процесс предварительной отбраковки. Поскольку проверка малых значений выполняется весьма быстро, о времени ее выполнения можно не беспокоиться. Конечно, вполне законным является и вопрос: не слишком ли рано мы запускаем отбраковку? Но мы видели, что даже отбраковка 100 000 чисел с помощью огромной таблицы, построенной с использованием первых 25 000 простых чисел, занимает не более 24 минут. Так что на предварительную отбраковку одного числа требуется совсем немного времени. Ну а выигрыш отбраковка дает уже при весьма малых значениях п. Поэтому не стоит беспокоиться из-за того, что на некотором, весьма малом отрезке отбраковка может оказаться менее эффективной, чем прямая проверка. Для таблицы, построенной с использованием первых 25 000 простых чисел, перерасход времени, например, не превысит десятка-двух секунд. Так что беспокойство по поводу того, что временные затраты на отбраковку малых значений n могут быть больше, чем на прямую проверку, необоснованно. Но перед запуском отбраковки нужно сгенерировать таблицу t. Даже если таблица строится с использованием первых 10 000 простых чисел, генерация таблицы занимает немного более полторы минуты. Кроме того, в данной программе есть еще один источник непредсказуемых временных затрат. Это показатели n, кратные 4. Давайте сначала посмотрим, что это будут за показатели. Вот нужная нам программа.
prog401 : = Block[{nn=8000,Yes=0,No=0,j=4), While[j<=nn, If[fx[j,t],Print[j,","];Yes++,No++]; j +=4]; Print["Yes=",Yes,";   No=",No]]
В пределах первых восьми тысяч будет найдено 401 число, кратное 4, такое, что его не удастся отсеять с помощью таблицы, для генерации которой использовались 25 000 простых чисел. (Отсеяно будет 1599 чисел.) Далее приведен список тех чисел, которые кратны 4 и выдерживают отсев.

Если показатель я равен одному из этих чисел, то наименьший делитель числа 3*2n+ 1, больший единицы, больше наибольшего простого модуля, использованного для составления таблицы. Вот примеры.
Factor-Integer [3*2АЗб+1] {(206158430209,1}} Factorlnteger[3*2Л92+1] {{1132314641089,1},{13119392730926401,1}} Factorlnteger[3*2л96+1] {{392840481939253,1},{605040718740253,1}} Factorlnteger [3*2лЮ8+1] {{805213,1},{1781311693519,1},{678750386188027,1}}
Таким образом, нетривиальные делители чисел 3*2n +1, показатели которых выдержали отсев, довольно велики. (Ведь если число 3*2n +1 делится на простой модуль, использованный при построении таблицы, то показатель не выдержал отсев.) Поэтому при применении функции PrimeQ могут произойти временные задержки. Давайте сосчитаем количество таких показателей в пределах первых 300 000. Вот нужная нам программа.
prog4:= Block[{nn=30000,Yes=0,No=0, j=4}, While[j<=nn,
If[fx[j,t],Yes++,No++]; j+=4]; Print["Yes=",Yes,"; No=",No]]
Вот результаты счета.

prog4//Timing Yes= 14451 ; N0= 60549

{1893.61 Second,Null}

Как видите, в пределах первых 300000 придется проверить только 14451 число, кратное 4. Правда, отсев остальных чисел при этом у меня занял примерно 1893 секунды, т.е. чуть больше получаса.

Хотя, как и в случае чисел 5*2n +1, время на построение таблицы окупается, поведение программы может показаться несколько неожиданным из-за временных затрат на построение таблицы и возможных временных затрат на вычисление функции PrimeQ в случае показателей, кратных 4. Действительно, сначала очень быстро печатаются первые значения п, а затем после n = 534 (на моем компьютере) происходит весьма существенная задержка. Она связана как с отсутствием искомых значений и вплоть до п = 2 208, так и с тем, что для части чисел (чуть более 11 %) приходится выполнять полную проверку, которая теперь уже требует весьма ощутимых временных затрат, причем зачастую весьма непредсказуемых для показателей п, кратных 4. Затем выдача замедляется, но поведение программы выглядит вполне стабильно, пока не будет напечатан показатель n = 3 912. После этого идет большой перерыв, вплоть до w = 20 909. В приведенной ниже программе учтено все, о чем говорилось выше.
М3n[n_]:= Block[{j=l,   nnprime=25000,primen=Prime[nnprime],nn=Min[n,primen]}, While   [MMkn[3,j]<= primen   &&   j<=n,    {If[MM3nPrime01[j],Print[j]], J++H; If[j]<n, t=Sort[DeleteCases[Array[f3,nnprime, 3] , Null]] ;  While[j]=nn, If[fx[j,t],   If[MM3nPrime01[j],Print[j]]] ; j++]]]
Эта программа существенно быстрее первоначальной, но если хотите по ней найти n = 20 909, понадобится время и ангельское терпение.

Как видим, предварительная отбраковка и в данном случае значительно повышает эффективность программы. Впрочем, эффективность программы все еще сильно определяется неотсеянными показателями n, кратными 4. Хотелось бы и для них найти такое основание а, чтобы для установления простоты числа 3*2n +1 достаточно было проверить выполнение сравнения ат k =-1 (mod /a2n+l). (Ранее для показателей n, не кратных 4, мы полагали а = 5.)

Но все дело в том, что если показатель п делится на 4, то выбрать основание а несколько сложнее. Если показатель п делится на 4, в качестве основания а можно взять любой квадратичный невычет г по модулю 3 * 2n +1. Поэтому если п делится на 4, то сначала нужно найти (желательно небольшой) квадратичный невычет t по модулю 3*2n+1, а затем проверить, выполняется ли сравнение n2n-1=-1(mod k*2n+l). Но как найти квадратичный вычет? Для этого можно использовать квадратичный закон взаимности.

Если Р>1 n Q>1 — положительные нечетные взаимно простые числа, то







Начало  Назад  Вперед


Книжный магазин