Testy pierwszości. Liczby pseudopierwsze Fermata i liczby silnie pseudopierwsze

Z Henryk Dąbrowski
Przejdź do nawigacji Przejdź do wyszukiwania
11.11.2022



Potęgowanie modulo

Uwaga L1
Z twierdzenia Fermata (zobacz J21) wynika, że jeżeli liczby [math]\displaystyle{ a }[/math] i [math]\displaystyle{ m }[/math] są względnie pierwsze oraz [math]\displaystyle{ m }[/math] nie dzieli liczby [math]\displaystyle{ a^{m - 1} - 1 }[/math], to [math]\displaystyle{ m }[/math] jest liczbą złożoną. Każde twierdzenie pozwalające wykryć złożoność liczby może być wykorzystane do badania pierwszości liczb. Twierdzenia takie nie dają całkowitej pewności, że badana liczba jest pierwsza. Mamy na przykład [math]\displaystyle{ 341 = 11 \cdot 31 }[/math], ale [math]\displaystyle{ 341 \mid (2^{340} - 1) }[/math], bo

[math]\displaystyle{ 2^{340} - 1 = 2239744742177804210557442280568444278121645497234649534899989100963791871180160945380877493271607115775 }[/math]
[math]\displaystyle{ \;\;\; = 341 \cdot 6568166399348399444449977362370804334667582103327417990909058947107894050381703652143335757394742275 }[/math]

Widzimy, że nawet dla niewielkiej liczby [math]\displaystyle{ 341 }[/math], potęga [math]\displaystyle{ 2^{340} - 1 }[/math] jest liczbą ogromną. Jeśli ta metoda ma mieć jakiekolwiek zastosowanie, to musimy znaleźć inny sposób obliczania reszty z dzielenia [math]\displaystyle{ a^n }[/math] przez [math]\displaystyle{ m }[/math], czyli potęgowania modulo [math]\displaystyle{ m }[/math].


Uwaga L2
Wykorzystując wzór rekurencyjny

