Симметричные кортежи из последовательных простых чисел

[Yadryara]

Программа, по которой это посчитано:

(PARI)

{print();
t0=getwalltime();
fpr = 71; n = vecprod(primes([2,fpr]));
pp=vector(97); Li=vector(4); summc=vector(4);
Li[2] = 1.0451637801174927848445888891946131365;
Li[4] = 1246.1372158993884596927711075290597925;
spis = [ 1, 2, 3, 5, 6, 7,10,11,13,14,15,17,19,21,22,23,26,29,30,31,33,34,35,37,38,39,41,42,43,46,47,51,53,55,57,58,59,61,62,65,66,67,69,70,71,73,74,77,78,79,82,83,85,86,87,89,91,93,94,95 ];
mu   = [ 1,-1,-1,-1, 1,-1, 1,-1,-1, 1, 1,-1,-1, 1, 1,-1, 1,-1,-1,-1, 1, 1, 1,-1, 1, 1,-1,-1,-1, 1,-1, 1,-1, 1, 1, 1,-1,-1, 1, 1,-1,-1, 1,-1,-1,-1, 1, 1,-1,-1, 1,-1, 1, 1, 1,-1, 1, 1, 1, 1 ];
print(fpr,"#");
print();
for(i=1, #spis,
po=spis[i];
x=n^(1/po);
if (x < 2, next,
if ( x >= 10^4,
c1 = Li[4] + intnum(t = 10^4, x, 1/log(t)),
c1 = Li[2] + intnum(t = 2, x, 1/log(t)));
c3 = -log(2);
c4 = intnum(t = x, oo, 1/t/(t^2-1)/log(t));
mc1 = mu[i]/po*c1; summc[1]+=mc1;
mc3 = mu[i]/po*c3; summc[3]+=mc3;
mc4 = mu[i]/po*c4; summc[4]+=mc4;
printf("%u                     %0.5f                            %0.5f       %0.5f\n",po,mc1,mc3,mc4);
));
print();
printf("%0.5f     %0.5f      %0.5f       %0.5f\n",summc[1],pp[fpr]-(summc[1]+summc[3]+summc[4]),summc[3],summc[4]);
print();
print(summc[1]+summc[3]+summc[4]);
print();
}quit;

Поскольку нами уже было замечено, что PARI весьма неточно считает интегралы для маленьких границ, пришлось вплоть до 1e4 прибегнуть к ухищрениям, в том числе считать на Вольфрам Альфе. Здесь Li[2] это и есть Li(2), но вот Li[4] это Li(1e4). При счёте от 1e4 уже было совпадение всех 38 знаков.

И все знаки с книжными совпали:

https://www.phantastike.com/math/prime_obsession/djvu/view/

Страница примерно 400-я. Примерно — потому что разная нумерация.

[Dmitriy40]

Э, не забываем что интегральный логарифм завышает количество простых, как раз где-то на 330 миллиардов около 5e26.

Насчёт формулы для паттернов я сомневаюсь - тогда бы получив её для паттерна 0,2 доказали бы бесконечность близнецов - но может из формулы не получить асимптотику и она полезна лишь для конечных интервалов ... Пробуйте, очень интересно посмотреть. Но вообще это задача для Дт(М) (даже не ПРР).

Книжку я уже утром скачал и стал читать, правда подряд, хоть пока почти всё и известно.

-- 04.04.2025, 12:33 --

