Eksponencijalna funkcija matrice i linearni dinamički sustavi

Darko Žubrinić
Zavod za primijenjenu matematiku, FER, Zagreb, 2020.

Nastavno pomagalo uz kolegije Linearna algebra i Diferencijalne jednadžbe za studente 1. godine diplomskog studija FER-a.
Svrha: vizualizacija faznih portreta linearnih dinamičkih sustava
Ovo je nastavak prikaza Matrične transformacije ravnine vizualizirane

U Vašem pregledniku treba biti omogućen JavaScript (preporučam da rabite Chrome). Ako nije omogućen, potražite upute putem Googlea, s upitom "enable javascript". U pripremi članka korišten je MathJax (mrežna inačica LaTeX-a) i Processing. Prikaz je u još razvoju.

Ovo je nastavak prikaza

Matrične transformacije ravnine vizualizirane

Linearni dinamički sustavi

Matrice i matrične funkcije su važne u teoriji dinamičkih sustava u ravnini. Najjednostavniji primjer dinamičkog sustava je linearni homogoeni sustav u ravnini, zadan s pomoću ove diferencijalne jednadžbe i početnog uvjeta (tzv. Cauchyjev problem):
$$ \mathbf x'(t) = \mathbf A\mathbf x(t),\quad \mathbf x(0) = \mathbf x_0, $$ gdje je $\mathbf x_0=(x_1,x_2)^{\top}$ zadani početni vektor u $\mathbb R^2$. Matrica $\mathbf A$ je zadana kvadratna matrica reda 2 s konstantnim koeficijentima. S lijeve strane je derivacija vektorske funkcije $$ \mathbf x(t)=\left[\begin{matrix}x_1(t)\\x_2(t) \end{matrix}\right], $$ koja se dobije tako da deriviramo komponente $x_1$ i $x_2$ vektorske funkcije $\mathbf x(t)$, tj. $$ \mathbf x'(t)=\left[\begin{matrix}x_1'(t)\\x_2'(t) \end{matrix}\right]. $$ Za zadani $\mathbf x_0$ (tj. za zadani početni uvjet $\mathbf x_0$ u $\mathbb R^2$), cilj je izračunati vektorsku funkciju $\mathbf x(t)$.

Diferencijalna jednadžba $\mathbf x'(t) = \mathbf A\mathbf x(t)$ ekvivalentna sustavu dviju diferencijalnih jednačaba s dvije nepoznate funkcije $x_1=x_1(t)$ i $x_2=x_2(t)$: \begin{aligned} x_1'(t)&=ax_1(t)+bx_2(t)\\ x_2'(t)&=cx_1(t)+dx_2(t). \end{aligned} Pritom je, naravno, \begin{equation} {\mathbf A}= \left[ \begin{matrix} a&b\\ c&d \end{matrix} \right] \end{equation} matrica tog sustava.

Eksponencijalna funkcija matrice

Može se pokazati da je rješenje Cauchyjeva problema moguće zapisati u obliku: $$ \mathbf x(t) =e^{\mathbf At}\, \mathbf x_0,\tag{1} $$ gdje je matrica $e^{\mathbf At}$ t. zv. eksponencijalna funkcija matrice $\mathbf At$: $$ e^{{\mathbf A}t}=\sum_{k=1}^{\infty}\frac{({\mathbf A}t)^k}{k!}=1 + {\mathbf A}t + \frac{{\mathbf A}^2t^2}{2!}+\dots+\frac{{\mathbf A}^kt^k}{k!}+\dots \tag{2} $$ Pokazuje se da ovaj red matrica konvergira za svaku zadanu kvadratnu matricu $\mathbf A$ i bilo koji realan broj $t$. (Pogledajte [Aglić Aljinović, Elezović, Žubrinić].) Matrična eksponencijalna funkcija $e^{\mathbf At}$ se često označuje i sa $\exp(\mathbf At)$.

Matrična funkcija $e^{\mathbf At}$ u (2), kao funkcija realnog parametra (vremena) $t$, zove se evolucijska matrica početnog Cauchyjeva problema (ponekad i fundamentalna matrica). Razlog za taj naziv je što je za svaki početni položaj $\mathbf x_0$ jednoznačno određena evolucija $\mathbf x(t)$ za sve $t\in\mathbb R$, kako je opisano u jednakosti (1). Evolucijska matrica $e^{\mathbf At}$ kao da "razmazuje" početnu točku $\mathbf x_0$.

Čvor (stabilan i nestabilan)

Za matricu oblika $\mathbf A=a\mathbf I$, gdje je $a$ realan broj, dobivamo da je pripadajuća evolucijska matrica $e^{a\mathbf It}=e^{at}\mathbf I$, pa je $e^{a\mathbf It}\mathbf x_0=e^{at}\mathbf x_0$, za bilo koju početnu točku ravnine $\mathbb R^2$ s radius-vektorom $\mathbf x_0$. Ako je $a<0$, trajektorije su zrake iz ishodišta orijentirane prema ishodištu, te samo ishodište (koje je stacionarna točka, ili ravnotežna točka); za $a=0$, tj. za nulmatricu, sve točke ravnine ostaju na miru (stacionarne su); za $a>0$ trajektorije su zrake iz ishodišta orijentirane suprotno od ishodišta, te samo ishodište.

Ovdje je prikazan fazni portret sustava $\mathbf x'(t)= a\mathbf x(t)$, gdje je $\mathbf A=a\mathbf I$, sa $a=-1/2$. Pritisnite strjelicom miša u odabranu točku i držite pritisnuto u toj točki, da biste vidjeli odgovarajuću trajektoriju. Svaka trajektorija je dio zrake, gdje početna točka $\mathbf x_0$ konvergira prema ishodištu kada $t\to+\infty$, ali ju nikada ne dosegne. Brzina konvergencije prema ishodištu je eksponencijalna, jer je $\|\mathbf x(t)\|=e^{-t/2}\|\mathbf x_0\|$. Ishodište je stabilna stacionarna (mirna) točka: koju god početnu točku $\mathbf x_0$ odabrali, nakon nekog vremena ona će biti vrlo blizu ishodišta. Preciznu definiciju stabilne stacionarne točke vidi u [Korkut, Županović].