[math]\displaystyle{ a^n = \left\{ \begin{array}{cll} a & & \text{gdy } n = 1 \\ (a^2)^{{\large\frac{n}{2}}} & & \text{gdy } n \text{ jest parzyste} \\ a \cdot (a^2)^{{\large\frac{n - 1}{2}}} & & \text{gdy } n \text{ jest nieparzyste} \\ \end{array} \right. }[/math]


możemy napisać w PARI/GP prosty program do potęgowania modulo:

modPower(a, n, m) = 
\\ a - podstawa, n - wykładnik, m - moduł
{
local(w);
if( m == 1, return(0) );
a = a % m;
w = 1;
while( n > 0,
       if( n % 2 == 1, w = (w * a) % m; n = n - 1); \\ gdy n jest nieparzyste, wyłączamy a i zmniejszamy n o jeden
       a = (a*a) % m; \\ wyliczamy nową podstawę modulo m
       n = n/2; \\ dla nowej podstawy wykładnik jest dwa razy mniejszy
     );
return(w);
}


Czytelnik łatwo sprawdzi, że w funkcji modPower() nie występują wyrażenia o wartości większej od [math]\displaystyle{ m^2 }[/math].

Zauważmy jeszcze, że PARI/GP umożliwia szybkie potęgowanie modulo i nie musimy korzystać z funkcji modPower(). Wystarczy napisać

lift( Mod(a, m)^d )

Co ważniejsze, powyższe polecenie jest wykonywane znacznie szybciej niż nasza funkcja modPower(). Podaliśmy kod funkcji dlatego, że jest ona bardzo ważna i Czytelnik powinien wiedzieć, jak jest w praktyce realizowana.


Uwaga L3
Wykorzystując wzór rekurencyjny

[math]\displaystyle{ a \cdot b = \left\{ \begin{array}{cll} a & & \text{gdy } b = 1 \\ 2 a \cdot \frac{b}{2} & & \text{gdy } b \text{ jest parzyste} \\ a + 2 a \cdot \frac{b - 1}{2} & & \text{gdy } b \text{ jest nieparzyste} \\ \end{array} \right. }[/math]


możemy napisać w PARI/GP prosty program do mnożenia modulo:

modMult(a, b, m) = 
\\ a, b - czynniki, m - moduł
{
local(w);
if( m == 1, return(0) );
a = a % m;
b = b % m;
w = 0;
while( b > 0,
       if( b % 2 == 1, w = (w + a) % m; b = b - 1 );  \\ gdy b jest nieparzysty, wydzielamy a i zmniejszamy b o jeden
       a = (2 * a) % m;  \\ wyliczamy nowy czynnik a modulo m
       b = b / 2;  \\ dla nowego czynnika a czynnik b jest dwa razy mniejszy
     );
return(w);
}


Czytelnik może zapytać, po co nam program do obliczania iloczynu modulo. Istotnie, jeśli piszemy programy w PARI/GP, to liczby całkowite mogą być ogromne i nie mamy powodu do zmartwienia (między innymi dlatego podajemy przykłady programów w PARI/GP). Jeżeli jednak będziemy potrzebowali napisać program w innym języku – powiedzmy w C – to ten problem stanie się nagle bardzo ważny. W C możemy przeprowadzać obliczenia dla bardzo dużych liczb całkowitych. Zmienne całkowite zadeklarowane jako uint32_t mogą przyjmować wartości z przedziału [math]\displaystyle{ [0, 2^{32} - 1] }[/math], a zmienne całkowite zadeklarowane jako uint64_t mogą przyjmować wartości z przedziału [math]\displaystyle{ [0, 2^{64} - 1] }[/math]. Liczba [math]\displaystyle{ 2^{64} \approx 1.84 \cdot 10^{19} }[/math] jest na tyle duża, że możemy wiele problemów liczyć, pisząc programy w C, co zapewnia większą szybkość obliczeń. W takich przypadkach funkcja modMult() może być bardzo użyteczna.

Zauważmy, że wykonując potęgowanie modulo, obliczamy iloczyny (w * a) % m i (a * a) % m. Jeżeli [math]\displaystyle{ m \lt 2^{32} }[/math], to nie napotkamy problemu: obydwa iloczyny są mniejsze od [math]\displaystyle{ 2^{64} }[/math] i będziemy mogli je wyliczyć. Ale w przypadku większych modułów już tak nie będzie i jeżeli chcemy zwiększyć zakres obliczeń, to musimy mnożenie wykonywać przy użyciu funkcji modMult(). Wystarczy założenie, że moduł [math]\displaystyle{ m \lt 2^{63} }[/math], aby suma (w + a) % m i iloczyn (2 * a) % m mogły zostać wyliczone.



Liczby pseudopierwsze Fermata

Liczby złożone nieparzyste spełniające równanie Fermata, otrzymały własną nazwę.

Definicja L4
Jeżeli [math]\displaystyle{ m }[/math] jest liczbą złożoną nieparzystą i dla pewnego [math]\displaystyle{ a \in \mathbb{Z} }[/math] prawdziwa jest kongruencja

[math]\displaystyle{ a^{m - 1} \equiv 1 \pmod m }[/math]

to powiemy, że [math]\displaystyle{ m }[/math] jest liczbą pseudopierwszą Fermata przy podstawie [math]\displaystyle{ a }[/math] lub krótko: [math]\displaystyle{ m }[/math] jest PSP([math]\displaystyle{ a }[/math]).


Uwaga L5
Zauważmy, że w definicji liczb pseudopierwszych Fermata nie musimy dodatkowo dołączać warunku [math]\displaystyle{ \gcd (a, m) = 1 }[/math], bo wynika on z przyjętej definicji. Mamy

[math]\displaystyle{ \gcd (a^{m - 1}, m) = \gcd (1, m) = 1 }[/math]

Zatem [math]\displaystyle{ \gcd (a, m) = 1 }[/math].

Możemy też łatwo pokazać, że jeżeli [math]\displaystyle{ \gcd (a, m) = d \gt 1 }[/math], to liczba [math]\displaystyle{ m }[/math] nie może być liczbą pseudopierwszą Fermata przy podstawie [math]\displaystyle{ a }[/math]. Istotnie, gdyby tak było, to mielibyśmy

[math]\displaystyle{ a^{m - 1} \equiv 1 \pmod{m} }[/math]

Ponieważ [math]\displaystyle{ d \mid m }[/math], to jest również

[math]\displaystyle{ a^{m - 1} \equiv 1 \pmod{d} }[/math]

Ale modulo [math]\displaystyle{ d }[/math] otrzymujemy natychmiast

[math]\displaystyle{ 0 \equiv 1 \pmod{d} }[/math]

Co jest niemożliwe, czyli [math]\displaystyle{ m }[/math] nie jest PSP([math]\displaystyle{ a }[/math]).


Twierdzenie L6
Dla każdej podstawy [math]\displaystyle{ a \geqslant 2 }[/math] istnieje nieskończenie wiele liczb pseudopierwszych Fermata.

Dowód

Niech [math]\displaystyle{ a \in \mathbb{Z} }[/math] i [math]\displaystyle{ a \geqslant 2 }[/math]. Jeżeli [math]\displaystyle{ p }[/math] jest liczbą pierwszą nieparzystą, to

[math]\displaystyle{ a^p - 1 = (a - 1) (a^{p - 1} + a^{p - 2} + \ldots + a^2 + a + 1) }[/math]

oraz

[math]\displaystyle{ a^p + 1 = (a + 1) (a^{p - 1} - a^{p - 2} + \ldots + a^2 - a + 1) }[/math]

Czyli [math]\displaystyle{ a - 1 \mid a^p - 1 }[/math] oraz [math]\displaystyle{ a + 1 \mid a^p + 1 }[/math].


Jeżeli przez [math]\displaystyle{ R_2 (a) }[/math] oznaczymy resztę z dzielenia liczby [math]\displaystyle{ a }[/math] przez [math]\displaystyle{ 2 }[/math] równą [math]\displaystyle{ 0 }[/math] lub [math]\displaystyle{ 1 }[/math], to prawdziwe są kongruencje

[math]\displaystyle{ a \equiv R_2 (a) \pmod 2 }[/math]

oraz

[math]\displaystyle{ a^n \equiv R_2 (a) \pmod 2 }[/math]

dla dowolnej liczby całkowitej dodatniej [math]\displaystyle{ n }[/math]. Zatem modulo [math]\displaystyle{ 2 }[/math] jest

[math]\displaystyle{ {\small\frac{a^p - 1}{a - 1}} \equiv R_2 (a) \cdot (p - 1) + 1 \equiv 1 \pmod 2 }[/math]
[math]\displaystyle{ {\small\frac{a^p + 1}{a + 1}} \equiv 1 \pmod 2 }[/math]

Co oznacza, że

[math]\displaystyle{ m = {\small\frac{a^p - 1}{a - 1}} \cdot {\small\frac{a^p + 1}{a + 1}} \equiv 1 \pmod 2 }[/math]

Czyli [math]\displaystyle{ m }[/math] jest złożoną liczbą nieparzystą. Pozostaje pokazać, że [math]\displaystyle{ a^{m - 1} \equiv 1 \pmod m }[/math].


Z twierdzenia Fermata wiemy, że

[math]\displaystyle{ a^p - 1 \equiv a - 1 \pmod p }[/math]

Ponieważ [math]\displaystyle{ (a - 1) \mid (a^p - 1) }[/math], to możemy napisać

[math]\displaystyle{ (a - 1) \cdot \left( {\small\frac{a^p - 1}{a - 1}} - 1 \right) \equiv 0 \pmod p }[/math]

Z założenia [math]\displaystyle{ p \nmid (a - 1) }[/math], zatem liczba pierwsza [math]\displaystyle{ p }[/math] musi dzielić [math]\displaystyle{ {\small\frac{a^p - 1}{a - 1}} - 1 }[/math] i otrzymujemy

[math]\displaystyle{ {\small\frac{a^p - 1}{a - 1}} \equiv 1 \pmod p }[/math]

Postępując analogicznie jak wyżej, dostajemy

[math]\displaystyle{ a^p + 1 \equiv a + 1 \pmod p }[/math]
[math]\displaystyle{ (a + 1) \cdot \left( {\small\frac{a^p + 1}{a + 1}} - 1 \right) \equiv 0 \pmod p }[/math]
[math]\displaystyle{ {\small\frac{a^p + 1}{a + 1}} \equiv 1 \pmod p }[/math]

Wynika stąd, że [math]\displaystyle{ m \equiv 1 \pmod p }[/math].

Zbierając mamy [math]\displaystyle{ 2 \mid (m - 1) }[/math] i [math]\displaystyle{ p \mid (m - 1) }[/math], zatem [math]\displaystyle{ 2 p \mid (m - 1) }[/math], czyli

[math]\displaystyle{ m = {\small\frac{a^{2 p} - 1}{a^2 - 1}} \equiv 1 \pmod{2 p} }[/math]

Oznacza to, że [math]\displaystyle{ m = 1 + 2 k p }[/math] dla pewnej liczby całkowitej [math]\displaystyle{ k \gt 0 }[/math].


Zauważmy teraz, że z definicji liczby [math]\displaystyle{ m }[/math] mamy [math]\displaystyle{ (a^2 - 1) m = a^{2 p} - 1 }[/math]. Rozpatrując to równanie modulo [math]\displaystyle{ m }[/math], otrzymujemy

[math]\displaystyle{ a^{2 p} \equiv 1 \pmod m }[/math]

Zatem modulo [math]\displaystyle{ m }[/math] jest

[math]\displaystyle{ a^{m - 1} = a^{2 k p} = (a^{2 p})^k \equiv 1^k \equiv 1 \pmod m }[/math]

Ponieważ dowolna liczba pierwsza [math]\displaystyle{ p \gt a^2 - 1 }[/math] nie dzieli [math]\displaystyle{ a^2 - 1 }[/math], to dla każdego [math]\displaystyle{ a \geqslant 2 }[/math] istnieje nieskończenie wiele liczb, które są PSP([math]\displaystyle{ a }[/math]). Co należało pokazać.


Przykład L7
Z dowodu twierdzenia L6 wynika, że jeżeli liczba [math]\displaystyle{ p }[/math] jest liczbą pierwszą nieparzystą i [math]\displaystyle{ p \nmid (a^2 - 1) }[/math], to liczba [math]\displaystyle{ m = {\small\frac{a^{2 p} - 1}{a^2 - 1}} }[/math] jest PSP([math]\displaystyle{ a }[/math]). Poniżej przedstawiamy przykłady takich liczb, dla kolejnych liczb pierwszych nieparzystych [math]\displaystyle{ p }[/math] takich, że [math]\displaystyle{ p \nmid (a^2 - 1) }[/math].

for(a=2, 5, s=1; d=a^2-1; forprime(p=3, 50, if( d%p == 0, next() ); m=(a^(2*p)-1)/d; print("a= ", a, "   m= ", m); s++; if( s>6, break() ) )) 


Uwaga L8
Wykorzystując funkcję potęgowania modulo, możemy napisać prosty program do testowania pierwszości liczb w oparciu o twierdzenie Fermata.

isPrimeOrPSP(m, a) = 
{
if( modPower(a, m-1, m) == 1, return(1), return(0) );
}


Przykład L9
Poniższa tabela zawiera najmniejsze liczby pseudopierwsze Fermata dla podstaw [math]\displaystyle{ a }[/math] od [math]\displaystyle{ 2 }[/math] do [math]\displaystyle{ 15 }[/math]

for(a=2, 15, s=1; forstep(m=1, 2000, 2, if( isPrimeOrPSP(m, a)  &&  !isprime(m), print("a=", a, "  m=", m); s++ ); if( s>5, break() ) ))


Przykład L10
Tabela pokazuje ilość liczb pseudopierwszych Fermata dla podstaw [math]\displaystyle{ a }[/math] od [math]\displaystyle{ 2 }[/math] do [math]\displaystyle{ 15 }[/math]

for(a=2, 15, s=0; forstep(k=1, 10^6, 2, if( isPrimeOrPSP(k, a)  &&  !isprime(k), s++ )); print("a= ", a, "   ", s))


Uwaga L11
Można pokazać, że jeżeli [math]\displaystyle{ m }[/math] jest liczbą nieparzystą złożoną i istnieje przynajmniej jedna liczba [math]\displaystyle{ a }[/math] względnie pierwsza z [math]\displaystyle{ m }[/math], taka że

[math]\displaystyle{ a^{m - 1} \not\equiv 1 \pmod m }[/math]

to dla co najmniej połowy liczb [math]\displaystyle{ b \in [1, m - 1] }[/math] względnie pierwszych z [math]\displaystyle{ m }[/math] jest

[math]\displaystyle{ b^{m - 1} \not\equiv 1 \pmod m }[/math]

Zatem przeprowadzając test Fermata, możemy z prawdopodobieństwem nie mniejszym niż [math]\displaystyle{ \tfrac{1}{2} }[/math] twierdzić, że liczba, która przeszła test, jest liczbą pierwszą. Wykonując test [math]\displaystyle{ k }[/math] razy dla [math]\displaystyle{ k }[/math] różnych podstaw z przedziału [math]\displaystyle{ [1, m - 1] }[/math] możemy z prawdopodobieństwem większym niż [math]\displaystyle{ 1 - \left( \tfrac{1}{2} \right)^k }[/math] twierdzić, że badana liczba [math]\displaystyle{ m }[/math] jest pierwsza.

Niestety, istnieją liczby złożone [math]\displaystyle{ m }[/math] takie, że

[math]\displaystyle{ a^{m - 1} \equiv 1 \pmod m }[/math]

dla każdego [math]\displaystyle{ a }[/math] względnie pierwszego z [math]\displaystyle{ m }[/math]. Liczby te nazywamy liczbami Carmichaela i jest ich nieskończenie wiele. Pokazano, że dla dostatecznie dużych [math]\displaystyle{ x }[/math] ilość liczb Carmichaela mniejszych od [math]\displaystyle{ x }[/math] przekracza [math]\displaystyle{ x^{1/3} }[/math][1][2][3]. Test Fermata jest zatem zbyt zawodny, aby można było go stosować.


Przykład L12
Oto wszystkie liczby Carmichaela mniejsze od [math]\displaystyle{ 100 000 }[/math]



Test Millera-Rabina

Rozpoczniemy od udowodnienia prostego twierdzenia

Twierdzenie L13
Jeśli [math]\displaystyle{ m }[/math] jest liczbą pierwszą nieparzystą i [math]\displaystyle{ x^2 \equiv 1 \pmod m }[/math], to albo [math]\displaystyle{ x \equiv - 1 \pmod m }[/math], albo [math]\displaystyle{ x \equiv 1 \pmod m }[/math].

Dowód

Z założenia

[math]\displaystyle{ x^2 \equiv 1 \pmod m }[/math]

Zatem [math]\displaystyle{ m \mid (x^2 - 1) }[/math], czyli [math]\displaystyle{ m \mid (x - 1) (x + 1) }[/math]. Ponieważ z założenia [math]\displaystyle{ m }[/math] jest liczbą pierwszą nieparzystą, to [math]\displaystyle{ m }[/math] dzieli dokładnie jedną z liczb [math]\displaystyle{ x - 1 }[/math] i [math]\displaystyle{ x + 1 }[/math]. Istotnie, gdyby [math]\displaystyle{ m \mid (x - 1) }[/math] [math]\displaystyle{ \text{i} \;\, m \mid (x + 1) }[/math], to [math]\displaystyle{ m }[/math] dzieliłaby również ich różnicę równą [math]\displaystyle{ 2 }[/math], co jest niemożliwe w przypadku gdy [math]\displaystyle{ m }[/math] jest liczbą pierwszą nieparzystą.


Prace Gary'ego Millera[4] i Michaela Rabina[5] pozwoliły sformułować znacznie silniejszy test. Podstawą tego testu jest następujące twierdzenie

Twierdzenie L14
Jeżeli [math]\displaystyle{ m }[/math] jest liczbą pierwszą nieparzystą i [math]\displaystyle{ m - 1 = 2^r d }[/math], gdzie [math]\displaystyle{ d }[/math] jest liczbą nieparzystą, to dla dowolnego [math]\displaystyle{ a \in [1, m - 1] }[/math] jest albo

[math]\displaystyle{ a^d \equiv 1 \pmod m }[/math]

albo

[math]\displaystyle{ a^{2^k d} \equiv - 1 \pmod m }[/math]

dla pewnego [math]\displaystyle{ k \in [0, r - 1] }[/math].

Dowód

Rozważmy ciąg [math]\displaystyle{ r + 1 }[/math] liczb zdefiniowanych następująco

[math]\displaystyle{ \begin{array}{l} u_0 = a^d \\ u_1 = a^{2 d} = (a^d)^2 \\ u_2 = a^{2^2 d} = (a^{2 d})^2 \\ \ldots \ldots \ldots \ldots \ldots \ldots \ldots \ldots \ldots \ldots \\ u_{r - 1} = a^{2^{r - 1} d} = (a^{2^{r - 2}})^2 \\ u_r = a^{2^r d} = (a^{2^{r - 1} d})^2 = a^{m - 1} \\ \end{array} }[/math]

Wyrazy ciągu [math]\displaystyle{ (u_i) }[/math] są dane wzorem ogólnym

[math]\displaystyle{ u_i = a^{2^i d} }[/math]

gdzie [math]\displaystyle{ i = 0, 1, \ldots, r }[/math]

Zauważmy, że mogą zdarzyć się następujące sytuacje

a) żaden z wyrazów ciągu [math]\displaystyle{ (u_i) }[/math] nie przystaje do [math]\displaystyle{ 1 }[/math] modulo [math]\displaystyle{ m }[/math]
b) wszystkie wyrazy ciągu [math]\displaystyle{ (u_i) }[/math] przystają do [math]\displaystyle{ 1 }[/math] modulo [math]\displaystyle{ m }[/math]
c) [math]\displaystyle{ u_k }[/math] jest pierwszym wyrazem ciągu [math]\displaystyle{ (u_i) }[/math], który przystaje do [math]\displaystyle{ 1 }[/math] modulo [math]\displaystyle{ m }[/math]


