Prędkość ucieczki do nieskończoności - czyli zabawy z rekurencją (część 2)

Dziś ciekawostka w nawiązaniu do wpisu z dnia 20 października 2015 roku "Liczba PI ukryta w zbiorze Mandelbrota" - przeanalizujmy równanie rekurencyjne dla liczb rzeczywistych:

x_n=\begin{cases}x_{n-1}^2+\frac{1}{4}+\epsilon,&\text{dla}\quad n>0\\0,&\text{dla}\quad n=0\end{cases}

Zbliżanie się do "ostrza" zbioru Mandelbrota

Powyższe wyrażenie powstaje na bazie równania zbioru Mandelbrota z_n=z_{n-1}^2+c (w liczbach zespolonych), jeśli ograniczymy się do prostej rzeczywistej (dlatego użyłem zapisu x_n) i będzie nas interesowało zbliżanie się elementu x_1=\frac{1}{4}+\epsilon do "ostrza" (ang. "cusp") zbioru o współrzędnych (\frac{1}{4},0).

Mandelbrot - Ostrze

Szybkość ucieczki do nieskończoności

Ustalając odpowiednio małe \epsilon>0 decydujemy jak bardzo chcemy się zbliżyć do "ostrza", teraz zadanie polega na znalezieniu pierwszego n, dla którego x_n>=2. Takie minimalne n jest dobrą miarą prędkości ucieczki x_n do nieskończoności w zależności od wybranego \epsilon. Na marginesie dodam, że zbiór Juli dla równania Mandelbrota jest na powyższym obrazku oznaczony kolorem czarnym, który reprezentuje punkty "nieuciekające" do nieskończoności w trakcie nieskończonej iteracji (ta tematyka jest sama w sobie bardzo ciekawa i zapewne kiedyś coś napiszę o atraktorach).

Rekurencja na rekurencji

Zapiszmy nasze zadanie wykorzystując rekurencję w celu poszukiwania rozwiązania:

f(n,\epsilon)=\begin{cases}f(n+1,\epsilon),&\text{dla}\quad x_n<2\\n,&\text{dla}\quad x_n>=2\end{cases}

Niestety zdefiniowaliśmy rekurencję na rekurencji - to zły znak dla wydajności - cóż - sprawdźmy to używając mXparsera.

/* Definicja funkcji rekurencyjnej */
Function x = new Function("x(n,c) = if( n>0, x(n-1,c)^2+0.25+c, 0 )");
Function f = new Function("f(n,c) = if( x(n,c) >= 2, n, f(n+1, c) )", x);

/* Obliczenia i wyświetlenie wyniku */
mXparser.consolePrintln( "c = 0.01" + ", f(0,c) = " + f.calculate(0, 0.01) + ", czas = " + f.getComputingTime() + " s" );
mXparser.consolePrintln( "c = 0.0001" + ", f(0,c) = " + f.calculate(0, 0.0001) + ", czas = " + f.getComputingTime() + " s" );
mXparser.consolePrintln( "c = 0.000001" + ", f(0,c) = " + f.calculate(0,0.000001) + ", czas = " + f.getComputingTime() + " s" );
mXparser.consolePrintln( "c = 0.00000001" + ", f(0,c) = " + f.calculate(0,0.00000001) + ", czas = " + f.getComputingTime() + " s" );

+ wyczekiwany wynik


c = 0.01, f(0,c) = 30.0, czas = 0.224 s
c = 0.0001, f(0,c) = 312.0, czas = 1.532 s
c = 0.000001, f(0,c) = 3140.0, czas = 37.343 s
c = 0.00000001, f(0,c) = 31414.0, czas = 4068.338 s


Wzorzec prędkości ucieczki

Wow - jaki przedziwny wzorzec liczby iteracji wymaganych, aby przekroczyć 2!! Dostajemy coś, co przypomina \pi, jednak wymaga postawienia przecinka w odpowiednim miejscu! Można również zauważyć, że 100-krotne zmniejszenie \epsilon zwiększa niezbędną liczbę iteracji około 10-krotnie. Zmniejszając \epsilon otrzymamy liczbę coraz bardziej przypominającą \pi 🙂

Pozdrowienia,

Mariusz Gromada

Zobacz również:

  1. Polowanie na czarownice - czyli zabawy z rekurencją (część 1)
  2. Naiwny test pierwszości - czyli zabawy z rekurencją (część 3)
  3. Rekurencja pośrednia - czyli zabawy z rekurencją (część 4)
  4. Liczba PI ukryta w zbiorze Mandelbrota
Views All Time
Views All Time
534
Views Today
Views Today
1

Dodaj komentarz

Twój adres email nie zostanie opublikowany. Pola, których wypełnienie jest wymagane, są oznaczone symbolem *