Privlačno radijalno vektorsko polje $\mathbf f$ koje svakoj točki s radius-vektorom $\mathbf x$ pridružuje vektor brzine $\mathbf f(\mathbf x)=a \mathbf x$ za $a<0$ (na platnu je $a=-1/2$), orijentiran prema ishodištu. Prikazani su vektori polja brzina samo u točkama koncentričnih kružnica oko ishodišta.

Ovdje je prikazan fazni portret sustava $\mathbf x'(t)= a\mathbf x(t)$, gdje je $\mathbf A=a\mathbf I$, sa $a=1/2$. Ishodište je nestabilna stacionarna točka: koliko god odabrali početnu točku $\mathbf x_0$ blizu ishodišta (ali različitu od nje), nakon nekog vremena ona će se udaljiti od ishodišta.

Odbojno radijalno vektorsko polje $\mathbf f$ koje svakoj točki s radius-vektorom $\mathbf x$ pridružuje vektor $\mathbf f(\mathbf x)=a\mathbf x$, gdje je $a>0$ (na platnu je $a=1/2$), orijentiran suprotno od ishodišta. Prikazani su vektori polja samo u točkama koncentričnih kružnica oko ishodišta.

Pogledajmo slučajeve kada je $\mathbf A$ dijagonalna matrica, koja na dijagonali ima brojeve istog predznaka (ne nužno jednake, kao što je bio slučaj gore). Pogledajemo najprije slučaj kada su obje vlastite vrijednosti matrice negativne. Pogledajmo slučajeve kada je $\mathbf A$ dijagonalna matrica, koja na dijagonali ima brojeve istog predznaka (ne nužno jednake, kao što je bio slučaj gore). Pogledajemo najprije slučaj kada su obje vlastite vrijednosti matrice negativne.

U ovom primjeru je matrica $\mathbf A$ dijagonalna, s negativnim brojevima $a$ i $b$ na dijagonali (ovdje $a=-1$ i $b=-2$). Ovdje je $\mathbf x(t)=(e^{at}x_{10},e^{bt}x_{20})^{\top}$, gdje je $\mathbf x(0)=\mathbf x_0=(x_{10},x_{20})^{\top}$ početna točka zadana vrhom strjelice miša (trajektorija se prikazuje sve dok miša držimo pritisnutim u toj točki). Budući da su vlastite vrijednosti matrice negativne, ishodište je stabilna točka. Za konkretne vrijednosti $a=-1$ i $b=-2$, trajektorije su "polovice" parabola oblika $x_2=Cx_1^2$ orijentirane prema ishodištu, gdje je $C$ realna konstanta (za $x_1<0$ ili $x_1>0$), i samo ishodište.

Pogledajemo sada slučaj kada su obje vlastite vrijednosti matrice negativne.

U ovom primjeru je matrica $\mathbf A$ dijagonalna, s pozitivnim brojevima $a$ i $b$ na dijagonali (ovdje $a=1$ i $b=2$). Ovdje je $\mathbf x(t)=(e^{at}x_{10},e^{bt}x_{20})^{\top}$, gdje je $\mathbf x(0)=\mathbf x_0=(x_{10},x_{20})^{\top}$ početna točka zadana vrhom strjelice miša. Budući da su vlastite vrijednosti matrice pozitivne, ishodište je nestabilna točka (jer ako je poečetna točka blizu ihodišta, onda kroz neko vrijeme točka $\mathbf x(t)$ biti daleko od njega). Za konkretne vrijednosti $a=1$ i $b=2$, trajektorije su "polovice" parabola oblika $x_2=Cx_1^2$ orijentirane suprotno od ishodišta, gdje je $C$ realna konstanta (za $x_1<0$ ili $x_1>0$), i ishodište.

U svim prethodnim slučajevima (kada matrica $\mathbf A$ ima vlastite vrijednosti istoga predznaka) kažemo da je ishodište čvor (engl. "node"). Čvor je stabilan ako su vlastite vrijednosti matrice $\mathbf A$ negativne, a nestabilan ako su pozitivne.

Hiperboličko vektorsko polje

Ovdje promatramo sustav $\mathbf x'(t)=\mathbf A\mathbf x(t)$ u kojem je matrica $\mathbf A$ dijagonalna, s vrijednostima suprotnih predznaka na dijagonali: $$ \mathbf A= \left[ \begin{matrix} a&0\\ 0&b \end{matrix} \right],\quad\mbox{gdje je $a>0$ i $b<0$ (ili obratno).} $$ Pogledajmo vektorsko polje koje točki s radius-vektorom $\mathbf x=(x_1,x_2)^{\top}$ pridružuje vektor $\mathbf A\mathbf x=(ax_1,bx_2)^{\top}$.

U ovom primjeru je $a>0$ i $b<0$, čime dobivamo hiperboličko vektorsko polje brzina oko ishodišta. Duž vodoravne $x_1$-osi ishodište je nestabilna (odbojna) točka. Duž vertikalne $x_2$-osi ishodište je stabilna (privlačna) točka. U ovom primjeru je $a=2$ i $b=-1/2$. Za $x_1$-os kažemo da je stabilna mnogostrukost, dok je $x_2$-os nestabilna mnogostrukost (primijetite da su to upravo vlastiti podprostori matrice $\mathbf A$). Trajektorije će općenito biti hiperbole.
Indeks vektorskog polja u odnosu na ishodište jednak je $-1$: ako pritisnutim mišem napravimo jedan okret oko ishodišta u pozitivnom smjeru (tj. u smjeru suprotnom od smjera kazaljke na satu), strjelica napravi jedan okret u negativnom smjeru (tj. u smjeru kazaljke na satu). Općenitu definiciju vidi u [Korkut, Županović].

Pogledajmo linearni dinamički sustav $\mathbf x'(t)= \mathbf A\mathbf x(t)$, gdje je matrica $\mathbf A=\mbox{diag}\,(a,b)$ dijagonalna, s brojevima $a>0$ i $b<0$ na dijagonali (primijetite da su ti brojevi ujedno i vlastite vrijednosti matrice $\mathbf A$). Evolucijska matrica $e^{\mathbf At}$ je onda također dijagonalna, $e^{\mathbf At}=\mbox{diag}\,(e^{at}, e^{bt})$, pa je rješenje dinamičkog sustava uz zadani početni uvjet $\mathbf x_0=(x_{01},x_{02})^{\top}$ jednako $\mathbf x(t)=e^{\mathbf At}\mathbf x_0$, tj. $$ x_1(t)=e^{at}x_{01},\quad x_2(t)=e^{bt}x_{02}. $$ Iz $a>0$ i $b<0$ slijedi da za bilo koju odabranu početnu točku $\mathbf x_0\ne0$ vrijedi da, ako $t\to+\infty$, onda