Co możemy zapisać jako

a) [math]\displaystyle{ u_i \not\equiv 1 \pmod m \quad }[/math] dla każdego [math]\displaystyle{ i \in [0, r] }[/math]
b) [math]\displaystyle{ u_i \equiv 1 \pmod m \quad }[/math] dla każdego [math]\displaystyle{ i \in [0, r] }[/math]
c) [math]\displaystyle{ u_k \equiv 1 \pmod m \quad }[/math] dla pewnego [math]\displaystyle{ k \in [0, r] }[/math]


Z definicji każdy wyraz ciągu [math]\displaystyle{ (u_i) }[/math] jest kwadratem poprzedniego. W szczególności oznacza to, że jeżeli dla pewnego [math]\displaystyle{ k \in [0, r] }[/math] jest [math]\displaystyle{ u_k \equiv 1 \pmod m }[/math], to [math]\displaystyle{ u_i \equiv 1 \pmod m }[/math] dla wszystkich [math]\displaystyle{ i \geqslant k }[/math]. Ten fakt pozwala doprecyzować zapis poszczególnych przypadków.

a) [math]\displaystyle{ u_i \not\equiv 1 \pmod m \quad }[/math] dla każdego [math]\displaystyle{ i \in [0, r] }[/math]
b) [math]\displaystyle{ u_0 \equiv 1 \pmod m }[/math]
c) [math]\displaystyle{ u_i \not\equiv 1 \pmod m \quad }[/math] dla każdego [math]\displaystyle{ i \in [0, k - 1] \quad }[/math] oraz [math]\displaystyle{ \quad u_i \equiv 1 \pmod m \quad }[/math] dla każdego [math]\displaystyle{ i \in [k, r] , \quad }[/math] gdzie [math]\displaystyle{ k \geqslant 1 }[/math]


