Generowanie paproci Barnsleya i smoka Heighwaya

Choć niektóre fraktale wydają się bardzo skomplikowane, ich konstrukcję można opisać w dość prosty sposób dzięki funkcjom IFSIFSIFS. Przyjrzyjmy się bliżej dwóm przykładom takich fraktali.

Smok Heighwaya

Z historią tego fraktala związanych jest wiele nazwisk. Pracowali nad nim naukowcy z NASA: John Heighway, Bruce Banks i William Harter, a opisany i spopularyzowany zastał przez Martina Gardnera w roku 1967. Był także przedmiotem badań Chandlera Davisa oraz Donalda Knutha (znanego z algorytmu Knutha-Morrisa‑Pratta), którzy odkryli wiele jego cech.

Generowanie tego fraktala opiera się na wykorzystaniu dwóch par funkcji, zgodnie z którymi obliczamy współrzędne kolejnych punktów. Te pary funkcji to:

\begin{array}{l} f_{1} (x) = 0, 5 x - 0, 5 y \\ f_{1} (y) = 0, 5 x + 0, 5 y \end{array}

\begin{array}{l} f_{2} (x) = - 0, 5 x - 0, 5 y + 1 \\ f_{2} (y) = 0, 5 x - 0, 5 y \end{array}

Funkcji tych należy użyć w sposób losowy w równych proporcjach. Jest to metoda Monte Carlometoda Monte CarloMonte Carlo, polegająca na losowaniu dostatecznej liczby danych w taki sposób, że po ich zestawieniu otrzymujemy pewne przybliżenie modelowanego zjawiska. Punktem startowym jest punkt $(0, 0)$ .

Każdy wygenerowany w ten sposób punkt rysujemy w układzie współrzędnych – kształt, który otrzymamy, jest nazywany smokiem Heighwaya. Im więcej punktów wygenerujemy, tym dokładniejszy będzie uzyskany kształt.

Wszystkie przedstawione kroki możemy zapisać za pomocą następującego pseudokodu:

// punkt startowy
x = 0,0
y = 0,0

// liczba generowanych punktów
n = 10000

dla i = 0, 1, 2 ... n - 1 wykonuj:
	// kopia zmiennej x
	x0 = x
    
    g = losuj liczbę z [0, 1]
    
    jeżeli g == 0 wykonaj:
        x = 0,5 * x - 0,5 * y
        y = 0,5 * x0 + 0,5 * y
    w przeciwnym razie wykonaj:
        x = - 0,5 * x - 0,5 * y + 1
        y = 0,5 * x0 - 0,5 * y
    
    rysuj punkt (x, y)

Zwróćmy uwagę na fakt, że potrzebna jest nam kopia zmiennych x, aby podczas generowania zmiennej y (w linijce 16. oraz 19.) nie użyć nowej wartości nadanej zmiennej x.

R10tWEZwGAvSA

Ilustracja przedstawia zakrzywioną linię, która tworzona jest na podstawie kwadratu, ale nigdy go nie zamyka, jedynie przechodzi przez 4 rogi, tak aby następnie przejść do następnego fragmentu. W ten sposób tworzona wypełnia przestrzeń. — Smok Heighwaya
Źródło: William Sikkema, licencja: CC BY-SA 4.0.

Polecenie 1

Napisz program w języku C++, który wygeneruje fraktal smok Heighwaya.

Pokaż podpowiedź

Pseudokod zapisany w C++:

Linia 1. kratka include otwórz nawias ostrokątny iostream zamknij nawias ostrokątny. Linia 2. kratka include otwórz nawias ostrokątny cstdlib zamknij nawias ostrokątny. Linia 3. kratka include otwórz nawias ostrokątny ctime zamknij nawias ostrokątny. Linia 4. using namespace std średnik. Linia 6. int main otwórz nawias okrągły zamknij nawias okrągły otwórz nawias klamrowy. Linia 7. prawy ukośnik prawy ukośnik punkt startowy. Linia 8. double x znak równości 0 kropka 0 średnik. Linia 9. double y znak równości 0 kropka 0 średnik. Linia 11. prawy ukośnik prawy ukośnik liczba generowanych punktów. Linia 12. int n znak równości 10000 średnik. Linia 14. prawy ukośnik prawy ukośnik inicjalizacja generatora liczb losowych. Linia 15. srand otwórz nawias okrągły time otwórz nawias okrągły NULL zamknij nawias okrągły zamknij nawias okrągły średnik. Linia 17. for otwórz nawias okrągły int i znak równości 0 średnik i otwórz nawias ostrokątny n średnik i plus plus zamknij nawias okrągły otwórz nawias klamrowy. Linia 19. double x0 znak równości x średnik prawy ukośnik prawy ukośnik kopia x. Linia 20. double g znak równości rand otwórz nawias okrągły zamknij nawias okrągły procent 2 średnik prawy ukośnik prawy ukośnik losowanie 0 lub 1. Linia 22. if otwórz nawias okrągły g znak równości znak równości 0 zamknij nawias okrągły otwórz nawias klamrowy. Linia 23. x znak równości 0 kropka 5 asterysk x minus 0 kropka 5 asterysk y średnik. Linia 24. y znak równości 0 kropka 5 asterysk x0 plus 0 kropka 5 asterysk y średnik. Linia 25. zamknij nawias klamrowy else otwórz nawias klamrowy. Linia 26. x znak równości minus 0 kropka 5 asterysk x minus 0 kropka 5 asterysk y plus 1 kropka 0 średnik. Linia 27. y znak równości 0 kropka 5 asterysk x0 minus 0 kropka 5 asterysk y średnik. Linia 28. zamknij nawias klamrowy. Linia 30. prawy ukośnik prawy ukośnik wypisanie punktu. Linia 31. cout otwórz nawias ostrokątny otwórz nawias ostrokątny x otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów cudzysłów otwórz nawias ostrokątny otwórz nawias ostrokątny y otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów lewy ukośnik n cudzysłów średnik. Linia 32. zamknij nawias klamrowy. Linia 34. return 0 średnik. Linia 35. zamknij nawias klamrowy.

#include <iostream>
#include <cstdlib>
#include <ctime>
using namespace std;

int main() {
    // punkt startowy
    double x = 0.0;
    double y = 0.0;

    // liczba generowanych punktów
    int n = 10000;

    // inicjalizacja generatora liczb losowych
    srand(time(NULL));

    for (int i = 0; i < n; i++) {

        double x0 = x;   // kopia x
        double g = rand() % 2;  // losowanie 0 lub 1

        if (g == 0) {
            x = 0.5 * x - 0.5 * y;
            y = 0.5 * x0 + 0.5 * y;
        } else {
            x = -0.5 * x - 0.5 * y + 1.0;
            y = 0.5 * x0 - 0.5 * y;
        }

        // wypisanie punktu
        cout << x << " " << y << "\n";
    }

    return 0;
}

Pokaż odpowiedź

Aby zwizualizować obliczone punkty, należy ich współrzędne zapisać do pliku, a potem zaimportować do Excela.

Uruchom program w C++, powstanie plik smok.csv.
Otwórz Excel.
Wybierz Dane → Z pliku tekstowego/CSV.
Wskaż plik smok.csv.
Kliknij Wstaw → Wykres punktowy (XY).
Wybierz wersję bez linii, same punkty.

I gotowe - zobaczysz pięknego, gęstego smoka Heighwaya.

#include <iostream>
#include <fstream>
#include <cstdlib>
#include <ctime>
using namespace std;

