книги / Mathematica 5. ╨б╨░╨╝╨╛╤Г╤З╨╕╤В╨╡╨╗╤М
.pdf23564163
241077871
252063689
263957809
277603553
2814630843
2928192750
3054400028
31105097565
32203280221
33393615806
34762939111
351480206279
362874398515
375586502348
3810866266172
3921151907950
4041203088796
4180316571436
42156661034233
43305761713237
44597116381732
451166746786182
462280998753949
474461632979717
PrimePi::largp |
Argument 281474976710656 |
||
|
in PrimePi[281474976710656] |
is too |
|
48 |
large for this implementation. More... |
||
PrimePi[281474976710656] |
|
||
PrimePi::largp |
Argument 562949953421312 |
||
|
in PrimePi[562949953421312] |
is too |
|
49 |
large for this implementation. More... |
||
PrimePi[562949953421312] |
|
||
PrimePi::largp |
Argument 1125899906842624 |
||
|
in PrimePi[1125899906842624] |
is too |
|
|
large for this implementation. More... |
||
General::stop |
|
|
|
Further output of PrimePi::largp will be |
|||
|
suppressed |
during this calculation. More... |
50PrimePi[1125899906842624]
51PrimePi[2251799813685248]
Тем не менее не следует воспринимать ситуацию уж очень пессимистически.
£
Давайте, например, построим графики функций и(дг) и — на промежутке (0, 100).
1пх
P lo t f f P r im e P i[х] , |
----------- 1 , { х , 1 , 1 0 0 } , |
||
|
|
|
L o g [x ] ‘ |
P lotR an g e |
-> |
( 0 , |
2 6 } , |
P lo t s t y le |
-> |
|
|
{{}, ( D a s h in g [ { 0 .0 2 , 0 . 0 2 } ] } } ] ;
Арифметика: простые числа |
161 |
видно, что она довольно неплохо приближает TT(JC) . Впрочем, вот график с интеграль
ным логарифмом.
Plot[{PrimePi[х],Loglntegral[х]},{х,1,100},PlotRange->{0,26},
PlotStyle->{{},{Dashing[{0.02,0.02}]}}];
Теперь даже в нуле неприятностей нет. Но вот вопрос: будет ли заметна ступенча тость на больших интервалах? На самом деле ступенчатость незаметна уже на интер вале (0, 10000).
P lo t [ { P r im e P i[ х ] , |
Xr |
} , {х, 1 , 1 0 0 0 0 ), |
|
|
|
L o g [x ] |
J |
P lo t R a n g e -> |
{0, |
1 20 0}, |
|
P l o t s t y l e -> |
|
|
|
{{}, { D a s h in g [ { 0 .0 2 , 0 . 0 2 } ) } } ] ;
162 |
Гпава 5 |
1 2 0 0
1000
800
600
400
200
2000 |
4000 |
6000 |
8000 |
10000 |
х
Но какая же функция приближает TC(JC) лучше: — или Li(x)? Конечно, Li(x).
InJC
Plot [{PrimePi[x],Loglntegral[x]},{x, 1,10000},
PlotRange->{0, 1200}, PlotStyle->{{},{Dashing[{0.02,0.02}]}}];
Если же построить график на интервале (0, 10000), то графики тг(;с) и Li(jc) прак тически сольются.
Plot[{PrimePi[х] ,Loglntegral[х]},{х,1,50000},
PlotRange->{0, 6000}, Plotstyle->{{},{Dashing[{0.02,0.02}]}}];
Арифметика: простые числа |
163 |
6 0 0 0
Раз уж речь зашла о приближении TT(JC) , давайте построим графики разности
функции тг(дг) и известных ее аппроксимаций. В качестве первой (самой простой) ап-
проксимирующей функции выберем — (приближение Чебышева и Лежандра), в ка- lnJC
честве второй — Li(jc) (приближение Гаусса), а в качестве третьей — R(x) (приближение Римана). Сначала определим функцию R(x).
R iem a n n [ х _ ] :=
L a s t [ B l o c k [ { R = L o g I n t e g r a l[ x ] , y = L o g I n t e g r a l [ S q r t [ x ] ] / 2 , n = 2 },
{ W h ile [ y > 0 . 0 0 0 0 0 0 0 0 1 , {R=R-y; n = n + l ; y = L o g I n t e g r a l [ x A ( 1 / n ) ] /n }
],R }
]]
Теперь можем построить нужные нам графики.
P l o t [ { P rim e P i [х] -R iem ann [ х ] , P r i m e P i[ х ] - х / (L o g [ х ] - 1 ) ,
L o g l n t e g r a l [ х ] - P r i m e P i [ х ] },
{ х , 1 , 1 0 л7 } , P lo t R a n g e - > { - 1 0 0 , 3 0 0 } ,
P l o t S t y l e - > { { } , { D a s h in g [ { 0 . 0 2 , 0 . 0 2 } ] } } ] ;
Как видим, из трех наилучшим является приближение Римана. |
|
|||
И ещ е одно |
замечание. Как |
мы видели, график функции |
TT(JC) |
выглядит вполне |
гладко, хотя на |
самом деле он |
является ступенчатым. Более |
того, |
существуют сколь |
164 |
Глава 5 |
угодно длинные интервалы, на которых он горизонтален. Однако длина таких интер валов незначительна по сравнению с их расстоянием до начала координат. Это интер валы, на которых функция TT(JC) не изменяется. Иными словами, это интервалы, на
которых нет простых чисел. Кстати, а как сосчитать количество простых чисел на ин тервале (а , Ь]1
К о л и ч ест в о п д о с т ы х ч и с е л н а о т к р ы т о м с л е в а о т р е з к е (а , Ь]
С помощью функции PrimePi несложно подсчитать и количество простых чисел на открытом слева отрезке (а , Ь\. Помните только, что если вы пользуетесь выражени ем к(Ь)-п(а), т.е. выражением PrimePisJb] -PrimePi [а], то в случае простоты простое
число Ь будет посчитано, а простое число а — нет. Иными словами, подсчет осущест вляется на открытом слева отрезке (а , Ь]. А что, если нужно посчитать простые числа на замкнутом с обоих концов отрезке [а, Ь]1 Ничего сложного: в качестве аргумента можно взять не а , а несколько меньшее число, ведь аргументом функции PrimePi
может быть любое вещественное положительное число. Давайте посчитаем, например, количество простых чисел в 1-й, 2-й, , 20-й сотне миллионов натуральных чисел. Вот как это можно сделать.
delta=l00000000; PrimePiAB[delta_Integer?Positive,xk_Integer?Positive]:= Block[{k=0, x=delta, kt=0 },
While[x<=xk,
{kt=PrimePi[x]; Print[x,":"rktf":",kt-k]; k=kt ;
x=x+delta}]]
PrimePiAB[delta,10 delta]
Полученные результаты удобно представить в виде таблицы (табл. Б.23).
В приведенной выше программе мы воспользовались тем, что интервалы примы кают друг к другу. Однако так бывает не всегда. Иногда нужно подсчитать количество простых чисел в интервалах заданной длины, причем начала интервалов расположены на числовой оси так, что конец очередного интервала не совпадает с началом сле дующего. Подсчитаем, например, количество простых чисел в интервалах длиной
150 000, причем пусть начинаются эти интервалы в точках х = 106, 107, |
1014 |
Для этого случая пригодится следующая функция.
beltaPi [х_, delta_J :=
Block [{xk=x+delta, k=PrimePi [xk] },
Pr i n t [ x , ,xk,":",к-PrimePi[x]]]
Спомощью этой функции нетрудно написать и нужную нам программу.
delta=150000;
Do[DeltaPi [10Лп, delta] ,{n, 6,14 }]
Результаты отформатированы в виде таблицы (табл. Б.24).
Резюме
Хотя каждый математик знает, ох, как это непросто, нужно признать, что система Mathematica дружит с простыми числами. В ней есть функции, проверяющие, являет ся ли число простым; функция, генерирующая л-е простое число; функция, генери рующая следующее простое число; функция, генерирующая предыдущее простое число;
Арифметика: простые числа |
165 |
а также знаменитая функция л(х), подсчитывающая количество простых чисел, не превосходящих заданного числа х . Конечно, реализовать все эти функции весьма не просто, и применяя их, нужно учитывать ограничения на их аргументы. Но как бы то ни было, в системе Mathematica они реализованы на современном уровне, с учетом совсем недавно доказанных теорем и новых методов.
Задачи /
Задача 5.1. Среди квадратов простых чисел нет. А вот часто ли среди чисел вида п2 +1 встречаются простые? Найдите все такие натуральные л, не превышающие 7V, что
п2 +1 — простое. Иными словами, нужно составить список всех натуральных n£N , для которых число п2 4-1 является простым. Рассмотрите случай N = 1 000.
Задача 5.2. Как вы видели, для многих п из первой тысячи числа вида п2 + 1 явля ются простыми. А как подсчитать, сколько таких п в первом миллионе?
Задача 5.3 (числа Каллена). Числа вида п2л + \ называются числами Каллена. Долгое время было известно только два значения л, для которых соответствующие числа Каллена были простыми: п = 1 и п = 141. Найдите все л<10 000, для которых соответствующее число Каллена является простым.
Задача 5.4. Для каждого ли к найдется такое л, что число /:2я+1 будет простым? Хотя есть много примеров, подтверждающих это высказывание, оно ложно! Оказыва ется, существуют такие к , что число /:2я +1 является составным при всех натуральных п. Составьте программу, которая для каждого £<1000 находит все такие л<3000, что число /:2я+ 1 является простым.
Задача 5.5 (простые числа, в двоичной записи которых ровно три единицы). Для каж дого ли т найдется такое простое двоичное т -значное число, в двоичной записи кото рого ровно три единицы? Иными словами, для каждого ли т найдется такое и, \<п<т,
что число 2т+ 2я +1 будет простым? (Составных чисел такого вида, как легко сообра зить, бесконечно много: например, составными являются все числа вида 22я+ 2я+| + 1 = (2я+1)2.) Составьте программу, которая для каждого /и<1000 находит
все такие я, \<п<т , что число 2т+ 2Я+1 является простым.
Задача 5.6. Для каких т<2 200 не существует простого двоичного m-значного чис ла, в двоичной записи которого ровно три единицы? Иными словами, для каких т
при всех п (1<п<т) числа 2т + 2яН-1 являются составными?
Задача 5.7. Постройте таблицу простых чисел, близких к степеням двойки.
Задача 5.8. Постройте таблицу простых чисел, близких к степеням тройки, однако в таблице укажите не сами простые числа, а разности между соответствующей степе нью тройки и самим простым числом.
166 |
Гпава 5 |
Глава 6
Арифметика:
наибольший общий делитель и наименьшее
общее кратное
Вэтой главе...
♦Наибольший общий делитель
♦Наименьшее общее кратное — функция LCM
♦Резюме
♦Задачи
Во множестве всех общих делителей данных чисел имеется такой, который делится на всякий другой общий делитель этих чисел: это — общий наибольший делитель данных чисел.
Л. К. Сушкевич. Теория чисел. Глава 1. О делимости чисел. Теорема 9
Имея каноническое разложение чисел на простые множители, легко найти их наи больший общий делитель (НОД) и наименьшее общее кратное (НОК). Однако еще со времен Евклида хорошо известно, что для вычисления наибольшего общего делителя существуют гораздо более эффективные алгоритмы, чем разложение на простые мно жители. И система Mathematica их применяет.
Наибольший общий делитель
Для нахождения наибольшего общего делителя чисел (целых, рациональных или гауссовых) в системе Mathematica предусмотрено две функции: GCD и ExtendedGCD.
Наибольший общий делитель — функция GCD
Функция GCD находит наибольший общий делитель в области целых, рациональ ных и гауссовых чисел.
Наибольший общий делитель в кольце целых чисел
Чтобы найти наибольший общий делитель чисел л,, th, ...» можно использовать
функцию GCD [ л,, п2, ...]. Вот примеры ее применения для нахождения наибольшего
общего делителя двух чисел.
GCD[36,45]
9
GCD[2200 + 3, З300 + 80]
349 а=177Л5+30621*173Л3-173А5
177309584821 Ь=173л5+30621*177л3-177л5 151037867129 с=173А4+30621А2+177л4 2814896923
GCD[a,b]
30637
GCD [а,с]
30637
GCD [Ь,с]
30637
Пример 6.1. Графики функции GCD.
Давайте теперь построим несколько графиков функции GCD. Поскольку это функ ция двух аргументов, построим изображения поверхности z = GCD[lntegerPart [х], IntegerPart [у] ]. Для этого используем функцию Plot3D.
Plot3D[GCD[IntegerPart[x] ,IntegerPart[у]],{x,1,50}, {у, 1,50}, PlotPoints->100, PlotRange->All, Mesh->False] ;
А вот вид вблизи.
Plot3D[GCD[IntegerPart[x ],IntegerPart[y]],{x,1,10},{у,1,10},
PlotPoints->100,PlotRange->AllfMesh->False] ;
168 |
Гпава 6 |
10
Вот пример нахождения наибольшего общего делителя нескольких чисел.
GCD [ F ib o n a c c i[9 4 5 ],F ib o n a c c i[ 9 0 ] ,F ib o n a c c i[450]]
1134903170
Раз уж мы заговорили о числах Фибоначчи, значит, мы не можем обойти случаи, которые традиционно считаются “наихудшими” для алгоритмов нахождения наиболь шего общего делителя.
“Наихудшие” случаи для алгоритмов нахождения наибольшего общего делителя
“Наихудшим” случаем для классического алгоритма Евклида нахождения наи большего общего делителя считается случай соседних чисел Фибоначчи. Чтобы оце нить быстродействие, попытаемся, например, определить, при каком п будет заметна задержка при вычислении GCD[Fibonacci [n+2], Fibonacci [п+1] ]. Возьмем, на пример, следующую программу.
Do[Print[n," GCD[Fibonacci[n+2],Fibonacci[n+1]]],{n, 10000}]
Увы! Алгоритм, реализованный в системе Mathematica, столь эффективен, что даже на весьма посредственном ПК эта программа выполняется без каких бы то ни было задержек. Если же программу модифицировать так:
Do[Print[п,":",GCD[Fibonacci[2An+2],Fibonacci[2Ап+1]]],{п,10000}]
то задержка будет |
весьма |
ощутимой уже при |
п = 25. Но не |
забывайте, |
что |
2^+2 = 33554434, |
а число |
Fibonacci [?Л25+2] |
имеет 7 012 462 |
цифры! Так |
что |
алгоритм, реализованный в системе Mathematica, справляется с классическим “наихудшим” случаем не просто вполне удовлетворительно, а блестяще! Правда, известны алгоритмы, для которых наихудший случай представляют числа и = ал +ап_{
и v= а где (1 + л/2)"-(1->/2)" . Давайте исследуем и этот случай. Нам
2-Jl
понадобятся следующие определения.
a[nj :=FullSimplify [((1+Sqrt[2])An-U-Sqrt[2])An) / (2 Sqrt[2])] u : = a [ n ] + a [ n - l ]
v : = a [ n - l]
Арифметика: наибольший общий делитель... |
169 |
Попытаемся выполнить программу.
D o [P rint[n," : " , GCD[u,v]],{n,100}]
Окажется, что уже при п = 98 выражения “не хотят” упрощаться! Поэтому опреде ление а[п_] лучше изменить.
а [п_] : = (Expand[ (1+Sqrt [2] ) Лп] -Expand[ (1-Sqrt [2] ) лп] ) / (2 Sqrt [2] )
Однако при таком подходе понятно, что значительное время тратится также на раскрытие скобок. Однако даже при и = 70 000 вычисления не занимают слишком много времени. При этом значении п числа и и v являются 26 794-значными, и ос новное время уходит на их вычисление, а не на вычисление их наибольшего общего делителя, который вычисляется в мгновенье ока! Отсюда можем сделать вывод, что функция GCD реализована весьма эффективно, и если при ее использовании возника ют неприятности, то, скорее всего, они связаны с вычислением не самой функции, а ее аргументов!
Пример 6.2. Взаимная простота чисел Ферма. Числа Ферма возрастают довольно быстро и являются взаимно простыми. Поэтому их наибольший общий делитель дол жен быть равен 1. Выясним, до какого номера п система Mathematica сможет в при емлемое время проверить взаимную простоту этих чисел.
Вот нужная нам функция.
ferm at[п_] :=2Л(2Ап)+1
Теперь можно написать и программу.
D o[P rint[n," : " , GCD[fermat[n+1], ferm at[n]] ] , {n,100}]
Задержка начинает ощущаться при п = 28. Однако вычисления заканчиваются (и притом довольно быстро, если учитывать количество цифр в числах) и при п = 29. Напомню, что уже 22-е число Ферма имеет более миллиона (1 262 612) знаков, а 23-е число Ферма является 2 525 223-значным! Шутки ради 24-е число Ферма, являющееся 5 050 446-значным, я перетащил в Word. Текстовый редактор с записью этого числа работает медленнее, чем система Mathematica! Вот уж воистину супервычислитель: решает задачу быстрее, чем иные успевают произнести ее формулировку! 25-е число Ферма является 10 100 891-значным, и Word его не выносит: замедляется работа ме ню, даже курсор движется так медленно, что работать становится практически невоз можно.
Пример 6.3. Взаимная простота чисел Т -1 и 2м-1 при взаимно простых п и т. Как
известно, числа 2я-1 и 2т-1 взаимно просты, если взаимно просты п и т . Выяс ним, до какого номера п (предполагая, что т<п) система Mathematica сможет в при емлемое время проверить взаимную простоту этих чисел.
Вот нужная нам функция. f n [п_] : =2Лп-1
Теперь можно написать и программу.
Do [
Do [
If[GCD[n, m]==1,If[GCD[fn[n], fn[m]]!=l, P r i n t [ n , m ] ]]
/ {m,n-1}] {n,10000}]//Timing
Для первой тысячи чисел все проверки на моем весьма слабеньком ПК выполня ются всего лишь за 8,125 секунды! Правда, для проверки этого утверждения в преде лах первых десяти тысяч чисел понадобится уже 3398,78 секунды. Здесь, очевидна
ПО |
Глава { |