новое решение сравнения 2^n = 3 (mod n)

[Руст]

Это необходимое условие. Дальше я пользуюсь достаточным условием допустимости, которое определяется по уровням.

[maxal]

Руст писал(а):

Далее решение сводится взятию радикала от $3^{T_b}=1\mod \prod_i p_i$, где $p_i|3^{T_b}-1$ допустимые простые делители этого числа.

Здесь и далее видимо подразумевается $b=3$? Или откуда берется тройка?

Что есть "радикал от $3^{T_b}=1\mod \prod_i p_i$"? Имеется ли в виду радикал $3^{T_b}-1$, за исключением уже известных допустимых простых, или что-то иное?

[Руст]

Да далее более конкретная ситуация рассматривается b=3, хотя метод от конкретного значения не зависит.
Нет, именно так. Это похоже вычислению именно радикала. Например для вычисления квадратного корня от a (являющегося квадратичным вычетом) по модулю p=3(mod 4) возводим в (р+1)/4 ую степень, т.е. $\sqrt a (mod p) =a^{(p+1)/4}\mod p.$

[dmd]

Руст писал(а):

Алгоритм поиска решений

$$(1) \ \ a^n=b\mod n.$$

Здесь а натуральное число (случай а=1 тривиальный и можно не рассматривать), b - целое.
Обозначим периоды числа а и b по модулю простого числа р через T(a,p),T(b,p). Для полного решения нужны так же $T(a,p^k),T(b,p^k)$. При этом простое число р делит некоторое решение n тогда и только тогда, когда $T(b,p)|T(a,p)$. Для степеней появляется допольнительное условие $p\not |T(b,p^k)$. Таких чисел очень мало, поэтому в дальнейшем рассмотрим только случаи, когда n свободно от квадратов.
Разделим простые числа на допустимые. Число р назовём допустимым нулевого порядка, если Т(а,р)=Т(b,p). Число р назовём допустимым k порядка (k>0), если T(b,p) делит T(a,p) и частное $\frac{T(a,p)}{T(b,p)}=p_1p_2...p_m$ разлагается на произведение различных (мы исключили квадратосодержащиу решения) допустимых простых чисел уровня меньше k и хотя бы один из них имеет уровень k-1. Это соответствует тому, что при рассмотрении по модулю p получается сравнение $n=x_p\mod T(a,p), \ (x_p,T(a,p))=p_1p_2...p_m$. Соответственно каждому допустимому простому числу высокого порядка соответствует целое дерево подчинённых ему простых чисел, состоящее из простых чисел в разложении T(a,p)/T(b,p) и их подчинённых.
Для каждого решения исходного сравнения (1) имеется общий модуль $T_b=lcm(T(b,p)| \ p|n)$ по всем простым делителям n. Для эффективного поиска решений будем искать решения пробегая по удобным модулям $T_b$ например $T_b=2^k3^l5^m7^l$ и рассмотрим мультипликативную группу G чисел по модулю $T_b$ с числом элементов $N=\phi(T_b)$. Эта группа есть прямая сумма (в аддитивной записи) циклических групп

$$Z_{2^{k-2}},Z_2,Z_{2*3^{l-1}},Z_{4*5^{m-1}},Z_{6*7^{l-1}}$$

(понятно каких для более общей группы). Далее решение сводится взятию радикала от $3^{T_b}=1\mod \prod_i p_i$, где $p_i|3^{T_b}-1$ допустимые простые делители этого числа. Определим для каждого допустимого простого числа p его координату в группе G в аддитивной записи, предварительно фиксировав некоторые образующие в примарных компонентах $T_b$. Например, если $q^k|T_b, (q,r_q)=1$, где $r_q=\frac{T_b}{q^k}$ и g мультипликативная образующая по модулю $q^k$ то координата $z_q$ допустимого простого числа p определяется по следующей схеме. Пусть $c_p=\frac{T(a,p)}{T(b,p)}$ (из определения допустимых простых $(c_p,T_b)=1$), тогда $g^{z_q}=y\mod q^k , \ 2^{pc_p}=3^y\mod p.$
Тогда произведение множества допустимых простых есть решение тогда и только тогда, когда вместо с каждым простым числом р это множество содержит так же все подчинённые и сумма координат равна 0 в каждой компоненте (сумма считается по модулю длины цикла).
Для эффективного поиска таких допустимых множеств упорядочим подгруппы G согласованным образом с включением (самая нижняя нулевая подгруппа меньше всех) и находим все комбинации множеств допустимых простых переводящих в более нижнюю подгруппу.
Дальнейшие детали дам в ходе обсуждения.

