O (chyba nieudanym) klonowaniu baroków - Part I
Mirek
Ale Włodek mnie obudził swoim blogiem.
Wprawdzie nie spałem, cały czas czegoś poszukiwałem, czytałem
literaturę. Ale mimo iż chciałem napisać wcześniej, cały czas
wstrzymywałem się, gdyż czekałem na jakieś spektakularne rezultaty z
mojej strony :)
Ale ponieważ (na razie) nie nadchodzą, to podzielę się moimi
przemyśleniami - może kogoś zapłodnią. Zgóry zaznaczam, że nie piszę o
SA.
===
Troszkę teorii.
Jak wiemy funkcja sd(x) jest częściowo multiplikatywna. Wspominana
dotychczas na grupie multiplikatywność dotyczyła jedynie dodawania
nowych czynników pierwszych do liczby x, tzn. gdy gdc(p,x)=1, to
sd(x*p^n)=sd(x)*sd(p^n).
To jednak niejedyny obszar, w którym można znaleźć coś podobnego do
multiplikatywności. Dla danych p,a istnieją pewne liczby b,c, dla
zachodzą dwie ciekawe (być może;) zależności:
(1) sd(p^(a*b))=sd(p^a)*sd(p^b)*XXX.
(2) sd(p^(a+c))=sd(p^a)*YYY.
Gdyby dopuścić wymierność XXX czy YYY, to zależności stają się
trywialne. Mnie chodzi jednak o przypadek całkowitoliczbowych XXX i YYY
- daje to niezanikanie (jak to nazwać?) dzielników pierwszych sd(p^a)
czy sd(p^b).
Dla jakich b lub c to się dzieje?
Odpowiedź na to pytanie daje struktura cyklotomiczna funkcji sd(p^n).
Jak wiemy
sd(p^n) = (p^(n+1)-1) / (p-1)
Natomiast p^m-1 możemy wyrazić za pomocą iloczynu wszystkich
wielomianów cyklotomicznych o rzędach d będących dzielnikami n
(symbol Phi_d(p) oznacza wielomian cyklotomiczny rzędu d obliczony dla
argumentu p):
Pn: p^n-1 = Prod_(d|n) Phi_d(p)
Jeżeli m=k*n dla k>0, to m zawiera wszystkie dzielniki n, więc zachodzi
Pm: p^m-1 = (p^n-1) * ZZZ
czyli Pm zawiera wszystkie dzielniki Pn plus jakieś dodatkowe.
Na podstawie powyższej własności możemy odpowiedzieć, jakie wartości
a,b,c dają nam dodatkową (nazwijmy ją cyklotomiczną) multiplikatywność
funkcji sd():
Dla (1) warunek brzmi:
gcd( a+1, b+1 ) == 1
Dla (2) sytuacja wydaje się zawiła, ale tak nie jest:
(a+1) | (a+c+1)
Gdyby w (2) operację dodawania a+c zamienić mnożeniem a*c, to warunek
całkowitoliczbowości YYY przyjmie postać:
(a+1) | (a*c+1)
wydaje mi się to jednak mało użyteczne.
Trochę się rozpisałem, więc na tym na razie kończę. W zamyśle to miał
być jeden list (a do tej pory nic nie wyjaśniłem tematu - zaczyna to
wyglądać jak kryminał rozwiązanie zagadki na końcu serialu), ale być
może będą jeszcze dwie części:
- technikalia: techniki programowe - walka z ograniczonością zasobów
obliczeniowych
- wyjaśnienie tematu i może jakieś rezultaty.
Mirek
Wlodzimierz Holsztynski
> Wprawdzie nie spałem, cały czas czegoś poszukiwałem, czytałem
> literaturę. [...]
Mogłeś też poczytać na dany temat
coś niecoś pod psem.
> Jak wiemy funkcja sd(x) jest częściowo multiplikatywna.
W teorii liczb po prostu mówimy, że sd
jest multyplikatywna (bez "częściowo").
Taka jest klasyczna, standardowa,
od dawna ustalona terminologia.
> Wspominana dotychczas na grupie multiplikatywność
> dotyczyła jedynie dodawania nowych czynników
> pierwszych do liczby x, tzn. gdy gdc(p,x)=1, to
>
> sd(x*p^n)=sd(x)*sd(p^n).
Powiem ostrożnie, że gdy dodasz nasze forum,
+ pod psem + alt.pl.matematyka, to na pewno
było o wiele więcej, a nie "jedynie".
> Dla danych p,a istnieją pewne liczby b,c, dla
> zachodzą dwie ciekawe (być może;) zależności:
>
> (1) sd(p^(a*b))=sd(p^a)*sd(p^b)*XXX.
> (2) sd(p^(a+c))=sd(p^a)*YYY.
Chyba miałeś na myśli:
(1) sd(p^((a*b)-1))=sd(p^(a-1))*sd(p^(b-1))*XXX.
(2) sd(p^(a+c-1))=sd(p^(a-1))*YYY.
Pisałem o takich tożsamościach w ostatnich
miesiącach z okazji baroku, a dawniej
z okazji uogólnionych liczb Mersenna,
jako że M(n) := 2^n-1 = sd(2^(n-1)).
Kiedyś pod psem rozpatrywałem
ogólniej M_a(n) := sd(a^(n-1)). W obu
okresach pisałem dosyć systematycznie,
ale raz napiszę wciąż pełniej. Potrafię
też chyba odnaleźć swój printout notki,
w której te sprawy rozpatrywałem
na poziomie odpowiadającym de facto
twierdzeniu o wzajemności kwadratowej
(ale wszystko robiłem "ręcznie", bez
powoływania się na tw. Gaussa).
Najlepiej takie rzeczy opisać ogólnie,
w terminach waluacji, ale mniejsza.
Cieszę się Mirku, że zająłeś się tymi
podzielnościami. Lwią część dostaje
się za pomocą szkolnej tożsamości
algebraicznej. Dalsze, detaliczniejsze
podzielności wymagają nieco więcej
delikatności.
Pozdrawiam,
Włodek
Mirek
On Mon, Jul 09, 2007 at 06:34:10PM +0200, Mirek wrote:
> Trochę się rozpisałem, więc na tym na razie kończę. W zamyśle to miał
> być jeden list (a do tej pory nic nie wyjaśniłem tematu - zaczyna to
> wyglądać jak kryminał rozwiązanie zagadki na końcu serialu), ale być
> może będą jeszcze dwie części:
Trzeciej części na razie nie planuję (planowaną tematykę zawarłem
tutaj), gdyż na razie mój pomysł powoli umiera (ze względu na słabą
efektywność). Ale opiszę co zrobiłem.
=============
Pomysł był taki:
1) Dla n1<n2 (n1,n2 to numery mpn-ów w tabeli 5000+) sprawdzam "pełną
cyklotomiczną obcość" zdefiniowaną jako
x=mpn[n1], y=mpn[n2] - wartości mpn-ów
P={p, p|x OR p|y} - unia dzielników pierwszych x i y
gcd( log_p(x)+1, log_p(y)+1 ) == 1, forall p in P
Sprawdziłem, że dla n1 >= 1618, n2>n1 wszystkie mpn-y są
cyklotomicznie spokrewnione!
2) Dla pary obcych sobie x i y zdefiniowałem operację krzyżowania
(wykorzystanie zależności (1) z poprzedniego listu):
z = Prod_(p in P) p^( (log_p(x)+1) * (log_p(y)+1) - 1 )
Takie z ma następujące własności:
x|z, y|z, sd(x)|sd(z), sd(y)|sd(z)
Udało mi się zaobserwować (ale nie przejrzałem całej tabeli, w zasadzie
tylko częściowo dla n1={2,3}) tylko jedną parę (n1=2, n2=3), której
skojarzenie dało baroka:
"2^2*7", "2*3" -> "2^5*3*7"
3) Czerpiąc z zależności (2) z poprzedniego listu zdefiniowałem wzrost
potomka. Mając P={p, p|z OR p|sd(z)} dla każdego p poszukiwałem
możliwości zwiększenia potęgi p^a do p^c w liczbie "z" bez zaburzenia
istniejącej struktury cyklotomicznej:
Mając:
P={p, p|z OR p|sd(z)}
a=log_p(z)
b=log_p(sd(z))
Dla każdego p poszukiwałem c spełniającego:
(a+1)|(c+1)
a<=c<=a+b
4) I tu zaczęły się schody:(
Czasami liczba kobinacji c jest zabójcza (z uwzględnieniem rekurencji
wzrostu z). Plus problemy z faktoryzacją dużych wartości sd(z), ale o
tym w technikaliach.
Tak naprawdę właśnie uaktualniłem oprogramowanie oraz usunąłem
subiektywne ograniczenia na kombinacje "c". Więc obliczenia rozpoczęły
się na nowo i znacznie się wydłużyły:(
Znałazłem tą metodą tylko (czyli nic nowego):
c 1.8733791 2^5*3*7
c 2.5778519 2^9*3*11*31
c 3.0498799 2^13*3*11*43*127
c 3.2049844 2^14*5*7*19*31*151
d 3.2485794 2^14*3*5*7*19*31*151
d 3.9828826 2^25*3^3*5^2*19*31*683*2731*8191
f 4.7729739 2^50*3^4*5^4*7^3*11^3*13*17*31*61*67*71*103*139*2143*11119*131071
=============
TECHNIKALIA
Program napisałem w pari-gp i aby przyspieszyć obliczenia użyłem
kompilatora gp2c.
1) Analizowane wartości przechowywane są w postaci sfaktoryzowanej. Ale
nadal zachodzi potrzeba rozkładu wartości sd():
- wykorzystuję multiplikatywność i cyklotomiczność sd()
- badam podzielność czynnika przez dzielniki "z"
- dopiero potem co zostaje leci na faktoryzację
Można by się też zastanowić nad:
- implementacją szczególnych czynników rozkładu wielomianów
cyklotomicznych: Gaussa i Aurifeuille'a.
- wykorzystaniem tablic typu Cunningham Project itp.
2) Faktoryzacja - pari ma uniwersalny engine do tego, wg helpa robi
takie kroki
(0) dzielenie przez małe liczby pierwsze i takie w tablicy addprimes
(1) MPQS
(2) first-stage ECM
(3) Pollard-Brent Rho and Shanks SQUFOF
(4) final ECM
... ale ... to jest za wolne. Stwierdziłem, że funkcja factor(x), jest
użyteczna dla max c40 (c40 oznacza liczbę złożoną 40-to cyfrową), no w
porywach do około c50-c55. Więc zrobiłem heurystyczną faktoryzację
wielostopniową:
(a) Jeżeli x<c40 , to wołam
normalnie factor(x), jak nie to tylko krok (0) z powyższego
schematu (w pari to: factorint(x,15) )
(b) x<c40 -> factor(x), x<c65 -> factorE(x) (to zewnętrzny program, o
tym dalej), w przeciwnym wypadku trial division do 10^6 (??? nie
jestem tego pewien, wołam po prostu factor(x,10^6) - help opisuje
jako partial factorization up to limit)
(c) x<c40 -> factor(x), x<c65 -> factorE(x) otherwise factor(x,10^8)
(d) factorE(x)
Kazałem też funkcji factor automatycznie dodawać do tablicy addprimes
(stwierdziłem, że wiele przypadków się powtarza, więc widać zysk).
3) Co to jest factorE - to moja funkcja, która woła zewnętrzny program
(via skrypt konwertujący wyniki), którym jest obecnie msieve 1.25.
Mała dygresja: Jason Papadopoulos umieścił mnie w pliku Changes :)
A mój udział to zlokalizownie błędu (podczas pracy nad barokami).
Dlaczego msieve - wydaje mi się, że jest to obecnie (z pośród ogólnie
dostępnych) najszybszy program klasy SIQS (np. autor qsieve twierdzi,
że msieve jest max 2x szybszy, ale to pisał dla v1.20, a mnie wychodzi,
że 1.25 ma współczynnik 8x). Próbowałem też ppmpqs28, ppsiqs11 oraz
kilka innych, ale odpadły w przedbiegach.
Msieve nie jest też idealny:
- nie jest to uniwersalny engine
- ma wprawdzie też ECM, ale sam Jason pisze o tym "early stage", a i
tak można go włączyć dopiero od c96
- do c65 działa błyskawicznie, do c80 (taki mam ustawiony próg
działania) to zależy od szczęścia
Ponieważ jak wcześniej spominałem zdarzają się powtórzenia, to mój
wraper msieve.sh też ma bazę wyników, a w dodatku funkcja factorE
dodaje do tablicy największy czynnik p otrzymany z msieve - efekty
widać gołym okiem!
Znalazłem też aplet javie na www.alpertron.com.ar - wprawdzie w nazwie
ECM, ale jest to full engine. Jest niesamowicie szybki i ma dodatku
wsparcie dla (cyklotomicznych) Cunningham Project tables. Kilka liczb
c75-c77, z którymi męczył się msieve, zrobił w kilka minut. Ale
przyznaję też, że trafiłem na taką c72, przy ktorej oszacował mi czas
>24h :(
Ma dwie podstawowe wady: jest apletem webowym, a w dodatku nawet
cut&paste wyników (ze względu na format) jest koszmarem.
Nie marzy mi się nawet przeniesienie tego na C, ale gdy ktoś podjął się
zrobienie z tego kodu konsolowego (nie znan się na javie), tak aby
czytał i pisał stdin/stdout - to skakał bym z radości.
Mirek
Mirek
>
> On 9 Lip, 09:34, Mirek <m...@zind.ikem.pwr.wroc.pl> wrote:
>
> > Wprawdzie nie spałem, cały czas czegoś poszukiwałem, czytałem
> > literaturę. [...]
>
>
> Mogłeś też poczytać na dany temat
> coś niecoś pod psem.
Tak wiem, na podstawie (ulotnych) wspomnień czyniłem dalsze poszukiwania.
A ciągle się zbieram, aby posty z pod psa skompletować w moim mailboksie
liczby-barokowe.
> > Jak wiemy funkcja sd(x) jest częściowo multiplikatywna.
>
> W teorii liczb po prostu mówimy, że sd
> jest multyplikatywna (bez "częściowo").
> Taka jest klasyczna, standardowa,
> od dawna ustalona terminologia.
Wybacz nieświadomemu, ale termin "partially multiplicative" spotkałem w
jakimś artykule - tylko nie pomnę czy to był jakiś tuz czy skromny adept
numerologii.
> Powiem ostrożnie, że gdy dodasz nasze forum,
> + pod psem + alt.pl.matematyka, to na pewno
> było o wiele więcej, a nie "jedynie".
Dzięki - poszukam jeszcze raz, ale twoje posty (to chyba były z onetu?)
są trudne do wuszukania, ze względu na zniekształcony adres.
> > Dla danych p,a istnieją pewne liczby b,c, dla
> > zachodzą dwie ciekawe (być może;) zależności:
> >
> > (1) sd(p^(a*b))=sd(p^a)*sd(p^b)*XXX.
> > (2) sd(p^(a+c))=sd(p^a)*YYY.
>
> Chyba miałeś na myśli:
>
> (1) sd(p^((a*b)-1))=sd(p^(a-1))*sd(p^(b-1))*XXX.
Sorry, tak to bywa, jak z programu robi się opis :(
U mnie naprawdę jest:
sd( p^((a+1)*(b+1)-1) )=sd(p^a)*sd(p^b)*XXX
co jest równoważne twojej postaci po zmniejszeniu o -1 wartości a,b .
> (2) sd(p^(a+c-1))=sd(p^(a-1))*YYY.
Tu nie ma błędu - postacie są równoważne, ale warunek podałem dla mojej
wersji.
> Pisałem o takich tożsamościach w ostatnich
> miesiącach z okazji baroku, a dawniej
> z okazji uogólnionych liczb Mersenna,
Dzięki za dodatkową wskazówkę, postaram się odkopać i skopiować posty do
skrzynki (już parę razy chciałem coś zapytać).
A zawartość mailboxa i swoich programów udostępnię na www (jak wreszcie
uruchomię z powrotem serwer).
Wlodzimierz Holsztynski
> > > Jak wiemy funkcja sd(x) jest częściowo multiplikatywna.
>
> > W teorii liczb po prostu mówimy, że sd
> > jest multyplikatywna (bez "częściowo").
> > Taka jest klasyczna, standardowa,
> > od dawna ustalona terminologia.
>
> Wybacz nieświadomemu, ale termin "partially multiplicative"
> spotkałem w jakimś artykule - tylko nie pomnę czy to był
> jakiś tuz czy skromny adept numerologii.
W teoriik liczb, funkcje f : N --> G ze zbioru liczb
naturalnych w dowolny monoid przemienny (G * 1)
-- powiedzmy w ciało liczb zespolonych G = C --
nazywamy multyplikatywną <=:=>
gcd(x y) = 1 ==> f(x*y) = f(x) * f*y)
dla dowolnych argumentów x y. Funkcja
multyplikatywna f jest jednoznacznie wyznaczona
przez swoje wartości dla potęg p^k liczb pierwszych.
Spróbuję się "starczymi oczamy WGRYŹĆ się"
w Twoje teksty,
Pozdrawiam,
Włodek
Mirek
>
> Trzeciej części na razie nie planuję (planowaną tematykę zawarłem
> tutaj), gdyż na razie mój pomysł powoli umiera (ze względu na słabą
> efektywność). Ale opiszę co zrobiłem.
Włodek zarzucał mi m.in. brak znajomości różnych wcześniejszych (nawet
z poprzenich lat) wątków dotyczących naszej tematyki. Niestety
odkrywanie koła zdarza się ciągle i wszędzie :)
Co zaś do merytorycznej wartości mojego pomysłu - to sam sądzę, że
jest to ślepa uliczka, choć kto wie trochę jescze policzę.
A sama technologia może się przydać do SA.
Ale przecież podobieństwo MPN-ów legło też u podstaw Włodkowych
TWR-ów. Ja zaś zaobserwowałem nie tylko podobieństwo czynników
pierwszych, ale także podobieństwo (cyklotomiczne) ich potęg.
Większość liczb w tabeli 5000+ jest wg mojego kryterium ze sobą
spokrewniona.
Z czego to wynika, nie mnie rozsądzać, ale odrobiłem trochę lekcji i
wiem, że w numerologii obowiązują Prawa Małych Liczb: każda skończona
sekwencja jest za krótka, aby wyjaśnić ogólne reguły nią rządzące.
Czy to jest ogólna prawidłowość???
Czy może za mała baza - przypuszczam, że to drugie może to być przyczyną.
Większość aktywnych w ostanich latach na polu baroków badaczy (sądzę po
googlu, a nie po publikacjach):
Fred Helenius,
Jason Moxham,
George Woltman,
Ron Sorli (+Graeme Cohen)
Achim Flammenkamp
Joseph Chein
(kolejność mniej więcej historyczna) bazuje, jeżeli nawet nie na tym
samym kodzie obliczeniowym (prawdobodobnie pochodzącym od Moxhama), to
na koncepcji tzw. sigma-chains (muszę to jeszce przestudiować).
====
Wracają do głównego tematu wątku.
Mam zaszczyt oznajmić, że od wczoraj mój program doznał około 13-to
krotnego przyspieszenia (uwzględniając wcześniejsze zagadnienie
faktoryzacji szacuję, że osiągnąłem znacznie ponad stukrony wzrost
względem wersji pierwotnej).
Osiągnałem to dokunując dwóch zmian
1) Pierwsza jest stosunkowo mało istotna. Zlikwidowałem sprawdzanie,
czy w drzewie rekurencji nie pojawiło się już dane rozwiązanie. To
wg moich szacunków ma mały wpływ (w obie strony +- ) na predkość
działania, ale kilkukrotnie redukuje zajętość pamięci.
2) Jak wspominałem wcześniej, liczbę n podawaną wzrostowi
cyklotomicznemu przechowywałem w postaci sfaktoryzowanej - oznaczmy
to fn. Natomiast wspólczynnik barokowy q=brq(n) (a dokładniej
q=brq(fn) ) obliczałem na bieżąco. A ponieważ w moim zagadnieniu i
tak potrzebuję czynników pierwszych sd(n), to obliczenia trochę
trwały.
Zdecydowałem, że sd(n) nie jest mi w ogóle potrzebna. Używam oprócz
fn tylko sfaktoryzowanej postaci q: fq. Obróbce podlegają tylko
tablice fn i fq (w pari-gp strukturą faktoryzacyjną jest tablica
dwukolumnowa: [p,n] - p czynniki pierwsze, n ich potęgi).
Tak, tak wiem, że q w ogólnym przypadku jest liczbą wymierną,
w takim przypadku czynniki pierwsze mianownika pojawiają się
w postaci ujemnych wartości potęg.
Zmieniając w fn potęgę danej liczby pierwszej p z n=a na n=b
nie obliczam na nowo całej fq tylko modyfikuję ją. Podstawą tej
operacji jest wzór:
X = p^(a-b) * ( p^(b+1) - 1 ) / ( p^(a+1) - 1 )
Wartości potęg czynników pierwszych X dodaję do odpowiednich
pozycji w tablicy fq.
Mirek
Mirek
Jeszcze raz przypominam - PML (prawo małych liczb).
-------------
On Wed, Jul 11, 2007 at 10:54:32PM +0200, Mirek wrote:
> Większość liczb w tabeli 5000+ jest wg mojego kryterium ze sobą
> spokrewniona.
Nie wiem skąd mi się wziął tak mocny termin "spokrewniona".
U założeń mojego pomysłu był "brak pokrewieństwa".
Mimo podobieństwa semantycznego chciałbym jednak podkreślić (w moim
odczuciu) istotną różnicę.
Dla przypomnienia - brakiem podobieństwa cyklotomicznego dwóch liczb
m,n nazwałem "nieistnienie ani jednego wspólnego dzielnika" o postaci:
p^q (p,q liczby pierwsze).
Tabela 5000+ to mały zbiór. Wydaje mi się, że nie znamy żadnych
sensownych podstaw do budowy jego w sposób algorytmiczny oprócz pełnego
przeglądu zbioru liczb - nawet wspomiana metodyka sigma-chains, jest
chyba heurystyczna i mało efektywna. Brak tzw. Włodkowych "rakiet"
Jednak bez sprawadzenia, czy choćby w tej skończonej populacji istnieje
statystyczna nadreprezentatywność wspólnych dzielników p^q (w całym
zbiorze tych dzielników zaobserwowanych w tabeli 5000+) -
trudno mówić o cyklotomicznym pokrewieństwie, jako regule rządzącej tą
tabelą (a przynajmniej większością tej populacji).
Mirek