W przypadku a) mamy [math]\displaystyle{ u_r = a^{m - 1} \not\equiv 1 \pmod m }[/math], zatem liczba [math]\displaystyle{ m }[/math] byłaby liczbą złożoną, wbrew założeniu, że jest liczbą pierwszą.

Przypadek b) jest możliwy (np. dla [math]\displaystyle{ m = 41 }[/math] i [math]\displaystyle{ a = 10 }[/math]), ale nie pozwala powiedzieć nic więcej ani o liczbie [math]\displaystyle{ m }[/math], ani o wyrazach ciągu [math]\displaystyle{ (u_i) }[/math], które wszystkie przystają do [math]\displaystyle{ 1 }[/math] modulo [math]\displaystyle{ m }[/math].

W przypadku c) mamy [math]\displaystyle{ u_k \equiv 1 \pmod m }[/math], czyli [math]\displaystyle{ (u_{k - 1})^2 \equiv 1 \pmod m }[/math]. Z twierdzenia L13 wiemy, że musi być albo [math]\displaystyle{ u_{k - 1} \equiv - 1 \pmod m }[/math], albo [math]\displaystyle{ u_{k - 1} \equiv 1 \pmod m }[/math]. Ale drugi przypadek nie może zachodzić, bo założyliśmy, że [math]\displaystyle{ u_k }[/math] jest pierwszym wyrazem ciągu [math]\displaystyle{ (u_i) }[/math], który przystaje do [math]\displaystyle{ 1 }[/math] modulo [math]\displaystyle{ m }[/math]. Zatem musi być [math]\displaystyle{ u_{k - 1} \equiv - 1 \pmod m }[/math].

