Hlavní navigace

Rozšířené Markus-Ljapunovovy fraktály

Pavel Tišnovský 6. 3. 2007

V dnešním článku bude ukázán způsob vytváření takzvaných rozšířených Markus-Ljapunovových fraktálů, u kterých je možné aplikovat více funkcí růstu. Také si na příkladu předvedeme, jak je možné použít základní i rozšířené Markus-Ljapunovovy fraktály pro tvorbu trojrozměrných modelů terénu pomocí výškových polí.

Obsah

1. Rozšířené Markus-Ljapunovovy fraktály – přidání dalších funkcí do modelu růstu
2. Funkce používané v rozšířených Markus-Ljapunovových fraktálech
3. První demonstrační příklad – vykreslení rozšířených Markus-Ljapunovových fraktálů
4. Tvorba fraktálů s řídicím řetězcem „aab“ až „azb“ a „aaab“ až „azzb“
5. Použití Markus-Ljapunovových fraktálů pro vytváření 3D modelů terénů
6. Druhý demonstrační příklad – 3D terén vykreslený v POV-Rayi
7. Odkazy na další informační zdroje
8. Obsah dalšího pokračování tohoto seriálu

1. Rozšířené Markus-Ljapunovovy fraktály – přidání dalších funkcí do modelu růstu

V předchozí části tohoto seriálu jsme si ukázali způsob výpočtu a vykreslení dvojrozměrných Markus-Ljapunovových fraktálů. Přechod od jednorozměrného bifurkačního diagramu (ve kterém je použita pouze jedna nezávislá proměnná – míra růstu populace Gr) k dvojrozměrnému obrázku byl nápaditě vyřešen tak, že se místo jedné funkce růstu použily funkce dvě, které se v každé iteraci pravidelně střídaly. První funkce se podle konvence jmenuje a, druhá b. V každé funkci byla použita jiná míra růstu, kterou jsme označili symboly startx a starty. Plošný fraktální objekt vykreslený tímto způsobem je zobrazený na prvním obrázku.

fractals70_1
Obrázek 1: Klasický Markus-Ljapunovův fraktál vykreslený programem FractInt

Jedno z mnoha rozšíření (tj. vlastně přidání dalšího stupně volnosti při vytváření fraktálního objektu) tohoto způsobu výpočtu spočívá v tom, že se specifikuje vlastní posloupnost funkcí a a b, které jsou použity při iterativním výpočtu Ljapunovova exponentu. Specifikace vybrané posloupnosti funkcí může být zadána buď pomocí řetězce („ab“, „bba“, „abba“) nebo decimálním číslem, které vznikne po převodu binárních symbolů (0=a, 1=b) do sekvence bitů, které můžeme chápat jako binárně zapsané celé číslo. Toto rozšíření jako první použili Markus a Hess. Fraktální objekt vytvořený tímto způsobem je zobrazený na druhém obrázku.

fractals70_2
Obrázek 2: Markus-Ljapunovův fraktál vykreslený programem FractInt s odlišnou posloupností funkcí ab

Markus a Hess však navrhli ještě další způsob rozšíření iterativního výpočtu – místo pouhých dvou funkcí růstu a a b je možné při výpočtech vystřídat libovolné množství dalších funkcí, kterým bývají přiřazeny další znaky malé abecedy, tj. d–z. Tyto funkce provádějící výpočet xn+1=f(xn, Gr) musí splňovat pouze jedno kritérium: pro malé hodnoty xn by mělo platit, že xn+1>xn a naopak pro velké hodnoty xn by měla platit opačná nerovnost, tj. xn+1<xn. Tím je zajištěno, že nebude docházet k divergenci (může však docházet k oscilacím s rostoucí amplitudou) a také se tím splňuje podmínka pro chaotické chování celého výpočtu (není to však podmínka dostatečná). Funkcí splňujících toto kritérium je velké množství, některé z nich si vypíšeme v následující kapitole.

fractals70_3
Obrázek 3: Markus-Ljapunovův fraktál s řídicím řetězcem „xotogbg“

2. Funkce používané v rozšířených Markus-Ljapunovových fraktálech

