Lekce 39 - Úvod do fyzikálních simulací
V gravitačním poli se pokusíme rozpohybovat hmotný bod s konstantní rychlostí, hmotný bod připojený k pružině a hmotný bod, na který působí gravitační síla - vše podle fyzikálních zákonů. Kód je založen na nejnovějším NeHeGL kódu.
Pokud zvládáte fyziku a chcete používat kód pro fyzikální simulaci, tak Vám tento tutoriál může pomoci. Abyste ale mohli něco vytěžit, měli byste vědět něco o počítání s vektory v trojrozměrném prostoru a fyzikálních veličinách, jako je síla nebo rychlost. Tutoriál obsahuje popis velmi jednoduchého fyzikálního simulátoru.
Třída Vector3D
Návrh fyzikálního simulačního enginu není vždy jednoduchý. Ale je zde jednoduchá posloupnost závislostí - aplikace potřebuje simulační část a ta potřebuje matematické knihovny. Tady tuto závislost uplatníme. Naším cílem je získat zásobník na simulaci pohybu objektů v prostoru. Simulační část bude obsahovat třídy Mass a Simulation. Třída Simulation bude naším zásobníkem. Pokud vytvoříme třídu Simulation budeme schopni vyvíjet aplikace, které ji využívají. Ale předtím potřebujeme matematickou knihovnu. Knihovna obsahuje pouze jednu třídu Vector3D, která pro nás bude představovat body, vektory, pozice, rychlost a sílu ve 3D prostoru.
Vector3D tedy bude jediným členem naší matematické knihovny. Obsahuje souřadnice x, y, z v přesnosti float a zavádí operátory pro počítání s vektory ve 3D. Abychom byli konkrétní, přetížíme operátory sčítání, odčítání, násobení a dělení. Protože se tento tutoriál zaměřuje na fyziku a ne matematiku, nebudu podrobně vysvětlovat Vector3D. Podíváte-li se na jeho zdrojový kód, myslím si, že nebudete mít problémy porozumět.
Síla a pohyb
Abychom mohli implementovat fyzikální simulaci, měli bychom vědět, jak bude vypadat náš objekt. Bude mít polohu a rychlost. Pokud je umístěn na Zemi, Měsíci, Marsu nebo na jakémkoliv místě, kde je gravitace musí mít také hmotnost, která se liší podle velikosti působící gravitační síly. Vezměme si třeba knihu. Na Zemi váží 1 kg, ale na Měsíci pouze 0,17 kg, protože Měsíc na ni působí menší gravitační silou. My budeme uvažovat hmotnost na Zemi.
Poté, když jsme pochopili, co pro nás znamená hmotnost, měli bychom se přesunout k síle a pohybu. Objekt s nenulovou rychlostí se pohybuje ve směru rychlosti. Proto je jeden z důvodů změny polohy v prostoru rychlost. Ač se to nezdá, je další působící veličinou čas. Posunutí předmětu tedy závisí na tom, jak rychle se pohybuje, a na tom kolik času uplynulo od počátku pohybu. Pokud vám vztah mezi polohou, rychlostí a časem není jasný, tak asi nemá cenu pokračovat. Doporučuji si vzít učebnici fyziky a najít si kapitolu zabývající se Newtonovy zákony.
Rychlost objektu se mění, pokud na objekt působí nějaká síla. Její vektor je kombinací směru (počáteční a koncový bod) a velikosti. Velikost působení je přímo úměrná působící síle a nepřímo úměrná hmotnosti objektu. Změna rychlosti za jednotku času se nazývá zrychlení. Čím větší síla působí na objekt, tím více zrychluje. Čím má, ale větší hmotnost, tím je menší zrychlení.
zrychlení = síla / hmotnost
Odsud jednoduše vyjádříme sílu:
síla = hmotnost * zrychlení
Při přípravě prostředí simulace si musíte dávat pozor na to, jaké podmínky v tomto prostředí panují. Prostředí v tomto tutoriálu bude prázdný prostor čekající na zaplnění objekty, které vytvoříme. Nejdříve se rozhodneme, jaké jednotky použijeme pro hmotnost, čas a délku. Rozhodl jsem se použít kilogram pro hmotnost, sekundu pro čas a metr pro délku. Takže jednotky rychlosti budou m/s a jednotky zrychlení budou m/s^2 (metr za sekundu na druhou).
Abychom toto všechno využili v praxi, musíme napsat třídu, která bude reprezentovat objekt a bude obsahovat jeho hmotnost, polohu, rychlost a sílu, která na něho působí.
class Mass
{
public:
float m;// Hmotnost
Vector3D pos;// Pozice v prostoru
Vector3D vel;// Rychlosti a směr pohybu
Vector3D force;// Síla působící na objekt
V konstruktoru inicializujeme pouze hmotnost, která se jako jediná nebude měnit. Pozice, rychlost i působící síly se určitě měnit budou.
Mass(float m)// Konstruktor
{
this->m = m;
}
Aplikujeme silové působení. Objekt může současně ovlivňovat několik zdrojů. Vektor v parametru je součet všech sil působících na objekt. Před jeho aplikací bychom měli stávající sílu vynulovat. K tomu slouží druhá funkce.
void applyForce(Vector3D force)
{
this->force += force;// Vnější síla je přičtena
}
void init()
{
force.x = 0;
force.y = 0;
force.z = 0;
}
Zde je stručný seznam toho, co při simulaci musíme provést:
- Vynulovat sílu - metoda init()
- Vypočítat znovu působící sílu
- Přizpůsobit pohyb posunu v čase
Pro práci s časem použijeme Eulerovu metodu, kterou využívá většina her. Existují mnohem sofistikovanější metody, ale tahle postačí. Velmi jednoduše se vypočítá rychlost a poloha pro další časový úsek s ohledem na působící sílu a uplynulý čas. Ke stávající rychlosti přičteme její změnu, která je závislá na zrychlení (síla/m) a uplynulém čase (dt). V dalším kroku přizpůsobíme polohu - opět v závislosti na čase.
void simulate(float dt)
{
vel += (force / m) * dt;// Změna rychlosti je přičtena k aktuální rychlosti
pos += vel * dt;// Změna polohy je přičtena k aktuální poloze
}
};
Jak by měla simulace pracovat
Při fyzikální simulaci se během každého posunu opakuje totéž. Síly jsou vynulovány, potom znovu spočítány. V závislosti na nich se určují rychlosti a polohy předmětů. Tento postup se opakuje tolikrát, kolikrát chceme. Je zajišťován třídou Simulation. Jejím úkolem je vytvářet, ukládat a mazat objekty a starat se o běh simulace.
class Simulation
{
public:
int numOfMasses;// Počet objektů v zásobníku
Mass** masses;// Objekty jsou uchovávány v jednorozměrném poli ukazatelů na objekty
Simulation(int numOfMasses, float m)// Konstruktor vytvoří objekty s danou hmotností
{
this->numOfMasses = numOfMasses;// Inicializace počtu
masses = new Mass*[numOfMasses];// Alokace dynamické paměti pro pole ukazatelů
for (int a = 0; a < numOfMasses; ++a)// Projdeme všechny ukazatele na objekty
masses[a] = new Mass(m);// Vytvoříme objekt a umístíme ho na místo v poli
}
~Simulation()// Smaže vytvořené objekty
{
release();
}
virtual void release()// Uvolní dynamickou paměť
{
for (int a = 0; a < numOfMasses; ++a)// Smaže všechny vytvořené objekty
{
delete(masses[a]);// Uvolní dynamickou paměť objektů
masses[a] = NULL;// Nastaví ukazatele na NULL
}
delete(masses);// Uvolní dynamickou paměť ukazatelů na objekty
masses = NULL;// Nastaví ukazatel na NULL
}
Mass* getMass(int index)// Získání objektu s určitým indexem
{
if (index < 0 || index >= numOfMasses)// Pokud index není v rozsahu pole
return NULL;// Vrátí NULL
return masses[index];// Vrátí objekt s daným indexem
}
Proces simulace se skládá ze tří kroků:
- Init() nastaví síly na nulu
- Solve() znovu aplikuje síly
- Simulate(float dt) posune objekty v závislosti na čase
virtual void init()// Tato metoda zavolá init() metodu každého objektu
{
for (int a = 0; a < numOfMasses; ++a)// Prochází objekty
masses[a]->init();// Zavolání init() daného objektu
}
virtual void solve()
{
// Bez implementace, protože nechceme v základním zásobníku žádné síly
// Ve vylepšených zásobnících, bude tato metoda nahrazena, aby na objekty působila nějaká síla
}
virtual void simulate(float dt)// Výpočet v závislosti na čase
{
for (int a = 0; a < numOfMasses; ++a)// Projdeme všechny objekty
masses[a]->simulate(dt);// Výpočet nové polohy a rychlosti objektu
}
Všechny tyto metody jsou volány v následující funkci.
virtual void operate(float dt)// Kompletní simulační metoda
{
init();// Krok 1: vynulování sil
solve();// Krok 2: aplikace sil
simulate(dt);// Krok 3: vypočítání polohy a rychlosti objektů v závislosti na čase
}
};
Nyní máme jednoduchý simulační engine. Je založený na matematické knihovně. Obsahuje třídy Mass a Simulation. Používá běžnou Eulerovu metodu na výpočet simulace. Teď jsme připraveni na vývoj aplikací. Aplikace, kterou budeme vyvíjet využívá:
- Objekty s konstantní hmotností
- Objekty v gravitačním poli
- Objekty spojené pružinou s nějakým bodem
Ovládání simulace aplikací
Předtím než napíšeme nějakou simulaci, měli bychom vědět, jak se třídami zacházet. V tomto tutoriálu jsou simulační a aplikační části odděleny do dvou samostatných souborů. V souboru s aplikační částí je funkce Update(), která se volá opakovaně při každém novém framu.
void Update (DWORD milliseconds)// Aktualizace pohybu
{
// Ošetření vstupu z klávesnice
if (g_keys->keyDown [VK_ESCAPE] == TRUE)
TerminateApplication (g_window);
if (g_keys->keyDown [VK_F1] == TRUE)
ToggleFullscreen (g_window);
if (g_keys->keyDown [VK_F2] == TRUE)
slowMotionRatio = 1.0f;
if (g_keys->keyDown [VK_F3] == TRUE)
slowMotionRatio = 10.0f;
DWORD milliseconds je čas, který uplynul od předchozího volání funkce. Budeme počítat čas při simulacích na milisekundy. Pokud bude simulace sledovat tento čas, půjde stejně rychle jako v reálném čase. K provedení simulace jednoduše zavoláme funkci operate(float dt). Předtím než ji zavoláme musíme znát hodnotu dt. Protože ve třídě Simulation nepoužíváme milisekundy, ale sekundy, převedeme proměnnou milliseconds na sekundy. Potom použijeme proměnnou slowMotionRatio, která udává, jak má být simulace zpomalená vzhledem k reálnému času. Touto proměnnou dělíme dt a dostaneme nové dt. Přidáme dt k proměnné timeElapsed, která udává kolik času simulace už uběhlo (neudává tedy reálný čas).
float dt = milliseconds / 1000.0f;// Přepočítá milisekundy na sekundy
dt /= slowMotionRatio;// Dělení dt zpomalovací proměnnou
timeElapsed += dt;// Zvětšení uplynulého času
Teď už je dt skoro připraveno na použití v simulaci. Ale! je tu jedna důležitá věc, kterou bychom měli vědět: čím menší je dt, tím reálnější je simulace. Pokud nebude dt dostatečně malé, naše simulace se nebude chovat realisticky, protože pohyb nebude spočítán dostatečně precizně. Analýza stability se užívá při fyzikálních simulacích, aby zajistila maximální přijatelnou hodnotu dt. V tomto tutoriálu se nebudeme pouštět do detailů. Pokud vyvíjíte hru a ne specializovanou aplikaci, tato metoda bohatě stačí na to, abyste se vyhnuli chybám.
Například v automobilovém simulátoru je vhodné, aby se dt pohybovalo mezi 2 až 5 milisekundami pro běžné auto a mezi 1 a 3 milisekundami pro formuli. Při arkádovém simulátoru je možné použít dt v rozsahu od 10 do 200 milisekund. Čím nižší je dt, tím silnější procesor potřebujeme, abychom stíhali simulovat v reálném čase. To je důvod proč se u starších her nepoužívají fyzikální simulace.
V následujícím kódu nastavíme maximální hodnotu dt na 0.1 sekundy (100 milisekund). S touto hodnotou spočítáme kolikrát cyklus simulace při každém projití funkce zopakujeme. To řeší následující vzorec:
int numOfIterations = (int)(dt / maxPossible_dt) + 1;
NumOfIterations je počet cyklů, které při simulaci provedeme. Dejme tomu, že aplikace běží 20 framů za sekundu. Z toho plyne, že dt=0.05. numOfIterations tedy bude 1. Simulace se provede jednou po 0.05 sekundách. Pokud by dt bylo 0.12 sekund, pak numOfIterations bude 2. Pod v kódu uvedeným vzorcem můžete vidět, že dt počítáme ještě jednou. Podělíme ho počtem cyklů a bude dt = 0.12 / 2 = 0.06. dt bylo původně vyšší než maximální možná hodnota 0.1. Teď se tedy rovná 0.06. My ale provedeme dva cykly simulace, takže v simulaci uběhne čas 0.12 sekund. Prozkoumejte následující kód a ujistěte se, že všemu rozumíte.
// Abychom nepřekročili hranici kdy už se simulace nechová reálně
float maxPossible_dt = 0.1f;// Nastavení maximální hodnoty dt na 0.1 sekund
int numOfIterations = (int)(dt / maxPossible_dt) + 1;// Výpočet počtu opakování simulace v závislosti na dt a maximální možné hodnotě dt
if (numOfIterations != 0)// Vyhneme se dělení nulou
dt = dt / numOfIterations;// dt by se měla aktualizovat pomocí numOfIterations
for (int a = 0; a < numOfIterations; ++a)// Simulaci potřebujeme opakovat numOfIterations-krát
{
constantVelocity.operate(dt);// Provedení simulace konstantní rychlosti za dt sekund
motionUnderGravitation.operate(dt);// Provedení simulace pohybu v gravitaci za dt sekund
massConnectedWithSpring.operate(dt);// Provedení simulace pružiny za dt sekund
}
}
1. Objekt s konstantní rychlostí
Objekt s konstantní rychlostí nepotřebuje působení externí síly. Pouze vytvoříme objekt a nastavíme jeho rychlost na (1.0f, 0.0f, 0.0f), takže se bude pohybovat po ose x rychlostí 1 m/s. Třídu ConstantVelocity odvodíme od třídy Simulation.
class ConstantVelocity : public Simulation
{
public:
// Konstruktor nejdříve použije konstruktor nadřazené třídy, aby vytvořil objekt o hmotnosti 1 kg
ConstantVelocity() : Simulation(1, 1.0f)
{
masses[0]->pos = Vector3D(0.0f, 0.0f, 0.0f);// Nastavíme polohu objektu na počátek
masses[0]->vel = Vector3D(1.0f, 0.0f, 0.0f);// Nastavíme rychlost objektu na (1.0f, 0.0f, 0.0f) m/s
}
};
Když je volána metoda operate(float dt) třídy ConstantVelocity, vypočítá se nová polohu objektu. Tato metoda je volána hlavní aplikací před každým překreslením okna. Dejme tomu, že aplikace běží 10 framů za sekundu. To znamená, že při každém novém výpočtu bude dt 0.1 sekundy. Když se potom zavolá funkce simulate(float dt) daného objektu, k jeho pozici se přičte rychlost*dt, které se rovná:
Vector3D(1.0f, 0.0f, 0.0f) * 0.1 = Vector3D(0.1f, 0.0f, 0.0f)
Při každé frame se objekt pohne o 0.1 metru doprava. Po 10 framech to bude právě 1 metr. Rychlost byla 1.0 m/s. Takže to bude fungovat celkem slušně.
Když spustíte aplikaci, uvidíte objekt pohybující se konstantní rychlostí po ose x. Aplikace nabízí dvě rychlosti plynutí času. Stisknutím F2 poběží stejně rychle jako reálný čas. Stisknutím F3 poběží 10krát pomaleji. Na obrazovce uvidíte přímky znázorňující souřadnicovou plochu. Mezery mezi přimkami jsou 1 metr. Díky těmto přímkám uvidíte, že se objekt pohybuje 1 metr za sekundu v reálném čase a 1 metr za 10 sekund ve zpomaleném čase. Výše popsaná technika je způsob, jak udělat simulaci tak, aby běžela v reálném čase. Abyste ji mohli použít musíte se pevně rozhodnout, v jakých jednotkách simulace poběží.
Aplikace síly
Při simulacích s konstantní rychlostí jsme nepoužili sílu působící na objekt, protože víme, že pokud síla působí na objekt, tak mění jeho rychlost. Pokud chceme pohyb s proměnlivou rychlostí použijeme vnější sílu. Nejdříve musíme všechny působící síly sečíst, abychom dostali výslednou sílu, kterou v simulační fázi aplikujeme na objekt.
Dejme tomu, že chcete použít na objekt sílu 1 N ve směru x. Pak do solve() napíšete:
mass->applyForce(Vector3D(1.0f, 0.0f, 0.0f));
Pokud chcete navíc přidat sílu 2 N ve směru y, napíšete:
mass->applyForce(Vector3D(1.0f, 0.0f, 0.0f));
mass->applyForce(Vector3D(0.0f, 2.0f, 0.0f));
Na objekt můžete použít libovolné množství sil, libovolných směrů, abyste ovlivnili pohyb. V následující části použijeme jednoduchou sílu.
2. Pohyb v gravitaci
MotionUnderGravitation vytvoří objekt a nechá na něj působit sílu. Touto silou bude právě gravitace, která se vypočítá vynásobením hmotnosti objektu a gravitačního zrychlení:
F = m * g
Gravitační zrychlení na Zemi odpovídá 9.81 m/s^2. To znamená, že objekt při volném pádu zrychlí každou sekundu o 9.81 m/s dokud na něho nepůsobí žádná jiná síla než gravitace. Může jí být odpor vzduchu, který působí vždycky, ale to sem nepatří.
class MotionUnderGravitation : public Simulation
{
public:
Vector3D gravitation;// Gravitační zrychlení
Konstruktor přijímá Vector3D, který udává sílu a orientaci gravitace.
// Konstruktor nejdříve použije konstruktor nadřazené třídy, aby vytvořil 1 objekt o hmotnosti 1kg
MotionUnderGravitation(Vector3D gravitation) : Simulation(1, 1.0f)
{
this->gravitation = gravitation;// Nastavení gravitace
masses[0]->pos = Vector3D(-10.0f, 0.0f, 0.0f);// Nastavení polohy objektu
masses[0]->vel = Vector3D(10.0f, 15.0f, 0.0f);// Nastavení rychlosti objektu
}
virtual void solve()// Aplikace gravitace na všechny objekty, na které má působit
{
// Použijeme gravitaci na všechny objekty (zatím máme jenom jeden, ale to se může do budoucna změnit)
for (int a = 0; a < numOfMasses; ++a)
masses[a]->applyForce(gravitation * masses[a]->m);// Síla gravitace se spočítá F = m * g
}
V kódu nahoře si můžete všimnout vzorce F = m * g. Pro reálné působení gravitace byste měli předat konstruktoru Vectror3D(0.0f, -9.81f, 0.0f). -9.81 znamená, že má gravitace působit proti směru y, což způsobuje, že objekt padá směrem dolů. Můžete zkusit zadat kladné číslo a určitě poznáte rozdíl.
3. Objekt spojený pružinou s bodem
V tomto příkladě chceme spojit objekt se statickým bodem. Pružina by měla objekt přitahovat k bodu upevnění a tak způsobovat oscilaci objektu. V konstruktoru nastavíme bod upevnění a pozici objektu.
class MassConnectedWithSpring : public Simulation
{
public:
float springConstant;// Čím vyšší bude tato konstanta, tím tužší bude pružina
Vector3D connectionPos;// Bod ke kterému bude objekt připojen
// Konstruktor nejdříve použije konstruktor nadřazené třídy, aby vytvořil 1 objekt o hmotnosti 1kg
MassConnectedWithSpring(float springConstant) : Simulation(1, 1.0f)
{
this->springConstant = springConstant;// Nastavení tuhosti pružiny
connectionPos = Vector3D(0.0f, -5.0f, 0.0f);// Nastavení pozice upevňovacího bodu
masses[0]->pos = connectionPos + Vector3D(10.0f, 0.0f, 0.0f);// Nastavení pozice objektu na 10 metrů napravo od bodu, ke kterému je uchycen
masses[0]->vel = Vector3D(0.0f, 0.0f, 0.0f);// Nastavení rychlosti objektu na nulu
}
Rychlost objektu je nula a jeho pozice je 10 metrů napravo od úchytu, takže se bude pohybovat ze začátku směrem doleva. Síla pružiny se dá zapsat jako
F = -k * x
k je tuhost pružiny a x je vzdálenost od úchytu. Záporná hodnota u k značí, že jde o přitažlivou sílu. Kdyby bylo k kladné, tak by pružina objekt odpuzovala, což zcela jistě neodpovídá skutečnému chování.
virtual void solve()// Užití síly pružiny
{
// Použijeme sílu na všechny objekty (zatím máme jenom jeden, ale to se může do budoucna změnit)
for (int a = 0; a < numOfMasses; ++a)
{
Vector3D springVector = masses[a]->pos - connectionPos;// Nalezení vektoru od pozice objektu k úchytu
masses[a]->applyForce(-springVector * springConstant);// Použití síly podle uvedeného vzorce
}
}
};
Výpočet síly v kódu nahoře odpovídá vzorci, který jsme si uvedli (F = -k * x). Jenom je zde místo x trojrozměrný vektor a místo k je zde springConstant. Čím vyšší je springConstant, tím rychleji objekt osciluje.
V tomto tutoriálu jsem se snažil předvést základní prvky pro tvorbu fyzikálních simulací. Pokud vás zajímá fyzika, nebude pro vás těžké vytvořit vlastní simulace. Můžete zkoušet složitější interakce a vytvořit tak zajímavá dema a hry. Další v pořadí by měli být simulace pevných objektů, jednoduché mechaniky a pokročilé simulační metody.
napsal: Erkin Tunca <erkintunca (zavináč) icqmail.com>
přeložil: Václav Slováček - Wessan <horizont (zavináč) host.sk>
Zdrojové kódy
- Visual C++ kód této lekce.
- Borland C++ Builder 6 kód této lekce. ( Le Thanh Cong )
- Code Warrior 5.3 kód této lekce. ( Scott Lupton )
- Dev C++ kód této lekce. ( Warren Moore )
- Linux/GLut kód této lekce. ( Laks Raghupathi )
- Visual Studio .NET kód této lekce. ( Grant James )