$x_1(t)\to+\infty$ ako je $x_{10}>0$, a $x_1(t)\to-\infty$ ako je $x_{10}<0$.
$x_2(t)\to0$.

Ovdje je prikazan fazni portret sustava $\mathbf x'(t)= \mathbf A\mathbf x(t)$, gdje je matrica $\mathbf A$ dijagonalna, s brojevima $a=1$ i $b=-1/2$ na dijagonali. Sve trajektorije (putanje) su hiperbole, osim onih s početnom točkom na koordinatnim osima. Ishodište je hiperbolička točka: u odnosu na dinamiku na vodoravnoj osi ona je nestabilna (odbojna) točka, a u odnosu na dinamiku na vertikalnoj osi, ona je stabilna (privlačna) točka.

Degeneriran slučaj (kada matrica sustava nije regularna)

Pogledajmo slučaj kada je matrica $\mathbf A$ neinvertibilna, tj. barem jedna vlastita vrijednost joj je jednaka nula. Na primjer, $\mathbf A= \left[ \begin{matrix} a&0\\ 0&0 \end{matrix} \right]$.

U ovom primjeru je matrica $\mathbf A$ dijagonalna, s negativnim brojem $a$ i s $0$ na dijagonali (uzimamo $a=-1$ i $b=0$). Ovdje je $\mathbf x(t)=(e^{at}x_{10},x_{20})^{\top}$, gdje je $\mathbf x(0)=\mathbf x_0=(x_{10},x_{20})^{\top}$ početna točka zadana vrhom strjelice miša. Sve točke na $x_2$-osi su stacionarne točke. U slučaju kada je $a<0$ (kao na gornjem platnu), stacionarne točke su stabilne, a ako je $a>0$ (kao na dolnjem platnu) stacionarne točke su nestabilne. Trajektorije su dijelovi vodoravnih zraka, kao i svaka od stacionarnih točaka na vertikalnoj osi. Dinamika je u bitnom jednodimenzionalna, tj. odgovara rangu matrice $\mathbf A$ koji je jednak $1$.

U slučaju kada je $\mathbf A$ nul-matrica, sve točke ravnine su stacionarne. Dinamika je u bitnom nul-dimenzionalna, koja odgovara rangu matrice (rang je jednak $0$).

Matrična imaginarna jedinica $\mathbf J$

U slučaju kada je $\mathbf A$ jednak matričnoj imaginarnoj jedinici $\mathbf J$, tj. $$ \mathbf A=\mathbf J=\left[ \begin{matrix} 0&-1\\ 1&\phantom{-}0 \end{matrix} \right], $$ onda je $e^{{\mathbf A}t}$ matrica rotacije ravnine oko ishodišta za kut od $t$ radijana, tj. $$ e^{{\mathbf J}t}=\left[ \begin{matrix} \cos t&-\sin t\\ \sin t&\phantom{-}\cos t \end{matrix} \right]. $$

Primijetite da gornju jednakost možemo pisati na ovaj način: $$ e^{{\mathbf J}t}=\cos t\,\mathbf I+\sin t\,\mathbf J $$ gdje je $\mathbf I$ jedinična matrica, a $\mathbf J$ matrična imaginarna jedinica (jer je $\mathbf J^2=-\mathbf I$). Ovo podsjeća na Eulerovu formulu za kompleksne brojeve, $e^{it}=\cos t + i\sin t$ (gdje je $i$ imaginarna jedinica, $i^2=-1$), zar ne?

Drugim riječima, za bilo koji vektor $\mathbf x_0$ u ravnini, vektor $e^{{\mathbf J}t}\mathbf x_0$ nastaje iz $\mathbf x_0$ njegovom rotacijom za kut od $t$ radijana oko ishodišta.

Ovo je svojstvo u lijepoj analogiji s kompleksnim brojevima. Za bilo koji kompleksni broj $z_0$ u Gaussovoj ravnini, broj $e^{it}z_0$ nastaje iz $z_0$ rotacijom za kut od $t$ radijana. Doista, ako je $z_0=re^{i\varphi}$ trigonometrijski oblik broja $z_0$, onda je: $$ e^{it}z_0=e^{it}re^{i\varphi}=re^{i(t+\varphi)}. $$ Prema tome, matričnoj imaginarnoj jedinici $\mathbf A=\mathbf J$ odgovara imaginarna jedinica $i$, a matričnoj funkciji $e^{{\mathbf A}t}=\cos t\,\mathbf I+\sin t\,\mathbf J$ odgovara znamenita Eulerova formula: $$ e^{it}=\cos t+ i\sin t. $$

Kao što vidimo, u slučaju matrice $\mathbf A=\mathbf J$, skup svih radius-vektora oblika $e^{{\mathbf A}t}\mathbf x_0$, gdje je $t$ bilo koji realan broj, je kružnica polumjera $\|\mathbf x_0\|$, gdje je $\|\mathbf x_0\|$ Euklidska norma (tj. udaljenost odgovarajuće točke s tim radius-vektorom od ishodišta). Taj skup zovemo putanjom ili trajektorijom dinamičkog sustava $\mathbf x'=\mathbf A\mathbf x$. Skup svih trajektorija (u ovom slučaju svih kružnica) zove se fazni portret dinamičkog sustava.

Cirkularno (kružno) vektorsko polje; centar

U cijeloj ravnini ishodište je jedina točka koja se ne miče (sve ostale točke kruže po svojim pripadajućim kružnicama pozitivnog radiusa). Iz tog razloga točku $0$ koja miruje zovemo stacionarnom točkom (mirnom točkom). Još preciznije, tu stacionarnu točku zovemo centrom, jer su oko nje kružne trajektorije.