Co kończy dowód twierdzenia.


Definicja L15
Złożoną liczbę nieparzystą [math]\displaystyle{ m }[/math], która spełnia twierdzenie L14 dla pewnej liczby [math]\displaystyle{ a \in \mathbb{Z} }[/math], będziemy nazywali liczbą silnie pseudopierwszą przy podstawie [math]\displaystyle{ a }[/math] (w skrócie: SPSP([math]\displaystyle{ a }[/math])).


Uwaga L16
Niech [math]\displaystyle{ a }[/math] będzie liczbą całkowitą względnie pierwszą z [math]\displaystyle{ m }[/math] i [math]\displaystyle{ a \in [1, m - 1] }[/math]. Można pokazać, że jeżeli [math]\displaystyle{ m \neq 9 }[/math] jest liczbą nieparzystą złożoną, to co najwyżej [math]\displaystyle{ \tfrac{1}{4} }[/math] liczb [math]\displaystyle{ a }[/math] stanowią liczby silnie pseudopierwsze. Zatem w przypadku liczb silnie pseudopierwszych nie istnieje odpowiednik liczb Carmichaela. Czyli nie istnieją liczby złożone nieparzyste [math]\displaystyle{ m }[/math], dla których twierdzenie L14 byłoby prawdziwe dla wszystkich podstaw [math]\displaystyle{ a }[/math].

Wynika stąd, że jeżeli dla [math]\displaystyle{ k }[/math] różnych liczb [math]\displaystyle{ a }[/math] względnie pierwszych z [math]\displaystyle{ m }[/math] prawdziwe jest twierdzenie L14, to prawdopodobieństwo uznania liczby złożonej [math]\displaystyle{ m }[/math] za pierwszą wynosi [math]\displaystyle{ \left( \tfrac{1}{4} \right)^k }[/math].