int main() {
    ofstream out("smok.csv");
    if (!out) {
        cout << "Nie moge utworzyc pliku.\n";
        return 1;
    }

    // nagłówki do Excela
    out << "x,y\n";

    double x = 0.0;
    double y = 0.0;

    int n = 50000;

    srand(time(NULL));

    for (int i = 0; i < n; i++) {
        double x0 = x;
        int g = rand() % 2;

        if (g == 0) {
            x = 0.5 * x - 0.5 * y;
            y = 0.5 * x0 + 0.5 * y;
        }
        else {
            x = -0.5 * x - 0.5 * y + 1.0;
            y = 0.5 * x0 - 0.5 * y;
        }

        out << x << "," << y << "\n";
    }

    out.close();
    cout << "Zapisano punkty do pliku smok.csv\n";
    return 0;
}

Paproć Barnsleya

Odkrywcą tego fraktala jest Michael Barnsley, który opisał go w książce Fractals Everywhere. Publikacja ta stanowi obszerne kompendium wiedzy na temat fraktali, szczególnie jednak traktuje o generowaniu fraktali z wykorzystaniem funkcji IFS. Autor przedstawił w niej też smoka Heighwaya.

Stworzona przez Barnsleya definicja paproci opiera się na czterech przekształceniach afinicznych. Ponieważ jednak posługiwał się on skomplikowaną terminologią matematyczną (związaną z rachunkiem macierzy, funkcjami wielu zmiennych, liczbami zespolonymi oraz krzywymi hiperbolicznymi), przedstawiony przez nas sposób generowania paproci Barnsleya nie jest oryginalnym zapisem stosowanym przez autora.

Podobnie jak miało to miejsce w przypadku smoka Heighwaya, generowanie tego fraktala opiera się na wykorzystaniu pewnej liczby par funkcji, które należy stosować losowo w określonych proporcjach. W tym wypadku mamy cztery pary funkcji:

\begin{array}{l} f_{1} (x) = 0, 85 x + 0, 04 y \\ f_{1} (y) = - 0, 04 x + 0, 85 y + 1, 6 \end{array}

\begin{array}{l} f_{2} (x) = - 0, 15 x + 0, 28 y \\ f_{2} (y) = 0, 26 x + 0, 24 y + 0, 44 \end{array}

\begin{array}{l} f_{3} (x) = 0, 20 x - 0, 26 y \\ f_{3} (y) = 0, 23 x + 0, 22 y + 1, 6 \end{array}

\begin{array}{l} f_{4} (x) = 0 \\ f_{4} (y) = 0, 16 y \end{array}

Funkcji tych należy użyć losowo według następujących proporcji:

85 : 7 : 7 : 1

Punktem startowym jest punkt $(0, 0)$ .

Pseudokod, którym posłużymy się podczas pisania programu generującego paproć Barnsleya, może wyglądać następująco:

// punkt startowy
x = 0.0
y = 0.0

// liczba generowanych punktów
n = 10000

dla i = 0, 1, 2 ... n-1 wykonuj:
	// kopia zmiennej x
	x0 = x
    
    g = losuj liczbę z [0, 1, ..., 99]
    
    jeżeli g < 1 wykonaj:
    	x = 0;
    	y = 0.16 * y;
    w przeciwnym razie jeżeli g < 86 wykonaj:
    	x = 0.85 * x + 0.04 * y;
    	y = -0.04 * x0 + 0.85 * y + 1.6;
    w przeciwnym razie jeżeli g < 93 wykonaj:
    	x = 0.2 * x - 0.26 * y;
    	y = 0.23 * x0 + 0.22 * y + 1.6;
    w przeciwnym razie wykonaj:
    	x = -0.15 * x + 0.28 * y;
    	y = 0.26 * x + 0.24 * y + 0.44;
    
    rysuj punkt (x, y)

W wyniku tworzenia paproci Barnsleya według pseudokodu możemy uzyskać na przykład taki kształt:

R18GjGseEJBIv

Ilustracja przedstwia liść paproci, którego każdy mniejszy listek to kolejny liść paproci, jedynie mniejszy. Liść zakrzywia się w prawo. — Źródło: Laug, licencja: CC BY-SA 4.0.

Polecenie 2

Napisz program w języku C++, który wygeneruje fraktal paproć Barnsleya.

Pokaż podpowiedź

Pseudokod zapisany w C++:

Linia 1. kratka include otwórz nawias ostrokątny iostream zamknij nawias ostrokątny. Linia 2. kratka include otwórz nawias ostrokątny cstdlib zamknij nawias ostrokątny. Linia 3. kratka include otwórz nawias ostrokątny ctime zamknij nawias ostrokątny. Linia 5. using namespace std średnik. Linia 7. int main otwórz nawias okrągły zamknij nawias okrągły otwórz nawias klamrowy. Linia 8. srand otwórz nawias okrągły time otwórz nawias okrągły 0 zamknij nawias okrągły zamknij nawias okrągły średnik. Linia 10. prawy ukośnik prawy ukośnik punkt startowy. Linia 11. double x znak równości 0 kropka 0 średnik. Linia 12. double y znak równości 0 kropka 0 średnik. Linia 14. const int n znak równości 10000 średnik prawy ukośnik prawy ukośnik liczba generowanych punktów. Linia 16. for otwórz nawias okrągły int i znak równości 0 średnik i otwórz nawias ostrokątny n średnik i plus plus zamknij nawias okrągły otwórz nawias klamrowy. Linia 18. double x0 znak równości x średnik prawy ukośnik prawy ukośnik kopia x. Linia 19. int g znak równości rand otwórz nawias okrągły zamknij nawias okrągły procent 100 średnik prawy ukośnik prawy ukośnik losuj liczbę 0 kropka kropka 99. Linia 21. if otwórz nawias okrągły g otwórz nawias ostrokątny 1 zamknij nawias okrągły otwórz nawias klamrowy. Linia 22. x znak równości 0 średnik. Linia 23. y znak równości 0 kropka 16 asterysk y średnik. Linia 24. zamknij nawias klamrowy. Linia 25. else if otwórz nawias okrągły g otwórz nawias ostrokątny 86 zamknij nawias okrągły otwórz nawias klamrowy. Linia 26. x znak równości 0 kropka 85 asterysk x plus 0 kropka 04 asterysk y średnik. Linia 27. y znak równości minus 0 kropka 04 asterysk x0 plus 0 kropka 85 asterysk y plus 1 kropka 6 średnik. Linia 28. zamknij nawias klamrowy. Linia 29. else if otwórz nawias okrągły g otwórz nawias ostrokątny 93 zamknij nawias okrągły otwórz nawias klamrowy. Linia 30. x znak równości 0 kropka 2 asterysk x minus 0 kropka 26 asterysk y średnik. Linia 31. y znak równości 0 kropka 23 asterysk x0 plus 0 kropka 22 asterysk y plus 1 kropka 6 średnik. Linia 32. zamknij nawias klamrowy. Linia 33. else otwórz nawias klamrowy. Linia 34. x znak równości minus 0 kropka 15 asterysk x plus 0 kropka 28 asterysk y średnik. Linia 35. y znak równości 0 kropka 26 asterysk x plus 0 kropka 24 asterysk y plus 0 kropka 44 średnik. Linia 36. zamknij nawias klamrowy. Linia 38. prawy ukośnik prawy ukośnik wypisanie punktu. Linia 39. cout otwórz nawias ostrokątny otwórz nawias ostrokątny x otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów cudzysłów otwórz nawias ostrokątny otwórz nawias ostrokątny y otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów lewy ukośnik n cudzysłów średnik. Linia 40. zamknij nawias klamrowy. Linia 42. return 0 średnik. Linia 43. zamknij nawias klamrowy.

#include <iostream>
#include <cstdlib>
#include <ctime>

using namespace std;

int main() {
    srand(time(0));

    // punkt startowy
    double x = 0.0;
    double y = 0.0;

    const int n = 10000;   // liczba generowanych punktów

    for (int i = 0; i < n; i++) {

        double x0 = x;          // kopia x
        int g = rand() % 100;   // losuj liczbę 0..99

        if (g < 1) {
            x = 0;
            y = 0.16 * y;
        }
        else if (g < 86) {
            x = 0.85 * x + 0.04 * y;
            y = -0.04 * x0 + 0.85 * y + 1.6;
        }
        else if (g < 93) {
            x = 0.2 * x - 0.26 * y;
            y = 0.23 * x0 + 0.22 * y + 1.6;
        }
        else {
            x = -0.15 * x + 0.28 * y;
            y = 0.26 * x + 0.24 * y + 0.44;
        }

        // wypisanie punktu
        cout << x << " " << y << "\n";
    }

    return 0;
}

