KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB doc. dr. Tit Turnšek KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Naslov publikacije: Kvalitativna analiza diferencialnih enačb Nosilec avtorskih pravic: Visoka šola za upravljanje podeželja Grm Novo mesto, Ljubljanska cesta 28, 8000 Novo mesto 1. Izdaja Avtor publikacije: doc. dr. Tit Turnšek Lektoriranje: doc. dr. Alenka Divjak Recenzija: doc. dr. Franc Brcar, ddr. Janez Usenik Ime in sedež založnika: Visoka šola za upravljanje podeželja Grm Novo mesto, Ljubljanska cesta 28, 8000 Novo mesto Leto izida: 2025 Leto izdelave: 2025 Število izvodov: Izdaja v E-obliki CIP podatki Kataložni zapis o publikaciji (CIP) pripravili v Narodni in univerzitetni knjižnici v Ljubljani COBISS.SI-ID 235622915 ISBN 978-961-94364-7-9 (PDF) KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB KAZALO UVOD 1 1. SMERNO POLJE 3 2. KONVEKSNOST IN KONKAVNOST FUNKCIJ – TOČKA PREVOJA 16 3. KVALITATIVNA ANALIZA NAVADNIH DIFERENCIALNIH ENAČB 23 4. KVALITATIVNA ANALIZA SISTEMA LINEARNIH DIFERENCUIALNIH ENAČB PRVEGA REDA 56 5. KVALITATIVNA ANALIZA SISTEMA DVEH NELINEARNIH DIFERENCIALNIH ENAČB 105 6. GRAFIČNA ANALIZA STABILNOSTI – DIAGRAM PAJČEVINE 119 ZAKLJUČEK 125 LITERATURA 126 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Recenzija – doc. dr. Franc Brcar Dr. Tit Turnšek je v preteklosti že izdal dve knjigi, tako da je knjiga »Kvalitativna analiza diferencialnih enačb« nadaljevanje njegovega dela. Knjiga je odlično dodatno branje ali dodatno študijsko gradivo predvsem za študente na višjih stopnjah študija s poudarkom na področju upravljanja podeželja in kmetijstva. Vsebina: Uvod 1) Smerno polje 2) Konkavnost in konveksnost, točka prevoja 3) Kvalitativna analiza navadnih diferencialnih enačb 4) Kvalitativna analiza sistema linearnih diferencialnih enačb prvega reda 5) Kvalitativna analiza sistema dveh nelinearnih diferencialnih enačb 6) Grafična analiza stabilnosti – diagram pajčevine Zaključek Literatura Vsebina knjige je zahtevna in terja od bralca poglobljeno znanje in razmišljanje. Teorija je ustrezna in dodatno razložena s številnimi primeri za lažje razumevanje. Strokovni izrazi so ustrezno razloženi, tako da knjiga predstavlja pomemben prispevek k slovenskemu strokovnemu izrazoslovju. V prvem poglavju je prikazana grafična analiza diferencialnih enačb in uporaba programske opreme Maxima. Drugo poglavje je namenjeno konkavnosti, konveksnosti in točki prevoja. Diferencialne enačbe so grafično prestavljene s programom Maxima, ki je uporabljen tudi v naslednjih poglavjih. V tretjem poglavju je razložen namen kvalitativne analize navadnih diferencialnih enačb. Navadne diferencialne enačbe včasih nimajo eksplicitne oblike in v takšnih primerih je smiselna izvedba kvalitativne analize. Tak pristop pogosto da boljši vpogled v delovanje sistema. Poglavje je obogateno z dodatnimi primeri in programom Maxima. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Četrto in peto poglavje obravnavata kvalitativno analizo sistema linearnih diferencialnih enačb prvega reda in kvalitativno analizo sistema dveh nelinearnih diferencialnih enačb. Tudi v teh dveh poglavjih je obilo rešenih grafičnih primerov. V šestem poglavju je prikazana grafična analiza stabilnosti kot dodatna metoda grafične analize diferencialnih enačb. V zaključku avtor povzame vsebino knjige. Knjiga obravnava zahtevno tematiko in poglobljeno znanje in je pomemben prispevek k razumevanju kvalitativne analize diferencialnih enačb. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Recenzija – ddr. Janez Usenik Knjigo Kvalitativna analiza diferencialnih enačb lahko razumemo kot matematični dodatek knjigama Sistemska dinamika in Populacijska dinamika, ki ju je avtor dr. Tit Turnšek izdal leta 2020 in 2024. Pri tem pa je avtor poudaril, da svoje novo delo predlaga v branje kot dodatno gradivo zlasti magistrskim in doktorskim študentom s področja upravljanja podeželja, ekologije, biologije in kmetijstva. K temu bi jaz dodal, da je to koristno branje tudi za študente tehniških fakultet. Vsebino knjige je avtor razdelil na 6 poglavij: Uvod 1 Smerno polje 2 Konkavnost in konveksnost, točka prevoja 3 Kvalitativna analiza navadnih diferencialnih enačb 4 Kvalitativna analiza sistema linearnih diferencialnih enačb prvega reda 5 Kvalitativna analiza sistema dveh nelinearnih diferencialnih enačb 6 Grafična analiza stabilnosti – diagram pajčevine Zaključek Literatura Vsa poglavja so podprta s številnimi primeri, s katerimi avtor ilustrira teoretično razlago. Primeri so dobrodošli pri branju in razumevanja teksta, ki ni enostaven in od bralca terja zbranost in osredotočenost na vsebino. Pohvalno je še dejstvo, da so v besedilu navedeni tudi angleški nazivi uporabljenih pojmov, pri čemer je lahko marsikateri slovenski izraz avtorjev predlog za slovensko besedišče. Ob tem bi bilo morda smiselno omeniti možnost kreiranja večjezikovnega slovarja uporabljenih strokovnih terminov, npr. poleg slovenskega in angleškega morda še nemški, francoski, hrvaški … izraz. V prvem poglavju (Smerno polje) avtor izčrpno opiše grafično analizo diferencialnih enačb. V štirih dodatkih prikaže možnost uporabe kratkih računalniških programov z aplikacijo Maxima, s katerimi uspešno rešujemo diferencialne enačbe. V drugem poglavju avtor vpelje pojme konkavnost, konveksnost, točka prevoja. Tudi tu vse nazorno ilustrira z grafičnimi predstavitvami diferencialnih enačb, v KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB dodatku poglavja pa predstavi reševanje konkretnih primerov, spet z uporabo Maxime. Tretje poglavje Kvalitativna analiza navadnih diferencialnih enačb obravnava možnost, da diferencialna enačba nima eksplicitne oblike ali pa morda takšne rešitve niti ne potrebujemo. V tem primeru si lahko pomagamo s kvalitativno analizo, ki marsikdaj lahko daje boljši vpogled v obnašanje sistema kot pa eksplicitna rešitev v obliki neke formule. Z veliko primeri, ki so analizirani tudi grafično, postane razumljiv avtorjev namen bralcu pojasniti smisel kvalitativne analize. Poglavje zaključuje sedem ilustrativnih in koristnih dodatkov, kjer si spet lahko pomagamo z uporabo Maxime. Naslednje poglavje obravnava kvalitativno analizo sistema linearnih diferencialnih enačb prvega reda. Avtor vpelje nekaj novih pojmov, ki jih potrebuje za pojasnitev celotnega besedila tega poglavja. Tudi tu najdemo veliko rešenih primerov in grafov, pojasnjenih tudi v kar trinajstih dodatkih. Program Maxima je uporabljen tudi tu. V petem poglavju z naslovom Kvalitativna analiza sistema dveh nelinearnih diferencialnih enačb avtor preide na splošno obravnavo sistema diferencialnih enačb. S tremi dodatki dodatno ilustrira uporabo aplikacije Maxima. V zadnjem poglavju z naslovom Grafična analiza stabilnosti – diagram pajčevine pa avtor bralca spozna še z eno grafično metodo analize diferencialnih enačb. V prvem dodatku je za ilustracijo prikazan še diskretni populacijski model, ki ga opišemo z diferenčnimi enačbami. Ob koncu dodatka sledi še kratka obrazložitev Maxima ukaza »staircase«. V zaključku avtor na kratko (in jedrnato) povzame snov svoje knjige. Uporabljena in bralcu priporočena literatura našteje šest naslovov, vse v angleškem jeziku. Celoten tekst je precej zahteven in terja poglobljeno branje. Knjige ni mogoče prebirati brez pomoči svinčnika in papirja, kar velja za matematiko nasploh. Zaželena je tudi uporaba računalnika, težko gre brez njega. Poudariti je potrebno, da je dr. Turnšek eden od pionirjev računalništva v Sloveniji, saj je bil »računalnikar« že pred šestdesetimi leti. Prav to daje knjigi še posebno dodatno kvaliteto. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB UVOD Kvalitativna analiza diferencialnih` enačb je dodatno gradivo višje matematike, namenjeno študentom magistrskega in doktorskega študija s področja upravljanja podeželja, ekologije, biologije in kmetijstva. Z matematičnega vidika je to področje diferencialnih enačb. Velikokrat ni mogoče najti eksplicitne rešitve diferencialne enačbe. Zelo pogosto pa nam povedo kvalitativni aspekti diferencialne enačbe mnogo več kot pa eksplicitna rešitev. Kaj pa so kvalitativni aspekti diferencialne enačbe? Grafična podoba rešitve diferencialne enačbe - integralna krivulja v vektorskem polju - v ravnini, ki jo tvorita odvisna in neodvisna spremenljivka ali pa tudi dve med seboj odvisni spremenljivki. Če pa gre za sistem diferencialnih enačb, govorimo o fazni ravnini, ki jo tvorita obe spremenljivki. Njuno medsebojno odvisnost imenujemo parametrične krivulje ali trajektorije. Integralna krivulja, prav tako tudi parametrična krivulja, pove bistveno več kot pa eksplicitna rešitev diferencialne enačbe oziroma sistema diferencialnih enačb. Ali je na nekem danem območju integralna krivulja ali pa trajektorija – ali narašča ali upada? Ali je na nekem območju integralna krivulja konkavna ali konveksna? Zanima nas, kje je točka prevoja. Zanima nas, ali je rešitev enačbe stabilna ali ne. Zanimajo nas stacionarne rešitve ali kritične točke. Zanima nas odvisnost integralnih krivulj ali trajektorij od začetne vrednosti ali pa od vrednosti posameznih parametrov. Vsemu temu je namenjena `Kvalitativna analiza diferencialnih enačb.` Kvalitativno analizo obravnavamo v naslednjih poglavjih: 1. Uvodno poglavje, imenovali smo ga `Smerno polje`, pojasnjuje grafično analizo diferencialnih enačb. 2. V drugem poglavju `Konkavnost konveksnost - točka prevoja` se pojasnjujeta dve grafični diferencialni enačbi. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 3. V tretjem poglavju `Navadne diferencialne enačbe` obravnavamo kvalitativno analizo navadnih diferencialnih enačb. 4. V četrtem poglavju `Kvalitativna analiza sistema linearnih diferencialnih enačb` preidemo na obravnavo sistema dveh, sicer linearnih diferencialnih enačb. 5. V petem poglavju `Kvalitativna analiza sistema dveh nelinearnih diferencialnih enačb` preidemo na splošno obravnavo sistema diferencialnih enačb. 6. V šestem poglavju `Grafična analiza stabilnosti - diagram pajčevine` je prikazana še ena od grafičnih metod analize, ki je večkrat uporabljena. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 1. SMERNO POLJE Naj bo y F x y , y  y x funkcija neodvisne spremenljivke x  , torej   ali   0 za I a b x :     dy na intervalu x y . Geometrijski pomen odvoda funkcije dx v dani točki x y , predstavlja naklon (angl. Slope) oziroma smerni koeficient   tangente na graf funkcije v tej točki, to je tangens kota, ki ga oklepa tangenta z dy  tan    osjo X (slika 1.1). dx  dy  tg( )  Slika 1.1: Ilustracija h geometrijskemu pomenu odvoda dx Primer 1 (Dawkins, 2007, str. 9). Vzemimo preprosto diferencialno enačbo 1: dy  9.8 0.196   y . (En. 1) dx y  dy  30 3, 92 Pri neki vrednosti x naj bo , tedaj je , odvod je torej pozitiven, dx to pomeni, da je tukaj funkcija y naraščajoča. Pri različnih vrednostih y lahko določimo naklone tangent. Ker je desna stran diferencialne enačbe 1 zvezna funkcija, se predznak naklona tangente spremeni, ko gre le ta skozi vrednost 0, torej pri vrednosti y  50 (glej sliko 1.2): 0  9.8 0.196   y  y  50 . Rešiti diferencialno enačbo   y f x y  ,  na intervalu a   x b pomeni najti družino krivulj, ki predstavljajo rešitve enačbe. Oziroma gledano geometrično, KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB najti družino krivulj, od katerih vsaka predstavlja rešitev in ima v vsaki točki naklon y . Te krivulje imenujemo integralne krivulje (angl. Integral curves ali Solution curves). Če že ne moremo najti analitične rešitve diferencialne enačbe y  f x y  ,  , pa lahko narišemo majhne usmerjene daljice v vsaki točki ravnine X-Y, kjer je x znotraj intervala I , ki kažejo nagib integralne krivulje v tej točki. Na primer: če je v neki točki x , y y    2 a a  , potem v tej točki narišemo kratko usmerjeno daljico z naklonom 2 . Te majhne usmerjene daljice imenujemo linijski elementi (angl. Line elements ali Lineal elements). Vsi ti elementi pa predstavljajo smerno polje (angl. Direction field). V primeru diferencialne enačbe 1 bi lahko sestavili tabelo: y : 70, y   3.92 ; y : 60, y  1.96 ; y : 50, y  0 ; y : 40, y  1.96 ; y : 30, y   3.92 … in nato skicirali smerno polje, kot je prikazano na sliki 1.2. Slika 1.2: Risanje smernega polja Smerno polje in integralne krivulje lahko narišemo s pomočjo različnih programskih paketov. V naših primerih je uporabljen programski paket Maxima. Glej sliki 1.3 in 1.4 (dodatek 1.1). 100 80 60 40 20 0 0 1 2 3 4 5 Slika 1.3: Dif. enačba 1 - smerno polje KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 1.4: Dif. enačba 1 - integralne krivulje Primer 2 (Dawkins, 2007, str. 16). Vzemimo diferencialno enačbo 2: dy   y x (En. 2) dx Če postavimo za odvod konstantne vrednosti (enačba 3), C   y x (En. 3) potem enačba 2 pri različnih vrednostih konstante C predstavlja družino krivulj, ki jih imenujemo izokline (angl. Isoclines). Izokline povezujejo točke, v katerih imajo integralne krivulje enak naklon, kot je prikazano na sliki 1.5. Izokline, kjer je , imenujemo ničte izokline dy  0 (angl. dx Null clines). Celotno smerno polje je prikazano na sliki 1.6 in integralne krivulje na sliki 1.7 (dodatek 1.2). 3 2 1 0 -1 -2 -3 -3 -2 -1 0 1 2 3 Slika 1.5: Dif. enačba dy   y x dx - izokline KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 2 1.5 1 0.5 0 -0.5 -1 -1.5 -2 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 Slika 1.6: Dif. enačba dy   y x dx - smerno polje Slika 1.7: Dif. enačba dy   y x dx - integralne krivulje Primer 3 (Dawkins, 2007, str. 13). Analizirajmo in določimo smerno polje in integralne krivulje diferencialne enačbe 4: dy 2 2     f y  y     y 2    1 y (En. 4) dx Ugotovimo, kje so ničte izokline, to se pravi točke, kjer je odvod 0. f y     2 y     y 2   1 y 2  0   y   2   y    1   y2 1  0 Ničte izokline so: y  2 y 1 y ,   1 in . Sedaj že lahko določimo linijske elemente v teh točkah (slika 1.8). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 1.8: Linijski elementi za y=2, y=1 in y=1 Za y  dy  2 0 je - linijski elementi so usmerjeni navzgor. Za 1   y 2 je dx dy dy  0    1 y 1  - linijski elementi so usmerjeni navzdol. Za 0 je - dx dx linijski elementi so usmerjeni navzdol. Za y   dy 1 0  je - linijski elementi dx so usmerjeni navzgor. Celotno smerno polje je prikazano na sliki 1.9 in integralne krivulje na sliki 1.10 (dodatek 1.3). 3 2 1 0 -1 -2 -3 0 1 2 3 4 5 Slika 1.9: Dif. enačba 2 2 ( y     y 2) (1 y )- smerno polje Slika 1.10: Dif. enačba 2     2 ( y y 2) (1 y )- integralne krivulje KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Vzemimo avtonomno diferencialno enačbo dy  f y   . Rešitve enačbe dx dy  f y    0 so neodvisne od spremenljivke x , so torej konstante. Označevali dx jih bomo z y r in imenovali stacionarne rešitve (angl. Stationary solutions) ali ravnotežne rešitve (angl. Equilibrium solutions). Različni avtorji uporabljajo tudi sledeče sinonime: fiksne točke (angl. Fixed points) ali kritične točke (angl. Critical points) ali točke ravnotežja (angl. Equilibrium points) (Nagy, 2015, str. 242). Primer 4 Analizirajmo diferencialno enačbo 5: dy   a y (En. 5) dx Rešitev enačbe je:   a x  y C e in je določena z začetno vrednostjo  x , y x 0 0   0 ter z vrednostjo parametra a . Pri tem je C poljubna konstanta. Glede na predznak parametra a imamo tri možnosti: 1. a  a x  lim C e    0  , če je C  0 in je  , če je C  x 0  2.    C a  a x 0  C e 3. a  lim a x  C e   0 0  x  Če je a  0 in se samo malo spremeni, se posledično radikalno spremeni tudi obnašanje rešitve enačbe. Pravimo, da ima diferencialna enačba y    a x bifurkacijsko točko pri a  0 . Integralne krivulje za a  0 so podane na sliki 1.11 in za a  0 na sliki 1.12. (dodatek 1.4). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 1.11: Dif. enačba 5 a  0 - integralne krivulje Slika 1.12: Dif. enačba 5 a  0- integralne krivulje KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 1.1 Primer 1 Diferencialna enačba: dy  9.8 0.196   y dx Smerno polje: slika 1.3 (Maxima): load(draw); load(drawdf); set_draw_defaults(xrange = [-1.5,1.5],yrange = [-1.5,1.5], grid = true, line_width = 2); drawdf([9.8-0.196*y],[x,y] ,[x,0,5],[y,0,100],color=black, field_grid=[10,10]); Integralne krivulje: slika 1.4 (Maxima): diff(y,x)=9.8-0.196*y; load("plotdf"); plotdf([9.8-0.196*y],[x,y],[y,0,80],[x,0,8]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 1.2 Primer 2 Diferencialna enačba: dy   y x dx Poiščimo rešitev diferencialne enačbe: y    y x y y x      Splošna oblika diferencialne enačbe   y p x tega tipa je:     y g x   in njena rešitev je: y x         x  g x dx C   1     p x dx C    2 x e    , pri tem je   x V obravnavnem primeru je p x   g x   x   1 in     dx C 2        x C  x e 2  x  e        x x  C e 3            x x  g x dx  C e    3   x dx  a x  Formula: a x  e x e   dx     2 (a x 1) a Za: a  1  x dx   e    3  C 3     x x  x C e      x 1    C e   x1 3            x x g x dx   C e  x 1 ; y x   C e   x    1 3  x  C4   x C e  3 Rešitev diferencialne enačbe je: y x     C e x 1 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Določimo izokline za različne vrednosti k   y x . k je konstanta po izbiri. k  3  y   3 x   x  0, y  3; x  3, y  0  , k  2  y   2 x   x  0, y  2; x  2, y  0  , k  1  y   1 x   x  0, y  1; x  1, y  0  , k  0  y  x , k  1  y    1 x   x  0, y  1; x  1, y  0  , k  2  y    2 x   x  0, y  2; x  2, y  0  , k  3  y    3 x   x  0, y  3; x  3, y  0  . Izokline: slika 1.5 (Maxima) load(draw); load(drawdf); set_draw_defaults( xrange=[-3, 3], yrange=[-3, 3], grid=true, line_width=1, drawdf([y-x], [x,-3,3],[y,-3,3],field_grid=[8,8], color=blue, line_width = 3,explicit(3+x, x, -3,3), explicit(2+x, x, -3,3),explicit(1+x, x, -3,3), explicit(x, x, -3,3),explicit(-1+x, x, -3,3), explicit(-2+x, x, -3,3),explicit(-3+x, x, -3,3)); Smerno polje: slika 1.6 (Maxima) load(draw); load(drawdf); set_draw_defaults( KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB xrange=[-2, 2], yrange=[-2, 2], grid= true, line_width=1, color=black); multiplot_mode(screen)$ drawdf([y-x], [x,-2,2],[y,-2,2],field_grid=[15,15]); Integralne krivulje: slika 1.7 (Maxima) 'diff(y,x)=y-x; load("plotdf"); plotdf([y-x],[x,y],[y,-5,5],[x,-5,5]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 1.3 Primer 3 Diferencialna enačba: dy 2   y     y 2   1 2 y dx Določimo ničte izokline :  2 y     y 2   1 y 2  0   y 2 1  0  y  1 ;  2 y    y 2  0  y  1 in y  1 Smerno polje: slika 1.9 (Maxima) load(draw); load(drawdf); set_draw_defaults( xrange=[0,1], yrange=[0,2], grid=true, line_width=2); drawdf([(y^2-y-2)*(1-y)^2], [x,0,5], [y,-3,3],field_grid=[15,15], line_width=3, color=red, polygon([[0,2],[3,2],[5,2]]), polygon([[0,1],[3,1],[5,1]]), polygon([[0,-1],[3,-1],[5,-1]])); Integralne krivulje: slika 1.10 (Maxima) 'diff(y,x)=(y^2-y-2)*(1-y)^2; load("plotdf"); plotdf([(y^2-y-2)*(1-y)^2],[x,y],[y,-2,3],[x,0,5]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 1.4 Primer 4 Diferencialna enačba: dy   a y dx Rešitev diferencialne enačbe: dy   a dx y   dy       a dx C   a x C a x  1  ln y    a x C y  e 1 y   C e 1   y Rešitev diferencialne enačbe je torej: y   a x  C e . Pri tem je C poljubna konstanta in določena z začetno vrednostjo  x , y x 0 0   0 . Integralne krivulje za a  0 : slika 1.11 (Maxima) eq: 'diff(y,x)=a*y; load("plotdf"); plotdf([0.5*y],[x,-10,10], [y,-10,10] ); Integralne krivulje za a  0 a  0 : slika 1.12 (Maxima) eq: 'diff(y,x)=a*y; load("plotdf"); plotdf([-0.5*y],[x,-10,10], [y,-10,10] ); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 2. KONVEKSNOST IN KONKAVNOST FUNKCIJ – TOČKA PREVOJA Konveksnost - definicija 1 Funkcija f x x  ,   a b je konveksna ( navzgor ukrivljena ) na intervalu  , če je izpolnjen pogoj: f x    f a   f b    f a    (En. 1) x a  b a  Enačba sekante: f x f b    f a    f a          x a (En. 2) b a  Slika 2.1: Konveksna funkcija Graf funkcije f x leži pod premico (sekanto), ki gre skozi točki: a f a , in       b f b ,   (slika 2.1). Enačbo 1 lahko zapišemo še malo drugače. Izberimo neko pozitivno število h , tako da velja: 0   h b  a    a h b . Pogoj za konveksnost funkcije lahko zapišemo: f a f b    f a     h   f a     h (En. 3) b a  Z    x  a označimo razmerje: , ki je manjše a1i kvečjemu enako: 1 b  a KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB  x  a         b  a x a  x   .  1      a  b Označimo b  a    1  1    in.    x      a  b 2 1 2 Prav tako velja tudi  f x    f a     f b    f a   f x    f a      f b      f a    f x       1     f a    f b    f x      f a      f b 1 2   Enačbo 1, ki je pogoj za konveksnost funkcije f x , lahko sedaj zapišemo tudi   v sledeči obliki (enačba 4): f        a  b  f a    f b 1 2  1   2   (En. 4) Konkavnost - definicija 1 Funkcija f x je konkavna (navzdol ukrivljena) na intervalu x     a b , , če je izpolnjen pogoj f x f b    f a    f a          x a (En. 5) b a  Slika 2.2: Konkavna funkcija Graf funkcije f x leži nad premico (sekanto), ki gre skozi točki: a f a ,      in b f b , (slika 2.2).     KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Na enak način, kot smo izpeljali še drugi dve obliki pogoja za konveksnost, ju lahko izpeljemo tudi za konkavnost. Pogoj za konkavnost funkcije lahko zapišemo, podobno kot za konveksnost, (glej enačbo 3) v obliki enačbe 6: f a f b    f a     h   f a     h (En. 6) b a  Pri tem je h neko pozitivno število h , tako da velja 0   h b  a    a h b Lahko pa pogoj za konkavnost zapišemo tudi v obliki enačbe 7 (glej tudi enačbo 4). f        a  b  f a    f b 1 2  1   2   (En. 7) Pri tem je: 1           1 1 , 2 in  1 2 ter: x      a  b f x    f a    f b 1 2 in   1   2   Preden nadaljujemo, napravimo predpostavko, da je f x funkcija dvakrat   odvedljiva. Konveksnost - definicija 2 Funkcija f x je na intervalu a b , konveksna, če v vsaki točki tega     intervala tangente na to funkcijo ležijo pod krivuljo f x (slika 2.3).   Slika 2.3: Konveksna funkcija KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Prvi odvod funkcije  f   x z naraščajočim x narašča od negativnih vrednosti k pozitivnim vrednostim. To pomeni, da je funkcija na intervalu a b ,   konveksna natanko takrat, ko je drugi odvod funkcije pozitiven f    x  0 . Konkavnost - definicija 2 Funkcija f x je na intervalu a b , konkavna, če v vsaki točki intervala     a b ,  tangente te funkcije ležijo nad krivuljo f x   (slika 2.4). Slika 2.4: Konkavna funkcija Prvi odvod funkcije  f   x z naraščajočim x upada od pozitivnih vrednosti k negativnim vrednostim. To pomeni, da je funkcija f x na intervalu a b ,     konkavna natanko takrat, ko je drugi odvod funkcije negativen f    x  0 . Zaključek Ko analiziramo ukrivljenost funkcije f x , je ta lahko konveksna ali konkavna.   Pri tem obstajajo štiri možnosti: 1. f  x  0 f  x  0   ,   - krivulja je konveksna in padajoča; 2. f  x    0 , f   x  0 - krivulja je konveksna in naraščajoča; 3. f  f  x 0   x  0  ,   - krivulja je konkavna in padajoča; 4. f  x  0 f x    ,    0 - krivulja je konkavna in naraščajoča. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Točka prevoja Če ima funkcija f x točko, v kateri preide funkcija iz konveksne v konkavno   ali obratno, iz konkavne v konveksno, potem jo imenujem točka prevoja (angl. Inflection point). Glej sliko 2.5. Večkrat se uporabljajo tudi izrazi: prevojna točka, prevoj, obračaj ali točka infleksije. Označimo točko prevoja z a . V tej točki je  f   a  0 . Na eni strani točke prevoja je drugi odvod pozitiven in na drugi negativen. Slika 2.5: Točka prevoja Brez dokazovanja: če je f    a  0 , pa to še ni zadosten pogoj, da je v točki a res prevoj  . Pogoj je tudi, da je tretji odvod f   a  0 . Če je tudi tretji odvod v točki a enak 0 , potem odvajamo funkcijo f x toliko časa, da dobimo odvod,   ki je različen od nič. Če je red višjega odvoda v točki a liho število, potem je v točki a prevoj, če pa je sodo število, v točki a ni prevoja. Za ilustracijo si poglejmo sledeči primer. Vzemimo funkcija   2 x f x   e (na sliki 2.6: modra krivulja). Prvi odvod funkcije (na sliki 2.6 : rdeča krivulja): f x  2 x          f x 2 x e f    x  0    x 0 f  x  0 0 x za in   za     Drugi odvod f x  f funkcije       x :  (na sliki 2.6: črna krivulja). f      2 2 2 x   x     x 2 x e     x f  x    2 e       2 x  2 x e   KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB f       2    x x 2 2 x 4 e f  1 1       x  x     0 x za in 2 2 f    1 1    x  0 x za 2 2 Pogoj za točko prevoja je določen z drugim odvodom, ki je v točki prevoja nič. Točka prevoja:  1  2    4 2  2 x 2 x e  0     4 x 2   x   0  1, 2 2 Točki prevoja 1 1 (označimo ju z a in a ) sta:   in   (dodatek 2.1). 1 2 1 1 a a 2 2 Slika 2.6: Primer točke prevoja: f x   (modra), f    x (rdeča), f   x črna Funkcija je 1 1 konveksna na intervalu od  do  in od do  . Na 2 2 intervalu od 1 1  do pa je funkcija konkavna. 2 2 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 2.1 Funkcija : f x  2 x    e , prvi odvod  2 x : f    x     2 x e , drugi odvod f  2  x x  4  x   2 e :    2 . Konkavnost, konveksnost, primer - slika 2.6 (Maxima) load(draw); set_draw_defaults( xrange=[−3, 3], yrange=[−2, 1.5], grid=true, line_width=3, color=black); draw2d ( color=blue, yrange=[−2, 1.5], line_width=2, explicit (%e^−(x^2), x, −3, 3), color=red, explicit (−2·x·%e^−(x^2), x, −3, 3), color=black, explicit ((4·x^2−2)·%e^−(x^2), x, −3, 3), line_width=1, polygon ([[−3,0], [0,0], [3,0]]), polygon ([[0,−2], [0, 0], [0, 1.5]]), polygon ([[−0.707, 0], [−0.707,1], [−0.707, −1]]), polygon ([[0.707, 0],[0.707, 1], [0.707, −1]]))$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 3. KVALITATIVNA ANALIZA NAVADNIH DIFERENCIALNIH ENAČB Uvod Niso vse navadne diferencialne enačbe take, da bi lahko našli njihovo eksplicitno rešitev. Vendar vedno ne potrebujemo rešitve, pač pa pogosto zadostujejo tudi nekateri kvalitativni aspekti rešitve. Tudi če poznamo natančno rešitev, pa nam lahko kvalitativna analiza da boljši vpogled v obnašanje kot rešitev v obliki formule. Najpogosteje nas zanimajo ravnotežne rešitve in ali so te stabilne. Kvalitativna analiza avtonomnih diferencialnih enačb V nadaljevanju tega poglavja bomo obravnavali navadne nelinearne diferencialne enačbe prvega reda z eno neznano funkcijo. Pretežno bomo obravnavali avtonomno diferencialno enačbo. To je diferencialna enačba, kjer je neodvisna spremenljivka le v odvodu. Posplošeno je njena oblika sledeča (enačba 1): dy  y   f y ( ) (En. 1) dx Včasih seveda lahko dano diferencialno enačbo rešimo eksplicitno. Kvalitativno analizo diferencialne enačbe bomo prikazali na nekaj konkretnih primerih. Primer 1 Diferencialna enačba y   f y ( ) naj bo: y    a y b (En. 2) Naj bosta oba parametra pozitivna a b ,   0 . Eksplicitno rešitev podaja enačba 3 (glej dodatke ob koncu poglavja od 3.1 do 3.7): y   a x  b C e  (En. 3) a Integralne krivulje (diferencialne enačbe 2 ) so prikazane na sliki 3.1. Pri začetni vrednosti C x  y   b 0 je rešitev enačbe 3: . Če vzamemo, da je a konstanta C  0 , potem je rešitev diferencialne enačbe neodvisna od x y   b , to pa je tudi ničta izoklina, določena z y  0 (glej enačbo 2). a KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Če je konstanta C  0 , bo funkcija eksponentno naraščala, če pa je C  0 , bo eksponentno padala, kar je razvidno tudi s slike 3.1. Slika 3.1: Smerno polje in integralne krivulje dif. enačbe 2; a=2, b=1 Primer 2 Avtonomna diferencialna enačba y   f y ( ) naj bo: y   sin y (En. 4) Implicitno rešitev diferencialne enačbe 4 podaja enačba 5 (glej dodatek 3.2): 1 cos  y   x C e (En. 5) sin y Graf funkcije f y    sin y je podan na sliki 3.2. Identificirajmo ničle funkcije f y . Z drugimi besedami, iščemo točke, kjer je y    0 , oziroma točke, kjer je sin y  0 . Te točke so y   n  n  n ,  ... 2, 1, 0, 1, 2,.....    . To pa so istočasno tudi rešitve diferencialne enačbe 4, neodvisne od x . Te od x neodvisne rešitve imenujemo tudi stacionarne rešitve, pa tudi ravnotežne rešitve (angl. Equilibrium solutions), fiksne točke (angl. Fixed points) ali tudi kritične točke (angl. Critical points). Ravnotežne točke bomo ' klasificirali'. Kaj to pomeni? Če je rešitev 'blizu' ravnotežne rešitve in se bo od nje oddaljevala (angl. Move away), potem bomo tako ravnotežno rešitev imenovali nestabilna ravnotežna rešitev ali točka (angl. Unstable equilibrium solution point). Če pa bo rešitev 'blizu' ravnotežne rešitve in se ji bo približevala, pa jo bomo imenovali asimptotično stabilna ravnotežna rešitev (angl. Asymptotical stable equilibrium solution/point). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Uporabili smo izraz 'blizu'. Kaj to pomeni? Za trenutek pozabimo na naš primer, h kateremu se bomo vrnili. Vzemimo, da imamo neko avtonomno diferencialno enačbo, ki ima dve ravnotežni rešitvi (točki), in sicer y  0 10 r 0 in y  r 1 . Ti dve točki razdelita številčno premico na tri področja:     y 0 , 0   y 10 in 10    y . 'Blizu' ravnotežne točke y  10 0   y 10 10 y   r 1 pomeni, da je ali . 'Blizu' ravnotežne točke y  0     y 0 0   y 10 r 0 pa pomeni, da je ali . Naslednje vprašanje pa je, kako spoznati, ali je ravnotežna točka (asimptotično) stabilna ali nestabilna. Ena možnost je, da narišemo smerno polje in vanj nekaj integralnih krivulj ter pogledamo, h katerim točkam y rn se pri naraščajoči neodvisni spremenljivki integralne krivulje približujejo in od katerih se oddaljujejo (glej sliko 3.3). Druga možnost je bolj formalna. Vrnimo se h grafu funkcije f y    sin y (slika 3.2 f y  0 ) in identificirajmo intervale, kjer je   , in intervale, kjer je f y    0 f y  0 . Na intervalih, kjer je   , narišemo na številčni premici puščice v desno f y in na intervalih, kjer je    0 , puščice v levo. Na intervalih, kjer je f y  0 f y  0   , bo rešitev upadala in na intervalih, kjer je   , bo rešitev naraščala. To je logično, saj je y   f y   . Slika 3.2: Funkcija f y    sin y Ravnotežna rešitev (točka) je stabilna, če puščice z leve in desne kažejo proti tej točki, in je nestabilna, če puščice iz leve in desne kažejo proč od te točke. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB S slike 3.2 je razvidno, da sta točki y   y  r 4 in   r 4 stabilni in točke y    2  , 0, 2 r 1, 2, 3 r r nestabilne. To potrjujejo tudi smerno polje in integralne krivulje, prikazane na sliki 3.3 (dodatek 3.2). Slika 3.3: Primer 2 Smerno polje in integralne krivulje dif. enačbe y  sin y Primer 3 Naslednji primer avtonomne diferencialne enačbe je enačba tako imenovane logistične rasti (enačba 6), ki je eden od osnovnih modelov populacijske dinamike: dy  y      r y  1  (En. 6) dt  k  Eksplicitno rešitev enačbe 6 podaja enačba 7 (glej dodatek 3.3): y  k 1     r t , (En. 7) A e pri čemer je A je poljubna konstanta. Z biološkega vidika r  0 in k  0 nimata smisla. Naredimo graf funkcije dy  y   f       r y  1  za k  10 r ,  0, 5 . Prikazan je na sliki 3.4 (dodatek dt  k  3.3). Iz enačbe 6 in tudi s slike 3.4 lahko razberemo, da ima rešitev dve ravnotežni točki, to sta: y  0 y  10 r 0 in r 1 . V graf funkcije vnesemo puščice, ki kažejo v desno  , kjer je f y   0 (rešitev enačbe narašča - integralne krivulje naraščajo), KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB in v levo, kjer je f y    0 (rešitev enačbe upada - integralne krivulje upadajo). Prva ravnotežna točka je nestabilna ravnotežna rešitev, ker sta puščici z leve in desne usmerjeni od točke yr0 . Druga je asimptotično stabilna ravnotežna rešitev, ker sta puščici z leve in desne usmerjeni k točki yr1 . Slika 3.4: Graf funkcije y f y ( )     r y (1 ) ; parametri: k=10, r=0.5 k Določimo še področja, kjer so integralne krivulje konkavne (vbočene) in kjer so konveksne (izbočene) ter padajoče ali naraščajoče. Določimo torej še ukrivljenost integralnih krivulj. To lahko ugotovimo, ko izračunamo še drugi odvod funkcije y y , torej : y   y   r    f y              2 r y  1  y f y     r y y    k   k  y    2  r r  y (En. 8) k Točke, ki jih določa y   0 , so točke, kjer se spremeni ukrivljenost in ki jih imenujemo točke prevoja (angl. Inflection points). Premica, ki določa točke prevoja, je: 0 2  r 2  r k   r  y   y r y    k k 2 Za k  10 y  5 y je in razdeli ravnino X-Y na področja, kjer je pozitiven – integralne krivulje so konveksne  in kjer je y negativen – integralne krivulje so konkavne. Obstajajo štiri možnosti (Davis, 1962, str. 25-30): KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 1. y  0 y ,   0- integralne krivulje so konkavne in padajoče; 2. y  0   , y 0 - integralne krivulje so konkavne in naraščajoče; 3. y  0 y  0 , - integralne krivulje so konveksne in padajoče; 4. y 0 y   ,   0 - integralne krivulje so konveksne in naraščajoče. Grafa obeh funkcij 𝑓′(𝑦) in 𝑓"(𝑦) prikazuje slika 3.5. Slika 3.5: Grafa funkcij 𝒇′(𝒚) 𝒊𝒏 𝒇"(𝒚) S ′ ″ slike 3.5 se da hitro identificirati področja, kjer sta 𝒇 (𝒚) in 𝒇(𝒚) večja odnosno manjša od 0: za y : 10 do 15 y   0 , y  0 za y : 5 do 10 y 0 y 0   ,   za y : 0 do 5 y  0 , y  0 za y : ˗5 do 0 y   0 y 0 ,   Tako lahko ugotovimo ukrivljenost integralnih krivulj. To potrjuje tudi slika 3.6 (glej dodatek 3.3), ki prikazuje smerno polje in integralne krivulje logistične diferencialne enačbe. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 3.6: Primer 3 Smerno polje in integralne krivulje logistične dif. enačbe y f y ( )     r y (1 ) ; parametri: r=0.5, k=10 k Transformacija logistične diferencialne enačbe v brez-dimenzionalno obliko Na splošno je numerična rešitev diferencialne enačbe odvisna od vrednosti parametrov. Kvantitativno obnašanje enačbe pa je odvisno od manjšega števila parametrov in njihovih medsebojnih kombinacij (Olinick, 2011, str. 68). Za potrebe kvalitativne analize, posebno kadar enačba modelira kakšen pojav (biološki, družbeni, itd.), pogosto izluščimo parametre in enačbo (model) transformiramo v brez-dimenzionalno obliko. V logistični enačbi 6 redefinirajmo najprej spremenljivko stanja y , ki ima pri populacijski dinamiki dimenzijo 'populacija'. Prav tako ima tudi parameter k , ki predstavlja zgornjo mejo rasti, dimenzijo 'populacija'. Vpeljemo novo brez- x  y dimenzionalno spremenljivko stanja . Enačbo 6 delimo s k : k 1 dy 1  y  dx       r y 1         r x 1 x in dobimo: . k dt k  k  dt V enačbi 6 ima parameter r dimenzijo '1/čas'. Vpeljimo še spremenljivko    r t in upoštevajmo: d   r dt in tako dobimo: dx    x   d 1 x . (En. 9) Enačbo logistične rasti (enačba 6) smo transformirali v enačbo 9. Primer 4 Poglejmo si sledečo avtonomno diferencialno enačbo: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y       y 12   2 f y y  4 (En. 10) Graf te funkcije je na sliki 3.7. Ni težko videti, da ima enačba tri ravnotežne točke: y  2 y   y  2 r 0 , 1 y r 1 in r 2 . S slike razberemo, da je r 0 nestabilna ravnotežna rešitev y y in r 2 asimptotično stabilna rešitev . Kaj pa r1 ? Ta ravnotežna rešitev je z ene strani asimptotično stabilna, rešitve se ji asimptotično približujejo, z druge strani pa je asimptotično nestabilna, rešitve se od nje oddaljujejo. Tako ravnotežno rešitev imenujemo semi-stabilna ravnotežna rešitev (angl. Semi-stable equilibrium solution). Vse povedano pa lahko razberemo tudi iz smernega polja in integralnih krivulj (slika 3.8), (glej dodatek 3.4). Slika 3.7: Graf funkcije f y ( ) Slika 3.8: Smerno polje in integralne krivulje dif. enačbe        2 y 1  2 y f y  y  4 Primer 5 Avtonomna diferencialna enačba naj bo: y  2 y  y      r y  1   (En. 11)   2 k 1  y KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Poiščimo ravnotežne rešitve: r y   2 y y     1    0 (En. 12)    2 k 1 y Ena ravnotežna rešitev je y  0 r 0 , ki je nestabilna. Druge ravnotežne rešitve so odvisne od vrednosti parametrov   r in k . Graf funkcije y f y   za k  10 in r  0, 48 je prikazan na sliki 3.9. S slike 3.9 je razvidno (glej dodatek 3.5), da sta ravnotežni rešitvi y y y r 1 in r 3 stabilni , medtem ko je ravnotežna rešitev r2 nestabilna. Slika 3.9: Funkcija y  f y ( ) ; vrednosti parametrov : k = 10 in r =0,48 Ostale ravnotežne rešitve (poleg y  0 ) y  0 lahko dobimo tudi z rešitvijo enačbe 12: r y  2 y  y  y  y    1       0 r 1     2    2 k 1  y k 1   y Na f y  y   f    r sliki 3.10  1 sta prikazani funkciji 1 2 in  1  2 z vrednostjo y  k  parametra k  10 in različnih vrednosti parametra r (0.38, 0.45, 0.56, 0.7). S slike 3.10 je razvidno, da dobimo poleg ravnotežne rešitve yr 0 še eno, dve ali tri ravnotežne rešitve (kjer se sekata krivulji f1 in f 2 ). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 3.10: Funkciji y  y  in r    1 pri različnih vrednostih parametra r   2 1 y  k  Smerno polje z integralnimi krivuljami je prikazano na sliki 3.11 (glej dodatek 3.5). Slika 3.11: Smerno polje in integralne krivulje; Dif. enačba:   2 y y y      r y 1  ; r=0,48 in k=10      2 k 1 y Bifurkacija ali razcep Razcep ali bifurkacija (angl. Bifurcation ) je pojav, ko se število ravnotežnih točk spremeni kot posledica spremembe vrednosti katerega od parametrov diferencialne enačbe. Če se pri določeni vrednosti katerega od parametrov ta samo malo spremeni, bo obnašanje rešitve diferencialne enačbe bistveno drugačno. Odvisnost ravnotežnih točk od parametrov določimo s pomočjo bifurkacijskega diagrama (angl. Bifurcation diagram) (dodatek 3.5). Ta nam pomaga najti tudi tako imenovane razcepne ali bifurkacijske točke, v katerih se kvalitativna narava rešitve bistveno spremeni. Vrnimo se k enačbi 12 , ki določa ravnotežne točke: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB r y  2 y  y     1    0    2 k 1 y Predpostavimo, da je konstanta: k  10 in se ne spreminja. Iz enačbe 12 dobimo odvisnost kritične točke od parametra r . Kritična točka y  0 je trivialna in nas ne zanima. r y  2 y  y  y  y    1     0  r   1      2    2 k 1 y  k 1  y r   y  1 0.1 y    1  2 0 y Razcepni oziroma bifurkacijski diagram je prikazan na sliki 3.12. Slika 3.12: Bifurkacijski diagram – enačba 12; parameter k=10 Tudi iz razcepnega diagrama je razvidno, da ima enačba 12 eno, dve ali tri stabilne točke. Prav tako lahko razberemo, da ima enačba dve razcepni točki, in sicer r 1 in r 2 . Primer 6 Diferencialna enačba: dy  f y      y 1 y    h (En. 6) dx Določimo točke ravnotežja: y    2 1 y    2 h 0   y    y h 0  y    y h 0  KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y , y  r 1  1 4   h 1 r 2 . 2  1 Točki ravnotežja sta dve: y y y r 1 in r 2 ali ena sama : r3 , če je diskriminanta 2 kvadratne enačbe 0 . Mogoče pa je, da ravnotežne točke sploh ni, če je diskriminanta manjša od 0 . Za h  2 je y  1 y  r 1 = in 2 r 2 . Stabilnost ravnotežnih točk, pri dani vrednosti parametra h , ugotovimo iz funkcije f y      y  1 y   h , ki je prikazana na sliki 3.13, pri tem je parameter. h  2 y y . S slike 3.13 je razvidno, da je točka ravnotežja r 2 stabilna in točka r1 nestabilna (glej dodatek 3.6). Slika 3.13: Dif. enačba f y      y  1 y   h ; parameter h=-2 Bifurkacijski diagram  za funkcijo y  g h 1  1 4   h   r je prikazan na 2 sliki 3.14. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 3.14: Bifurkacijski diagram. Dif. enačba f y      y  1 y   h Smerno polje in integralne krivulje diferencialne enačbe 13 pri različnih vrednostih parametra h  0.125. 0.25, 0.8 so prikazane na sliki 3.15 - dve točki ravnotežja, 3.16 - ena semi stabilna točka in 3.17 - ni točke ravnotežja. Slika 3.15: Smerno polje in integralne krivulje - Dif. enačba f y      y  1 y  h ; parameter: h=0,125 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 3.16: Smerno polje in integralne krivulje - Dif. enačba f y      y  1 y   h ; parameter; parameter: h=0, 25 Slika 3.17: Smerno polje in integralne krivulje - Dif. enačba f y      y  1 y   h ; parameter: h=0,8 Integralne krivulje na vseh treh slikah 3.15, 3.16 in 3.17 kažejo popolnoma različna obnašanja rešitev diferencialne enačbe 13. Primer 7 (Davis, 1962, str. 25-29). S primerom želimo prikazati še eno metodo določanja oblik integralnih krivulj v ravnini (X-Y). V primeru 7 bomo analizirali ne-avtonomno diferencialno dy  f x y ,  enačbo prvega reda, katere splošna oblika je: . Konkretno bomo dx analizirali enačbo 14 : dy    x y  y   2 (En. 14) dx KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Splošno rešitev enačbe 14 podaja enačba 15 (dodatek 3.7): y  2 1   2 (En. 15) x C e Oziroma pri danih začetnih pogojih: x0 , y0 y 2  y0    (En. 16) y    y  2 2 x  x  0 2 e 0 0 Določimo področja, na katerih je y  0 , in področja, na katerih je y  0 . Postavimo y   0 y   0 , x y    y  2  0  x  0 , y  0 , y  2 Dobili smo tri premice, ki razdelijo ravnino (X-Y) na 6 področij (slika 3.18), na katerih je y  0 y 0 ali   : 1. Področje: ( x  0 in y  0 )  y  0 2. Področje: ( x  0 y 0 y 0 in  )    3. Področje :( x  0 in 0   y 2 )  y   0 4. Področje :( x  0 0   y 2 y  0 in )  5. Področje: ( x  0 in y  0 )  y  0 6. Področje: ( x  0 n y  0 )  y  0 Slika 3.18: Področja, kjer je y  pozitiven in kjer je negativen KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Izračunajmo drugi odvod Drugi odvod je podan z enačbo 17 (dodatek 3.7): y   d y 2   y            y 2 2 2 1 2 x y 2 x 2 (En. 7) dx Postavimo  y 2 2 0 y  y        x y x 0    2 1 2 2      y  0 y  2 , 1 2      2 1 2 y   1 x y 2 x  0   2 2 x Dobili smo dve premici in krivuljo, ki razdelijo ravnino (X-Y) na 7 področij, ki določajo, kje je y  0 in kje je y  0 . Glej sliko 3.19, slika 3.20 pa dodatno prikazuje tudi smerno polje. V točkah prevoja ( y  0 ) se spremeni ukrivljenost integralnih krivulj. . Slika 3.19: Področja, kjer sta y in y pozitivna in kjer sta negativna KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 3.20: Smerno polje. Diferencialna enačba y     x y ( y 2) Smerno polje z integralnimi krivuljami je prikazano na sliki 3.21, ki dodatno potrjuje zgornjo analizo poteka integralnih krivulj. Slika 3.21: Smerno polje in integralne krivulje diferencialne enačbe 14 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.1 Primer 1 Diferencialna enačba y    a y b  dy 1 Formula:   ln ax b  a y b   a Eksplicitna rešitev  dy 1    a x C    dx C   ln a y   b   C a y b    e 0   0    x 0  a y b a     a x   a C   a x      a y b e e 0  a x C e C e  b 1  a y  y   a x b  C e  a Integralne krivulje: slika 3.1 (Maxima) 'diff (y,x) = a*y+b; a: 2; b: 1; load("plotdf"); plotdf([a·y+b],[x,y],[y, −3, 3],[x, 0, 8]); ode2(%i1,y,x); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.2 Primer 2 Diferencialna enačba: y   sin y (En. 1, D3.2) Poiščimo implicitno rešitev diferencialne enačbe: dy dy   y  sin y   dx dx sin y  dy  csc y dy   dx C    x C (En. 2, D3.2)   sin y csc y pomnožimo in delimo s csc y  cot y   csc csc y  csc y  cot y  y dy     dy  (En. 3, D3.2) csc y  cot y  Označimo: csc y  cot y  u (En. 4, D3.2) Izračunajmo odvod od u , torej u u   csc y  cot y csc   1   cos y 1 cos y   y      csc y    cot y 2  sin y  sin y sin y sin y cot  cos y   sin y  sin y  cos y  cos y 1    y      2 2  sin y  sin y sin y  2 u    csc y  cot y  csc y (En. 5, D3.2) Vstavimo v enačbo 3, D3.2 označbi za u (enačba 4, D3.2) in u  (enačba 5, D3.2).  csc csc y   csc y  cot y  y    dy  csc y  cot y   ln u y cot x C    du  ln csc   y    (En. 6, D3.2) u KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB  x csc y  cot y   C e Implicitna rešitev diferencialne enačbe je: 1 cos  y x   C e (En. 7, D3.2) sin y Diferencialna enačba: slika 3.2 (Maxima) load(draw); set_draw_defaults( xrange=[−2·%pi, 2·%pi], yrange=[−1.5,1.5], grid=true, line_width=2, color=black); draw2d (yrange=[-1.5,1.5],color=red, explicit(sin(x), x,-2*%pi,2*%pi))$ Integralne krivulje: slika 3.3 (Maxima) 'diff(y,x) = sin(y); load(draw); plotdf( [sin(y)], [x, −2·%pi, 2·%pi], [y, −2·%pi, 2·%pi])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.3 Primer 3 Diferencialna enačba: dy  y      r y  1  dt  k   dy  y  r     2 2 r r y  1      r y y     r y a y a   pri tem je dt  k  k k   r dt dy dy   r dt      C    r t C 1 y  y    1 y    1  y    1   k   k  .-.-.-.-.-.  dy k dy  k   y y dy y      dy       y  y k    y  y k    y  y  k  y  y    1   k   dy y  ln   y  ln  k  y  y  y  k  y    ln 1  k   .-.-.-.-.-.- ln y y r t r t     r t C   C e 1   y    C e    k y   k  y k  y y     r t y 1 C e    r t C e  r t  k  k    C e  k y     r t   1   1 C e  r t e  1 C Označimo . Eksplicitna rešitev 1  A diferencialne enačbe je: C y  k 1     r t A e Vrnimo se k diferencialni enačbi in izračunajmo še drugi odvod: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB dy  y  2 r     2 d y r r y 1          r y y r 2 y  2 dt  k  k dt k Diferencialna enačba : slika 3.4 (Maxima) r=0.5; k=10; load(draw); set_draw_defaults (xrange=[−5, 15], yrange=[−0.5, 2], grid=true, line_width=4, color=black); draw2d( color=blue, explicit(0.5·y·(1−y/10), y,−5,15)); Diferencialna enačba in drugi odvod: slika 3.5 (Maxima) r=0.5; k=10; load (draw); set_draw_defaults( xrange=[−2,12], yrange=[−1,2], grid=true, line_width=4, color =black); draw2d(color=blue, explicit(r·y·(1−y/k), y, −2, 12), color=red, explicit(r−2·r·y/k, y,−1,12)); Integralne krivulje: slika 3.7 (Maxima) 'diff(y,t)=r*y*(1-y/k); r=0.5; KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB k=10; load("plotdf"); plotdf([0.5·y·(1−y/10)],[t,y], [y,0,15], [t,0,15]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.4 Primer 4 Diferencialna enačba. y    y  12   2 y  4 Diferencialna enačba: slika 3.7 (Maxima) 'diff(y,x)=(y^2-4)*(y+1)^2; load(draw); set_draw_defaults( xrange=[−3,3], yrange=[−15,8], grid=true, line_width=3, color=black); multi_plot_mode (screen); draw2d (color=blue, yrange=[−15,8], explicit ((x^2−4)·(x+1)^2, x,−3,3)); Integralne krivulje: slika 3.8 (Maxima) 'diff (y, x) = (y^2-4)*(y+1)^2; load("plotdf"); plotdf ([(y^2−4)·(y+1)^2], [x,0,5], [y,−3,3]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.5 Primer 5 Diferencialna enačba: y  2 y  y      r y  1      2 k 1 y Graf enačbe f y  2 y  y     r y    1      2 (slika 3.9) k 1 y Diferencialna enačba: slika 3.9 (Maxima) Opomba: oznaka je spremenjena: y  x k: 10; r: 0.48; set_draw_defaults( xrange=[0, 10], yrange=[−5, 5], grid=true, line_width=3, color=black); draw2d(color=blue, yrange=[−0.5,0.6], explicit(r·x·(1−x/k)−x^2/(1+x^2), x,0,10), color=black, line_width=2, polygon([[0,0],[5,0],[10,0]])); Graf funkcij f y  y  y 1    f    2 in   y r 1 za k  10 in različne 1  2   y  k  vrednosti r  0.38, 0.45, 0.56, 0.7 Funkciji f y y 1 in f 2 : slika 3.10 (Maxima)     Opomba: oznaka je spremenjena y  x k: 10; r: 0.48; KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB load(draw); set_draw_defaults( xrange=[0,10], yrange=[−5,5], grid=true, line_width=4, color=black); draw2d( color=blue, yrange=[0,0.6], explicit (x/(1+x^2), x,0.1,10), color=red, explicit (0.38·(1−x/k), x,0.1,10), explicit (0.45·(1−x/k), x,0.1,10), explicit (0.56·(1−x/k), x,0.1,10), explicit (0.7·(1−x/k), x,0.1,10)); Integralne krivulje: slika 3.11 (Maxima) Opomba: vrednosti parametrov k  10 in r  0.48 'diff(y,x)=r·y·(1−y/k)−y^2/(1+y^2); k: 10; r: 0.48; load("plotdf"); plotdf ([r·y·(1−y/k)−y^2/(1+y^2)], [x,y], [x,0,100],[y,0,10]); Bifurkacijski diagram: slika 3.12 (Maxima) Izhodišče: 0   2 y y  y  y     r y 1    r   1  2    ; k  10  k 1    2 y  k 1  y Opomba: sprememba oznake: r  x load(draw); set_draw_defaults( xrange=[0,0.8], KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB yrange=[0,10], grid=true, line_width=4,color=black); draw2d (grid=true, color=blue,line_type=solid, implicit(x·(1−0.1·y)=y/(1+y^2), x,0,0.8, y,0,10)); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.6 Primer 6 Diferencialna enačba: dy  f y      y 1 y    h dx Točke ravnotežja y 2 1  1 4   h ,    1 y    h 0  y    y y  y h 0  r 1 r 2 . 2 Točka ravnotežja je ena sama: r3 , ko je diskriminanta D     y  1 1 4 h 0 , to 2 h  1 pa je ko je . 4 Funkcija f y    y  h   y  1  , pri h  2 je prikazana na sliki 3.13 (Maxima). Opomba: sprememba oznake: y  x h: −2; load(draw); set_draw_defaults( xrange=[−2,3], yrange=[−0.5,2.5], grid=true, line_width=3, color=black); draw2d (color=blue, explicit(x·(1−x)−h, x,−2,3), color=black, line_width= 2, polygon( [[−2,0],[0,0],[3,0]]), polygon( [[−1,−0.5], [−1, 0], [−1, 2.5]]), polygon( [ [2, −0.5], [2, 0], [2, 2.5]]), KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB polygon( [ [0, −0.5], [0, 0], [0, 2.5]])); Točke ravnotežja so funkcije parametra y  h 1  1 4   h : r 2 Bifurkacijski diagram: slika 3.14 (Maxima) Opomba: sprememba oznake: y  x load(draw); set_draw_defaults( xrange=[−6,2], yrange=[−4,4], grid=true, line_width=4, color=black); draw2d (color=blue, explicit(0.5·(1+sqrt(1−4·x)),x,−6,2), explicit(0.5·(1−sqrt(1−4·x)),x,−6,2))$ Integralne krivulje: slike: 3.15, 3.16 in 3.17 (Maxima) Opomba: oznaka spremenljivke je spremenjena: y  x ; spreminja se vrednost parametra h  0.125, 0.25, 0.8 . 'diff(x,t)=x·(1−x)−h; h: 0.125; load ("plotdf"); plotdf ([x·(1−x)−h],[t, x],[x,−0.5,1.5], [t,0,15]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 3.7 Primer 7 Diferencialna enačba: dy    x y  y   2 dx Poiščimo rešitev diferencialne enačbe: dy dy    x y y   2     x dx  y dx   y  2 1 A B    A   y      y  2 B y 1   y  2 y y   2  A B       y 2 A 1  Postavimo: A B   1 1 B  0 in    A   2 A 1  in 2 2  x dx y 1 1 dy 1 dy  dy          y  2  2 y 2  y  2    1 dy 1     ln y ; 2 y 2 1 dy 1     ln  y  2  2 y  2 2   1 dy 1 dy 1  y  2  1  2      y  2  2 ln    x  C ln  x    C  2 y 2 y  2 2  y  2  y  y  2 2 2 x   2 x C e x y    1 C e  2  y     2 y C e    y Eksplicitna splošna rešitev je: 2 x   y  2 1 C e Pri začetnih pogojih: x y 0 in 0 določimo konstanto C KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y 2 y  2 2 x   0 0 0   2  y    y C e  C 2 2 x 0 0  x 1 C e 0 y e  0 0 Eksplicitna rešitev pri danih začetnih pogojih je: y  2 2 2   y e  x 0 1   x 0   x 0 x y e  y 2 e 0 0 2 e  y  y 0  2 2 2 2   y e  x0 0 y 2  y0      2 2 x  x y  y  0 2 e 0 0 d y 2 Poiščimo drugi odvod 2 dx dy 2 d y f f dy  f x y  , 2 )    x y    2 x y     dx 2 dx x y dx f f 2  y   2 y      2 x y 2 x , ,  x  y f dy 2          2 x y 2 x   x y    2 x y    y dx f dy 2          2 2   x y x x y    2 x y   y dx d y 2   y y         2 x y y 2       2 2 x y 2 x  dx d y 2              y 2 2 y 2 1 2 x y 2 x dx 2 Razdelitev ravnine X-Y po kriteriju  y x    0 , y x     0 , slika 3.18 (Maxima) lord(draw); set_draw_defaults ( xrange=[−2,2], KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB yrange=[−2,4], grid=true, line_width = 4, color=black); draw2d (color=black,line_width=3, line_width=3, polygon([[−2,2],[0,2],[2,2]]), polygon([[−2,0],[0,0],[2,0]]), polygon([[0,−2],[0,0],[0,4]]))$ Razdelitev ravnine X-Y po kriteriju y x  ( )  0 , y x  ( )  0 : slika 3.19 (Maxima) load(draw); set_draw_defaults( xrange=[−3,3], yrange=[−3,3], grid=true, line_width=4, color=black); draw2d (color=black, line_width=3, line_width=3, polygon([[−3,2],[0,2],[3,2]]), polygon([[−3,0],[0,0],[3,0]]), polygon([[0,−3],[0,0],[0,3]]),color=blue, polygon([[−3,1],[0,1],[3,1]]), explicit(1−1/(2·x^2),x,−3,3))$ draw2d (color = black, line_width = 1, polygon([ [0,-3], [0, 0], [0,3]]), line_width = 3, polygon([ [-3, 2], [0, 2], [3,2]]), polygon([ [-3, 0], [0, 0], [3,0]]), color = blue, explicit( 1-1/(2*x^2), x, -3, 3))$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Smerno polje: slika 3.20 (Maxima) load(draw); load(drawdf); set_draw_defaults( xrange=[−3,3], yrange=[−3,3], grid=true, line_width=1, color =black); drawdf( [x·y·(y−2)],[x,−3,3],[y,−3, 3], field_grid=[20,20], line_width=3,color=black, polygon([[−3,2],[0,2],[3,2]]), polygon([[−3,0],[0,0],[3,0]]), polygon([[0,−3],[0,0],[0,3]]),color=blue, polygon([[−3,1],[0,1],[3,1]]), explicit(1−1/(2·x^2),x,−3,3))$ Integralne krivulje: slika 3.21 (Maxima) 'diff(y,x)= x·y·(y−2); load("plotdf"); plotdf([ x·y·(y−2)],[x,y],[y,−3,3],[x,−2,2]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 4. KVALITATIVNA ANALIZA SISTEMA LINEARNIH DIFERENCIALNIH ENAČB PRVEGA REDA Kvalitativno analizo sistema diferencialnih enačb pričnemo s homogenim sistemom linearnih diferencialnih enačb prvega reda s konstantnimi koeficienti. Primer takega sistema prikazujeta enačbi 1 in 2. To je povezan sistem (angl. y y y Coupled ), kjer v odvodu funkcije 1 nastopata funkciji 1 in 2 in v odvodu funkcije y y y 2 prav tako nastopata funkciji 1 in 2 . V tem poglavju bomo obravnavali povezan sistem samo dveh funkcij. dy 1  y     a y b y 1 1 2 (En. 1) dt dy2  y      c y d y 2 1 dt 2 (En. 2) Enačbi 1 in 2 zapišimo v matrični obliki - enačbe 3. y   A y (En. 3) Pri tem je y   y 1   y  1a b y     ,    A in    y    2   y  c d 2  a , b , c in d so poljubne konstante. y 1 kot tudi y2 sta funkciji neodvisne spremenljivke t . Grafično funkciji y 1 in y2 predstavljata integralni krivulji v ravnini (Y-T). Velikokrat pa nam več pove medsebojna odvisnost obeh funkcij, ki jo grafično predstavimo v tako imenovani fazni ravnini (angl. Phase plane) (Y1-Y2) kot parametrično krivuljo , ki jo včasih imenujemo tudi trajektorija. Niz parametričnih krivulj imenujemo fazni portret (angl. Phase portrait). Definicija Fazni portret je reprezentativen niz rešitev, ki so grafično predstavljene kot niz parametričnih krivulj v kartezični ravnini (Y1-Y2), imenovani fazna ravnina. Podobno kot smerno polje je tudi fazni portret primerno orodje za predvidevanje obnašanja rešitev sistema diferencialnih enačb. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Splošno rešitev sistema diferencialnih enačb 1 in 2 podaja enačba 4 (dodatek 4.1): y t    t   t C e  1  C e  2 1 2 (En. 4) Pri tem sta 1 in 2 lastni vrednosti in in lastna vektorja matrike A     1 2 (dodatek 4.1). Iz enačb 1 in 2 lahko v vsaki točki fazne ravnine P(Y1-Y2) izračunamo smer dy2 tangente na parametrično krivuljo, ki gre skozi to točko (enačba 5). dy 1 dy c y    d y 2  12 (En. 5) dy a y    b y 1 1 2 Za risanje faznih portretov smo uporabili programski paket Maxima. Kritična točka  P y  y 0 je točka, kjer sta 1 in 2 enaka 0 . V tej točki je smerni koeficient tangente na parametrično krivuljo nedoločen: dy 0 2  (Kreyszig, dy 0 1 2011, str. 142-149). Iz enačbe 5 hitro razberemo, da je to točka P 0 0, 0 .   V   y     A y dodatku 4.1 smo, izhajajoč iz enačbe 3  , določili lastni vrednosti     det A     I 0 1 in 2 :   . Pri tem sta 1 in 2 lastni vrednosti in in lastna vektorja matrike A     . Za 1 2 vsako lastno vrednost i izračunamo lastni vektor :     i    i A     i I   =0 Pri tem ni enak 0 .    i Obstajajo različni tipi kritičnih točk, ki jih imenujemo: 1. Pravilno vozlišče (angl. Proper node ali Star point), 2. Nepravilno vozlišče (angl. Improper node), 3. Degenerirano vozlišče (angl. Degenerate point), 4. Sedlo (angl. Saddle point), KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 5. Center (angl. Center), 6. Spiralna točka (angl. Spiral point); njim pripadajoče fazne portrete bomo predstavili v naslednjih primerih. 1. Pravilno vozlišče Pravilno vozlišče je asimptotično stabilno, če sta lastni vrednosti enaki in manjši od 0 -    0 1 2 , in je nestabilno, če sta lastni vrednosti enaki in večji od 0 :     0 1 2 . Pogoj, da je vozlišče pravilno, je, da je matrika A večkratnik identitetne matrike, torej: A   0    , pri tem je  katera koli od 0 različna konstanta (primer 1a in 0    primer 1b). To pa pomeni, da diferencialni enačbi nista med seboj povezani . Lastna vektorja sta:    1   0        1 2 in (En. 6)   0   1 Primer 1a Nestabilno pravilno vozlišče. Sistem diferencialnih enačb: y   1 0     y y   y y  y oziroma 1 1 , 22 (En. 7) 0 1  Obe lastni vrednosti sta   1 realni, enaki in večji od 0 :   1 2 (dodatek 4.2). Lastna vektorja sta:    1   0        1 2 in   0   1 Splošno rešitev prikazuje enačba 8: y    1   0 t  t C     e C   e 1 2   (En. 8)   0   1 S faznega portreta na sliki 4.1 vidimo, da je kritična točka nestabilna. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 4.1: Fazni portret - Nestabilno pravilno vozlišče Primer 1b Stabilno pravilno vozlišče Sistem diferencialnih enačb: y         1 0  0  y oziroma y y y   y 2 2 , (En. 9)   1 1  1 Obe lastni vrednosti sta   realni , enaki in negativni    0 1 2 (dodatek 4.3). Lastna vektorja sta:    1   0        1 2 in .   0   1 Splošno rešitev prikazuje enačba 10: y   1   0  t t  C     e  C     e 1 2 (En. 10)   0   1 S faznega portreta na sliki 4.2 vidimo, da je kritična točka asimptotično stabilna. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 4.2: Fazni portret - Stabilno pravilno vozlišče 2. Nepravilno vozlišče Nepravilno vozlišče je asimptotično stabilno, če sta lastni vrednosti realni, različni in manjši od 0 :     0 1 2 (primer 2a). Včasih uporabimo izraz ponor (angl. Sink). Nepravilno vozlišče je nestabilno, če sta lastni vrednosti realni, ralični in večji od 0 :     0 1 2 (primer 2b). Včasih uporabimo izraz izvor (angl. Source). Primer 2a Asimptotično stabilno nepravilno vozlišče Sistem diferencialnih enačb: y      3 1   1  y (En. 11)  3  Lastni vrednosti sta različni, realni in obe manjši od 0 .    2,     4 0 1 2 (dodatek 4.4) Lastna vektorja sta:    1   1        1 2 in .    1   1 Splošno rešitev prikazuje enačba 12: y    1   1     2 t     4t C e C  e 1   2   (En. 12)   1    1 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB S faznega portreta na sliki 4.3 vidimo, da je kritična točka asimptotično stabilna. Slika 4.3: Fazni portret - Asimptotično stabilno nepravilno vozlišče Primer 2b Nestabilno nepravilno vozlišče Sistem diferencialnih enačb: y   4 2     1 3  y (En. 13) Lastni vrednosti sta različni, realni in večji od 0 :   2,    5 0 1 2 (dodatek 4.5). Lastna vektorja sta:    1   2        1 2 in .    1   1 Splošno rešitev prikazuje enačba 14: y    1   2 2  t 5t C     e  C     e 1  2 (En. 14)   1   1 S faznega portreta na sliki 4.4 (dodatek 4.5) je razvidno, da je kritična točka nestabilna. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 4.4: Fazni portret - Nestabilno pravilno vozlišče 3. Degenerirano vozlišče Degenerirano vozlišče je nestabilno, če sta lastni vrednosti realni, enaki (ponavljajoči se) in pozitivni: 0     1 2 (primer 3a). Degenerirano vozlišče je stabilno, če sta lastni vrednosti realni, enaki in manjši od 0: 0     1 2 (primer 3b). Primer 3a Nestabilno degenerirano vozlišče Sistem diferencialnih enačb: y   7 1     y (En. 15) 4 3  Lastni vrednosti sta    5 0 realni , enaki in pozitivni   1 2 (dodatek 4.6). Splošno rešitev, ko sta lastni vrednosti enaki, prikazuje enačba 16 (dodatka 4.1 in 4.6): y  t t t 5     C e  C   t e   e   1 2   (En. 16) Lastni vektor je:    1   0           in vektor 2   1 Rešitev prikazuje enačba 17: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y    1    1   0  5    t  5   t 5  t C e  C  e     e 1 2      (En. 17)    2     2   1  S faznega portreta na sliki 4.5 vidimo, da je kritična točka nestabilna. Slika 4.5: Fazni portret - Nestabilno degenerirano vozlišče Primer 3b Stabilno degenerirano vozlišče Sistem diferencialnih enačb: y  3   1   2      y (En. 18)  1    2    6  Lastni vrednosti sta  3        realni , enaki in manjši od  0 0 1 2   2  (dodatek 4.7). Lastni vektor je:         3    6 in vektor (dodatek 4.7)      1   1 Rešitev prikazuje enačba 19: y    3 3 3 3   t     3   t    6  t  2 2 2  C   e  C    t    e   e 1 2      (En. 19)   1    1   0  KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB S faznega portreta na sliki 4.6 vidimo, da je kritična točka asimptotično stabilna. Slika 4.6: Fazni portret - Stabilno degenerirano vozlišče 4. Sedlo Sedlo je nestabilna kritična točka, če sta lastni vrednosti realni in imata nasproten predznak:     0  12  (primer 4a in primer 4b). Primer 4a (Kreyszig, 2011, str. 143). Sedlo Sistem diferencialnih enačb (dodatek 4.8): y    1 0    0  y (En. 20)  1  Lastni vrednosti sta realni, po absolutni vrednosti enaki in imata nasprotna predznaka 1,       1 1 2 . Lastna vektorja sta:    0   1        1 2 in   1   0 Rešitev prikazuje enačba 21: y    1   0 t   t  C     e C e 1 2   (En. 21)   0   1 S faznega portreta na sliki 4.7 vidimo, da je kritična točka nestabilna. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 4.7: Fazni portret – Sedlo Primer 4b (Dawkins, 2007, str. 280-281) Sedlo Sistem diferencialnih enačb (dodatek 4.9): y   1 2     3 2  y (En. 22) Lastni vrednosti sta realni, različni, ena je manjša in druga večja od 0   1,   4 1 2 Lastna vektorja sta:     1   2  1    in 2      1   3 Rešitev podaja enačba 23: y     1   2  t 4t C   e  C   e 1   2   (En. 23)   0   3 S faznega portreta na sliki 4.8 vidimo, da je kritična točka nestabilna. Slika 4.8: Fazni portret – Sedlo KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 5. Center Center je stabilna kritična točka, toda ne asimptotično stabilna, včasih se uporabi tudi izraz: nevtralno stabilna (angl. Neutrally stable). Lastni vrednosti sta imaginarni: i         1 in i 1 . Primer 5a Center Sistem diferencialnih enačb (dodatek 4.10): y  0 1      y . (En. 24) 4 0  Lastni vrednosti sta imaginarni števili:    i 2     1 in i 2 2 . Lastna vektorja sta:   1   1        1 2 in  i     2  i 2 Splošno rešitev podaja enačba 25: y  cos 2    t   sin 2    t   C     C  1     2 (En. 25)   2 sin 2  t 2 cos 2        t      S faznega portreta na sliki 4.9 vidimo, da je kritična točka center, ki je stabilna. Trajektorije obkrožajo kritično točko (dodatek 4.10). Slika 4.9: Fazni portret – Center Primer 5b Center Sistem diferencialnih enačb (dodatek 4.11): KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y     3  9    4  y (En. 26) 4  Lastni vrednosti sta imaginarni:     i 3 3      i 1 in 3 3 2 Lastna vektorja sta:         1 in   3   3   1  3  2 i  1  3  i         Splošno rešitev podaja enačba 27: y   3 cos 3    3  t     C   1   cos 3   3   t  3 sin 3    3  t          3 sin 3   3  (En. 27) t   C  2       sin 3  3   t 3 cos 3   3  t     S faznega portreta na sliki 4.10 vidimo, da je kritična točka center, ki je stabilna. Trajektorije obkrožajo kritično točko (dodatek 4.11). Slika 4.10: Fazni portret - Center 6. Spiralna točka Spiralna točka je nestabilna kritična točka, če je realni del kompleksne in konjugirano kompleksne lastne vrednosti   ,     i  1 2 pozitiven:   0 (primer 6a). Včasih se uporabi izraz spiralni izvor (angl. Spiral source) in je stabilna, če je realni del kompleksne in konjugirano kompleksne lastne vrednosti KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB negativen 0   (primer 6b). Včasih se uporabi izraz spiralni ponor (angl. Spiral sink). Primer 6a Spiralna točka Sistem diferencialnih enačb (dodatek 4.12): y   2 3      y (En. 28) 3 2  Lastni vrednosti sta kompleksno in konjugirano kompleksno število, realni del je pozitiven: 2     i 3     2 i 3 1 in 2 . Lastna vektorja sta:     i   i  1    in 2      1   1 Splošno rešitev podaja enačba 29: y  sin 3    t   cos 3    t 2  t  C  e C 1        (En. 29) cos 3    2 t sin 2           t  S faznega portreta na sliki 4.11 (dodatek 4.12) je razvidno, da je kritična točka - spiralni izvor -nestabilna. Slika 4.11: Fazni portret – Spiralni izvor Primer 6b Spiralna točka Sistem diferencialnih enačb (dodatek 4.13): KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y  1   1     2  y   (En. 30) 1   1  Lastni vrednosti sta kompleksno in konjugirano kompleksno število, realni del je negativen.  1 7 1 7     i      1 i in 2 4 4 4 4 Lastni vektor za lastno vrednost 1 je        1 3 7    0 1    i        4    4  Splošno rešitev podaja enačba 31: y   1    1  0   t      7 7       C e  4   cos   t  sin  1 3  7  t              4   4         4  4   C e  1    1  0   (En. 31) t      7 7 4         sin   t  cos  2 3   7  t               4   4      4  4   S faznega portreta na sliki 4.12 (dodatek 4.13) je razvidno, da je kritična točka - spiralni ponor - asimptotično stabilna. Slika 4.12: Fazni portret – Spiralni ponor KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.1 Homogen sistem linearnih diferencialnih enačb s konstantnimi koeficienti Homogen sistem diferencialnih enačb s konstantnimi koeficienti zapišimo v matrični obliki: y    A y (En. 1, D4.1), pri tem je y   dy 1 1 A je matrika dimenzije n n  , katere elementi niso odvisni dt od spremenljivke y a y t . Vzemimo, da je n  1 , in pridemo do enačbe:    , katere rešitev je:   a t  y C e . Po analogiji lahko predpostavimo, da je rešitev sistema y    A y y  t     e (En. 2, D4.1) To je za sedaj predpostavka, pri tem moramo in šele določiti.   y t    t tt      e  y            e  e A e       A (En. 3, D4.1) Tako pridemo do enačbe 4, D4.1:  A I        0 (En. 4, D4.1) Predpostavimo, da 0   , in pridemo do enačbe 5, D4.1:  A I      0 (En. 5, D4.1) Da bo imela A I 0 enačba 5, D4.1 neničelno rešitev, mora biti matrika       singularna. Matrika je singularna, če je njena determinanta enaka nič, sicer je nesingularna. Tako pridemo do enačbe 6, D4.1: det A I       0 (En. 6, D4.1) Od tod izračunamo lastne vrednosti  i matrike A ( n n  ). V tem poglavju obravnavamo samo matrike 2 2  . Za vsako lastno vrednost potem izračunamo lastni vektor  ( ) i KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB        ( ) i A I  0 i (En. 7, D4.1) Lastne vrednosti so lahko vse med seboj različne, nekatere so med seboj različne ali pa so vse enake. Predpostavimo, da ima matrika A splošno n realnih in med seboj različnih lastnih vrednosti, 2 , …   n , katerim pripada n linearno neodvisnih lastnih vektorjev (1) (2) n ( ) :  , , ….   V tem primeru je splošna rešitev sistema diferencialnih enačb: y   1 1  t 2  t n  t 1   C   e  C    e 2  .........  C    e n 1 2 n (En. 8, D4.1) V primeru enakih lastnih vrednosti (angl. Repeated Eigenvalues) enačba 8, D4.1 ni splošna rešitev. V nadaljevanju povzemamo po Dawkinsu (2018, str. 310). Naj bo A matrika 2 2  in obe lastni vrednosti sta enaki. Matrika A ima torej dvojno lastno vrednost 1      2 (glej tudi primera 4a in 4b) in en sam lastni vektor . Torej dobimo tudi eno samo rešitev sistema diferencialnih enačb:  y   t   e Želimo imeti dve linearno neodvisni rešitvi, da lahko oblikujemo splošno rešitev. Dobiti moramo še drugo rešitev. K prvotni rešitvi dodajmo še spremenljivko t , tako bi bila druga rešitev: y     t t e . Poglejmo, ali smo sedaj res dobili splošno rešitev. Izhajamo iz enačbe 2, D4.1: y          t A y  y A t e Odvajajmo drugo rešitev: y  t      t t       t e      y e     t   e  Odvod vstavimo v enačbo 3, D4.1 in dobimo enačbo 9, D4.1:   t t    t e     t   e     A t  e (En. 9, D4.1) KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Poglejmo, ali je druga rešitev  t y    t  e res prava rešitev. Na levi strani enačbe t t      9, D4.1 imamo dve funkciji:       e in t e , ki morata zadostiti enačbi 3, D4.1 oziroma 4, D4.1.      A A I       0 ( 3, D4.1 )    (4, D4.1), pri tem pa   0 Izenačimo koeficiente na levi in desni strani enačbe in dobimo: 1. A            A I    0 2. 0    Pri enačbi 1. je  lastna vrednost in  lastni vektor. Pri enačbi 2. je   lastni vektor, ki ne more biti 0 . Če zaključimo,     t y t e ne more biti rešitev diferencialne enačbe. Predpostavimo, da je druga rešitev: y  t     t  e   (En. 10, D4.1), pri tem je še neznan vektor, ki ga moramo šele določiti.  Odvajajmo enačbo 10, D4.1 in dobimo: y            t      t   t e t e e y in y vstavimo v enačbo 1, D4.1: y   A y y                   t     t     t   e t e e    t t A t e e         t   t t t e     t e     A  t e    A  e (En. 11, D4.1) Ponovno izenačimo koeficiente leve in desne strani enačbe. Koeficienta pri t e t       A A      I  0    (En. 12, D4.1) KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Koeficient pri et        A    A      I    (En. 13, D4.1) Iz enačbe 12, D4.1, ki je enaka enačbi 4, D4.1, izračunamo vektor   A      I   0 Iz enačbe13, D4.1 izračunamo vektor   A      I    (En. 14, D4.1) Splošna rešitev sistema diferencialnih enačb je:   t  t t           C e C t e e   y 1 2 (En. 15, D4.1) Poglejmo še primer, ko sta lastni vrednosti kompleksni števili. Če je lastna vrednost sistema diferencialnih enačb  kompleksno število, potem je tudi konjugirano kompleksno število lastna vrednost. Torej, če je lastna  vrednost i     , potem je tudi    i  lastna vrednost. V primeru, da je lastna vrednost i     1 imaginarno število, potem je druga lastna vrednost      i 2 (glej primera 6a in 6b). -.-.-.- Eulerjeve formule ei  cos   i sin  e   i  cos    i sin  e e e     i  t t i  t t   e  cos       t isin     t e e e e     i  t t    i t t    cos      t  i sin    t  -.-.-.- Vzemimo, da je i           a i b 1 in pripadajoči lastni vektor 1 . Ena od rešitev je: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y   a t 1 t  y  e     e        1 1 1 1  cos   t i sin     t    a t  y e     a i b cos          t i sin     t 1  (En. 16, D4.1):     y a t   a t e                   t a cos         1 t b sin t i e a sin t b cos Enačbo lahko zapišemo tudi: y  u t ( )   i v t ( ) u t 1 . Pri tem sta   in v t ( ) dve linearno neodvisni rešitvi . Splošno rešitev lahko torej zapišemo: y  C u t  ( )  C v t  ( ) 1 2 (En. 17, D4.1) (En. 18, D4.1): y  C   a cos 1 1        t b sin     t  C   2 a sin        t b cos     t KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.2 Primer 1a Nestabilno pravilno vozlišče Sistem dveh diferencialnih enačb: y   1 0 A     0 1  y Korak 1 Izračun lastnih vrednosti 1   0 2  0   1     0  karakteristični polinom 0 1      2 2  4 4  1  2 0        2  , 1 0    1 2  2     1 1 2 Obe lastni vrednosti sta enaki in pozitivni 0      1 2 . Korak 2 Poiščimo lastna vektorja. Diferencialni enačbi nista povezani. y   t y y   C e 1 1  1 1 y  t y y  C e  2 2  2 2 Iz sistema diferencialnih enačb je razvidno, da sta lastna vektorja:    1   1   2   0          in 0   1 Korak 3 Splošna rešitev sistema je: y    1   0 t t C      e C   e 1 2     0   1 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Izračun lastnih vrednosti, lastnih vektorjev in fazni portret: slika 4.1 (Maxima) en1: 'diff (y1, t) = y1; en2: 'diff (y2, t) = y2; load ("plotdf"); A: matrix ([1,0],[0,1]); expand (charpoly (A, λ)); solve (λ^2−2 λ); v: eigenvectors (A); plotdf([y1,y2], [y1,y2]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.3 Primer 1b Stabilno pravilno vozlišče Sistem diferencialnih enačb: y      1 0   0  y  1  Korak 1 Izračun lastnih vrednosti   1  0 2  0     1    0  karakteristični polinom 0   1     2   2 4 4  1  2 0          2  1 0  1, 2  2     1 2 Obe lastni vrednosti sta enaki in negativni 0      2 . Korak 2 Poiščimo lastna vektorja Diferencialni enačbi nista povezani: y t    y  y   C e 1 1 1 1 y    t y  y  C e  2 2 2 2 Iz sistema diferencialnih enačb je razvidno, da sta lastna vektorja:    1   1   2   0          in 0   1 Korak 3 Splošna rešitev sistema je: y   1   0   t t  C    e  C   e 1 2     0   1 Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.2 (Maxima) KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB en1: 'diff (y1, t)=−y1; en2: 'diff (y2, t)=−y2; load ("plotdf"); A: matrix ( [−1, 0], [0, −1] ); expand (charpoly (A, λ)); solve (λ^2+2·λ+1); v: eigenvectors (A); plotdf ([−y1,−y2], [y1, y2]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.4 Primer 2a Asimptotično stabilno nepravilno vozlišče Sistem diferencialnih enačb: y      3 1   1  y  3  Korak 1 Izračun lastnih vrednosti   3  1 2  0     3     1 0  karakteristični polinom 1   3       2   6 36 4 8   3   2 1 0         6  ,  8 0  1 2 2    2    1 in 4 2 Obe lastni vrednosti sta realni, različni in manjši od 0 ,   ,  0 1 2 . Korak 2 Poiščimo lastna vektorja. Za 2    1    3 2 1     0 1           vzemimo prvo vrstico   1   3 2       0 2      0 1     1 2   1  1 2     . Izberimo  1      1   1 1 in 1    1    dobimo: 1 .   1 Za 4    2 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB    3 4 1       0  1        vzemimo prvo vrstico  1   3 4        0 2       1      0 1 2      1 2 2 1         . Izberimo   1 in    2 1    1  1     1    dobimo 2 .    1    1   1        Lastna vektorja sta torej: 1 in 2    1   1 Korak 3 Splošna rešitev sistema je: y   1 1   2 t     4t  C   e  C   1   2   e 1      1 Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.3 (Maxima) en1: 'diff(y1,t)=-3*y1+1*y2; en2: 'diff(y2,t)=1*y1-3*y2; load("plotdf"); A: matrix([-3,1],[1,-3]); charpoly (A, λ), expand solve (λ^2+6·λ+8); plotdf([-3*y1+1*y2,1*y1-3*y2],[y1,y2]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.5 Primer 2b Nestabilno nepravilno vozlišče Sistem dveh diferencialnih enačb: y  4 2     1 3  y Korak 1 Izračun lastnih vrednosti 4   2  0  4       3     1 3  2 0  karakteristični polinom     7  49 4 10   2    ,  7  10  0  1 2  2   2   5 1 in 2 Obe lastni vrednosti sta realni in večji od 0    ,  0 1 2 . Korak 2 Poiščimo lastna vektorja. Za 2   1 4 2  2       0  1        vzemimo prvo vrstico  1 3 2       0 2   2     2  0 1      1   1 2 2 1     . Izberimo     1    1  1    1    in dobimo 1 .    1 Za 5   2 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 4 5  2       0  1        vzemimo drugo vrstico   1 3 5        0 2    2  0 1 2     2 2 1 2     .   2 2   2  Izberimo 1   2 2 in dobimo       2   . 1    1   2        Lastna vektorja sta 1 in . 2    1   1 Korak 3 Splošna rešitev sistema je: y    1 2 2   t   5  t C  1   e C  2   e    1   1 Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.4 (Maxima) en2: 'diff(y2,t)=1·y1+3·y2; load ("plotdf"); A: matrix ([4,2],[1,3]); charpoly (A, λ), expand; en1: 'diff(y1,t)=4·y1+2·y2; solve (λ^2−7·λ+10); v: eigenvectors (A); plotdf ( [4·y1 + 2·y2, 1·y1 + 3·y2], [y1, y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.6 Primer 3a Nestabilno degenerirano vozlišče Sistem dveh diferencialnih enačb y    7 1      y 4 3  Korak 1 Izračun lastnih vrednosti 7   1  0   7       3      4 0  karakteristični polinom 4 3    2 10  100 4 25         10  , 25  0  1 2  2      5 1 2 Obe lastni vrednosti sta realni, enaki in večji od 0      0 1 2 . Korak 2 Poiščimo lastni vektor. Za 5   7 5  1     1   0 2     0   1 2                vzemimo prvo vrstico   4 3 5            2   0  4 2 0 1 2     2      0      2   1 2 2 1     1 . Izberimo     1   2  1 1      1    in dobimo .    2 Korak 3 Splošna rešitev sistema: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y   1 5t  C     e 1 .    2 Ker imamo dvojno lastno vrednost, smo dobili eno samo rešitev, potrebujemo pa dve. V takem primeru je splošna rešitev: y   t   tt    C e   t e      e   1 2 C  dodatek 4.1 - enačba 14, D4.1 Korak 4 Določimo vektor  . 7 5  1     1 1                 2 1 , postavimo   0  1 2 1 in 4 3 5          2 2    0    dobimo vektor .   1 Korak 5 Splošna rešitev sistema diferencialnih enačb je: y    1    1 0  5  t 5    t 5         t C   e C t 1 2   e     e   2       2   1  Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.5 (Maxima) en1: 'diff(y1,t)=7·y1+1·y2; en2: 'diff(y2, t)=−4·y1+3·y2; load ("plotdf"); A: matrix ( [7, 1], [−4, 3] ); charpoly(A, λ), expand; solve (λ^2−10·λ+25); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB v: eigenvectors (A); plotdf ( [7·y1 + 1·y2, −4·y1 + 3·y2], [y1, y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.7 Primer 3b Stabilno degenerirano vozlišče. Sistem dveh diferencialnih enačb: y  3   1   2      y  1    2    6  Korak 1 Izračun lastnih vrednosti   3 1   1 2 1  0   1        2     0   4 2  6  karakteristični polinom   9 3 9 4    2 9 4     3  0   ,   1 2  4 2  3      1 2 2 Obe lastni vrednosti sta realni, enaki in manjši od 0      0 1 . Korak 2 Poiščimo lastni vektor.  3   Za 2  1 2 2       0      0 2 2  1          1  3           1 3    3 3       1 3         0   0     1 1   2 1 1       6 2   6 2  KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB   1 3       1 2  2 2     0 3        vzemimo prvo vrstico  0   1 1          1 1 6 2  1 3 3          0  1 2      3  2   1 1 2     . Izberimo 2 2   2 2  in dobimo     . 1    3   Korak 3 Splošna rešitev sistema: y     3 3  t 2     C e 1   1 Ker imamo dvojno lastno vrednost, smo dobili eno samo rešitev. Potrebujemo pa dve. V takem primeru je splošna rešitev: y t  t  t   C e      C     1 2 t e      e  dodatek 4.1 – enačba 14, D4.1 Korak 4 Določimo vektor  . y          3 2 2  1 1 3         1 3   1  1        3    1 1 2 , postavimo   0 2   2   2 2    6 2  potem dobimo vektor:        6   . 0 Lastni vektor     3    6   je:    in vektor  je:    .   1   0 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Korak 5 Splošna rešitev sistema diferencialnih enačb: y    3   3   3 3   t   3   t  6  t  C    e 2  C   t  e 2    e 2 1 2         1     1   0  Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.6 (Maxima) en1: 'diff(y1,t)=−1·y1 + 3/2·y2; en2: 'diff (y2, t) = −1/6·y1 − 2·y2; load ("plotdf"); A: matrix ([−1,3/2], [−1/6, −2] ); charpoly (A, λ)),expand; solve (λ^2+3·λ+9/4); eigenvectors (A); plotdf ([−1·y1+3/2·y2,−1/6·y1−2·y2],[y1,y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.8 Primer 4a Sedlo Sistem diferencialnih enačb: y   1 0   0  y  1  Korak 1 Izračun lastnih vrednosti: 1   0  0   1        1    0   0  karakteristični polinom 1   2   1 0    ,   1 1 2    1 1 in   1 2 Lastni vrednosti sta realni, po absolutni vrednosti enaki, ena je manjša in druga večja od 0    ,   1 2 1  . Korak 2 Poišči lastni vektor. Za   1 1  1 1  0     0 1           vzemimo drugo vrstico   0   1 1     o 2   0     2  0  0     2  0   1  1 2 1 2 . Izberimo za 1 in 2 mora biti 0  0 :  2 in dobimo:    1    1   0 Za 1    2  1 1  0       0  1        vzemimo prvo vrstico  0   1 1      o 2   2     0  0 0     2   1 2  2 1 . Izberimo za 1 2 in 1 mora biti 0 .   0 1 in dobimo KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB    0    2   1 Korak 3 Splošna rešitev sistema y    1   0 t t  C   e C     e 1 2     0   1 Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.7 (Maxima) (%i1) en1: 'diff (y1, t) =y1; (%i2) en2: 'diff (y2, t) =−y2; (%i3) load ("plotdf"); (%i4) A: matrix ( [1, 0], [0, −1] ); (charpoly(A, λ)), expand; solve(λ^2−1=0); v: eigenvectors (A); plotdf ( [y1, y2], [y1, y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.9 Primer 4b Sedlo Sistem diferencialnih enačb: y    1 2     2 3  y Korak 1 Izračun lastnih vrednosti 1   2  0  1       2     3 2  6 0  karakteristični polinom  1      2     6 0   2 3  9   4 4     3  ,  4 0    1 2 2    1   4   4 1 in 2 2 Lastni vrednosti sta realni, različni, pri tem ena je manjša in druga večja od 0    0,   0 1 2  . Korak 2 Poiščimo lastna vektorja. Za   1 1  1 1  2     0 1           vzemimo prvo vrstico   3 2 1      0 2   2     2  0     1 2  12     izberimo za   1     1    2   2  2 in dobimo       1    1   . 1 Za   4 2  1 4  2       0  1        vzemimo prvo vrstico   3 2 4      0 2   KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB     2    2   2      2 3  2   0 3        3 1 2  2   1 2  1 2    3    2   Izberimo za   3 2    2   2    in dobimo .   3 Korak 3 Splošna rešitev sistema: y     1   2  t 4t C   e  C   e 1   2     1   3 Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.8 (Maxima) en1: 'diff (y1,t)=y1+2·y2; en2: 'diff (y2,t)=3·y1+2·y2; load ("plotdf"); A: matrix([1,2],[3,2]); (charpoly(A,λ)), expand; solve (λ^2−3·λ−4=0); v: eigenvectors (A); plotdf([y1+2·y2,3·y1+2·y2],[y1, y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.10 Primer 5a Center Sistem diferencialnih enačb: y   0 1       y 4 0  Korak 1 Izračun lastnih vrednosti 0   1  0  0       0      4 0  karakteristični polinom 4 0    2   2 4 0    4    ,   4 1 2     2 i     2 i 1 in 2 Lastni vrednosti sta imaginarni. Korak 2 Poiščimo lastna vektorja. Za    2 i 1    2 i 1     0 1           vzemimo prvo vrstico   4   2 i      0 2       2 i    0      2 i  1 2 2 1    . Izberimo   1 2    1   1     i1  1    1  1     in dobimo . 2   i  Za     2 i 2 dobimo drugi lastni vektor       2  1    . 2 i  Korak 3 Splošna rešitev sistema: Rešitev, ki jo dobimo od prve lastne vrednosti    2 i 1 , je: y   1  1    1  1  2     i t y            e   t t   cos 2 sin 2  2    i   2 i  KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB .-.-.- Formula et  cos      i sin    e   t  cos      i sin    .-.-.-  y         1 cos 2   t sin 2   t        i   1  y    u i v   2 sin 2    t 2 cos 2        t      Pri tem je: u  cos 2    t   sin 2    t  v       in   2 sin 2     t  2 cos 2      t      Funkciji u t v t in sta dve linearno neodvisni rešitvi     sistema diferencialnih enačb (Dawkins, 2007, str. 291). To pomeni, da ju lahko uporabimo kot realni rešitvi, da oblikujemo splošno rešitev. Splošna rešitev sistema diferencialnih enačb z imaginarnimi oziroma kompleksnimi lastnimi vrednostmi je: y  C u t    v t 1   C 2   Splošna rešitev obravnavanega sistema diferencialnih enačb je: y  cos 2    t   sin 2    t   C     C  1 2      2 sin 2    t  2 cos 2      t      Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.9 (Maxima) en1: 'diff(y1,t)=y2; en2: 'diff (y2,t)=−4·y1; load ("plotdf"); A: matrix ([0,1],[−4,0]); (charpoly (A, λ)),expand; KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB solve(λ^2+4=0); v: eigenvectors (A); plotdf([y2,−4·y1],[y1,y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.11 Primer 5b (Dawkins, 2007, str. 290-293). Center Sistem diferencialnih enačb: y 3  9      y 4   3  Korak 1 Izračun lastnih vrednosti: 3    9  0  3        3    36  4   0 3   karakteristični polinom  2     9 4 9 0   2  27  0    1 ,2   27   1   3 3  i in 2    3 3  i Lastni vrednosti sta imaginarni. Korak 2 Poiščimo lastna vektorja. Za 1    3 3  i   3 3   3  i  9      1   0         vzemimo prvo vrstico            0  4 3 3 3 i  2 3 1   1  3     i  9  0       1    2 1 2  1 3 i   3      1  1    1    3   izberimo :    1 3   i  1 in dobimo  1   3     3  1       1  3 i     Za 3     3  i 2 dobimo drugi lastni vektor. Ta vektor je: KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB    3  2       1  3  i     Rešitev , ki jo dobimo od prve lastne vrednosti    3 3  i 1 y    3  1           3   i  3 3 i t e 1       3  1        y     cos 3  3    t i sin 3  3  t  1  3  i    y    3 cos(3 3 ) 3     t isin(3 3 )  t  1         cos(3                  3 ) t i 3 cos(3 3 ) t i sin(3 3 ) t 3 sin(3 3 ) t  y   3 cos(3   3 )  t   3 sin(3   3 )  t  1      i     cos(3  3 )   t 3 sin(3   3 )  t   sin(3  3 )   t 3 cos(3   3 )  t      y   1 u t  ( )   i v t ( )  3 cos(3  3 )  t   3 sin(3  3 )  t  Pri tem je u    , v    cos(3  3 )   t 3 sin(3     3 ) t     sin(3  3 )   t 3 cos(3  3 ) t   Splošna rešitev sistema je: y  C u t  ( )  C v t  ( ) 1 2 . Za obravnavani primer je splošna rešitev:  3 cos(3 3 ) t   3 sin(3 3 ) t        y  C     C  1 2   cos(3  3 )   t 3 sin(3   3 )  t sin(3  3 )   t 3 cos(3   3 )     t  Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.10 (Maxima) en1:'diff(y1,t)=3·y1−9·y2; en2:'diff(y2,t)=−4·y1−3·y2; load("plotdf"); A: matrix([3,-9],[4,−3]); (charpoly(A, λ)),expand; KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB solve(λ^2+27=0); v: eigenvectors (A); plotdf([3·y1−9·y2,4·y1−3·y2],[y1,y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.12 Primer 6a Spiralni izvor Sistem diferencialnih enačb: y      2 3     y 3 2  Korak 1 Izračun tipične vrednosti: 2   3  0  2       2      9 0  karakteristični polinom 3 2       4  13  0    2 4  16   4 13 4  36 ,  1 2    ,  1 2  2 2     2 3 i     2 3 1 in i 2 Lastni vrednosti sta konjugirano kompleksni števili. Korak 2 Poiščimo lastna vektorja. Za     2 3 i 1 2 2 3    i 3       0 1   3i 3       0    1                3 2 2 3    i  2 0            3 3 i    2   0 vzemimo drugo vrstico       3  3 i  0 1      i 2 12     1    i2     izberimo:   1 in dobimo lastni vektor   2  2       1    i   1 Za     2 3i 2    2   i    je lastni vektor   1 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Korak 3 Splošna rešitev sistema Rešitev, ki jo dobimo od prve lastne vrednosti     2 3 i 1 y   1    i   1    i    2    t 3 i  2 t   e  y e   cos 3        t i sin 3    t    1   1 y   1 sin 3    t   cos 3    t  2  t      2  t e i e    cos 3    t sin 3       t      Splošna rešitev je: y   sin 3    t   cos 3    t 2  t 2         t C e C e 1    2  cos 3 t sin 3        t    Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.11 (Maxima) en1:'diff(y1,t)=2·y1+3·y2; en2:'diff(y2,t)=−3·y1+2·y2; load("plotdf"); A: matrix ( [2, 3], [−3, 2] ); (charpoly (A,λ)),expand; solve(λ^2−4·λ+13=0); v: eigenvectors (A); plotdf([2·y1+3·y2,−3·y1+2·y2],[y1,y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 4.13 Primer 6b Spiralni ponor Sistem diferencialnih enačb: y   1   1   2  y   1   1  Korak 1 Izračun lastnih vrednosti: 2  1 0          1   1      1  1 0  1 2    1   karakteristični polinom 1 1 7   4  2  1 1 1 0      4 2 4     ,    ,     2 1 2 1 2 2 4 2 4 2      1 i  7 1 i  7     1 2 in 4 2 4 2 Lastni vrednosti sta kompleksno in konjugirano kompleksno število (vedno nastopata v paru). Korak 2 Poiščimo lastna vektorja. Za    1  1 i  7 4 2 y     1 1 i  7     1 1       0  2 4 4  1               0 1 i 7   1    2  1   4 4  KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y     3 i  7    1 1       0  4 4    1         vzemimo prvo vrstico     0 3 i 7 2    1     4 4      3 i  7 3 i  7             0  1     2  2  1   4 4   4 4    1     1      3  i  7       1   izberimo za in dobimo   1  1 4   4    1     1      3  i  7       4 4     Lastni vektor zapišimo v obliki:     1 0  1          1 i          . i 3 7       4    4  Korak 3 Splošna rešitev sistema Rešitev, ki jo dobimo od prve lastne vrednosti    1  1 i  7 4 2 y    1 i 7  1    0       t 1    4 4         e    i  3 7            4     4    y   1    1 0   i  7 1   t      e 4  e 4      i  3 7            4     4    y   1  1 0       1 t             7 7  4 e    i   cos    t i sin     t    3 7              4 4           4     4   KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y   1  0     1   1 t      7 7       e 4   cos   t  sin  3   t        7       4  4         4     4       1 0   1   t      7   7 4    i e    sin   t  cos     t 3      7       4  4         4     4   Splošna rešitev sistema diferencialnih enačb je: y   1    1 0     1 t      7 7  4      C e    cos   t  sin  t  1 3         7       4  4         4     4   C   1    0 1   t      7 7       e 4   sin   t  cos     t 2 3     7        4  4         4     4   Izračun lastnih vrednosti, lastnih vektorjev, fazni portret: slika 4.12 (Maxima) en1:'diff(y1,t)=1/2·y1−1·y2; en2: 'diff(y2,t)=−1·y1−1·y2; load("plotdf"); A: matrix([1/2,−1],[1,−1]); (charpoly(A,λ)),expand; solve (λ^2−λ/2+1/2=0); v: eigenvectors (A); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB plotdf([1/2·y1−1·y2,1·y1−1·y2], [y1, y2])$ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 5. KVALITATIVNA ANALIZA SISTEMA DVEH NELINEARNIH DIFERENCIALNIH ENAČB Odstavek je povzet po Panfilovu (2010, str. 83-86). Sistem dveh nelinearnih diferencialnih enačb lahko splošno zapišemo v sledeči obliki: dx  f x y   , dt dy (En. 1)    g x y , dt Kritične točke ali točke ravnotežja so točke, kjer je sistem stacionaren, ko pride sistem v stacionarno stanje, ostane v tem stanju in se ne spreminja več. Točke ravnotežja so določene z enačbama f x y ,    0 in g x y  ,   0 . Točka  x , y , 0 g x 0 r r  je točka ravnotežja, če sta f x  y  r r  in  , y  r r  . Sistem ima lahko več točk ravnotežja, ki vplivajo na parametrične krivulje. Ravnino X-Y imenujemo fazna ravnina, večkrat uporabljamo tudi izraz vektorsko polje (angl. Vector field). V fazni ravnini imajo vektorji, ki kažejo smeri tangent na parametrične ozirom fazne krivulje, v posameznih točkah dve komponenti: V  dx dy    v , v y  , x    .  dt dt  Vsaka komponenta je lahko pozitivna ali pa negativna. Tako dobimo v vektorskem polju štiri področja, označimo jih z I, II, III in IV (slika 5.1). Na prehodu med poljema I in II ter poljema III in IV, ki se razlikujeta po predznaku. dx dx  0 , je x-ničta izoklina (angl. Null-cline), kjer je . dt dt Na x-ničti izoklini med področjema I in II so vektorji usmerjeni navpično navzgor in na ničti izoklini med področjema III in IV so vektorji usmerjeni navpično navzdol. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 5.1: Smeri vektorjev (Panfilov, 2005, str. 83) Na prehodu med poljema I in III ter poljema II in IV, ki se razlikujeta po dy dy  0 predznaku , je y-ničta izoklina, kjer je Na y-ničti izoklini med dt dt področjema I in III so vektorji usmerjeni vodoravno desno in na ničti izoklini med področjema II in IV so vektorji usmerjeni v nasprotni smeri. Definicija (Panfilov, 2005, str. 84): x-ničta izoklina je množica točk, ki zadostijo pogoju f x y ,    0 . Geometrično so to krivulje v vektorskem polju, kjer je horizontalna komponenta vektorja v  0 x , in razmejujejo področja, kjer vertikalna komponenta vektorja vy spremeni predznak. y-ničta izoklina je množica točk, ki zadostijo pogoju g x y ,    0 . Geometrično so to krivulje v vektorskem polju, kjer je vertikalna komponenta vektorja v  0 y , in razmejujejo področja, kjer horizontalna komponenta vektorja v x spremeni predznak. Kritične ali ravnotežne točke so na presečiščih x-ničtih izoklin in y-ničtih izoklin. Primer 1 (Panfilov, 2010, str. 85-86) Sistem diferencialnih enačb podajata enačbi 2 in 3: dx  x      r x 1      b x y (En. 2) dt  k  dy      c x y d y (En. 3) dt To je tako imenovan Lotka-Volterra model , dopolnjen z upoštevanjem logistične rasti plen populacij. V nadaljnjem nas ne bo zanimala biološka KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB interpretacija modela. Vzemimo, da imajo konstante sledeče vrednosti: r  3, k  1 , b  1.5 , c  0.5 in d  0.25 . Enačbi 2 in 3 sedaj zapišimo kot enačbi 4 in 5: dx     3 x  1 x  1.5   x y (En. 4) dt dy  0.5    x y 0.25  y (En. 5) dt Poiščimo točke ravnotežja. Iščemo torej rešitve sledečega sistema enačb, enačbi 6 in 7: 3     x  1 x  1.5   x y 0 (En. 6) 0.5    x y 0.25   y 0 (En. 7) x-ničti izoklini določimo iz enačbe 6. Prva x 0 x-ničta izoklina je takoj razvidna iz enačbe 6 :  ( y os). Določimo še drugo x-ničto izoklino: 3     x  1 x  1.5   x y 0  3 3    x1.5   y 0  y    2 2 x y-ničti izoklini določimo iz enačbe 7. Prva y-ničta izoklina je prav tako takoj razvidna iz enačbe 7: y  0 ( x os). Določimo še drugo y-ničto izoklino: 0.5    x y 0.25   y 0  0.5  x 0.25  0  x  0.5 (vertikala) Imamo torej: x-ničti izoklini:: x  0 in y    2 2 x y-ničti izoklini: x  0.5 y  0 in Določimo točke ravnotežja, te so na presečiščih x-ničtih in y-ničtih izoklin.   0, 0  y  0 presečišče izoklin: x  0 in   1, 0  presečišče izoklin y    2 2 x in y  0 0.5, 1  presečišče izoklin y    2 2 x in x  0.5 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Vektorsko polje z ničtimi izoklinami (dodatek 5.1) je prikazano na sliki 5.2. Z biološkega vidika je model interakcij med populacijami smiseln le v prvem kvadrantu. Slika 5.2: Vektorsko polje z ničtimi izoklinami – enačbi 4 in 5 ( x-ničti izoklini – modra, y-ničti izoklini - rdeča) S faznega portreta na sliki 5.3 razberemo, da je točka 0.5,1 stabilna.   Slika 5.3: Fazni portret; Dif. enačbi: dx      3 x (1 x) 1.5   x y in dt dy  0.5    x y 0.25  y dt Primer 2 Sistem diferencialnih enačb podajata enačbi 8 in 9: dx    x  1 x    dt x y (En. 8) dy  y      2 y  1     3 x y (En. 9) dt  2  KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Omejimo se na področje x  0 in y  0 . dx dy  0  Ničte-izokline 0 so opredeljene z in , (enačbi 10 in 11): dt dt x       1 x  x y 0 (En. 10) 2 y 1         y    3 x y 0 (En. 11)  2  x-ničti izoklini določimo iz enačbe 10. Prva x-ničta izoklina je takoj razvidna iz enačbe 10: x  0 ( y os). Določimo še drugo x-ničto izoklino: x       1 x  x y 0  1    x y 0  y   1 x y-ničti izoklini določimo iz enačbe 11. Prva y 0 y-ničta izoklina je razvidna iz enačbe 11  ( x os). Določimo še drugo y-ničto izoklino: 2  y     y 1       3 x y 0  2     y 3 x 0  y    2 3 x  2  Imamo torej: x-ničti izoklini: x  0 in y   1 x . y-ničti izoklini y 0 y x :     2 3 in . Določimo točke ravnotežja, te so na presečiščih x-ničtih in y-ničtih izoklin:   0, 0  presečišče izoklin x  0 in y  0 .   1, 0  presečišče izoklin y   1 x in y  0 . 0, 2  presečišče izoklin x  0 in y    2 3 x . 0.5, 0.5  presečišče izoklin y   1 x in y    2 3 x . Ničte izokline razdelijo vektorsko polje na 4 področja ( dodatek 5.2). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 5.4: Vektorsko polje in ničte izokline sistema; Dif. enačbi dx dy y      x (1 x ) x y in     2 y (1 ) 3    x y dt dt 2 Fazni portret je prikazan na sliki 5.5. Slika 5.5: Fazni portret; Dif. enačbi dx dy y      x (1 x ) x y in     2 y (1 ) 3    x y dt dt 2 Primer 3 Sistem diferencialnih enačb: dx    x  1 x    dt x y (En. 12) dy  2 y      2 2 y 1      3 x y (En. 13) dt  2  dx dy  0  Ničte izokline 0 so opredeljene z in (enačbi 14 in 15): dt dt x       1 x  x y 0 (En. 14) KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 2     2 y  2 y 1       3 x y 0 (En. 15)  2  x-ničti izoklini določimo iz enačbe 14. Prva x-ničta izoklina je x  0 ( y os). Določimo še drugo x-ničto izoklino: 1   x  y 0  y   1 x y-ničti izoklini določimo iz enačbe 15. Prva y  0 y-ničta i zoklina je ( x os). Določimo še drugo y-ničto izoklino: 2  2 y  2        y  1  3 x y 0 2 2 2  2  y   3 x  2 0  3  x  y  2  2   druga y-ničta izoklina je elipsa. Imamo torej: x-ničti izoklini y x : x  0   1 in y-ničti izoklini 2 2 : y  0 in 3  x  y  2 Določimo točke ravnotežja, te so na presečiščih x-ničtih in y-ničtih izoklin:   0, 0  presečišče izoklin x  0 in y  0   1, 0  presečišče izoklin y   1 x in y  0 0, 1.41 2  presečišče izoklin x  2 0 in 3  x  y  2 0.80, 0.20 in  0.31, 1.31  presečišče izoklin y   1 x in 3  2 2 x  y  2 (dodatek 5.3). Točke ravnotežja so:   0, 0  ,  0, 1.41  1, 00.31, 1.31 ,   0.80, 0.20 ,   in  . Vektorsko polje z ničtimi izoklinami je prikazano na sliki 5.6, fazni portret pa je prikazan na sliki 5.7 (dodatek 5.3). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 5.6: Ničte izokline in smerno polje; Dif. enačbi 2 dx dy y      2 x (1 x ) x y in     2 y (1 ) 3    x y dt dt 2 Slika 5.7: Fazni portret; Dif. enačbi: dx 2      x (1 x) x y     in dy y 2 2 y (1 ) 3    x y dt dt 2 KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 5.1 Primer 1 Sistem diferencialnih enačb: dx  x      r x  1     b x y dt  k  dy       c x y d y dt Vrednosti konstant so r  3 , k  1 , b  1.5 , c  0.5 in d  0.25 . x-ničti izoklini  dx  0 x  0 (vertikala, ki sovpada z y-osjo) in dt y    2 2 x - premica, ki gre skozi točke 0, 2 ,   1, 0 , 1.5, 1  . y-ničti izoklini  dy  0 0.5    x y 0.25   y 0 enačba ima dve rešitvi dt y  0 (horizontala, ki sovpada z x-osjo) in x  0.5 (vertikalna linija). Točke ravnotežja so na presečiščih x-ničtih izoklin z y-ničtimi izoklinami. Točke ravnotežja so torej 0, 0 , 0.5, 1 1, 0 in       1, 0   . Opomba: V poglavju 3 smo za ponazarjanje faznega portreta uporabili ukaz: 'plotdf’, v tem poglavju bomo uporabili tudi ukaz: ’drawdf’. Vektorsko polje in ničte izokline: slika 5.2 (Maxima) load (drawdf); drawdf ([3·x·(1−x)−1.5·x·y, 0.5·x·y−0.25·y],[x,−1,2.2],[y,−1, 2.2], field_degree=1, line_width=2, field_grid=[15, 15] , color=blue, line_width=3, explicit (2−2·x, x,−1, 2), polygon ([[0,−1],[0,0],[0,2.2]]),color = red, polygon ([[0.5,−1],[0.5,0],[0.5,2.2]]), polygon ([[−1,0],[1,0],[2.2,0]]), color = black, point_at(0,0),point_at(0.5,1),point_at(1,0)); Fazni portret: slika 5.3 (Maxima) load (draw); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB load (drawdf); set_draw_defaults ( xrange=[−3, 3], yrange=[−3, 3], grid=true, line_width=1, color=black); drawdf ([3·x·(1−x)−1.5·x·y, 0.5·x·y−0.25·y], [x,0,2.2], [y,0,2.2], color=black, field_degree=1, field_grid=[20, 20], line_width=3, point_at(0, 0), color=red, solns_at ([0.5, 1.5],[0.8, 1.2],[1, 1.5],[2, 1.7],[1, 0.5], [1.5,0.8],[0.5, 0.3],[0.3, 0.7] )); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 5.2 Primer 2 Sistem diferencialnih enačb: dx    x  1 x    x y dt dy  y      2 y  1     3 x y dt  3  x-ničti izoklini x  0 (vertikala, ki sovpada z y-osjo) in y   1 x (premica, ki gre skozi točki 0, 1 in 1, 0 )     y-ničti izoklini y  0 (horizontala, ki sovpada z x-osjo in y    2 3 x (premica, ki gre skozi točki  2  , 0   0, 2  in   3  Točke ravnotežja so 0, 0 , 0, 2 , 1, 0 in 0.5, 0.5 .         Vektorsko polje. Slika 5.4 (Maxima) lord (draw); lord (drawdf); set_draw_defaults( xrange=[0, 1], yrange=[0, 2], grid=true, line_width=1); drawdf([x·(1−x)−x·y,2·y·(1−y/2)−3·x·y],[x,0,1],[y,0,2], field_degree=1, line_width=3, field_grid=[15, 15],color=red explicit(1−x,x,0,1), polygon([[0,0],[0,1],[0,2]]), color=blue, polygon ([[0,0],[0.5,0],[1,0]]), explicit(2−3·x,x,0,1)); Fazni portret: slika 5.5 (Maxima lord (draw); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB lord (drawdf); set_draw_defaults( xrange=[0,1], yrange=[0,2], grid=true, line_width=1); drawdf ([x·(1−x )−x·y, 2·y·(1−y/2)−3·x·y], [x,0,1], [y,0,2], field_degree=1, line_width=4, field_grid=[15, 15], color=red, explicit (1−x, x,0,1), polygon ([[0,0],[0,1],[0,2]]), color=blue, polygon ([[0, 0],[0.5, 0],[1, 0]]), explicit(2−3·x, x, 0, 1), line_width=2,color=red, solns_at ([[0.2,1],0.3,1],[0.4,0.5],[0.4,1.5], [0.6,1],[0.6,0.3],[0.8,0.5], [0.8,1]] )); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 5.3 Primer 3 Sistem diferencialnih enačb: dx    x  1 x    x y dt dy  2 y  2     2 y 1       3 x y 0 dt  2  Imamo sledeče izokline: x-ničti izoklini x  0 y 1 x in   y-ničti izoklini y  2 2 0 in 3  x  y  2 Točke ravnotežja na presečiščih x-ničtih izoklin z y-ničtimi izoklinami. Določimo presečišča med y   1 x x-ničto izoklino in y-ničto izoklino 3 2 2  x  y  2 . Postavimo  2 2  2 y   1 x v 3 x y  2 2  3    x  1 x   2  4 2  x     2 x 1 0 . Izračunajmo x x ,  x 2  20 1  1 2   in x x 0.80 in x   0.31 2 1 2. 8 Iz enačbe y   1 x dobimo y  0.20  1 in y 1.31 2 , pripadajoči točki ravnotežja sta 0.80, 0.20    in 0.31, 1.31 . Ničte izokline in smerno polje: slika 5.6 (Maxima) load (draw); load (drawdf); set_draw_defaults ( xrange=[−1.5,1.5], yrange=[−1.5,1.5], grid=true, line_width=2); drawdf ([x·(1−x)−x·y,2·y·(1−y^2/2)−3·x^2·y],[x,−1.5,1.5], KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB [y,−1.5,1.5], field_degree=1, field_grid=[15, 15], point_at(0, 0), line_width=4, color=blue, explicit (sqrt(2−3·x^2), x,−1.5,1.5), explicit (−sqrt(2−3·x^2), x,−1.5,1.5), polygon ([[−1.5,0],[0,0],[1.5,0]]), color=red, explicit (1−x, x,−1.5,1.5), polygon ([[0,−1.5],[0,0],[0,1.5]])); Fazni portret: slika 5.7 (Maxima) en1:'diff (x,t)=x·(1−x)−x·y; en2:'diff (y,t)= 2·y·(1−y^2/2)−3·x^2·y; plotdf ([x·(1−x)−x·y, 2·y(1−y^2/2)−3·x^2·y], [x,y],[x,−1.5,1.5],[y,−1.5,1.5 ]); KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB 6. GRAFIČNA ANALIZA STABILNOSTI – DIAGRAM PAJČEVINE Stabilnost ali nestabilnost dinamičnega sistema lahko ugotavljamo tudi z grafično analizo. V ta namen se pogosto uporablja tako imenovanimi 'diagram pajčevine' (angl. Cobweb diagram). Diagram pajčevine prikaže grafično funkcijo prenove y  F x   , premico y  x ter izmenjujoče se zaporedje vertikalnih in horizontalnih segmentov, združujoč točke y , y y 0 0 , , y 01 ,      y , y y , y y  F x 1 1  ,  1 2  itd. Točke ravnotežja so na presečiščih krivulje   in premice y  x . Primer diagrama pajčevine prikazuje slika 6.1. Slika 6.1: Diagram pajčevine; cos(y); izhodiščna točka y0=2 Grafično analizo z uporabo diagrama pajčevine bomo prikazali na poenostavljenem diskretnem modelu logistične rasti (dodatek 6.1). Populacija N naj ima konstantno stopnjo rasti a in stopnjo umrljivosti proporcionalno z velikostjo populacij: b N  . V intervalu n je prirastek populacije  2 a N b N  n osebkov, odmre pa osebkov. V začetku intervala n 1 je velikost populacije: N  2 2 N   a N   b N N n  1 n n n     1  1 a  N   b N n  n n  N 1 n  1 1 a N       b     N n   n  (En. 1)  1 a  Vpeljimo novo spremenljivko y n : KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB y   N N   y n b 1  a 1  n n n   a b 1  1  a a   y    1 a    y y n  1    n  1n   b b y      1  1 a  y   1 y n nn  (En. 2) Vpeljimo še oznako c    1 a in pridemo do enačbe 3: y     c y 1  1 y n  nn  (En. 3) Diagrami pajčevine diskretnega logističnega modela (enačba 3), narisani z različnih izhodišč, so prikazani na slikah 6.2, 6.3 in 6.4 (dodatek 6.1). Točke ravnotežja so na presečišču krivulj, ki jih ponazarjata funkcija prenove y  F x   y  y in premica n  1n . Če je v bližini točke ravnotežja x ri naklon krivulje funkcije prenove manjši od naklona premice y  x , potem zaporedje konvergira proti točki ravnotežja. Torej je pogoj F x 1     ri zadosten za konvergentnost zaporedja pri tej točki ravnotežja. S slik je razvidno, da je točka ravnotežje y  2 0 y  r1 r nestabilna in 2 stabilna. 3 Slika 6.2 Diagram pajčevine; enačba 3; parameter c=3; izhodiščna točka 0.05. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Slika 6.3 Diagram pajčevine; enačba 3; parameter c=3; izhodiščna točka 0.2. Slika 6.4 Diagram pajčevine; enačba 3; parameter c=3; izhodiščna točka 0.5. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Dodatek 6.1 Diskretni populacijski modeli Definicije: 1. Generacije se ne prekrivajo. 2. Velikost populacije v nekem časovnem intervalu je odvisna od velikosti populacije v predhodnem časovnem intervalu (predhodne generacije). N  N  f N n  1 n   n N  F N n  1   n Funkcija f N n predstavlja »stopnjo rasti populacije« po osebku ali »neto   stopnjo reprodukcij e«. Uporabljali bomo izraz » fitnes funkcija« (angl. Fitness function). Funkcija F N n se imenuje »funkcija prenove«, ki je praviloma nelinearna in   ima smisel le, če je ne-negativna - N  0 . Diskreten model logistične rasti: N  Nn   N   a N   1 n  1 n n   (En. 1, D6.1)  K  konstantna stopnja rasti: a , nosilna kapaciteta okolja: K . b  a Označimo: K N  2 N   a N   b N n (En. 2, D6.1)  1 n n n Enačbi 1, D6,1 in 2, D6. nekoliko preoblikujmo v enačbo 3, D6.1: N  b     a N  n 1  1     1 N  n  n  (En. 3, D6.1) 1  a  Vpeljimo novo spremenljivko y n y   N N   y n b 1  a 1  n n n  a b in jo vstavimo v enačbo 3, D6.1 1  a 1  a  y   a n  1      y  1 1 y nn   y      1 y y n 1     b  a n 1 n b KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Vpeljimo še oznako: + c    1 a  Tako pridemo do enačbe 4, D6.1: y     c y 1 y n  1 n n  (En. 4, D6.1 4) Diferenčna enačba : y  y       y c y 1 y  y   y y    y n  1 n n  n  n   c  1n  1 n Funkcija prenove: f y       c y  1 y nn  Določimo točki ravnotežja: za c  3   y 0  y    r  c  1 y    1 r  0 Imamo dve točki ravnotežja: y  0 r 1 in c    1 y 2    y c  1 2  1 0 y  r  r 2  r 2 c 3 Diagrami pajčevine: slike 6.2, 6.3, 6.4 (Maxima) kill (y)$ load ("dynamics") $ a: logen: 3·y·(1−y); 3*(1-y)*y fixed: solve (logen − y); [y=2/3,y=0] kill (y)$ load ("dynamics") $ staircase (3·y·(1−y), 0.05, 8, [y, 0, 1]) $ staircase (3*y*(1-y), 0.2, 8, [y, 0, 1]) $ staircase (3*y*(1-y), 0.5, 8, [y, 0.5, 0.8]) $ KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB Kratka obrazložitev Maxima ukaza: Ukaz: n - Število vertikalnih segmentov [y, 0, 1] – Obseg grafa - velikost osi Ukaz nariše funkcijo f (v obravnavanem primeru je to funkcija prenove) in premico yn+1 = yn. Vertikalni segment se nariše od točke (y0, y0) (na premici) do točke, kjer se vertikala dotakne funkcije f. Od te točke se potegne horizontalni segment do točke (y1, y1) (na premici). Procedura se ponovi n-krat, do točke (yn, yn). KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB ZAKLJUČEK Kvalitativna analiza `diferencialnih enačb` obravnava grafično predstavitev navadnih diferencialnih enačb, to je integralnih krivulj, in njihovo odvisnost od vrednosti parametrov kot tudi od začetnih vrednosti. Poseben poudarek je na določitvi točk ravnotežja in njihovi stabilnosti. Poglavji 4 in 5 obravnavata sistem dveh diferencialnih enačb. V poglavju 4 so to linearne diferencialne enačbe prvega reda in v poglavju 5 nelinearne diferencialne enačbe. Večina modelov na področju ekologije, biologije in upravljanja podeželja obravnava interakcije med dvema spremenljivkama. Grafične rešitve so trajektorije v fazni ravnini. Za predstavitev grafičnih rešitev je uporabljena aplikacija Maxima. Uporaba te aplikacije je podana v dodatkih posameznih poglavij. Za konec naj opozorimo, da obstaja še ena metoda ugotavljanja stabilnosti točk ravnotežje in lastnosti teh točk. Ta metoda je tako imenovana 'linearna analiza stabilnost'. Bistvo linearne analize stabilnost je v tem, da diferencialno enačbo oziroma sistem diferencialnih enačb v okolici točk stabilnosti pretvorimo v linearno diferencialno enačbo, odnosno sistem linearnih diferencialnih enačb, in ga kot takega tudi obravnavamo. Linearizacija je zasnovana na Taylorjevem približku prvega reda. KVALITATIVNA ANALIZA DIFERENCIALNIH ENAČB LITERATURA 1. Davis, H., T. (1962). Introduction to Nonlinear Differential and Integral Equations. Dover Publications, Inc, New York. 2. Dawkins, P. (2007). Differential Equations. http://tutorial.math.lamar.edu/trms.aspx 3. Kreyszig, E., Kreyszig, H., & Norminton, E., J. (2011). Advanced Engineering Mathematics. Tenth edition, John Wiley & Sons. 4. Nagy, G. (2015). Ordinary Differential Equation. Mathematics Department. Michigan State University. East Lansing, MI, 48824. http://users.math.msu.edu/users/gnagy/teaching/ode.pdf 5. Panfilov, A. (2005). Qualitative Analysis of Differential Equation. Utrecht University, Utrecht. 6. Weckesser, W. (2005). Math 312 Lecture Notes Linear Two-dimensional Systems of Differential Equations. Department of Mathematics. Colgate University. http://math.colgate.edu/~wweckesser/math312Spring05/handouts/Linear2 x2.pdf