Поясните пожалуйста Ваш метод на примере одного из известных решений, например
$n=23333*205319*1746169$.
Периоды этих допустимых простых будут:
$T(b,23333)=2^2*19*307$
$T(b,205319)=251*409$
$T(b,1746169)=2^2*31*2347$
Их общий период $T_b=lcm(..)$ будет равен
$T_b=2^2*19*31*251*307*409*2347$
Как дальше вычисляется тот странный радикал и для чего?
Как затем определяется образующая $g$ для каждого из трех допустимых простых?
Затем вычисление координат $z_q$ из условия $g^{z_q}=y\mod q^k$, где так понимаю $y=1$ для этого примера и $q$ - простой делитель $T_b$, $k$ - его степень.
И в итоге, насколько понял, сумма вычисленных координат должна стать равной 0.

[dmd]

Такое свойство. Из $p|(2^h-1)$ и $p|(2^k-3)$ следует, что $p|(3^{h \backslash k} 2^{h \% k}-1)$.
И это дальше продолжаемо: обозначив $a=h \backslash k$ и $b=h \% k$, получаются $p|(3^{a(h \backslash b)}-2^{h \% b})$ и $p|(3^{a(k \backslash b)+1}-2^{k \% b})$ и т.д.
В связи с этим вопрос: есть ли в сети что-то наподобие Cunningham Project, где уже факторизуют числа вида $3^x-2^y$ и $3^x2^y-1$?

[maxal]

Нет. Но Joe Crump создал проект по факторизации чисел $2^k-3$:
http://www.immortaltheory.com/NumberTheory/2nm3_db.txt

[megamix62]

Исправил формулы в topic59242.html

[megamix62]

maxal в сообщении #75435 писал(а):

Цитата:

Ну, например, если 17 делит решение $2^n\equiv 3\pmod{n}$, то $\frac{n}{17}\equiv 14\pmod{16}.$ При этом, $n$ обязано делиться на gcd(14,16)=2, но этого быть не может. Поэтому 17 не может делить решения $2^n\equiv 3\pmod{n}$.

Можна рассказать откуда берется $\frac{n}{17}\equiv 14\pmod{16}.$

[megamix62]

<maxal >
Не дадите контакты - Крампа, Мантгомери, свой, R.K.Guy для личной переписки , если они у вас есть?
Как я вижу Вы преподователь , хотя по форуму не скажешь, не любите учить ?
Вас спрашивают, а Вы - молчок...
Или Вы не хотите , что решилась проблема решения сравнения $2^n = 3 \mod\ n$
(Department of Computer Science & Engineering University of South Carolina)

[maxal]

megamix62 в сообщении #581528 писал(а):

Можна рассказать откуда берется $\frac{n}{17}\equiv 14\pmod{16}.$

В том примере была описка - должно быть $2^n \equiv 13\pmod{n}$. Тогда из того, что 17 делит $n$ следует, что $2^n \equiv 13\pmod{17}$, что влечет $n\equiv 6\pmod{8}$ и т.п.

-- Fri Jun 15, 2012 07:26:34 --

megamix62 в сообщении #585293 писал(а):

Не дадите контакты - Крампа, Мантгомери, свой, R.K.Guy для личной переписки , если они у вас есть?

Контактов нету. Но подозреваю, что при желании их можно найти нагуглить. Мне можно писать в ЛС на форуме.

[megamix62]

maxal писал(а):