Pokaż odpowiedź

Aby zwizualizować obliczone punkty, należy ich współrzędne zapisać do pliku, a potem zaimportować do Excela.

Linia 1. kratka include otwórz nawias ostrokątny iostream zamknij nawias ostrokątny. Linia 2. kratka include otwórz nawias ostrokątny fstream zamknij nawias ostrokątny. Linia 3. kratka include otwórz nawias ostrokątny cstdlib zamknij nawias ostrokątny. Linia 4. kratka include otwórz nawias ostrokątny ctime zamknij nawias ostrokątny. Linia 6. using namespace std średnik. Linia 8. int main otwórz nawias okrągły zamknij nawias okrągły otwórz nawias klamrowy. Linia 9. const int N znak równości 100000 średnik prawy ukośnik prawy ukośnik liczba punktów. Linia 10. ofstream out otwórz nawias okrągły cudzysłów fern kropka csv cudzysłów zamknij nawias okrągły średnik. Linia 12. if otwórz nawias okrągły wykrzyknik out zamknij nawias okrągły otwórz nawias klamrowy. Linia 13. cerr otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów Nie mogę otworzyć pliku fern kropka csv lewy ukośnik n cudzysłów średnik. Linia 14. return 1 średnik. Linia 15. zamknij nawias klamrowy. Linia 17. out otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów x przecinek y lewy ukośnik n cudzysłów średnik. Linia 19. srand otwórz nawias okrągły time otwórz nawias okrągły 0 zamknij nawias okrągły zamknij nawias okrągły średnik prawy ukośnik prawy ukośnik prosta inicjalizacja generatora. Linia 21. double x znak równości 0 kropka 0 przecinek y znak równości 0 kropka 0 średnik. Linia 23. for otwórz nawias okrągły int i znak równości 0 średnik i otwórz nawias ostrokątny N średnik i plus plus zamknij nawias okrągły otwórz nawias klamrowy. Linia 25. prawy ukośnik prawy ukośnik losowa liczba z zakresu 0–1. Linia 26. double r znak równości rand otwórz nawias okrągły zamknij nawias okrągły prawy ukośnik otwórz nawias okrągły double zamknij nawias okrągły RAND podkreślnik MAX średnik. Linia 28. double xn przecinek yn średnik. Linia 30. if otwórz nawias okrągły r otwórz nawias ostrokątny 0 kropka 01 zamknij nawias okrągły otwórz nawias klamrowy. Linia 31. prawy ukośnik prawy ukośnik łodyga. Linia 32. xn znak równości 0 kropka 0 średnik. Linia 33. yn znak równości 0 kropka 16 asterysk y średnik. Linia 34. zamknij nawias klamrowy. Linia 35. else if otwórz nawias okrągły r otwórz nawias ostrokątny 0 kropka 86 zamknij nawias okrągły otwórz nawias klamrowy. Linia 36. prawy ukośnik prawy ukośnik główna część paproci. Linia 37. xn znak równości 0 kropka 85 asterysk x plus 0 kropka 04 asterysk y średnik. Linia 38. yn znak równości minus 0 kropka 04 asterysk x plus 0 kropka 85 asterysk y plus 1 kropka 6 średnik. Linia 39. zamknij nawias klamrowy. Linia 40. else if otwórz nawias okrągły r otwórz nawias ostrokątny 0 kropka 93 zamknij nawias okrągły otwórz nawias klamrowy. Linia 41. prawy ukośnik prawy ukośnik lewy liść. Linia 42. xn znak równości 0 kropka 20 asterysk x minus 0 kropka 26 asterysk y średnik. Linia 43. yn znak równości 0 kropka 23 asterysk x plus 0 kropka 22 asterysk y plus 1 kropka 6 średnik. Linia 44. zamknij nawias klamrowy. Linia 45. else otwórz nawias klamrowy. Linia 46. prawy ukośnik prawy ukośnik prawy liść. Linia 47. xn znak równości minus 0 kropka 15 asterysk x plus 0 kropka 28 asterysk y średnik. Linia 48. yn znak równości 0 kropka 26 asterysk x plus 0 kropka 24 asterysk y plus 0 kropka 44 średnik. Linia 49. zamknij nawias klamrowy. Linia 51. x znak równości xn średnik. Linia 52. y znak równości yn średnik. Linia 54. out otwórz nawias ostrokątny otwórz nawias ostrokątny x otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów przecinek cudzysłów otwórz nawias ostrokątny otwórz nawias ostrokątny y otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów lewy ukośnik n cudzysłów średnik. Linia 55. zamknij nawias klamrowy. Linia 57. out kropka close otwórz nawias okrągły zamknij nawias okrągły średnik. Linia 58. cout otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów Zapisano cudzysłów otwórz nawias ostrokątny otwórz nawias ostrokątny N otwórz nawias ostrokątny otwórz nawias ostrokątny cudzysłów punktów do fern kropka csv lewy ukośnik n cudzysłów średnik. Linia 60. return 0 średnik. Linia 61. zamknij nawias klamrowy.