V následující tabulce jsou vypsány některé vhodné rozšiřující funkce, které mohou být použity při tvorbě rozšířených Markus-Ljapunovových fraktálů. Přiřazení těchto funkcí znakům malé abecedy tak, aby bylo možné jednoduše zapsat řídicí řetězec, je převzato z diplomové práce Martina Fědora, mimochodem člena demogrupy Vectors, známé jak Amigistům, tak i PCčkářům (nazdar Defore!). Symbol startx odpovídá míře růstu Gr1 (ten se mění v horizontálním směru), symbol starty pak míře růstu Gr2 (měněný ve směru vertikálním). Zápis funkcí uvedených v tabulce přímo odpovídá zdrojovému kódu prvního demonstračního příkladu, který bude podrobněji popsán ve třetí kapitole.

Znak Funkce Nová hodnota r
‚a‘ z=(startyz)(1.0-z) r=starty
‚b‘ z=(startxz)(1.0-z) r=startx
‚c‘ z=(startyz)(1.0-zz) r=starty
‚d‘ z=(startxz)(1.0-zz) r=startx
‚e‘ z=(startyz)(1.0-sqrt(z)) r=starty
‚f‘ z=(startxz)(1.0-sqrt(z)) r=startx
‚g‘ z=(startyz)(1.0-sin(z)) r=starty
‚h‘ z=(startxz)(1.0-sin(z)) r=startx
‚i‘ z=(startyz)(1.0-cos(z)) r=starty
‚j‘ z=(startxz)(1.0-cos(z)) r=startx
‚k‘ z=(startyz)(1.0-zzz) r=starty
‚l‘ z=(startxz)(1.0-zzz) r=startx
‚m‘ z=(startyz)(1.0-sqrt(zzz)) r=starty
‚n‘ z=(startxz)(1.0-sqrt(zzz)) r=startx
‚o‘ z=(startyz)(1.0-sin(zz)) r=starty
‚p‘ z=(startxz)(1.0-sin(zz)) r=startx
‚q‘ z=(startyz)(1.0-cos(zz)) r=starty
‚r‘ z=(startxz)(1.0-cos(zz)) r=startx
‚s‘ z=(startyz)(1.0-sin(zzz)) r=starty
‚t‘ z=(startxz)(1.0-sin(zzz)) r=startx
‚u‘ z=(startyz)(1.0-cos(zzz)) r=starty
‚v‘ z=(startxz)(1.0-cos(zzz)) r=startx
‚w‘ z=(startyz)(1.0-sin(sqrt(zzz))) r=starty
‚x‘ z=(startxz)(1.0-sin(sqrt(zzz))) r=startx
‚y‘ z=(startyz)(1.0-cos(sqrt(zzz))) r=starty
‚z‘ z=(startxz)(1.0-cos(sqrt(zzz))) r=startx

fractals70_4
Obrázek 4: Markus-Ljapunovův fraktál s řídicím řetězcem „xxmnsn“

3. První demonstrační příklad – vykreslení rozšířených Markus-Ljapunovových fraktálů

Dnešní první demonstrační příklad, jehož zdrojový kód můžete získat zde, slouží k výpočtu a následnému vykreslení rozšířených Markus-Ljapunovových fraktálů. Céčková funkce, která výpočet těchto fraktálů provádí, je rozšířena o 24 nových funkcí růstu, které byly vypsány v předchozí kapitole. V řídicím řetězci, který je této funkci předáván, se mohou vyskytovat všechny znaky malé abecedy, tj. „a“„z“. Všechny ostatní znaky jsou ignorovány a výpočet žádným zásadním způsobem neovlivňují (až na zbytečně proběhlé iterace). Obrázky uvedené v tomto článku jsou vytvořeny právě tímto demonstračním příkladem. Ústředním bodem celého programu je funkce nazvaná lyapunov(), jejíž zdrojový kód má tvar:

//-----------------------------------------------------------------------------
// Vypocet Markus-Lyapunovych fraktalu
//-----------------------------------------------------------------------------
void recalcLyapunov(pixmap *pix,                // pracovni pixmapa
                    int    palette,             // cislo barvove palety
                    int    maxiter,             // maximalni pocet iteraci
                    double xmin,                // meze fraktalu v rovine
                    double ymin,
                    double xmax,
                    double ymax,
                    const char *order)          // poradi aplikace funkci
{
    int x, y, i, j;                             // pocitadla smycek
    double startx, starty;
    double z, r;
    double total;

    starty=ymin;
    // pro vsechny radky pixmapy
    for (y=0; y<pix->height; y++) {
        startx=xmin;
        printf("%d\t", y);
        // pro vsechny pixely na radku
        for (x=0; x<pix->width; x++) {
            startx+=(xmax-xmin)/pix->width;     // pocatecni hodnota prvni funkce
            z=0.5;
            total=0.0;
            j=0;                                // pocitadlo pro prisup k retezci
            // iteracni smycka
            for (i=0; i<maxiter; i++) {
                switch (order[j]) {             // vyber iteracni funkce
                    case 'a': r=starty; z=(starty*z)*(1.0-z);                break;
                    case 'b': r=startx; z=(startx*z)*(1.0-z);                break;
                    case 'c': r=starty; z=(starty*z)*(1.0-z*z);              break;
                    case 'd': r=startx; z=(startx*z)*(1.0-z*z);              break;
                    case 'e': r=starty; z=(starty*z)*(1.0-sqrt(z));          break;
                    case 'f': r=startx; z=(startx*z)*(1.0-sqrt(z));          break;
                    case 'g': r=starty; z=(starty*z)*(1.0-sin(z));           break;
                    case 'h': r=startx; z=(startx*z)*(1.0-sin(z));           break;
                    case 'i': r=starty; z=(starty*z)*(1.0-cos(z));           break;
                    case 'j': r=startx; z=(startx*z)*(1.0-cos(z));           break;
                    case 'k': r=starty; z=(starty*z)*(1.0-z*z*z);            break;
                    case 'l': r=startx; z=(startx*z)*(1.0-z*z*z);            break;
                    case 'm': r=starty; z=(starty*z)*(1.0-sqrt(z*z*z));      break;
                    case 'n': r=startx; z=(startx*z)*(1.0-sqrt(z*z*z));      break;
                    case 'o': r=starty; z=(starty*z)*(1.0-sin(z*z));         break;
                    case 'p': r=startx; z=(startx*z)*(1.0-sin(z*z));         break;
                    case 'q': r=starty; z=(starty*z)*(1.0-cos(z*z));         break;
                    case 'r': r=startx; z=(startx*z)*(1.0-cos(z*z));         break;
                    case 's': r=starty; z=(starty*z)*(1.0-sin(z*z*z));       break;
                    case 't': r=startx; z=(startx*z)*(1.0-sin(z*z*z));       break;
                    case 'u': r=starty; z=(starty*z)*(1.0-cos(z*z*z));       break;
                    case 'v': r=startx; z=(startx*z)*(1.0-cos(z*z*z));       break;
                    case 'w': r=starty; z=(starty*z)*(1.0-sin(sqrt(z*z*z))); break;
                    case 'x': r=startx; z=(startx*z)*(1.0-sin(sqrt(z*z*z))); break;
                    case 'y': r=starty; z=(starty*z)*(1.0-cos(sqrt(z*z*z))); break;
                    case 'z': r=startx; z=(startx*z)*(1.0-cos(sqrt(z*z*z))); break;
                    default:  r=0;      z=z;                                 break;
                }
                j++;
                if (j>strlen(order)) j=0;       // postupne projit celym retezcem order

                if (i>maxiter/2) {              // ve druhe polovine cyklu se pocita Lyapunuv exponent
                    if (fabs(r-2.0*r*z)>0.0)    // = takzvana filtrace
                        total+=log(fabs(r-2.0*r*z))/log(2.0);
                }
                if (z>500.0) break;
            }
            //total=0.7*total;
            total=255.0*total/maxiter;
            if (total<0) {                      // obarvit pixely se zapornym exponentem
                unsigned char r, g, b;
                if (total<-255) total=-255;
                mapPalette(palette, -total, &r, &g, &b); // barvu ziskat z vybrane palety
                putpixel(pix, x, y, r, g, b);
            }
            else                                // "chaoticke" pixely jsou cerne
                putpixel(pix, x, y, 0, 0, 0);
        }
        // prechod na dalsi obrazovy radek
        starty+=(ymax-ymin)/pix->height;
    }
} 

