Ahoj lidi!
Dnes jsem se rozhodl pro vás napsat návod na jeden zajímavý efekt,
který můžete s malou námahou a velkým efektem použít například na
pozadí titulků svého programu, jako to udělala např. Iguana v demu
Heartquake. Jedná se o kreslení pohybujících se vlnek na vodní
hladině. Začnu trochou matematiky, což se, jak brzy zjistíte, u mě
stane tradicí. Já jsem totiž šílenec a mám na to dokonce průkaz
(s nápisem MFF UK :-).
VLNOVÁ ROVNICE
Protože najít nějaký jiný výpočet, který by dával podobné výsledky je
velmi obtížné, tady je obecná vlnová rovnice, jak ji odvodili fyzici,
ve své plné kráse ;-) (pro vlnu na povrchu - 2D):
2 + 2 2 +
d p 2 | d p d p | dp
----- = V |----- + -----| - D ----
2 | 2 2 |
dt + dx dy + dt
Nebojte se, až to budete programovat, tak tam nebude ani jeden
součin, JENOM DVĚ ROTACE a budete to počítat v INTEGERECH! Nicméně,
pokud chcete pochopit, proč by to mělo vlastně fungovat, tak tady je
vysvětlení: Nejdříve popis symbolů: Symbol d znamená diferenciál,
p je okamžitá výchylka v určitém bodě a čase, x a y jsou prostorové
souřadnice a t souřadnice časová. Konstanta V je rychlostní konstanta
a D je tlumící konstanta. A nyní trochu polopatické vysvětlení, co
tahle rovnice znamená. Druhá derivace na levé staně představuje
zrychlení a rovnice o něm tvrdí, že je přímo úměrné zakřivení daného
povrchu a v případě tlumení se snižuje při vyšší rychlosti.
Pokud pořád ještě netušíte, co tahle rovnice vyjadřuje, nebo nevíte,
jak něco takového numericky počítat, uvádím podrobné vysvětlení
eulerovy metody řešení této rovnice.
Pro výpočet druhé derivace eulerovou metodou musíme rovnici rozdělit
na soustavu dvou rovnic s první časovou derivací. Asi takto:
+ 2 2 +
dv | d p d p |
---- = sqV |----- + -----| - D v
| 2 2 |
dt + dx dy +
dp
---- = v,
dt
kde v je rychlost v daném bodě a V^2 bylo nahrazeno symbolem sqV,
čímž chci říct, že tu konstantu přece nebudete pokaždé umocňovat.
Tyto dvě rovnice potom převedeme do diferenčního tvaru a obdržíme
následovný tvar:
+ +
v - v | p + p - 2p p + p - 2p |
n o | i+1,j i-1,j i,j i,j+1 i,j-1 i,j |
--------- = sqV |----------------------- + ---------------------- | -
| 2 2 |
⌂t | ⌂x ⌂y |
+ +
- D v
o
p - p
n o
--------- = v
n (<-- tady správně patří o, ale prakticky je to jedno)
⌂t
Pokud zvolíme vhodné jednotky, snadno se zbavíme jmenovatelů (budou
prostě 1), no a pak už stačí deklarovat příslušná pole a přepsat výše
uvedené rovnice do C syntaxe:
v[i][j] = sqV * (p[i+1][j] + p[i-1][j] + p[i][j+1] + p[i][j-1] -
p[i][j]) + jednaminusD * v[i][j];
p[i][j] += v[i];
Doufám, že jste pochopili, že změna konstanty D na jednaminusD je
způsobena tím, že se v rovnici objevuje dvojčlen v - D*v.
Ekvivalentně lze použít operátor += a původní konstantu. Konstanta
sqV musí být z důvodu přesnosti nejvýše 1/2, jinak to bude počítat
pěkný blbosti a ne vlny.
No a hurá programovat! Tak to zase prr! Ještě nám zbývá několik
problémů.
OKRAJOVÉ PODMÍNKY
Protože pole, ve kterém budete vlny počítat, někde končí, je nutné
ošetřit, co se bude dít na jeho okrajích, kde zjevně neexistují
všechna sousední políčka. Máte dvě základní možnosti:
Pevný okraj
Pevné okraje se vyskytují například u blány na povrchu bubnu a pro
simulaci jsou velmi jednoduché. Pevný okraj má tu vlastnost, že
výchylka p je zde vždy 0. Takže prostě na okraje pole uložíte nuly a
nebudete pro tyto uzly provádět ŽÁDNÝ výpočet.
Volný okraj
Na hladině vody se však vyskytují okraje volné. Jejich simulace
vypadá efektněji, ale staví nás před větší problémy. Volný konec je
obecně charakterizován tím, že zde je povrch v rovině. Tedy bod na
okraji má stejnou výchylku, jako bod těsně vedle něj (směrem
dovnitř). Takže pro okrajové uzly musíme provádět speciální operaci,
která jim přiřadí hodnoty jejich sousedů. Pochopitelně můžete mít
některé okraje volné a jiné pevné, nebo můžete na některý okraj
připojit generátor vln (budete tam postupně ukládat předem danou
posloupnost hodnot výchylky)
Průměrná výška a rychlost
Použití pouze volných okrajů nás staví ještě před jeden problém. Díky
nepřesnosti, které se přibližný numerický výpočet nevyhne, se velmi
rychle stane, že celý vlnící se povrch odcestuje nahoru, nebo dolu.
Proto musíte udržovat celý povrch v nulové průměrné výšce a
rychlosti, kterážto podmínka je odvozena z toho, že pod vodní
hladinou bývá konstantní objem nádoby a tento je nutno zachovat.
Proto musíme do výrazů pro výpočet v[i][j] a p[i][j] přidat člen vavg
a pavg respektive. Takže:
v[i][j] = ... - vavg;
p[i][j] += ... - pavg;
Hodnoty vavg a pavg potom vypočítáme tak, že sečteme hodnoty v,
respektive p ve všech bodech sítě a vydělíme je počtem uzlů
(idim*jdim).
POČÁTEČNÍ PODMÍNKY
Než můžete začít počítat, musíte pole naplnit nějakými hodnotami,
kterým se říká počáteční podmínky. Pokud tam dáte samé nuly, tak se
vám nic nerozvlní. Takže musíte buď na počátku naplnit pole
nestejnými hodnotami, nebo později vyvolat rozruch (dodatečně
pozměníte některé hodnoty)
No a teď se konečně můžete pustit do toho programování!
-Bulb-
PS: Na závěr jeden PŘÍKLAD. Přeji příjemnou zábavu!
výheň