Вот для затравки все тройки $(p,a,m)$ в пределах $p<1000$:

5 15 20
19 247 342
23 184 253
29 145 812
47 893 1081
53 901 2756
71 1136 2485
97 1843 4656
101 6969 10100
139 5699 19182
149 12963 22052
163 16463 26406
167 8016 13861
173 4671 29756
191 16999 18145
197 35657 38612
211 9073 44310
239 11472 28441
263 21829 34453
269 29321 72092
293 46001 85556
311 28612 48205
317 78933 100172
359 24771 64261
379 51923 143262
383 70472 73153
389 105419 150932
409 4499 83436
431 5603 18533
461 51171 212060
479 11975 114481
499 75349 82834
503 65390 126253
509 4581 258572
557 242295 309692
599 104825 179101
643 112525 137602
647 51760 208981
653 289279 425756
677 324283 457652
701 300729 490700
719 161056 258121
743 81730 275653
773 16233 596756

Но в Ваших тройках - нету
ни $p=11$
ни $p=17$
т.к.
$p=24-7(+7)$
$p=11(-11)$ тоже не должны входить ($F10- R.K.Guy$)

[maxal]

megamix62 в сообщении #586365 писал(а):

Но в Ваших тройках - нету
ни $p=11$
ни $p=17$
т.к.
$p=24-7(+7)$
$p=11(-11)$ тоже не должны входить ($F10- R.K.Guy$)

Не понял вашего вопроса. Легко убедиться, что ни 11, ни 17 не могут делить решения сравнения $2^n\equiv 3\pmod n$.

[dmd]

Простые в форме $T=\frac{3^p-1}{2}$ и $T=\frac{3^p+1}{4}$, где $p$ из последовательностей A028491 и A007658 соответственно, удовлетворяют сравнению $2^\frac{T-1}{2p}\equiv \pm 3^x \pmod T$. Гипотетически некоторые из таких $T$ могут оказаться подходящими и быть делителями решения сравнения $2^n \equiv 3 \pmod n$, при условии что $x=1$ либо $\frac{T-1}{2px}$ натурально. Тогда решением будет $n=\frac{T(T-1)}{2px}$. Есть ли смысл в этом, или можно сходу показать, что такие $T$ не могут являться подходящими простыми?

[dmd]

У меня тогда два варианта кода получилось для сбора-проверки подходящих чисел. Один на основе массива Set; поиск в нём быстр, но пополнение его очень не эффективно, т.к. каждый раз пересорт массива происходит.