fractals70_5
Obrázek 5: Markus-Ljapunovův fraktál s řídicím řetězcem „aadb“

4. Tvorba fraktálů s řídicím řetězcem „aab“ až „azb“ a „aaab“ až „azzb“

Vzhledem k tomu, že vyhledávání vhodného obsahu řídicího řetězce může být zdlouhavá činnost, můžeme si celou situaci zjednodušit tím, že si vytvoříme galerii fraktálních objektů pro vhodně zvolené řídicí řetězce. V praxi se ukazuje, že bývá vhodné, aby řídicí řetězec začínal symbolem a, končil symbolem b a obsahoval jeden až čtyři další znaky. Pokud tyto požadavky vezmeme v úvahu, můžeme vytvořit všechny fraktály definované řídicími řetězci „aab“„azb“ (těchto řetězců je 26), popř. všechny fraktály definované řídicími řetězci „aaab“„azzb“ (zde se již dostaneme k číslu 676 možných kombinací). Vygenerování prvních 26 řetězců zajistí následující fragment kódu (ten je možné použít v prvním demonstračním příkladu po odkomentování příslušných řádků ve funkci main()):

for (i='a'; i<='z'; i++) {
    int palette=RND(10);
    str[0]='a';
    str[1]=i;
    str[2]='b';
    // ukonceni retezce
    str[3]=0;
    puts(str);
    puts("*** lyapunov ***");
    puts("processing:");
    // vypocet leve horni ctvrtiny fraktalu
    recalcLyapunov(pix, palette, 500, 0.0, 0.0, 4.0, 3.0, str);
    // pridat k retezci koncovku a ulozit obrazek fraktalu do souboru
    strcat(str, ".tga");
    savePixmap(pix, str);
    puts("\ndone!");
} 

fractals70_6
Obrázek 6: Markus-Ljapunovův fraktál s řídicím řetězcem „aamb“

O vytvoření dalších 676 kombinací se postará následující fragment kódu, který je taktéž možné použít v prvním demonstračním příkladu po odkomentování příslušných řádků, které se nacházejí ve funkci main():

for (i='a'; i<='z'; i++) {
    for (j='a'; j<='z'; j++) {
        int palette=RND(10);
        str[0]='a';
        str[1]=i;
        str[2]=j;
        str[3]='b';
        // ukonceni retezce
        str[4]=0;
        puts(str);
        puts("*** lyapunov ***");
        puts("processing:");
        // vypocet leve horni ctvrtiny fraktalu
        recalcLyapunov(pix, palette, 500, 0.0, 0.0, 4.0, 3.0, str);
        // pridat k retezci koncovku a ulozit obrazek fraktalu do souboru
        strcat(str, ".tga");
        savePixmap(pix, str);
        puts("\ndone!");
    }
} 

fractals70_7
Obrázek 7: Markus-Ljapunovův fraktál s řídicím řetězcem „abfb“

5. Použití Markus-Ljapunovových fraktálů pro vytváření 3D modelů terénů

Markus-Ljapunovovy fraktály je možné mimo jejich přímého vykreslení použít i při vytváření trojrozměrných modelů terénu. Zde využijeme takzvaná výšková pole, která jsme si popsali v předchozích částech tohoto seriálu (v sekci, kde jsme se věnovali obrázkům plasmy). Výšková pole představovaná vypočtenými Markus-Ljapunovovými fraktály mohou sloužit k vytvoření hor s ostrými vrcholky, hřebeny a mnoha údolími. Takto vytvořené hory jsou zcela nepodobné horám vypočteným pomocí plasmy či Mandelbrotových množin, v praxi se však může provést kombinace například s obrázky plasmy, které do modelu vnesou jistou dávku náhodnosti. Konkrétní model složitější scény s výškovým polem bude ukázán v šesté kapitole.

6. Druhý demonstrační příklad – 3D terén vykreslený v POV-Rayi