Кстати из кода в OEIS узнал про real(-eint1(-log(n))=li(n) в PARI. Чтобы не выписывать сам интеграл intnum().

[Yadryara]

Dmitriy40 в сообщении #1681034 писал(а):

не забываем что интегральный логарифм завышает количество простых, как раз где-то на 330 миллиардов около 5e26.

Не понял к чему это сказано. Это ведь как раз и неизвестно насколько он завышает, потому что неизвестно точное количество простых от нуля до 5e26.

Dmitriy40 в сообщении #1681034 писал(а):

полезна лишь для конечных интервалов ... Пробуйте, очень интересно посмотреть. Но вообще это задача для Дт(М) (даже не ПРР).

Пока что речь именно о конечных интервалах. Всё-таки склоняюсь к ПР/Р.

Dmitriy40 в сообщении #1681034 писал(а):

Книжку я уже утром скачал и стал читать, правда подряд, хоть пока почти всё и известно.

Надо же, особо и не сомневался, что Вы её читали раньше. Тогда правильно, надо читать подряд, хотя разжевывание порой избыточное.

[Dmitriy40]

Yadryara в сообщении #1681045 писал(а):

Не понял к чему это сказано. Это ведь как раз и неизвестно насколько он завышает, потому что неизвестно точное количество простых от нуля до 5e26.

Вообще-то таблицы $\pi()$ есть до $10^{29}$. И насколько понимаю посчитаны они примерно как Вы и собираетесь, по неким формулам (я просто не разбирался), не интегральным логарифмом (или не только им). Соответственно и величины завышения тоже есть до $10^{29}$. Вот как пример какие алгоритмы использует primecount:

  -d, --deleglise-rivat  Count primes using the Deleglise-Rivat algorithm
  -g, --gourdon          Count primes using Xavier Gourdon's algorithm.
                         This is the default algorithm.
  -l, --legendre         Count primes using Legendre's formula
      --lehmer           Count primes using Lehmer's formula
      --lmo              Count primes using Lagarias-Miller-Odlyzko
  -m, --meissel          Count primes using Meissel's formula
      --Li               Approximate pi(x) using the logarithmic integral
      --Li-inverse       Approximate the nth prime using Li^-1(x)
  -n, --nth-prime        Calculate the nth prime
  -p, --primesieve       Count primes using the sieve of Eratosthenes
      --phi <X> <A>      phi(x, a) counts the numbers <= x that are not
                         divisible by any of the first a primes
      --Ri               Approximate pi(x) using Riemann R
      --Ri-inverse       Approximate the nth prime using Ri^-1(x)

Алгоритмом по умолчанию считает долго (даже многопоточно), $\pi(10^{20})$ порядка минут, другие алгоритмы не проверял.

Yadryara в сообщении #1681045 писал(а):

Надо же, особо и не сомневался, что Вы её читали раньше.

Похоже нет, я же даже близко не математик, да и проблематика вычисления $\pi()$ (и тем более гипотеза Римана) мне недоступна.

[Yadryara]

Dmitriy40 в сообщении #1681056 писал(а):

Вообще-то таблицы $\pi()$ есть до $10^{29}$.

Ну да, и в той же моей ссылке на Вики есть от нуля до всех степеней 10-ки до 29 включительно. Можно и в OEIS заглянуть. Если б gris заглянул, да повнимательнее, то и не стал бы писать, что $\pi(10^{29})$ нашли в этом году. Ведь чётко сказано что в начале 2022-го.

Но интересует-то не 1e29, а 5e26 и 6e26 в качестве соседей, а на самом деле 71#.

Dmitriy40 в сообщении #1681056 писал(а):

Вот как пример какие алгоритмы использует primecount:

И среди них нет того самого, не комбинаторного, а аналитического FKBJ17. (Давал ссылку на него в теме в ПР/Р)

Я хочу понять как считать по нему потому что понимаю как считать по формуле Римана, а он на ней основан.

Где брать нули вроде знаю.

[Dmitriy40]

Yadryara в сообщении #1681070 писал(а):

Но интересует-то не 1e29, а 5e26 и 6e26 в качестве соседей, а на самом деле 71#.

Интерполяцией. Кроме $\pi(10^n)$ есть таблицы $\pi(2^n)$, они покрывают и $71\#$ ($2^{89}$ всего на 11% больше $71\#$). А точность интерполяции оценить по известным значениям для $53\#, 59\#, 61\#, 67\#, 10^{20\ldots27}$. Это если задаться целью посчитать именно приближённое значение $\pi(71\#)$. Впрочем я так понимаю задача другая, эта лишь как тренировка.

-- 04.04.2025, 16:20 --

И кстати.

Yadryara в сообщении #1681025 писал(а):

Поскольку нами уже было замечено, что PARI весьма неточно считает интегралы для маленьких границ,

Это не так. Это не PARI криво считает интегралы, а лишь аппроксимация (асимптота, в смысле не прямой, а функции) по HL1 расходится с реальными значениями на начальном участке, что вполне нормально для асимптоты. Интегралы при этом правильные и PARI не виноват.

[Yadryara]

Ну вот стоило не написать ещё раз слово "точное" и пошло поехало...

Dmitriy40 в сообщении #1681034 писал(а):

не забываем что интегральный логарифм завышает количество простых, как раз где-то на 330 миллиардов около 5e26.

Мы не знаем (пока) точное значение $\pi(5\cdot10^{26})$. Стало быть не знаем насколько завышает.

Для $\pi(10^{6})$ — знаем: 78498. И знаем превышение:

$Li(10^{6})-\pi(10^{6})= 129.549159462181919862910747947261161...$

Тренироваться лучше на меньших числах. Только пока непонятно как.

[Yadryara]

Да, забыл сказать — вчера был наконец-то найден кортеж 15-228-2 и в низине. И список самых маленьких таких кортежей наконец-то разбавился нашим результатом:

Номер     Стартовое число       LGap RGap     Найден

1.    2079914861571286679         46   12       TBEG
2.    3665619319531504883         10    6       TBEG
----------------------------------------------------- 53#
3.    271541128585758431779       68    6      An+Dm
4.    1006882292528806742273       6    6      Jarek
----------------------------------------------------- 59#

Мною проверено почти 6% интервала $0-61\#$. И вот сравнение для общих количеств 15-228-2 :

Интервал         0-53#   0-59#   0-61#
Прогноз по HL1       2      49    1133
Найдено              2       4      89

Текущая статистика:

                 Посчитано      Найдено        Цепочек       Соотношение
                  в группе      цепочек        на юнит             между
                                 тысячи                         группами
1-й период 61# G19   100 %           98           2242            1.1112
1-й период 61# G20   100 %         2645           2025            1.1157
1-й период 61# G21   100 %        23366           1815            
1-й период 61# G22    25 %        22648           1598            1.1358
...
1-й период 61# G28   100 %        11231            813            1.1216
1-й период 61# G29   100 %         1566            731            1.1204
1-й период 61# G30   100 %           73            656            1.1198

Как видим кэф действительно растёт к центру — 1.1358 явно больше всех остальных в этом же столбце.

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

[Yadryara]

Кстати, решил сделать аналогичные таблички сравнения прогноза с фактом для всех трёх преемственных паттернов:

15-228-2

Интервал         0-53#   0-59#   0-61#
Обсчитано          35%      6%      6%
Прогноз по HL1       2      49    1133
Найдено              2       4      89


17-240-1

Интервал         0-59#   0-61#   0-67#   0-71#
Обсчитано         100%    0.4%*   0.4%   <0.1% 
Прогноз по HL1       0**    10     212    5273
Найдено              1       6      10      13

* — Jarek предположительно проверил не меньше 16%
** — 0.490


19-252

Интервал         0-67#   0-71#   0-73#
Обсчитано         100%    2.8%   <0.1%
Прогноз по HL1       1*     11     257
Найдено              0       1       1

* — 0.511

Напомню, что прогноз по HL1 для кортежей длиной 21 и больше мы пока толком сделать не смогли.

[Dmitriy40]

Yadryara в сообщении #1681025 писал(а):

Программа, по которой это посчитано:

Чуть оптимизировал:

s=0.0; n=1e6;\\Аргумент
for(i=1,oo,
m=moebius(i)/i; if(m==0, next);
y=sqrtn(n,i); if(y<2, break);
a=real(-eint1(-log(y)));\\Li(y)
b=0;\\Не сделано
c=-log(2);
d=intnum(t=y,oo, 1/t/(t^2-1)/log(t));
s+=m*(a+b+c+d);
printf("%3u:  %15.5f  %15.5f  %15.5f  %15.5f  = %15.5f\n",i,m*a,m*b,m*c,m*d,m*(a+b+c+d));
);
print(s);

Но самого интересного - вычисления второго члена - так и нет.

Если же говорить про паттерны, то в первом члене вместо Li(x) будет интеграл с теми константами и обратными степенями логарифма что и в HL1. Третий и четвёртый члены вообще без понятия нужны ли будут и чему будут равны. А вот что делать со вторым членом не очень ясно, то ли простая замена Li() на такой же интеграл с константами и обратными степенями логарифма, то ли нули дзета функции Римана вообще не в кассу (ведь даже кратных их вроде не обнаружено) и надо думать.
С другой стороны, первый член - это HL1 в чистом виде, остальные лишь малая поправка к ней (причём порядка квадратного корня на некую константу), по крайней мере для больших значений первого члена. Может на них тогда и плюнуть? Но тогда остается просто HL1.

-- 05.04.2025, 12:54 --

А, кстати, для паттернов достаточно лишь первой строки, все остальные должны быть нулевыми (конечно если нас не интересуют тысячи и миллионы кортежей) - ведь они берутся по корням из аргумента, а в такой низине мы точно знаем что кортежей нет и соответственно j()=0.
Вот и получается что всё сводится фактически к HL1 (первый член в j()) с малой поправкой к ней (остальные члены в j()). Порядок поправки неизвестен, но можно предположить что он будет ещё на порядки меньше чем для простых чисел, не корень из аргумента, а очень примерно в len-ую (длина паттерна) степень меньше. И пока нас интересуют малое количество кортежей на поправку можно и плюнуть. И остаётся HL1 в чистом виде.

-- 05.04.2025, 12:59 --

Э, собственно если мы точно знаем что ниже корня из аргумента кортежей нет, то $j(x)=\pi(x)$ (остальные члены в сумме зануляются) и соответственно никаких поправок к HL1 быть и не должно ... Кажется.

[Yadryara]

Сравнил две итоговых суммы для $10^6$, Вашу и мою :

78527.37378068578297055960454091 2354990
78527.37378068578297055960454091 6645019

Вроде пока некритично. И непонятно у кого точнее.

Dmitriy40 в сообщении #1681170 писал(а):

Но самого интересного - вычисления второго члена - так и нет.

Почитайте, там в районе 400-й страницы подробно рассказывается как считать по Риману. Я лет 12 назад делал на Бейсике.

Но сходится жутко медленно. А вдруг Вы придумаете как ускорить даже не читая никаких других работ.

Сам я перевёл Платта, всё-таки у него поменьше — 15 страниц, разбираюсь.

Dmitriy40 в сообщении #1681170 писал(а):

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

Я вообще пока настраиваюсь не плевать ни на что. Вот Вы почти полтора года в одиночку искали 19-252 в первом периоде 67#. А что если можно было бы посчитать аналитически и сказать, что кортежей в первом периоде 67# ровно 0, а во втором периоде 67# — ровно один. Вот это был бы высший класс. И за эти флуктуации, кстати, как раз нули дзета-функции и отвечают.

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

[Dmitriy40]

Yadryara в сообщении #1681182 писал(а):

И непонятно у кого точнее.

Да, я видел что eint вычисляется с погрешностью от intnum. Но думаю что встроенные функции точнее численного интегрирования. Можно сравнить с вольфрамом, например Li(2) в PARI:
? -eint1(-log(2))
%1 = 1.0451637801174927848445888891946131365 + 3.1415926535897932384626433832795028842*I
Вольфрам: 1.04516378011749278484458888919461313652261557815120157583290914407501320521...
Действительная часть полностью совпадает. Так что такому вычислению Li() верить можно.

Кстати primecount.exe --Ri моментально даёт именно эти же значения $\pi(x)$ (из таблицы и программы) без второго члена.

Стал разбираться как посчитать второй член, наткнулся на непонятку: $20^{\frac{1}{2}+14.134735i}=-0.302203-4.46191i$ на PARI примерно одинаково:
? 20^(1/2+14.134725*I)
%1 = -0.30230466031723746206585563990422198947 - 4.4619067552281367317176433520788076513*I
А вот $\operatorname{Li}(20^{\frac{1}{2}+14.134735i})=-0.105384+3.14749i$ в PARI даёт совсем другое значение:
? -eint1(-log(20^(1/2+14.1347251417346937904572519835625*I)))\\Li(20^(1/2+14.134735i))
%2 = 1.9979692368474822688136820603245564820 - 0.77224977227477630516373401329504663721*I
И вольфрам даёт это же значение, не то что в книге! Беда!

[Yadryara]

Dmitriy40 в сообщении #1681192 писал(а):

Действительная часть полностью совпадает. Так что такому вычислению Li() верить можно.

Отлично. Так не просто можно верить: мы же кучу времени потратили на то чтобы выбрать стартовую точку интегрирования для счёта по HL1 ! Ещё с vicvolfом спорили по этому поводу. Так небось можно было просто считать этой функцией прямо от двойки и не париться!

[Dmitriy40]

Yadryara в сообщении #1681196 писал(а):

Так небось можно было просто считать этой функцией прямо от двойки и не париться!

Нельзя: li(x) это интеграл от 1/ln(x), а нам нужно от того же в степени - и как её подставить в -eint1(-log(x)) непонятно. И считает она от нуля, не от двойки, это именно li(), не Li(), правда скорректировать легко.
Вот если бы по HL1 считали бы $\pi(x)$ как паттерна [0] - вот тогда степень логарифма была бы первой и тогда да, было бы можно.

Кстати, в PARI это выходит не Li(), а li(). А у меня в комментариях выше ошибка.
В книге же - каша (путаница), кажется тоже везде вместо Li() надо писать li() (раз уж совпадает с PARI), но не совсем уверен что прямо везде. Где-то в начале про это было примечание, но уже забыл. Нашёл, на 146 странице примечание что везде пишется Li(), но смысл у него li(), т.е. от нуля.

-- 05.04.2025, 17:54 --

Кажется $\int\frac{1}{\ln^k t}dt$ вполне берётся по частям ...

[Yadryara]

Dmitriy40 в сообщении #1681202 писал(а):

В книге же - каша (путаница), кажется тоже везде вместо Li() надо писать li()

В Альфе это одно и то же.

Dmitriy40 в сообщении #1681202 писал(а):

Нашёл, на 146 странице примечание что везде пишется Li(), но смысл у него li(), т.е. от нуля.

Так здесь же примечание есть, на 397-й (на 400-й бумажной)! То есть в Альфе надо писать Ei[(1/2+14.134725i) log[20]] и всё сходится.