#include <iostream>
#include <fstream>
#include <cstdlib>
#include <ctime>

using namespace std;

int main() {
    const int N = 100000;      // liczba punktów
    ofstream out("fern.csv");

    if (!out) {
        cerr << "Nie mogę otworzyć pliku fern.csv\n";
        return 1;
    }

    out << "x,y\n";

    srand(time(0));            // prosta inicjalizacja generatora

    double x = 0.0, y = 0.0;

    for (int i = 0; i < N; i++) {

        // losowa liczba z zakresu 0–1
        double r = rand() / (double)RAND_MAX;

        double xn, yn;

        if (r < 0.01) {
            // łodyga
            xn = 0.0;
            yn = 0.16 * y;
        }
        else if (r < 0.86) {
            // główna część paproci
            xn = 0.85 * x + 0.04 * y;
            yn = -0.04 * x + 0.85 * y + 1.6;
        }
        else if (r < 0.93) {
            // lewy liść
            xn = 0.20 * x - 0.26 * y;
            yn = 0.23 * x + 0.22 * y + 1.6;
        }
        else {
            // prawy liść
            xn = -0.15 * x + 0.28 * y;
            yn = 0.26 * x + 0.24 * y + 0.44;
        }

        x = xn;
        y = yn;

        out << x << "," << y << "\n";
    }

    out.close();
    cout << "Zapisano " << N << " punktów do fern.csv\n";

    return 0;
}

Excel / LibreOffice Calc

Otwórz fern.csv w Excelu.
Zaznacz kolumny x i y.
Wstaw wykres:
- Wykres punktowy (XY),
- typ: punktowy bez linii.
- Wykres punktowy (XY),
- typ: punktowy bez linii.
Ustaw:
- brak linii,
- bardzo małe punkty,
- kolor zielony,
- usuń osie.
- brak linii,
- bardzo małe punkty,
- kolor zielony,
- usuń osie.

Efekt będzie trochę pikselowy, ale paproć będzie widoczna.

Słownik

IFS

(z ang. iterated function system – system funkcji iterowanych) funkcje, które pozwalają generować konstrukcje samopodobne

metoda Monte Carlo

metoda modelowania zjawisk fizycznych lub matematycznych za pomocą losowania danych w określony sposób tak, aby ich zestawienie produkowało przybliżenie oczekiwanego rezultatu

Kurz Cantora w C++

Strefa wyzwań

Smok Heighwaya

Paproć Barnsleya

Słownik