Pogledajmo fazni portret harmonijskog oscilatora $x_2=x_2(t)$, čije je gibanje duž vertikalne $x_2$-osi zadano sa $$x_2''(t)=-x_2(t).$$ Označimo li trenutačnu brzinu sa $x_1(t)$, tj. $x_2'(t)=x_1(t)$, onda dobivamo dinamički sustav u koordinatnom sustavu $x_1,x_2$ (faznu ravninu čine brzina i položaj) opisan sa \begin{aligned} x_1'(t)&=-x_2(t)\\ x_2'(t)&=\phantom{-}x_1(t). \end{aligned} Drugim riječima, dobivamo sustav $$ \mathbf x'(t)=\mathbf J\mathbf x(t), $$ gdje je $\mathbf x=(x_1,x_2)^{\top}$, a $\mathbf J$ je matrična imaginarna jedinica: $\mathbf J=\left[ \begin{matrix} 0&-1\\ 1&\phantom{-}0 \end{matrix} \right]$.

Za matričnu imaginarnu jedinicu $\mathbf J$, trajektorije sustava $\mathbf x'(t)=\mathbf J\mathbf x(t)$, tj. $x_1'(t)=-x_2(t),\,\,x_2'(t)=x_1(t)$ su koncentrične kružnice s obzirom na ishodište. Početni uvjet se zadaje tako da miš postavite u bilo koju točku platna (canvas), te pritisnete i držite strjelicu miša stalno u toj početnoj točki.

Radi boljeg doživljaja traga trajektorije (putanje), generirani kružići biraju na slučajan način svaki svoju boju, među njih ukupno $256^3=16\,777\,216$. Naime, svaka od triju osnovnih boja RGB (red, green, blue - crvena, zelena, plava) dolazi u 256 nijansa, a svaka se boja dobiva njihovim miješanjem. Ovo ćete vidjeti uz predpostavku da niste daltonist. Ako ste daltonist, obratite pozornost samo na točku na platnu koja titra.

Kružna trajektorija (ciklus) $\mathbf x(t)=\exp(\mathbf Jt)\,\mathbf x_0$ će započeti u toj početnoj točki $\mathbf x_0$. Vektor $\mathbf x(t)=\exp(\mathbf Jt)\,\mathbf x_0$ s početkom u ishodištu dobiva se rotacijom radius-vektora $\mathbf x_0$ za kut od $t$ radijana. Po jedan kružni put se odvija kroz $2\pi$ (tj. oko 6.28) sekunda, tj. period ciklusa iznosi $T=2\pi$. Ciklus se prekida odpuštanjem miša.

Primijetite da zadavanjem početnog položaja miša zadajemo zapravo početnu brzinu (na vodoravnoj osi) i početni položaj (na vertikalnoj osi).

Da su trajektorije doista kružne i pozitivno orijentirane, može se vidjeti i prelaskom s Kartezijevih koordinata $(x_1,x_2)$ na polarne koordinate $r,\varphi$. U polarnim koordinatama, jednadžba promatranog sustava postaje iznenađujuće jednostavna (dovoljno je koristiti da je $r=\sqrt{x_1^2+x_2^2}$ i $\varphi=\arctan(x_2/x_1)$): \begin{aligned} r'(t)&=0\\ \varphi'(t)&=1. \end{aligned} tj. radijalni dio $r(t)$ je konstantan, a kutna brzina $\varphi'(t)$ jednaka je $1$. Drugim riječima, $$ r(t)=C_1,\quad\varphi(t)=t+C_2, $$ gdje je $C_1$ pozitivna konstanta ili nula (jer broj $r$, kao udaljenost točke od ishodišta, ne može biti negativan), a $C_2$ je bilo koja realna konstanta. Pripadajuća trajektorija je kružnica polumjera $C_1$, pozitivno orijentirana (tj. suprotno od smjera gibanja kazaljke na satu). Jedan krug svaka točka (osim ishodišta) obiđe kroz vrijeme od $2\pi\approx 6.28$ sekunda.

U ravnini $\mathbb R^2$, tj. u ravnini platna, iz svake njene točke s koordinatama $(x_1,x_2)$ polazi vektor $(-x_2,x_2)^{\top}$, određen desnom stranom sustava $\mathbf x'=\mathbf J\mathbf x$, gdje je $\mathbf x=\mathbf x(t)$. Na sljedećem platnu je ucrtano polje vektora, ili vektorsko polje, koje "vuče" svaku početnu točku. To polje vektora interpretiramo kao polje brzina u ravnini.

Cirkularno (kružno) vektorsko polje brzina definirano je u cijeloj ravnini, tj. iz svake točke izlazi po jedan vektor brzine, definiran desnom stranom sustava $\mathbf x'=\mathbf J\mathbf x$. Ovdje su radi jednostavnosti ucrtani samo neki od vektora s početcima na koordinatnim osima. Oko ishodišta se jasno nazire da je vektorsko polje "vrtložno". Vektor $\mathbf J\mathbf x$ je u ovom slučaju uvijek okomit radius-vektor $\mathbf x$.

Primijetite da su vektori $\mathbf J\mathbf x=(-x_2,x_1)^{\top}$ i $\mathbf x=(x_1,x_2)^{\top}$ istih duljina i međusobno okomiti, jer je $(\mathbf J\mathbf x\,|\,\mathbf x)=-x_2x_1+x_1x_2=0$. Vektor brzine "vuče" točku duž kružnice.

Pritisnite mišem odabranu početnu točku s radius-vektorom $\mathbf x_0$ i držite ju pritisnutom u istoj točki, da biste vidjeli gibanje točke, kao i pripadajuću (plavi) vektor brzine $\require{color}{\color{blue}\mathbf J\mathbf x}={\color{blue}(-x_2,x_1)}$ koji ju "vuče" duž kružne trajektorije (ciklusa). Položaj točke je određen njenim (zelenim) radius-vektorom $\require{color}{\color{green}\mathbf x}={\color{green}(x_1,x_2)^{\top}}$. Vektor $\require{color}{\color{blue}\mathbf J\mathbf x}$ je na njega uvijek okomit.

Indeks kružnog (pozitivno cirkularnog) vektorskog polja $\mathbf x\mapsto\mathbf{Jx}$ u odnosu na ishodište jednak je $1$: jednom punom okretu točke $\mathbf x$ oko ishodišta u pozitivnom smjeru, odgovara jedan puni okret (plavog) vektora $\mathbf J\mathbf x$ oko te točke u pozitivnom smjeru.

Duž $x_2$-osi imamo harmonijsko osciliranje odgovarajuće točke, koja se dobiva ortogonalnim projiciranjem točke koja se jednoliko giblje po kružnici.