Uwaga L17
Wykorzystując twierdzenie L14, możemy napisać w PARI/GP program wykonujący test Millera-Rabina dla ustalonej podstawy.

isPrimeOrSPSP(m, a) =
{
local(d, k, r, x);
r = valuation(m - 1, 2); \\ wykładnik, z jakim liczba 2 występuje w rozwinięciu na czynniki pierwsze liczby m - 1
d = (m - 1) / 2^r;
x = modPower(a, d, m);
if( x == 1 || x == m - 1, return(1) ); \\ x = m - 1 to przypadek k == 0
k = 0;
while( k++ < r,
       x = x^2 % m;
       if( x == m - 1, return(1) );
     );
return(0);
}


Zauważmy, że nie musimy sprawdzać, czy [math]\displaystyle{ \gcd (a, m) = 1 }[/math], bo jeśli tak nie jest, to dla takiej podstawy powyższy test i tak wykryje złożoność liczby [math]\displaystyle{ m }[/math]. Istotnie, jeżeli [math]\displaystyle{ \gcd (a, m) = g \gt 1 }[/math], to rozważając kongruencje z twierdzenia L14

[math]\displaystyle{ a^d \equiv 1 \!\! \pmod{m} }[/math]
[math]\displaystyle{ a^{2^k d} \equiv - 1 \!\! \pmod{m} }[/math]

modulo [math]\displaystyle{ g }[/math], otrzymujemy natychmiast

[math]\displaystyle{ 0 \equiv 1 \!\! \pmod{g} }[/math]
[math]\displaystyle{ 0 \equiv - 1 \!\! \pmod{g} }[/math]

Co jest niemożliwe. Jednak sprawdzenie, czy [math]\displaystyle{ \gcd (a, m) = 1 }[/math] jest wskazane, bo operacja ta jest wykonywana bardzo szybko. A jeśli mieliśmy tyle szczęścia, że [math]\displaystyle{ \gcd (a, m) = g \gt 1 }[/math], to jednocześnie znaleźliśmy dzielnik testowanej liczby [math]\displaystyle{ m }[/math], co w przypadku dużych liczb nie jest rzeczą prostą. Zatem program wykonujący [math]\displaystyle{ k }[/math] testów Millera-Rabina dla przypadkowych podstaw [math]\displaystyle{ a }[/math], gdzie [math]\displaystyle{ a \in [2, m - 2] }[/math], powinien wyglądać tak

PrimeTest(m, k) = 
{
local(a, d, j);
if( m < 2, return(0) );
if( m % 2 == 0, return(m == 2) );  \\ testowana liczba jest liczbą parzystą
setrand(getwalltime());  \\ ustawiamy ziarno (ang. seed) generatora liczb losowych
j = 0;
while( j++ <= k,
       a = random([2, m - 2]);  \\ a jest liczbą losową z przedziału domkniętego [2, m-2]
       d = gcd(a, m);
       if( d > 1, return(0) );  \\ testowana liczba jest liczbą złożoną podzielną przez d
       if( !isPrimeOrSPSP(m, a), return(0) );  \\ testowana liczba jest liczbą złożoną
     );
return(1);  \\ testowana liczba jest prawdopodobnie liczbą pierwszą
}

Testując dla pięciu przypadkowych podstaw, trudno znaleźć liczbę, dla której wartość funkcji PrimeTest() byłaby różna od isprime(). Nam się to nie udało.

forstep(j = 10^6+1, 10^7, 2, if( PrimeTest(j, 5) != isprime(j), print(j) ))


Zadanie L18
Pokazać, że jeżeli liczba [math]\displaystyle{ m }[/math] jest SPSP([math]\displaystyle{ a }[/math]), to jest PSP([math]\displaystyle{ a }[/math]).

Rozwiązanie

Z założenia [math]\displaystyle{ m }[/math] jest SPSP([math]\displaystyle{ a }[/math]), zatem spełniony jest dokładnie jeden z warunków

  • [math]\displaystyle{ a^d \equiv 1 \pmod m }[/math]
  • [math]\displaystyle{ a^{2^k \cdot d} \equiv - 1 \pmod m }[/math], dla pewnego [math]\displaystyle{ k \in [0, r - 1] }[/math]

gdzie [math]\displaystyle{ m - 1 = 2^r \cdot d }[/math], przy czym [math]\displaystyle{ d }[/math] jest liczbą nieparzystą.

Jeżeli spełniony jest pierwszy warunek, to

[math]\displaystyle{ (a^d)^{2^r} \equiv 1 \pmod m }[/math]
[math]\displaystyle{ a^{2^r \cdot d} \equiv 1 \pmod m }[/math]

Czyli [math]\displaystyle{ m }[/math] jest PSP([math]\displaystyle{ a }[/math]).

Jeżeli spełniony jest drugi warunek, to

[math]\displaystyle{ (a^{2^k \cdot d})^{2^{r - k}} \equiv (- 1)^{2^{r - k}} \pmod m }[/math]
[math]\displaystyle{ a^{2^r \cdot d} \equiv 1 \pmod m }[/math]

Czyli [math]\displaystyle{ m }[/math] jest PSP([math]\displaystyle{ a }[/math]).


Przykład L19
Pokażemy, że jeżeli [math]\displaystyle{ m }[/math] jest PSP([math]\displaystyle{ 2 }[/math]), to [math]\displaystyle{ 2^m - 1 }[/math] jest SPSP([math]\displaystyle{ 2 }[/math]).