pkh()=
{
local(ii, n, nl, p, pmax, k, h2, h3, g, l, z1, z2, o, dz, a, m, mm, a2, m2):int;
local(s):bool;
local(v, w, Y, Z):vec;
local(pr, st);

Y= read("pkh.dbt"):vec;
st= Set(Y[,1]);
Z= read("pkhset.dbt"):vec;
\\st= Set():vec;
\\Y= matrix(0,3):vec;

\\n= (5):int;
\\n= (52*10^5+1):int;
\\n= (10527*10^5+1):int;
\\n= (47907*10^5+1):int;
\\n= (52637*10^5+1):int;
\\n= (74121*10^5+1):int;
\\n= (407433*10^5+1):int;
\\n= (102116*10^6+1):int;

\\n= (125641*10^6+1):int;
n= (10^5+1):int;
\\pmax= 0:int;
pmax= vecmax(Y[,1]):int;

while(n < (10^6),
 if(n%3, if(n%7, if(n%11, if(n%13, if(n%17, if(n%31, if(n%37,
 if(issquarefree(n),
  s= 1:bool;
  w= factorint(n):vec;
  nl= #w[,1]:int;
  for(j=1, nl,
   if(w[j,1]<=pmax,
    if(!setsearch(st,w[j,1]), s= 0:bool; break)
   )
  );

  if(s,
   for(j=1, nl,
    p= w[j,1]:int;

    if(p<=pmax,
     ii= setsearch(st,p):int;
     a2= Z[ii,2]:int; m2= Z[ii,3]:int
     ,
     g= 0:int;
     h2= znorder(Mod(2,p)):int; h3= znorder(Mod(3,p)):int;
     if(h2<>h3,
      if(h2 % h3, s= 0:bool; break,
       g= (h2/h3):int;
       if(!issquarefree(g), s= 0:bool; break,
        v= factorint(g):vec;
        l= #v[,1]:int;
        for(i=1, l,
         if(!setsearch(st,v[i,1]), s= 0:bool; break(2))
        )
       )
      )
     );

     pr= znprimroot(p);
     z1= znlog(3:int, pr):int; z2= znlog(2:int, pr):int;
     o= eulerphi(p):int; dz= gcd(o,z2):int;
     if(z1 % dz, s= 0:bool; break(), k= ((z1/z2) % (o/dz)):int);
     if(!(h2%24),
      if(kronecker(6,lift(Mod(k,24)))==-1, s= 0:bool; print(n); break)
     );

     a2= (p*k):int; m2= (p*h2):int;

     if(g>1,
      for(i=1, l,
       ii= setsearch(st,v[i,1]):int;
       a= Z[ii,2]:int; m= Z[ii,3]:int;
       mm= gcd(m,m2);
       if(Mod(a,mm) == Mod(a2,mm),
        a2= lift(chinese(Mod(a,m),Mod(a2,m2))):int;
        m2= (m*m2/mm):int, s= 0:bool; break(2)
       )
      )
     );

     Y= concat(Y, [p, a2, m2]):vec;
     pmax= p:int;
    );

    if(j==1, a= a2:int; m= m2:int,
     mm= gcd(m,m2);
     if(Mod(a,mm) == Mod(a2,mm),
      a= lift(chinese(Mod(a,m),Mod(a2,m2))):int;
      m= (m*m2/mm):int, s= 0:bool; break
     )
    );

   )
  );

  if(s,
   if(a>10^5,
    if(a%2, if(a%3, if(a%7, if(a%11, if(a%13, if(a%17, if(a%31, if(a%37,
     if(lift(Mod(2,a)^a)==3,
      print("a= ",a,"  ",factorint(a),"  n= ",n,"  ",factorint(n));
      write("pkh.txt",a,"  ",factorint(a),"  ",n,"  ",factorint(n))
     )
    ))))))))
   )
  )

 )
 )))))));
 n+=2
);
write("pkhnew.dbt", Y);
}