Pozitivno orijentirana kružna trajektorija, dobivena iz sustava $\require{color}{\color{green}x_1'(t)=-x_2(t)}$, $\require{color}{\color{green}x_2'(t)=x_1(t)}$.
Na $x_2$-osi imamo gibanje harmonijskog oscilatora s jednadžbom $\require{color}{\color{blue}x_2''(t) = -x_2(t)}$. Plava točka na vertikalnoj osi dobiva kao $x_2$-koordinata zelene točke na kružnici. Njeno gibanje po $x_2$-osi je sinusoidalno: $x_2(t)=A\sin(t+\alpha)$, gdje se konstante $A$ (a to je ujedno i polumjer kružnice) i $\alpha$ lako određuju iz početnog položaja $x_2(0)$ i početne brzine $x_2'(0)$ plave točke.

Ako harmonijski oscilator zadamo duž vodoravne osi, tj. sa $x_1=x_1(t)$, onda je $$x_1''(t)=- x_1(t).$$ Označivši trenutačnu brzinu točke sa $x_2(t)$, tj. $x_1'(t)=x_2(t)$, dobivamo dinamički sustav u koordinatnom sustavu $x_1,x_2$ (fazni prostor čine položaj i brzina točke) opisan sa \begin{aligned} x_1'(t)&=\phantom{-}x_2(t)\\ x_2'(t)&=-x_1(t). \end{aligned} Drugim riječima, dobivamo sustav $$ \mathbf x'(t)=-\mathbf J\mathbf x(t). $$ Trajektorije opisane sa $\mathbf x(t)=\exp(-\mathbf Jt)\,\mathbf x_0$ su opet kružne, ali se točka giblje duž kružnice u negativnom smjeru (tj. u smjeru gibanja kazaljke na satu): kut zakreta kroz vrijeme $t$ iznosi $-t$ radijana oko ishodišta.

Odgovarajući sustav u polarnim koordinatama glasi\begin{aligned} r'(t)&=0\\ \varphi'(t)&=-1. \end{aligned} tj. radijalni dio $r(t)$ je konstantan, a kutna brzina $\varphi'(t)$ jednaka je $-1$. Drugim riječima, $$ r(t)=C_1,\quad\varphi(t)=-t+C_2, $$ pa se radi o kružnoj putanji polumjera $C_1$, orijentiranoj negativno (tj. u smjeru gibanja kazaljke na satu). Jedan krug svaka točka (osim ishodišta) obiđe kroz vrijeme od $2\pi\approx 6.28$ sekunda.

Ovo je vjerojatno najjednostavniji i najizravniji pristup. Ako definiramo kompleksnoznačnu funkciju $z(t)=x_1(t)+i\,x_2(t)$, onda je sustav $x_1'(t)=x_2(t)$, $x_2'(t)=-x_1(t)$ ekvivalentan s običnom diferencijalnom jednadžbom $$ z'(t)=-i\,z(t). $$ To se odmah vidi uspoređivanjem realnih i imaginarnih dijelova u jednakosti $(x_1(t)+i\,x_2(t))'=-i(x_1(t)+i\,x_2(t))$. Rješenje jednadžbe $z'=-iz$ glasi $$ z(t)=e^{-it}z_0, $$ gdje je $z_0=x_{10}+i\,x_{20}$ početna vrijednost, jer je $z'=e^{-it}(-i)z_0=-iz$. Vrijednost $z(t)=e^{-it}z_0$ dobiva se prema Eulerovoj formuli iz $z_0$ rotacijom tog kompleksnog broja za kut od $-t$ radijana. Doista, ako je $z_0=r_0e^{i\varphi_0}$ trigonometrijski oblik kompleksnog broja $z_0$, onda je $z(t)=r_0e^{i(\varphi_0-t)}$ trigonometrijski oblik kompleksnog broja $z(t)$.

Trajektorije sustava $\mathbf x'(t)=-\mathbf J\mathbf x(t)$, tj. $x_1'(t)=x_2(t),\,\,x_2'(t)=-x_1(t)$ su koncentrične kružnice s obzirom na ishodište.
Kružna trajektorija (ciklus) $\mathbf x(t)=\exp(-\mathbf Jt)\mathbf x_0$ će započeti u odabranoj početnoj točki $\mathbf x_0$, a orijentirana je negativno (tj. u smjeru gibanja kazaljke na satu). Ishodište je stacionarna točka, koja je stabilna: ako je početna točka blizu nje, onda cijela trajektorija ostaje blizu. (Preciznu definiciju stabilne stacionarne točke vidi u [Korkut, Županović].)

Indeks kružnog (negativno cirkularnog) vektorskog polja $\mathbf x\mapsto-\mathbf{Jx}$ u odnosu na ishodište jednak je također $1$: jednom punom okretu točke $\mathbf x$ oko ishodišta u pozitivnom smjeru, odgovara jedan puni okret (plavog) vektora $-\mathbf J\mathbf x$ oko te točke u pozitivnom smjeru.

Opći harmonijski oscilator

Općeniti harmonijski oscilator na vodoravnoj $x_1$-osi ima jednadžbu oblika $x_1''(t) = -kx_1(t)$, gdje je $k$ zadani pozitivan broj. U svjetlu trećeg Newtonovog aksioma (zakona akcije i reakcije), sila na jediničnu masu (određena lijevom stranom) jednaka je položaju pomnoženom s $-k$ (predznak minus nam znači da je riječ o privlačnoj sili). Označimo li sa $x_2$ brzinu, onda je gibanje opisano u faznoj ravnini $x_1,x_2$ (položaj i brzina) ovako: \begin{aligned} x_1'(t)&=\phantom{-}x_2(t)\\ x_2'(t)&=-kx_1(t). \end{aligned}

Trajektorije se sastoje od elipsa, za koje je omjer okomite i vodoravne poluosi uvijek jednak $\sqrt k$.