Z założenia [math]\displaystyle{ m }[/math] jest złożoną liczbą nieparzystą, zatem [math]\displaystyle{ N = 2^m - 1 }[/math] jest złożoną liczbą nieparzystą. Ponieważ [math]\displaystyle{ m }[/math] jest FPSP([math]\displaystyle{ 2 }[/math]), to

[math]\displaystyle{ 2^{m - 1} \equiv 1 \pmod m }[/math]

Wynika stąd natychmiast, że [math]\displaystyle{ 2^{m - 1} - 1 = k m }[/math], gdzie [math]\displaystyle{ k }[/math] jest liczbą nieparzystą. Zatem

[math]\displaystyle{ N - 1 = 2^m - 2 = 2 (2^{m - 1} - 1) = 2 k m }[/math]

Widzimy, że liczba [math]\displaystyle{ 2 }[/math] występuje w pierwszej potędze w rozwinięciu liczby [math]\displaystyle{ N - 1 }[/math] na czynniki pierwsze i łatwo otrzymujemy, że

[math]\displaystyle{ 2^{(N - 1) / 2} = 2^{k m} = (2^m)^k = (2^m - 1 + 1)^k = (N + 1)^k \equiv 1 \pmod N }[/math]

Czyli [math]\displaystyle{ N = 2^m - 1 }[/math] jest SPSP([math]\displaystyle{ 2 }[/math]).


Przykład L20
Tabela zawiera najmniejsze liczby silnie pseudopierwsze dla podstaw [math]\displaystyle{ a }[/math] od [math]\displaystyle{ 2 }[/math] do [math]\displaystyle{ 15 }[/math]

for(a=2, 15, s=1; forstep(m=3, 20000, 2, if( isPrimeOrSPSP(m, a)  &&  !isprime(m), print("a=", a, "  m=", m); s++ ); if( s>5, break() ) ))


Przykład L21
Tabela pokazuje ilość liczb silnie pseudopierwszych dla podstaw [math]\displaystyle{ a }[/math] od [math]\displaystyle{ 2 }[/math] do [math]\displaystyle{ 15 }[/math]

for(a=2, 15, s=0; forstep(k=3, 10^6, 2, if( isPrimeOrSPSP(k, a)  &&  !isprime(k), s++ )); print("a= ", a, "   ", s))

Widzimy, że liczb silnie pseudopierwszych jest znacznie mniej niż liczb pseudopierwszych Fermata.


Uwaga L22
Interesujące i pożyteczne będzie zbadanie najmniejszych liczb silnie pseudopierwszych dla wielu podstaw. Niech badanymi podstawami będą kolejne liczby pierwsze. Najmniejszą liczbę SPSP([math]\displaystyle{ 2 }[/math]) już znamy: [math]\displaystyle{ 2047 }[/math]. Najmniejszą liczbę, która jest jednocześnie SPSP([math]\displaystyle{ 2 }[/math]) i SPSP([math]\displaystyle{ 3 }[/math]) musimy poszukać. Prostym poleceniem

forstep(m=3, 10^7, 2, if( isprime(m), next() ); if( isPrimeOrSPSP(m, 2) && isPrimeOrSPSP(m, 3), print("m=", m) ))

znajdujemy, że [math]\displaystyle{ m = 1373653 }[/math]. Więcej czasu będzie wymagało znalezienie liczby jednocześnie silnie pseudopierwszej względem podstaw [math]\displaystyle{ 2, 3, 5 }[/math]. Poleceniem

forstep(m=3, 10^8, 2, if( isprime(m), next() ); if( isPrimeOrSPSP(m, 2) && isPrimeOrSPSP(m, 3) && isPrimeOrSPSP(m, 5), print("m=", m) ))

znajdujemy, że szukana liczba to [math]\displaystyle{ m = 25326001 }[/math].


Stosując bardziej wyrafinowane metody[6][7] znaleziono wartości liczb silnie pseudopierwszych względem wielu podstaw, które są kolejnymi liczbami pierwszymi, dla większej ilości liczb pierwszych[8]. Dla przykładu

Podane w prawej kolumnie liczby [math]\displaystyle{ m }[/math] są najmniejszymi liczbami jednocześnie silnie pseudopierwszymi względem podstaw [math]\displaystyle{ p_1, \ldots, p_n }[/math]. Zauważmy, że wyniki te mają bardzo praktyczne zastosowanie. Przykładowo, jeśli [math]\displaystyle{ m }[/math] przechodzi test Millera-Rabina dla siedmiu podstaw [math]\displaystyle{ a = 2, 3, 5, 7, 11, 13, 17 }[/math] i jest liczbą mniejszą od [math]\displaystyle{ 3.41 \cdot 10^{14} }[/math], to jest z pewnością liczbą pierwszą.


Uwaga L23
Pomysł przedstawiony w uwadze L22 ma proste uogólnienie. Niech [math]\displaystyle{ A_r = \{ a_1, \ldots, a_r \} }[/math] będzie zbiorem liczb naturalnych większych od [math]\displaystyle{ 1 }[/math]. Możemy teraz szukać takiego zbioru [math]\displaystyle{ A_r }[/math], dla którego najmniejsza liczba silnie pseudopierwsza jednocześnie względem podstaw [math]\displaystyle{ a_1, \ldots, a_r }[/math] będzie największa ze wszystkich rozpatrywanych przypadków.

Dla przykładu przyjmijmy, że

  •    [math]\displaystyle{ r = 3 }[/math]
  •    [math]\displaystyle{ a_k \lt 100 }[/math]
  •    [math]\displaystyle{ a_k }[/math] są liczbami pierwszymi