pkhset()=
{
local(l):int;
local(Y, Z, st):vec;

Y= read("pkh.dbt"):vec;
st= Set(Y[,1]):vec;
l= (#st):int;
Z= matrix(l,3):vec;
for(i=1, l, for(j=1, l, if(eval(st[i])==Y[j,1], Z[i,]= Y[j,]; print(i); break())));
write("pkhset.dbt", Z);
}

Другой на основе собственного бинарного поиска в массиве.

zsearch(p:int, Y:vec)=
{
local(j, i1, i2, it, t):int;
if(p<5, j= 0:int,
 i1= 1:int;
 i2= (#Y):int;
 t= -1:int;
 while(p!=t,
  it= ((i2+i1)\2):int;
  t= Y[it]:int;
  if(p<t, i2= (it-1):int,
   if(p>t, i1= (it+1):int,
    j= it:int; break
   )
  );
  if((i1>i2), j= 0:int; break)
 )
);
return(j:int)
};

pkh()=
{
local(Z, w):vec;
local(n, g, i, h, h3, z2, z3, o, dz, k, a, m, mm, a1, m1, a2, m2, d, L):int;
local(s):bool;
local(pr);

\\Z= matrix(0, 3):vec;
Z= read("mz78.dbt"):vec;
L= (#Z~):int;
if(L==0, n= (5):int, n= (Z[L,1]+2):int);

\\n= (526*10^4+1):int;
\\n= (105274*10^4+1):int;
\\n= (479070*10^4+1):int;
\\n= (741213*10^4+1):int;

while(n < 10^8,

 if(n%3, if(n%7, if(n%11, if(n%13, if(n%17, if(n%31, if(n%37,
 if(issquarefree(n),

  s= 1:bool;

  if(isprime(n),

   g= 0:int;
   h= znorder(Mod(2, n)):int; h3= znorder(Mod(3, n)):int;
   if(h<>h3,
    if(h%h3, s= 0:bool,
     g= (h/h3):int;
     i= zsearch(g, Z[,1]):int;
     if(!i, s= 0:bool)
    )
   );

   if(s,
    pr= znprimroot(n);
    z2= znlog(2:int, pr):int; z3= znlog(3:int, pr):int;
    o= eulerphi(n):int; dz= gcd(o, z2):int;
    if(z3%dz, s= 0:bool,
     k= ((z3/z2)%(o/dz)):int;
     if(!(h%24), if(kronecker(6, lift(Mod(k, 24)))==-1, s= 0:bool; print(n)))
    )
   );

   if(s,
    a= (n*k):int; m= (n*h):int;
    if(g,
     a2= (g*Z[i,2]):int; m2= (g*Z[i,3]):int;
     mm= gcd(m, m2):int;
     if(Mod(a, mm) == Mod(a2, mm),
      a= lift(chinese(Mod(a, m), Mod(a2, m2))):int;
      m= (m*m2/mm):int, s= 0:bool
     )
    )
   )

  ,

   w= divisors(n):vec;
   d= w[2]:int;
   i= zsearch(d, Z[,1]):int;
   if(i, a1= (d*Z[i,2]):int; m1= (d*Z[i,3]):int, s= 0:bool);

   if(s,
    L= (#w):int;
    d= w[L-1]:int;
    i= zsearch(d, Z[,1]):int;
    if(i, a2= (d*Z[i,2]):int; m2= (d*Z[i,3]):int, s= 0:bool)
   );

   if(s,
    mm= gcd(m1, m2):int;
    if(Mod(a1, mm) == Mod(a2, mm),
     a= lift(chinese(Mod(a1, m1), Mod(a2, m2))):int;
     m= (m1*m2/mm):int, s= 0:bool
    )
   );

  );

  if(s,

   if(a%2, if(a%3, if(a%7, if(a%11, if(a%13, if(a%17, if(a%31, if(a%37,
   if(lift(Mod(2,a)^a)==3,
    print("a= ",a,"  ",factorint(a),"  n= ",n,"  ",factorint(n));
    write("pkh.txt",a,"  ",factorint(a),"  ",n,"  ",factorint(n))
   )
   ))))))));

   Z= concat(Z, [n,a/n,m/n])
  );

 )
 )))))));

 n+=2
);
write("mz8.dbt", Z)
};

Мне хотелось создать код, который бы всю базу подходящих чисел держал в памяти, быстро её пополнял и быстро в ней искал. Нужно было выбрать, с каким массивом лучше работать - Set, вектор или матрицу; соответственно в каком формате сохранять. Нужно минимизировать дисковые операции, лучше однократно при входе читать из файла базы и при выходе писать данные в файл базы. Возможно ли сделать управляемый выход, т.е. чтоб в любой момент вручную остановить исполнение и корректно сохранить набранные данные? Какие слабые места наблюдаются в этом коде? Если кто понял алгоритм Руста, то не могли бы объяснить его на конкретном примере?

[dmd]

Немного поправил код, стало чуть быстрее. Для версии 2.3 уже не подходит, только для 2.4-2.6, т.к. изменился znlog.

zsearch(p:int, Y:vec)=
{
local(j, i1, i2, it, t):int;
if(p<5, j= 0:int,
 i1= 1:int;
 i2= (#Y):int;
 t= -1:int;
 while(p!=t,
  it= ((i2+i1)\2):int;
  t= Y[it]:int;
  if(p<t, i2= (it-1):int,
   if(p>t, i1= (it+1):int,
    j= it:int; break
   )
  );
  if((i1>i2), j= 0:int; break)
 )
);
return(j:int)
};

pkh()=
{
local(Z, w):vec;
local(n, g, i, j, h, h3, k, a, m, mm, a1, m1, a2, m2, d1, d2, L):int;
local(s):bool;

Z= matrix(0, 3):vec;
\\Z= read("pkh0.dbt"):vec;
L= (#Z~):int;
if(L==0, n= (5):int, n= (Z[L,1]+2):int);

\\n= (526*10^4+1):int;
\\n= (105274*10^4+1):int;
\\n= (479070*10^4+1):int;
\\n= (741213*10^4+1):int;

while(n < 10^6,

 if(n%3, if(n%7, if(n%11, if(n%13, if(n%17, if(n%31, if(n%37,
 if(issquarefree(n),

  s= 1:bool;

  if(isprime(n),

   g= 1:int;
   h= znorder(Mod(2, n)):int; h3= znorder(Mod(3, n)):int;
   if(!(h%h3),
    g= (h/h3):int;
    if(g-1,
     i= zsearch(g, Z[,1]):int;
     if(!i, s= 0:bool)
    )
    ,
    s= 0:bool
   );

   if(s,
    k= znlog(3, Mod(2,n), h);
    if(!(h%24), if(kronecker(6, lift(Mod(k, 24)))==-1, s= 0:bool; print(n)))
   );

   if(s,
    a= (n*k):int; m= (n*h):int;
    if(g-1,
     a2= (g*Z[i,2]):int; m2= (g*Z[i,3]):int;
     mm= gcd(m, m2):int;
     if(Mod(a, mm) == Mod(a2, mm),
      a= lift(chinese(Mod(a, m), Mod(a2, m2))):int;
      m= (m*m2/mm):int, s= 0:bool
     )
    )
   )

  ,

   w= divisors(n):vec;
   d1= w[2]:int;
   i= zsearch(d1, Z[,1]):int;
   if(i,
    d2= w[#w-1]:int;
    j= zsearch(d2, Z[,1]):int;
    if(j,
     a1= (d1*Z[i,2]):int; m1= (d1*Z[i,3]):int;
     a2= (d2*Z[j,2]):int; m2= (d2*Z[j,3]):int;
     mm= gcd(m1, m2):int; 
     if(Mod(a1, mm) == Mod(a2, mm),
      a= lift(chinese(Mod(a1, m1), Mod(a2, m2))):int;
      m= (m1*m2/mm):int
      ,
      s= 0:bool
     )
     ,
     s= 0:bool
    )
    ,
    s= 0:bool
   );

  );

  if(s,

   if(a%2, if(a%3, if(a%7, if(a%11, if(a%13, if(a%17, if(a%31, if(a%37,
   if(lift(Mod(2,a)^a)==3,
    print("a= ",a,"  ",factorint(a),"  n= ",n,"  ",factorint(n));
    write("pkh.txt",a,"  ",factorint(a),"  ",n,"  ",factorint(n))
   )
   ))))))));

   Z= concat(Z, [n,a/n,m/n])
  );

 )
 )))))));

 n+=2
);
write("pkh.dbt", Z)
};

После первого запуска файл pkh.dbt нужно переименовать в pkh0.dbt, закоментарить-раскоментарить соответственно строчки инициализации массива Z и увеличить верхнюю границу допустим до 10^7 в этом куске кода:

Z= matrix(0, 3):vec;
\\Z= read("pkh0.dbt"):vec;
L= (#Z~):int;
if(L==0, n= (5):int, n= (Z[L,1]+2):int);

\\n= (526*10^4+1):int;
\\n= (105274*10^4+1):int;
\\n= (479070*10^4+1):int;
\\n= (741213*10^4+1):int;

while(n < 10^6,

Чтоб проверить известные решения, нужно раскоментарить соответствующую строчку для n и изменить верхнюю границу до 10^10.

Мне по-прежнему интересно, возможно ли улучшить этот код, либо предложить другой более быстрый? Этот собирает базу подходящих до 10^6 примерно за 2 минуты в Винде, в Линуксе из-под gp2c будет примерно в 8 раз быстрее.