Da bismo pokazali da su trajektorije doista elipse, dovoljno je prvu jednadžbu sustava pomnožiti s $k x_1$, drugi s $x_2$, i zbrojiti. Dobivamo da je $kx_1'x_1+x_2'x_2=0$ za sve $t$, tj. $(kx_1^2+x_2^2)'=0$, tj. $kx_1^2+x_2^2=C$ (jasno je da konstanta $C$ ne može biti negativna). Drugim riječima, za $C>0$ vrijedi $$ \frac{x_1^2}{(C/k)}+\frac{x_2^2}{C}=1, $$ pa se radi o elipsi s vodoravnom poluosi $\sqrt{C/k}$ i vertikalnom $\sqrt{C}$. Znači, omjer vertikalne i vodoravne poluosi je za sve te elipse jednak $\sqrt k$. Nije teško uvjeriti se da je eliptička trajektorija negativno orijentirana. To se može lako vidjeti i crtanje odgovarajućeg vektorskog polja, po kvadrantima. Još preciznije, ako uvedemo eliptičke polarne koordinate $(r,\varphi)$ sa $$ x_1=r\cos\varphi,\quad x_2=\sqrt k\,r\sin\varphi, $$ onda deriviranjem izraza $kx_1^2+x_2^2=kr^2$ po $t$ odmah dobivamo da je $r'=0$, tj. $r=C_1$ (tj. trajektorija u početnom koordinatnom sustavu je elipsa $kx_1^2+x_2^2=kC_1^2$), dok iz $\sqrt k\,\tan\varphi=x_2/x_1$, tj. iz $\varphi=\frac1{\sqrt k}\arctan(x_2/x_1)$, deriviranjem po $t$ dobivamo nakon kratkog računa $$ \varphi'(t)=\frac1{1+\frac{x_2^2}{kx_1^2}}\cdot\frac{x_2'x_1-x_2x_1'}{\sqrt k\,x_1^2}= \frac{kx_1^2}{kx_1^2+x_2^2}\cdot\frac{-kx_1^2-x_2^2}{\sqrt k\,x_1^2} =-\sqrt k. $$ Budući da je generalizirana kutna brzina $\varphi'(t)$ negativna konstanta, eliptička putanja je negativno orijentirana (tj. suprotno od smjera gibanja kazaljke na satu). Jedan puni okret duž elipse ostvaruje se na primjer za $\varphi\in[0,2\pi]$, tj. kroz period od $$ T=\frac{2\pi}{\sqrt k} $$ sekunda. Naravno, što je koficijent $k$ veći, jača je privlačna sila (prema ishodištu), pa je period $T$ manji. 😊

Za one koji znaju rješavati linearne diferencijalne jednadžbe, moguć je i izravniji dokaz. Karakteristična jednadžba pridružena diferencijalnboj jednadžbi $x_1''(t)+kx_1(t)=0$ glasi $\lambda^2+k=0$, tj. $\lambda_{1,2}=\pm\sqrt k\,i$, gdje je $i$ imaginarna jedinica. Radi toga je $x_1(t)=C_1\cos(\sqrt k\,t)+C_2\sin(\sqrt k\,t)$, tj. $$ x_1(t)=A\sin(\sqrt k\,t+\alpha), $$ gdje su $A$ i $\alpha$ slobodne konstante ($A=\sqrt{C_1^2+C_2^2}\,$), pa je riječ o sinusoidalnom gibanju točke na $x_1$-osi, perioda $T=2\pi/\sqrt k$. Kao što vidimo, brzina iznosi $x_2=x_1'=A\sqrt k\,\cos(\sqrt k\,t+\alpha)$, pa je $x_1^2+x_2^2/k=A^2$, tj. dobivamo eliptičku trajektoriju u faznoj ravnini $(x_1,x_2)$ (položaj, brzina). Vektor brzine točke u faznoj ravnini je $$ \mathbf v=(x_1,x_2)'=(x_1',x_2')=(x_2,-kx_1), $$ gdje deriviramo po vremenu $t$. Odatle zbog pozitivnosti konstante $k$ lako vidimo je eliptička trajektorija negativno orijentirana. Na primjer u prvom kvadrantu (gdje su $x_1$ i $x_2$ pozitivni) je prva komponenta vektora brzine pozitivna, pa se točka giblje u desno. Ako smo u četvrtom kvadrantu, onda je prva komponenta vektora brzine negativna, pa se točka u ravnini giblje u lijevo. Itd.