V dnešním druhém demonstračním příkladu si ukážeme, jakým způsobem je možné použít obrázky s fraktály při renderingu složitější trojrozměrné scény s 3D povrchem. Jako zdroj pro výškové pole byl zvolen Markus-Ljapunovův fraktál zobrazený na obrázku číslo 8. V tomto obrázku odpovídá úroveň šedi jednotlivých pixelů přímo výšce terénu v daném místě. Takto vytvořený obrázek byl po nezbytných úpravách (rozmazání ostrých hran) předán programu POV-Ray, kde sloužil jako výškové pole při renderingu scény vyobrazené na devátém obrázku. Popis celé trojrozměrné scény v popisném jazyce POV-Raye je uveden níže. Ostré vrcholky jsou způsobeny nedokonalým vyhlazením obrázku, které mělo sloužit k odstranění aliasu. Na výsledném modelu se také negativně projevuje fakt, že výška je ve výškovém poli reprezentována pouze osmi bity; ideální by bylo použití šestnáctibitové hloubky

fractals70_8
Obrázek 8: Markus-Ljapunovův fraktál jako zdroj pro výškové pole

fractals70_9
Obrázek 9: Výsledná 3D scéna s výškovým polem
// Zapis sceny urcene pro rendering trojrozmerneho povrchu
// vytvoreneho pomoci Markus-Lyapunovova fraktalu

#include "shapes.inc"       // definice tvaru
#include "colors.inc"       // definice barev
#include "textures.inc"     // definice textur
#include "skies.inc"        // textury mraku apod.

// specifikace pohledu do sceny - kamery
camera  {
    location <-430,400,-630>
    up       <   0,  1,   0>        // orientace kamery
    look_at  <   0,100,   0>
}

//hlavni svetelny zdroj
light_source {
    <-500.0, 300.0, -60.0>
    color White
}

//vedlejsi svetelny zdroj
light_source {
    <-50.0, 800.0, -60.0>
    color red 0.5 blue 0.5 green 0.5
}

// deklarace 3D povrchu
height_field {
    png "lyap02.png" smooth             // soubor s obrazkem
    water_level 0.0
    pigment { color rgbf <0.0, 1.0, 0.2, 0.0> }
    finish  {
        //crand   0.025                   // dither
        ambient 0.2                     // tmave stiny
        diffuse 0.8
        phong   0.75                    // shiny
        phong_size 100
        specular 1.0
        roughness 0.005
    }
    scale <1400,300,1400>               // zmena meritka hory
    translate<-200, -30,-600>             // vycentrovat
}

// zobrazeni koule s nanesenymi mraky
// (textura obsahuje pruhledne casti)
sphere {<0,0,0> 1400
    texture {T_Cloud2
        scale 250
    }
}

// dalsi otexturovana koule, tentokrat s texturou mlhy
// (textura obsahuje pruhledne casti)
sphere {<0,0,0> 1250
    texture {T_Cloud1
        scale 250
    }
}

// mezi mraky se prosvecuje spinave modra
// obloha nanesená na kouli o vetsim polomeru
// nez koule s texturami mraku
sphere {<0,0,0> 1500
    pigment {
        color rgbf <0.2, 0, 0.8, 0>
    }
}

// hladina oceanu
plane { y, 40
    pigment {                           // barva hladiny
        Blue
    }
    normal  {                           // vlnky
        waves 0.05
        frequency 5000
        scale 15000
    }
    finish {                            // odrazy
        reflection 0.8
    }
}

//dodelavka hladiny oceanu
plane { y, 20
    pigment {                           // barevna vrstva pod hladinou
        Blue
    }
    finish  {
        ambient 1.0
        diffuse 0.0
    }
}

// finito 

