Fibonacci (ok. 1170-1250)
Aby rozwiązać ten problem doboru najbardziej złośliwych danych, należy spojrzeć na problem od drugiej strony: jak za pomocą najmniejszych liczb uzyskać z góry zadaną liczbę obrotów pętli? Widać, że aby pętla wykonała się raz, wystarczą dane \( (1, 1) \). Ale te dane nie mogą być wynikiem wcześniejszego obrotu pętli: pamiętajmy, że pierwszy z argumentów jest zawsze dzielnikiem z poprzedniego obrotu pętli, a skoro jest nim jedynka, to reszta nie mogła być równa jeden. Zatem przy więcej niż jednym obrocie pętli najmniejsze dane dla ostatniego obrotu pętli, to \( (2, 1) \) (po czym dostajemy od razu parę \( (1, 0) \) kończącą pętlę). Jakie były zatem dane przedostatnie? Widać, że dzieliliśmy przez \( 2 \) i dostaliśmy resztę \( 1 \). Zatem w poprzednim kroku mogliśmy mieć pary \( (3, 2) \) lub \( (5, 2) \) lub \( (7, 2),\ldots \). Z nich para \( (3, 2) \) jest najoszczędniejsza – są to liczby, które dają najmniejszą sumę, a liczba obrotów naszej pętli jest równa \( 2 \). Jak wyglądała zatem para dwa kroki wstecz? Dzielnik musiał być równy \( 3 \), a reszta równa \( 2 \), więc w grę wchodzą pary \( (5, 3), (8, 3), (11, 3),\ldots \). Z nich znów najoszczędniejsza jest para \( (5, 3) \). Dalej, cofając się i rozumując w analogiczny sposób, dostaniemy pary \( (8,5) \) i kolejno \( (13, 8), (21, 13),\ldots \). Widać zatem, że zawsze najoszczędniej będzie tak dobierać kolejną parę, aby wynik dzielenia był równy \( 1 \), czyli innymi słowy, jeśli w kolejnym obrocie pętli argumenty są równe \( a \) i \( b \), to w poprzednim powinny być równe \( a + b \) i \( a \). Wtedy bowiem \( b \) jest resztą z dzielenia \( a + b \) przez \( a \), natomiast iloraz równy jest 1.
Zbadajmy zatem, jakie najmniejsze argumenty dają zadaną liczbę obrotów pętli.
| Liczba obrotów | a | b |
|---|---|---|
| 0 1 2 3 4 5 6 ... |
1 1 3 5 8 13 21 |
0 1 2 3 5 8 13 |
Liczby, które występują w tabeli, są znane pod nazwą liczb Fibonacciego. Mają one w informatyce duże znaczenie i warto znać podstawowe fakty o nich. Definiuje się je następująco:
Liczby Fibonacciego
\( F_{0}=0 \) \( F_{1}=1 \) \( F_n = F_{n-1}+F_{n-2} \) dla \( n \ge 2 \)
Każda kolejna liczba Fibonacciego jest sumą dwóch swoich poprzedniczek. Parę początkowych wyrazów tego ciągu warto znać na pamięć.
| \( n \) | \( F_{n} \) |
|---|---|
| 0 1 2 3 4 5 6 7 8 ... |
0 1 1 2 3 5 8 13 21 |
Jacques Philippe Marie Binet (1786-1856)
Widzimy więc, że poczynając od drugiego wiersza tabeli (\( n\ge 1 \)), najbardziej złośliwych danych dla algorytmu Euklides2 mamy zawsze \( a = F_{n+2}, b= F_{n+1} \). Aby wymusić \( n \) obrotów pętli musimy wziąć zatem co najmniej \( n + 2 \)-gą i \( n + 1 \)-wszą liczbę Fibonacciego.
Leonhard Euler (1707-1783)
Liczby Fibonacciego rosną szybko. Konkretnie znany jest ogólny wzór na \( n \)-tą liczbę Fibonacciego przypisywany Binetowi, a znany jeszcze na pewno Eulerowi 100 lat przed Binetem.
\( F_n=\frac{1}{\sqrt{5}}((\frac{1+\sqrt{5}}{2})^n-(\frac{1-\sqrt{5}}{2})^n) \)
Dowód tego wzoru można znaleźć w #MatematykaDyskretna-xxx.
Wprowadzając oznaczenia \( \varphi = \frac{1+\sqrt{5}}{2}, \hat{\varphi} = \frac{1-\sqrt{5}}{2} \), otrzymujemy wzór w nieco bardziej czytelnej postaci
\( F_n=\frac{1}{\sqrt{5}}(\varphi^n-\hat{\varphi}^n) \).
Biorąc pod uwagę to, że \( \varphi = 1,618... \), zaś \( \hat{\varphi} = -0,318... \) i to, że w związku z tym składnik \( \hat{\varphi}^{n} \) w naszym wzorze bardzo szybko dąży do zera, możemy łatwo pokazać, że n-ta liczba Fibonacciego jest równa
\( F_n=[\frac{1}{\sqrt{5}}\varphi^n] \),
gdzie przez \( [x] \) oznaczamy zaokrąglenie \( x \), czyli liczbę całkowitą najbliższą danej liczby rzeczywistej \( x \). Zatem n-ta liczba Fibonacciego jest prawie dokładnie równa n-tej potędze \( \varphi \) podzielonej przez \( {\sqrt{5}} \). Zauważmy jeszcze, że skoro tak, to po zlogarytmowaniu obustronnie wzoru przy podstawie \( \varphi \) otrzymujemy wzór
\( \log_{\varphi} F_n = n - log_{\varphi}\sqrt{5} \).
Zatem indeks \( n \) zadanej liczby Fibonacciego \( F_{n} \) wynosi
\( n=\log_{\varphi}F_n + log_{\varphi}\sqrt{5} \).
Zapamiętajmy sobie niezwykle ważny wniosek:
Wniosek
Liczby Fibonacciego rosną wykładniczo szybko. Ich wzrost jest prawie identyczny, jak funkcji wykładniczej \( \frac{1}{\sqrt{5}} \varphi^n \).
Wprowadźmy oznaczenie \( \text{FIB}(a) \) na największą liczbę Fibonacciego mniejszą od \( a \), a przez \( [n]\text{FIB}(a) \) jej indeks. Ponieważ dla dowolnych argumentów \( (a, b) \) liczba obrotów pętli algorytmu Euklides2 nie przekracza liczby obrotów pętli dla argumentów najbardziej złośliwych, czyli \( \text{FIB}(a)=F_{[n]\text{FIB}(a)} \) oraz poprzedniej liczby Fibonacciego \( F_{[n]\text{FIB}(a)-1} \), a liczba obrotów pętli dla tych argumentów jest o \( 2 \) mniejsza, niż indeks większej z nich, to otrzymujemy szacowanie na liczbę obrotów pętli \( M(a, b) \) dla dowolnych argumentów \( a \) i \( b \):
\( M(a,b) \le M(\text{FIB}(a)), \text{FIB}(\text{FIB}(a)) = [n]\text{FIB}(a)-2, \)
skąd
\( M(a, b)+ 2 \le \log_{\varphi}F_n + log_{\varphi}\sqrt{5} \),
a biorąc pod uwagę, że
\( log_{\varphi}\sqrt{5}=1,67... \)
otrzymujemy
\( M(a,b)\le \log_{\varphi}a \log_{\varphi}a \le log_{\varphi}(\text{FIB}(a)) \).
Zapamiętajmy jeszcze jedną ważną prawidłowość.
Uwaga
Indeksy liczb Fibonacciego rosną logarytmicznie wolno w stosunku do wartości tych liczb.
Zatem funkcja [n]FIB(x) rośnie logarytmicznie ze względu na x.
Wracając do naszych danych trzydziestocyfrowych: możemy oszacować liczbę obrotów pętli przez \( log_{\varphi} 10^{30} = 148,33... \). Zatem wykonamy nie więcej niż 150 obrotów pętli, co oczywiście będzie w zasięgu nawet bardzo wolnego komputera. Pamiętajmy, że przy dużych liczbach możemy zapomnieć o wbudowanych w języki programowania procedurach arytmetycznych. O arytmetykę musimy zadbać sami. Własną arytmetyką dużych liczb zajmiemy się później.
Zauważmy pewną niedogodność. W algorytmie Euklides2 jeden krok pętli jest nieco trudniejszy. W poprzednim algorytmie Euklides1 mieliśmy tylko porównywanie liczb i ich odejmowanie. Tutaj musimy zaprogramować dzielenie z resztą. Jest to nie tylko trudniejsze, ale i ogólnie wolniejsze od odejmowania. Jeśli zdecydujemy się na algorytm szkolny dzielenia słupkowego, to trzeba będzie wykonać całą serię obliczeń realizujących kolejne kroki wyznaczania cyfr ilorazu. Każdy z tych kroków wymaga wyznaczenia stosownej cyfry wyniku, przemnożenia jej wartości przez dzielnik, a następnie odjęcia od fragmentu dzielnej tego wyniku. Robimy to z grubsza tyle razy, ile cyfr ma iloraz.
Jeżeli za miarę wielkości liczby przyjmiemy długość jej reprezentacji w systemie pozycyjnym (czyli liczbę cyfr), to o ile porównywanie oraz odejmowanie można zrobić w czasie proporcjonalnym do długości tej reprezentacji, to dzielenie może wymagać kwadratowego czasu. Niech długość dzielnej wynosi \( n \). Zauważmy, że jeśli dzielnik jest o połowę krótszy od dzielnej, (czyli mając długość \( n/2 \) jest mniej więcej równy pierwiastkowi kwadratowemu z dzielnej), to iloraz będzie miał długość podobną jak dzielnik, czyli \( n/2 \) i tyle razy będzie się musiała wykonać zasadnicza pętla algorytmu dzielącego, bo tyle cyfr trzeba wyznaczyć. Z kolei wyznaczenie każdej cyfry ilorazu wymaga odjęcia jakiejś niewielkiej wielokrotności dzielnika, a więc liczby również z grubsza \( n/2 \)-cyfrowej. A odejmowanie jest proporcjonalnie kosztowne do długości argumentów. Łącznie zatem \( n/2 \) cyfr ilorazu razy \( n/2 \) jednocyfrowych kroków przy odejmowaniu daje nam łącznie \( {n^2/4} \), więc kwadratowo wiele w stosunku do \( n \).
Liczba obrotów głównej pętli algorytmu jest też proporcjonalna do \( n \), bo w dowolnym systemie pozycyjnym liczba cyfr jest proporcjonalna do logarytmu z danej liczby przy podstawie będącej bazą systemu, czyli również proporcjonalna do logarytmu przy podstawie \( \varphi \), bo logarytmy o różnych podstawach różnią się od siebie tylko o czynnik stały.
Jeśli skupimy się na operacjach na pojedynczych cyfrach, to łączna liczba wszystkich operacji będzie rzędu co najwyżej \( n^3 \). W rzeczywistości możemy się pokusić o przypuszczenie, że będzie to nawet mniej. Zauważmy bowiem, że złośliwe dane dla głównej pętli, to kolejne liczby Fibonacciego, a te długością różnią się co najwyżej o 1, natomiast złośliwe dane dla algorytmu dzielenia z resztą, to dane różniące się długością dwukrotnie; liczby Fibonacciego można podzielić w czasie liniowym, a nie kwadratowym. I to liczby Fibonacciego będą się pojawiały cały czas, w każdym kroku algorytmu. Jeżeli zatem zaczniemy od pary kolejnych liczb Fibonacciego, to \( O(n) \)-krotnie wykonamy dzielenie kosztujące \( O(n) \), co nam da \( O(n^{2}) \). Jeśli natomiast będziemy się starali wywindować koszt dzielenia, to długości kolejnych par powinny być równe \( (n,n/2), (n/2,n/4 ), (n/4, n/8),\ldots \). Ale wtedy łączny koszt dzieleń będzie równy
\( \frac{n^2}{4}+\frac{n^2}{16}+ \frac{n^2}{64}+\ldots =O(n^2) \).
Zatem obie skrajności: złośliwe dane dla zewnętrznej i dla wewnętrznej pętli dają koszt kwadratowy. Ale czy nie można wypośrodkować złośliwości tak, aby uzyskać jednak złożoność sześcienną?
Czy można tak dobrać dane, żeby wymusić sześcienną złożoność algorytmu Euklides2?
Odpowiedź
NIE. Złożoność algorytmu Euklides2 jest jednak kwadratowa.
Powstaje pytanie, czy nie dałoby się znaleźć takiego algorytmu znajdowania największego wspólnego dzielnika tak, aby zachowując złożoność kwadratową korzystać z łatwiejszych operacji niż dzielenia z resztą. Rozwiązanie jest zaskakująco proste, jeśli zauważymy parę dość oczywistych faktów. Będziemy rozważać parzystość argumentów i redukować problem ze względu na tę właśnie własność. Oznaczmy zbiór liczb parzystych przez \( P \). Kluczem do algorytmu jest spostrzeżenie, że jeśli jedna liczba jest parzysta, a druga nie, to największy wspólny dzielnik nie zmieni się, jeśli parzysty argument podzielimy przez 2.
\( (a,b) = \begin{cases} a & \mbox{jeśli }b=0 \\ 2(\frac{a}{2},\frac{b}{2}) & \mbox{jeśli }a,b\in P \\ (\frac{a}{2},b) & \mbox{jeśli }a\in P, b\notin P \\ (a,\frac{b}{2}) & \mbox{jeśli }a\notin P, b\in P \\ (a-b,b) & \mbox{jeśli }a,b\notin P \end{cases} \)
Oczywiście, podobnie jak poprzednio, dbamy zawsze o to, żeby pierwszy argument nie był mniejszy od drugiego i w razie czego zamieniamy je miejscami.
Euklides 3
<b>Read</b>(a,b); //Wczytujemy a i b, zakładając że użytkownik wie, że a>=b, a+b>0 wynik:=1; <b>while</b> b > 0 <b>do</b> <b>begin</b> <b>if</b> a< b <b>then</b> zamień(a,b); //po wykonaniu tej instrukcji zawsze a>=b <b>if</b> parzyste(a) <b>and</b> parzyste(b) <b>then</b> <b>begin</b> wynik:=wynik*2 a:=a div 2; b:=b div 2 <b>end</b> <b>else</b> // w przeciwnym razie <b>if</b> parzyste(a) <b>and</b> <b>not</b> parzyste(b) <b>then</b> a:=a div 2 <b>else</b> <b>if</b> <b>not</b> parzyste(a) <b>and</b> parzyste(b) <b>then</b> b:=b div 2 <b>else</b> a:=a-b <b>end</b>; wynik:=wynik*a; Write(a)
W kodzie tym posługujemy się predykatem parzyste(x), który przyjmuje wartość true (prawda), jeśli x jest parzyste i false (fałsz) jeśli jest nieparzyste. Operacja div daje nam wynik dzielenia całkowitego (z ucięciem reszty).
Przyjrzyjmy się naszemu algorytmowi pod kątem liczby obrotów pętli. Zauważmy, że w każdym obrocie we wszystkich przypadkach, poza sytuacją obu argumentów nieparzystych, co najmniej jeden z argumentów jest połowiony. Natomiast w przypadku obu argumentów nieparzystych jeden z argumentów stanie się parzysty w wyniku odejmowania i w następnym obrocie pętli zostanie podzielony przez 2. Zatem co najmniej raz na dwa obroty pętli, co najmniej jeden z argumentów jest dzielony przez 2. Ale dzieleń przez 2 można wykonać tylko logarytmicznie dużo.
Uwaga
W dziedzinie całkowitoliczbowej możemy spojrzeć na logarytm jak na liczbę możliwych dzieleń przez 2, zanim dojdziemy do jedynki. Liczba ta jest równa podłodze z logarytmu rzeczywistego przy podstawie 2.
Zatem łączna liczba obrotów pętli nie przekracza \( 2(\log_2a + \log_2b) \). Co się dzieje w każdym obrocie pętli? Każda z operacji ma złożoność liniową. Sprawdzenie parzystości liczby wymaga zajrzenia do ostatniej cyfry. Sprawdzenie, czy liczba jest równa zero wymaga przejrzenia w najgorszym razie wszystkich jej cyfr jednokrotnie. Dzielenie przez 2 i mnożenie przez dwa, podobnie jak odejmowanie, mają złożoność liniową (czyli proporcjonalną do logarytmu z wartości a i b). Zatem złożoność algorytmu na poziomie operacji na cyfrach jest kwadratowa, czyli w tym przypadku proporcjonalna do kwadratu z logarytmu \( a \).