Za $k=1$, tj. kada je u faznoj ravnini $(x_1,x_2)$ riječ o kružnim trajektorijama, jasno je da je brzina točke duž kružnica konstantna (naime, kutna brzina je u tom slučaju, kao što smo već vidjeli, konstnatna, tj. $\varphi'(t)=0$). Pokažimo da to vrijedi i za eliptičke putanje, tj. za bilo koji $k>0)$, tj. da je $v=|\mathbf v|=\sqrt{(x_1')^2+(x_2')^2}$ konstanta. Doista, vrijedi $$(x_1')^2+(x_2')^2=(x_2)^2+(-kx_1)^2=k^2r^2\sin^2\varphi+k^2r^2\cos^2\varphi=k^2r^2=k^2C_1^2,$$ pa je $v=kC_1$.

Na sljedećem platnu možete vidjeti situaciju kada je $k=1/3$.

Negativno orijentirana eliptička trajektorija, dobivena iz sustava $\require{color}{\color{green}x_1'(t)=x_2(t)}$, $\require{color}{\color{green}x_2'(t)=-kx_1(t)}$ za $k=1/3$.
Na $x_1$-osi imamo gibanje harmonijskog oscilatora s jednadžbom $\require{color}{\color{blue}x_1''(t) = -kx_1(t)}$. Plava točka na vodoravnoj osi se dobiva kao $x_1$-koordinata točke na elipsi. Njeno gibanje je sinusoidalno u ovisnosti o vremenu: $x_1(t)=A\sin(\sqrt k\, t+\alpha)$. Konstante $A$ i $\alpha$ lako se određuju iz početnog položaja $x_2(0)$ i početne brzine $x_2'(0)$ plave točke. Elipsa ima vodoravnu i okomitu poluos jednaku $A$ i $A\sqrt k$.

Pitanjem kako u kvalitativnom smislu izgledaju fazni portreti linearnih sustava $\mathbf x'=\mathbf A\mathbf x$ (i općenitijih - nelinearnih sustava) bavi se teorija dinamičkih sustava. Na primjer, ako imamo bilo koju kvadratnu matricu $\mathbf A$ reda dva, onda će skup svih stacionarnih točaka činiti rješenja jednadžbe $\mathbf A\mathbf x$, a to je nul-podprostor (ili jezgra) matrice $\mathbf A$. Upravo te točke miruju, i niti koja druga.

Eksponencijalna spirala

U slučaju kada je matrica $\mathbf A$ nešto općenitijeg oblika $$ \mathbf A=\alpha \mathbf I+ \beta\mathbf J= \left[ \begin{matrix} \alpha&-\beta\\ \beta&\phantom{-}\alpha \end{matrix} \right]. $$ pri čemu su $\alpha$ i $\beta$ zadani realni brojevi različiti od $0$, onda je rješenje Cauchyjeva problema $\mathbf x'=\mathbf A\mathbf x$, s početnim uvjetom $\mathbf x(0)=\mathbf x_0$ oblika \begin{gathered} \mathbf x(t)=e^{{\mathbf A}t}\mathbf x_0=e^{\alpha t} \left[ \begin{matrix} \cos\beta t&-\sin\beta t\\ \sin \beta t&\phantom{-}\cos\beta t \end{matrix} \right]\mathbf x_0, \end{gathered} jer je $$ e^{{\mathbf A}t}=e^{\alpha t\,\mathbf I + \beta t\mathbf J}= e^{\alpha t\,\mathbf I}e^{\beta t\mathbf J}=e^{\alpha t}e^{\beta t\mathbf J}. $$ Pritom smo za drugu jednakost rabili i činjenicu da matrice $\mathbf I$ i $\mathbf J$ međusobno komutiraju: $\mathbf I\mathbf J=\mathbf J\mathbf I$, a u zadnjoj jednakosti smo koristili činjenicu da je $e^{\alpha t\,\mathbf I}=e^{\alpha t}\mathbf I$, koja se odmah dobije razvojem $e^{\alpha t\,\mathbf I}$ u red potencija.

Primijetite da su vlastite vrijednosti matrice $\mathbf A=\alpha \mathbf I+ \beta\mathbf J$ jednake $\lambda_{1,2}=\alpha\pm\beta i$. To vidimo odmah iz Teorema o preslikavanju spektra, jer matrica $\alpha \mathbf I$ ima vlastitu vrijednost $\alpha$, a matrica $\beta\mathbf J$ ima vlastite vrijednosti $\pm\beta i$ (pripadajuća funkcija za primjenu Teorema o preslikavanju spektra je $f(x)=\alpha + \beta x$, jer je onda $f(\mathbf J)=\alpha \mathbf I+ \beta\mathbf J$).

Također primijetite da je rješenje Cauchyjeva problema $\mathbf x'=\mathbf A\mathbf x$, $\mathbf x(0)=\mathbf x_0$, tj. sustava \begin{aligned} x_1'(t)&=\alpha x_1(t)-\beta x_2(t)\\ x_2'(t)&=\beta x_1(t) +\alpha x_2(t). \end{aligned} s početnim uvjetima $x_1(0)=x_{10}$ i $x_2(0)=x_{20}$ (gdje su brojevi $x_{10}$ i $x_{20}$ zadane početne vrijednosti) moguće napisati kao $$ \mathbf x(t)=e^{\alpha t}e^{\beta t\mathbf J}\mathbf x_0 $$ Kao što vidimo, rješenje je dobiveno kao umnožak eksponencijalne funkcije $e^{\alpha t}$ (koja predstavlja udaljenost odgovarajuće točke od ishodišta) i matrice rotacije ravnine za kut $\beta t$. Zanimljivo je primijetiti da je broj $\alpha$ (koji je realan dio vlastite vrijednosti $\lambda_1$) odgovoran za amplitudu (tj. udaljenost točke s radius-vektorom $\mathbf x(t)$ od ishodišta), a imaginarni dio $\beta$ (koji je imaginarni dio vlastite vrijednosti $\lambda_1$) je odgovoran za za kružne oscilacije, s periodom jednog okreta jednakom ${2\pi}/\beta$.

Kako u ovom slučaju izgledaju trajektorije sustava $\mathbf x'=\mathbf A\mathbf x$?

Ako su $\alpha>0$ i $\beta>0$, onda dobivamo spiralnu trajektoriju koja se eksponencijalnom brzinom udaljava od ishodišta kada vrijeme $t$ raste (u svakom trenutku $t$ je udaljenost od ishodišta jednaka točno $e^{\alpha t}$), a se istodobno okreće u pozitivnom smjeru oko ishodišta (tj. suprotno od smjera gibanja kazaljke na satu)), presjecajući na primjer $x$-os svakih ${2\pi}/\beta$ sekundi. To znači da je brzina sve veća i veća. Točnije, vektor brzine je $\mathbf x'=\mathbf A\mathbf x$ (tj. brzina je 'proporcionalna' s položajem $\mathbf x$), pa odatle vidimo da i brzina $\mathbf x'(t)$) raste eksponencijalnom brzinom kada $t$ teži u plus beskonačno, kao i $\mathbf x(t)$.

Ako je $\alpha<0$, onda dobivamo spiralnu trajektoriju koja eksponencijalnom brzinom konvergira prema ishodištu. Ako je pritom $\beta>0$, onda se spirala namata oko ishodišta 'rotirajući' u pozitivnom smjeru, a ako je $\beta<0$ onda se spirala namata oko ishodišta 'rotirajući u negativnom smjeru.

Pritisnite mišem odabranu početnu točku s radius-vektorom $\mathbf x_0$ i držite ju pritisnutom, da biste vidjeli gibanje točke, kao i pripadajuću trajektoriju (putanju). Ovdje je $\alpha =-1$ (pa se točka giblje prema ishodištu) i $\beta = 1$, pa se spirala namata oko ishodišta u pozitivnom smjeru. U polarnom sustavu $(r,\varphi)$, riječ je o eksponencijalnoj spirali zadanoj sa $r = r_0e^{-(\varphi-\varphi_0)}$, gdje su $(r_0,\varphi_0)$ polarne koordinate početne točke (ta se spirala često zove i logaritamskom spiralom, ako $\varphi$ izrazimo s pomoću $r$). Radi jednostavnosti, ovdje uvijek uzimamo istu početnu točku u prvom kvadrantu. Također, program crta samo petnaest prvih zaokreta (kako bi se štedjela memorija Vašeg računala).

Dvije eksponencijalne spirale, koje polaze istodobno iz dviju različitih početnih točaka, koje nisu na istoj spirali. Pripadajuće spiralne trajektorije su disjunktne, tj. nemaju niti jednu zajedničku točku. Ovo svojstvo imaju tzv. autonomni sustavi, koji su po definiciji oblika $\mathbf x'=\mathbf F(\mathbf x)$, tj. neprekinuto vektorsko polje na desnoj strani ne ovisi o vremenu, nego samo o položaju točke (neautonomni sustavi su oblika $\mathbf x'=\mathbf F(t,\mathbf x)$, tj. vektorsko polje na desnoj strani bitno ovisi i o vremenu). Ne zaboravimo da su trajektorije debljine NULA, ali s našim okom krivulje ne crtamo debljine nula (jer ih inače ne bi bilo moguće vidjeti), nego s nekom malom debljinom.