Okazuje się[7][9], że przy takich założeniach szukanym zbiorem jest [math]\displaystyle{ A_3 = \{ 2, 7, 61 \} }[/math]. Najmniejszą liczbą silnie pseudopierwszą jednocześnie względem podstaw [math]\displaystyle{ 2, 7, 61 }[/math] jest liczba [math]\displaystyle{ 4759123141 }[/math]. Łatwo znajdujemy, że dla [math]\displaystyle{ m \lt 3 \cdot 10^{10} }[/math] istnieje osiem liczb silnie pseudopierwszych jednocześnie względem podstaw [math]\displaystyle{ 2, 7, 61 }[/math]:

[math]\displaystyle{ 4759123141, 8411807377, 11207066041, 11711154457, 12015212653, 18074903681, 19632812033, 27913980641 }[/math]

Korzystając z tego rezultatu możemy napisać prosty program, który rozstrzyga w sposób pewny, czy badana liczba [math]\displaystyle{ m \lt 1.12 \cdot 10^{10} }[/math] jest liczbą pierwszą.

Ponieważ teraz podstawy [math]\displaystyle{ a }[/math] są ustalone, a testowana liczba [math]\displaystyle{ m }[/math] może być dowolna, to musimy wykluczyć sytuacje, że [math]\displaystyle{ \gcd (a, m) \gt 1 }[/math]. Musimy tak zrobić, bo pierwszość liczby [math]\displaystyle{ m }[/math] nie zostanie wykryta, gdy [math]\displaystyle{ \gcd (a, m) = g \gt 1 }[/math] (zobacz L17).

  •  Jeżeli podstawa [math]\displaystyle{ a }[/math] jest liczbą pierwszą, to wystarczy zbadać, czy [math]\displaystyle{ R_a (m) = 0 }[/math]. Jeśli tak, to tylko w przypadku, gdy [math]\displaystyle{ m = a }[/math] liczba [math]\displaystyle{ m }[/math] jest liczbą pierwszą.
  •  Jeżeli podstawa [math]\displaystyle{ a }[/math] jest liczbą złożoną, powiedzmy [math]\displaystyle{ a = p^{\alpha_1}_1 \cdot \ldots \cdot p^{\alpha_s}_s }[/math], to wystarczy zbadać, czy dla pewnego [math]\displaystyle{ i = 1, \ldots, s }[/math] jest [math]\displaystyle{ R_{p_i} (m) = 0 }[/math]. Jeśli tak, to tylko w przypadku, gdy [math]\displaystyle{ m = p_i }[/math] liczba [math]\displaystyle{ m }[/math] jest liczbą pierwszą.

Poniżej przedstawiamy odpowiedni kod w PARI/GP. Zauważmy, że wstępne sprawdzanie pierwszości nieprzypadkowo uwzględnia wszystkie liczby pierwsze [math]\displaystyle{ p \leqslant 61 }[/math]. Wybraliśmy taki zakres, aby zostały objęte podstawy [math]\displaystyle{ 2, 7, 61 }[/math].

MyIsPrime(m) = 
{
if( m < 2, return(0) );
forprime(p = 2, 61, if( m % p == 0, return(m == p) ));
if( m == 4759123141 || m == 8411807377, return(0) );
return( isPrimeOrSPSP(m, 2) && isPrimeOrSPSP(m, 7) && isPrimeOrSPSP(m, 61) );
}



Uzupełnienie

Uwaga L24
W funkcji isPrimeOrSPSP() wykorzystaliśmy zaimplementowane w PARI/GP funkcje:

  • gcd(a, b) – znajduje największy wspólny dzielnik liczb a i b
  • valuation(a, b) – znajduje największą wartość liczby [math]\displaystyle{ r }[/math] taką, że [math]\displaystyle{ b^r \mid a }[/math]

Wprowadzenie tych funkcji pozwoliło zwiększyć czytelność kodu, ale bez trudu możemy sami je napisać:

gcd2(a, b) =
{
local(r);
if( b == 0, return(abs(a)) );
r = a % b;
while( r > 0, a = b; b = r; r = a % b );
return(abs(b));
}


valuation2(a, b) =
{
local(s);
s = 0;
while(a % b == 0, s++; a = a / b);
return(s);
}








Przypisy

  1. W. R. Alford, Andrew Granville and Carl Pomerance, There are Infinitely Many Carmichael Numbers, Annals of Mathematics, 140, (1994), 703-722
  2. Glyn Harman, On the Number of Carmichael Numbers up to x, Bull. London Math. Soc. 37 (2005) 641–650
  3. Glyn Harman, Watt’s Mean Value Theorem and Carmichael Numbers, International Journal of Number Theory, Vol. 4, No. 2 (2008) 241–248
  4. Gary L. Miller, Riemann's Hypothesis and Tests for Primality, Journal of Computer and System Sciences 13, 300-317 (1976)
  5. Michael O. Rabin, Probabilistic Algorithm for Testing Primality, Journal of Number Theory 12, 128-138 (1980)
  6. Carl Pomerance, J. L. Selfridge and Samuel S. Wagstaff, Jr., The Pseudoprimes to 25*10^9, Mathematics of Computation, Vol. 35, No. 151 (1980), 1003-1026
  7. 7,0 7,1 Gerhard Jaeschke, On Strong Pseudoprimes to Several Bases, Mathematics of Computation, Vol. 61, No. 204 (Oct., 1993), 915-926
  8. On-Line Encyclopedia of Integer Sequences, Smallest odd number for which Miller-Rabin primality test on bases <= n-th prime does not reveal compositeness, (A014233)
  9. Wikipedia, Test Millera-Rabina, (Wiki-pl), (Wiki-en)