7. Odkazy na další informační zdroje

  1. Stránky raytraceru POV-Ray:
    http://www.po­vray.org
  2. Stránky programu FractInt:
    http://www.frac­tint.org
  3. Lyapunov exponent (Wikipedia):
    http://en.wiki­pedia.org/wiki/Ly­apunov_exponent
  4. Lyapunov fractal (Wikipedia):
    http://en.wiki­pedia.org/wiki/Ly­apunov_fractal
  5. The Lyapunov exponent:
    http://monet.u­nibas.ch/~elmer/pen­dulum/lyapexp­.htm
  6. Lyapunov Fractals:
    http://spanky­.triumf.ca/www/frac­tint/lyapunov_ty­pe.html
  7. Markus-Lyapunov Fractals:
    http://perso.o­range.fr/char­les.vassallo/en/ly­ap_art/lyapdoc­.html

fractals70_a
Obrázek 10: Markus-Ljapunovův fraktál s řídicím řetězcem „bhawdy“

fractals70_b
Obrázek 11: Markus-Ljapunovův fraktál s řídicím řetězcem „psfadp“

8. Obsah dalšího pokračování tohoto seriálu

V následujícím pokračování tohoto seriálu si ukážeme některé struktury, které sice na první pohled vypadají jako fraktály, ve skutečnosti se však o fraktály nejedná, protože není splněna základní vlastnost fraktálů – soběpodobnost, tj. nezávislost základních geometrických parametrů vytvářených objektů na změně měřítka.

fractals70_c
Obrázek 12: Markus-Ljapunovův fraktál s řídicím řetězcem „wgbusbmb“
Našli jste v článku chybu?

9. 3. 2007 14:45

No, ja ted pracovne take delam na DB nude :-) A dokonce i sefy mame podobne, ale jinak celkem dostavam pri praci volnost. Hint: pokud sef moc mluvi o UML/XML/Jave a "enterprise resenich", s velkou pravdepodobnosti do technologii nevidi a mate moznost ho ovlivnit (aby si myslel, ze to byl vlastne jeho napad) :-)
Měšec.cz: Zdravotní a sociální pojištění 2017: Připlatíte

Zdravotní a sociální pojištění 2017: Připlatíte

Lupa.cz: Co se dá měřit přes Internet věcí

Co se dá měřit přes Internet věcí

Lupa.cz: Propustili je z Avastu, už po nich sahá ESET

Propustili je z Avastu, už po nich sahá ESET

Podnikatel.cz: Přehledná titulka, průvodci, responzivita

Přehledná titulka, průvodci, responzivita

Root.cz: Pinebook: linuxový notebook za 89 dolarů

Pinebook: linuxový notebook za 89 dolarů

120na80.cz: Horní cesty dýchací. Zkuste fytofarmaka

Horní cesty dýchací. Zkuste fytofarmaka

DigiZone.cz: Sony KD-55XD8005 s Android 6.0

Sony KD-55XD8005 s Android 6.0

Měšec.cz: Kdy vám stát dá na stěhování 50 000 Kč?

Kdy vám stát dá na stěhování 50 000 Kč?

Měšec.cz: U levneELEKTRO.cz už reklamaci nevyřídíte

U levneELEKTRO.cz už reklamaci nevyřídíte

Lupa.cz: Teletext je „internetem hipsterů“

Teletext je „internetem hipsterů“

Vitalia.cz: Paštiky plné masa ho zatím neuživí

Paštiky plné masa ho zatím neuživí

Vitalia.cz: Jsou čajové sáčky toxické?

Jsou čajové sáčky toxické?

Měšec.cz: Jak levně odeslat balík přímo z domu?

Jak levně odeslat balík přímo z domu?

Podnikatel.cz: Babiše přesvědčila 89letá podnikatelka?!

Babiše přesvědčila 89letá podnikatelka?!

Podnikatel.cz: Udávání kvůli EET začalo

Udávání kvůli EET začalo

DigiZone.cz: Recenze Westworld: zavraždit a...

Recenze Westworld: zavraždit a...

Vitalia.cz: Říká amoleta - a myslí palačinka

Říká amoleta - a myslí palačinka

Vitalia.cz: Spor o mortadelu: podle Lidlu falšovaná nebyla

Spor o mortadelu: podle Lidlu falšovaná nebyla

Měšec.cz: Jak vymáhat výživné zadarmo?

Jak vymáhat výživné zadarmo?

Vitalia.cz: Znáte „černý detox“? Ani to nezkoušejte

Znáte „černý detox“? Ani to nezkoušejte