Eksponencijalna spirala može se prikazati u Gaussovoj ravnini (koju možemo tako]er gledati kao platno na gornjoj slici) kao $$ z(t)=e^{\alpha t}e^{\beta ti}z_0=r_0e^{\alpha t}e^{(\beta t+\varphi_0)i}, $$ pri čemu je $i$ imaginarna jedinica. Odavde bez problema odmah očitavamo parametarski zadanu krivulju u polarnom sustavu (s vremenskim parametrom $t$): $$ r=r_0e^{\alpha t}\quad\mbox{i}\quad\varphi=\beta t+\varphi_0. $$ Kao što vidimo, za jedan okret duž spirale treba vrijeme od $T = 2\pi/\beta$ sekunda. Eliminacijom vremena $t$ dobivamo da jednadžba spirale u polarnom sustavu glasi $$ r=r_0e^{\frac{\alpha}{\beta}(\varphi-\varphi_0)}, \quad\mbox{za $\varphi\in\mathbb R.$} $$ Eksponencijalna spirala $r = e^{-\varphi}$ približava se prema nuli kad $\varphi\to+\infty$ (a to će biti kad $t\to+\infty$ ako je $\beta>0$, ili pak kada $t\to-\infty$ ako je $\beta<0$.

Eksponencijalna spirala ima nekoliko iznimno zanimljivih svojstava.

Prvo, ona siječe bilo koju zraku iz ishodišta uvijek pod istim kutem (što se dobro vidi na zadnjem platnu). Taj kut ovisi samo o parametrima $\alpha$ i $\beta$.
Ako eksponencijalnu spiralu $r=e^{-\varphi}$ zoomirate prema ishodištu (tj. ako se približavate prema ishodištu), vidjeti ćete potpuno istu spiralu kao na početku, samo zarotiranu oko ishodišta za neki kut. To je lako razumjeti: zoomirati spiralu prema ishodištu (s nekim faktorom $k$ između nula i jedan) je isto što i promatrati spiralu $r=ke^{-\varphi}$. Ako napišemo $k=e^{c}$, gdje je $c=\ln k$, onda je $r=ke^{-\varphi}=r^{-(\varphi -c)}$, tj. početna je spirala zarotirana za kut od $c$ radijana. Isto vrijedi i ako se udaljavate od ishodišta. Prema tome, eksponencijalna spirala izgleda podjednako i iz bliza i iz daleka!
Eksponencijalna spirala ima beskonačno mnogo okretaja i duž nje točka putuje kroz beskonačno vrijeme. Međutim, ta je spirala (začudo) konačne duljine! Štoviše, njenu duljinu nije teško izračunati sredstvima integralnog računa (treba koristiti formulu za duljinu luka krivulje u polarnom sustavu, koju znamo iz Matematičke analize 2).
Kao funkcija kuta $\varphi$, udaljenost točke od ishodišta $r=e^{-\varphi}$ konvergira prema nuli eksponencijalnom brzinom kad $\varphi\to+\infty$, tj. iznimno brzo.

Ova svojstva ima ne samo spirala oblika $r=e^{-\varphi}$, nego i bilo koja druga eksponencijalna spirala oblika $r=a^{\varphi}$, gdje je baza $a$ bilo koji pozitivan broj različit od $1$. Na gornjem platnu je $a=e^{-1}$.

Fokus

U oba je slučaja ishodište $0$ jedina stacionarna točka (primijetite da je matrica $\mathbf A=\alpha\mathbf I+\beta\mathbf J$ regularna, jer joj je determinanta jednaka $\alpha^2+\beta^2>0$). Ta stacionarna točka je bitno različita od centra, jer je okružena sa spiralnim trajektorijama, a ne s kružnim. Za takvu stacionarnu točku kažemo da je fokus.

Stabilna matrica

Za slučaj matrice oblika $\mathbf A=\alpha \mathbf I+ \beta\mathbf J$, vektor $\mathbf x(t)=e^{\mathbf At}\mathbf x_0$ će konvergirati k 0 (tj. prema ishodištu) kad $t\to+\infty$ za bilo koji početni vektor $\mathbf x_0$, onda i samo onda ako je $\alpha<0$, tj. onda i samo onda ako vlastite vrijednosti matrice $\mathbf A$ imaju negativne realne dijelove. Razlog je taj što je u ovom slučaju $\mathbf x(t)=e^{\alpha t}e^{\beta t\mathbf J}\mathbf x_0$. Kažemo da je matrica $\mathbf A$ stabilna matrica.

Općenitije, za bilo koju kvadratnu matricu $\mathbf A$ reda $n$ kažemo da je stabilna matrica ako su sve njene vlastite vrijednosti strogo lijevo od imaginarne osi kompleksne ravnine (tj. realni dijelovi svih vlastitih vrijednosti od $\mathbf A$ su negativni). Pokazuje se da je matrica $\mathbf A$ stabilna onda i samo onda ako $e^{\mathbf At}$ konvergira prema nul-matrici kada $t\to+\infty$. (Dokaz vidi u [Aglić Aljinović, Elezović, Žubrinić].)

Ako je matrica $\mathbf A$ stabilna, onda u matrici $e^{\mathbf At}$ svi matrični koeficijenti konvergiraju (kao funkcije od $t$) prema nuli eksponencijalnom brzinom, kada $t$ teži u $+\infty$.

Nelinearni dinamički sustavi

Literatura

Neven Elezović: Kompleksni brojevi, Element, Zagreb
Neven Elezović, Andrea Aglić Aljinović, Darko Žubrinić: Linearna algebra, Element, Zagreb 2020.
Luka Korkut i Vesna Županović: Diferencijalne jednadžbe i teorija stabilnosti, Element, Zagreb 2009.
Steven Strogatz (istaknuti stručnjak za dinamičke sustave): Two dimensional linear systems, video snimak predavanja, traje 1 sat i 15 min.

Matematika za osnovne i srednje škole

Matematičke šale

History of Croatian science

Školovanje Nikole Tesle u Hrvatskoj