Obsah 1 Základní pojmy 1 2 Elementární metody řešení 5 2.1 Rovnice typu x = f(t) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.2 Exaktní rovnice . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.3 Rovnice se separovanými proměnnými x = f(t)g(x) . . . . . . . . . . . . . . . . . 5 2.4 Homogenní rovnice x = f x t . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2.5 Rovnice typu x = f at + bx + c t + x + . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2.6 Lineární rovnice x = a(t)x + b(t) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.7 Bernoulliova rovnice x = a(t)x + b(t)xr , r R . . . . . . . . . . . . . . . . . . . . 8 2.8 Rovnice nerozřešená vzhledem k derivaci F(t, x, x ) = 0 (implicitní rovnice) . . . . 8 2.9 Rovnice druhého řádu typu x = f(x) . . . . . . . . . . . . . . . . . . . . . . . . 10 2.10 Rovnice typu F(t, x(k) , x(k+1) , . . . , x(n) ) = 0, k {1, . . ., n - 1} . . . . . . . . . . . 10 2.11 Autonomní rovnice typu F x, x , x , . . . , x(n) = 0 . . . . . . . . . . . . . . . . . 11 2.12 Rovnice homogenní v x, x , x , . . . , x(n) . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.13 Ekvidimensionální rovnice . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.14 Cvičení . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 3 Obecné vlastnosti obyčejných diferenciálních rovnic 13 3.1 Vektorové a maticové funkce . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 3.2 Existence a jednoznačnost řešení systému ODR . . . . . . . . . . . . . . . . . . . . 14 3.3 Globální vlastnosti řešení systému ODR . . . . . . . . . . . . . . . . . . . . . . . . 17 3.4 Odhady řešení . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 4 Lineární rovnice 21 4.1 Systémy lineárních ODR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 4.2 Lineární diferenciální rovnice vyššího řádu . . . . . . . . . . . . . . . . . . . . . . . 27 4.3 Eulerova a Riccatiho diferenciální rovnice . . . . . . . . . . . . . . . . . . . . . . . 32 4.4 Cvičení . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 5 Autonomní systémy 35 5.1 Fázový prostor, trajektorie, stacionární body . . . . . . . . . . . . . . . . . . . . . 35 5.2 Autonomní systémy v rovině . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 5.3 Stabilita . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 6 Aplikace 49 6.1 Některé klasické elementární úlohy . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 6.2 Epidemiologický model Daniela Bernoulliho . . . . . . . . . . . . . . . . . . . . . . 62 6.3 Model populace produkující škodlivé odpady . . . . . . . . . . . . . . . . . . . . . 64 6.4 Lotkovy-Volterrovy systémy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 i Následující text je zápisem části přednášky předmětu M5858 Diferenciální rovnice a jejich užití. Má sloužit především k tomu, aby studentky/studenti nebyly/i nuceny/i si během přednášky dělat podrobné písemné poznámky, ale raději se soustředily/i na pochopení výkladu. Dále může být pomůckou k rychlému připomenutí toho, co člověk již zná. V žádném případě nemůže být považován za základní zdroj nahrazující standardní učební texty, z něhož by bylo možné se problematice diferenciálních rovnic naučit; sám o sobě bez komentářů během přednášky je málo srozumitelný až nesrozumitelný (aby byl s komentáři srozumitelný, je mým přáním a bude mou snahou). V současné chvíli se stále jedná o polotovar; určitě obsahuje i nějaké nedůslednosti, formulační nedostatky, překlepy nebo dokonce chyby. Budu vděčný každému, kdo mě na ně upozorní. Zdeněk Pospíšil prosinec 2009 Za základní učební texty k předmětu M5858 lze považovat: 1. J. Kalas, M. Ráb: Obyčejné diferenciální rovnice. MU, Brno 2001 (druhé vydání), 212 stran. Teorie obyčejných diferenciálních rovnic probraná důkladněji, než je v sylabu předmětu M5858. 2. J. Kalas, Z. Pospíšil: Spojité modely v biologii. MU, Brno 2001, 265 stran. Doplňky k teorii autonomních systémů, aplikace diferenciálních rovnic především v populační dynamice a teorii šíření epidemií. 3. M. Ráb: Metody řešení obyčejných diferenciálních rovnic. MU, Brno 1998, 96 stran. Popis základních elementárních metod řešení explicitních obyčejných diferenciálních rovnic. 4. R. Plch: Příklady z matematické analýzy. Diferenciální rovnice. MU, Brno 1995, 29 stran. Sbírka úloh z elementárních metod řešení explicitních i implicitních obyčejných diferenciálních rovnic. Je doplněna stručným popisem potřebných metod. Jako doplňující literaturu lze doporučit * P. Hartman: Ordinary Differential Equations. John Willey&Sons, New York-London-Sydney, 1964. Klasická monografie o teorii obyčejných diferenciálních rovnic. * E. Kamke: Differentialgleichungen Lösungsmethoden und Lösungen. Band I, Gewöhliche Differentialgleichungen. Akademische Verlagsgesellschaft Geest & Portig, Leipzig 1951. Důkladná příručka všech rovnic řešitelných elementárními metodami. * J. Kaucký: Elementární metody řešení obyčejných diferenciálních rovnic. Nakladatelství ČSAV, Praha 1952. Popis elementárních metod řešení obyčejných diferenciálních rovnic obsáhlejší než skripta 3 * N. F. Britton: Essential Mathematical Biology. Springer, London-Berlin-Heidelberg-New York-Hong Kong-Milan-Paris-Tokio, 2003 (second printing). Učebnice deterministických modelů v biologii; první tři kapitoly obsahují aplikace probírané v rámci předmětu M5858. * R. J. Barro, X. Sala-i-Martin: Economic growth. The MIT Press, Cambridge, MassachusettsLondon, England, 1999. Obsahuje aplikace obyčejných diferenciálních rovnic v ekonomii. ii Kapitola 1 Základní pojmy 1.0.1 Definice Buď G R2 množina s neprázdným vnitřkem, f : G R. Rovnice x = f(t, x) (1.1) se nazývá obyčejná diferenciální rovnice prvního řádu rozřešená vzhledem k derivaci. Řešením této rovnice se rozumí diferencovatelná funkce x : J R, kde J R je interval, která splňuje podmínky t, x(t) G, x (t) = f t, x(t) pro každé t J. Graf řešení rovnice (1.1) se nazývá integrální křivka. 1.0.2 Příklad G = {(t, x) R2 : t > 0}, f(t, x) = x t . Funkce x(t) = kt, kde k R je řešením rovnice x = x t na intervalu J = (0, ). Diferenciální rovnice může mít více řešení. 1.0.3 Definice Nechť G, f mají stejný význam jako v 1.0.1 a nechť (t0, x0) G je libovolný bod. Úloha najít řešení rovnice (1.1), které splňuje podmínku x(t0) = x0 (1.2) se nazývá počáteční nebo Cauchyova úloha, podmínka (1.2) se nazývá počáteční nebo Cauchyova podmínka. 1.0.4 Příklad Nechť t0 > 0, x0 R. Funkce x(t) = x0 t0 t je řešením úlohy x = x t , x(t0) = x0 na intervalu J = (0, ). 1 1.0.5 Definice Nechť x = x(t) je řešením úlohy (1.1), (1.2) na intervalu J a ~x = ~x(t) je řešením úlohy (1.1), (1.2) na intervalu ~J, t0 J ~J. Jestliže ~J J a pro každé t ~J je x(t) = ~x(t), řekneme, že řešení x = x(t) je prodloužením řešení ~x = ~x(t) a že řešení ~x = ~x(t) je zúžením řešení x = x(t). Jestliže řešení x = x(t) úlohy (1.1), (1.2) není zúžením žádného jiného řešení této úlohy, řekneme, že x = x(t) je úplným (neprodlužitelným) řešením úlohy (1.1), (1.2). V dalším budeme pod pojmem ,,řešení rozumět úplné řešení. 1.0.6 Příklad x(t) 0, x(t) = t3 , x(t) = 0, t < a (t - a)3 , t a , kde a 0 je libovolné číslo, jsou tři různá úplná řešení počáteční úlohy x = 3 3 x2, x(0) = 0. 1.0.7 Definice Buď g funkce dvou proměnných. Řekneme, že g je obecným řešením rovnice (1.1), jestliže ke každému (t0, x0) G existuje C0 R takové, že x = x(t) = g(t, C0) je řešením úlohy (1.1), (1.2). Řešení úlohy (1.1) (1.2) se nazývá partikulární řešení rovnice (1.1). 1.0.8 Příklad x = Ct je obecným řešením rovnice z příkladu 1.0.2. Rovnice z příkladu 1.0.2 nemá obecné řešení. 1.0.9 Geometrická interpretace Rovnice (1.1) přiřazuje každému bodu z G právě jednu hodnotu x = f(t, x), tedy každému bodu (t0, x0) G lze přiřadit směrový vektor tečny k integrální křivce v bodě (t0, x0), tj. přímky x - x0 = f(t0, x0)(t - t0). Tento vektor má souřadnice 1, f(t0, x0) . To znamená, že rovnice (1.1) definuje na G vektorové pole. Toto pole se nazývá směrové pole rovnice (1.1). Každá integrální křivka rovnice (1.1) je vektorovou čárou směrového pole. Směrové pole tedy poskytuje představu o průběhu řešení rovnice (1.1). Vrstevnice funkce f, (tj. křivky zadané rovnicí f(t, x) = c) se nazývají izokliny rovnice (1.1). Jsou to křivky, na nichž mají vektory ze směrového pole stejný směr. 1.0.10 Definice Buď G Rn+1 množina s neprázdným vnitřkem, f : G Rn . Rovnice x = f(t, x), kde x = (x 1, x 2, . . . , x n) (1.3) se nazývá systém n obyčejných diferenciálních rovnic prvního řádu nebo n-vektorová obyčejná diferenciální rovnice prvního řádu. Počáteční podmínku k rovnici (1.3) lze zadat: x(t0) = x0 = (x0 1, x0 2, . . . , x0 n). (1.4) Pojmy řešení, obecné řešení, partikulární řešení, úplné řešení rovnice (1.3) jsou analogiemi těchto pojmů z jednorozměrného případu. Obecné řešení závisí na n parametrech. 2 1.0.11 Definice Buď G Rn+1 množina s neprázdným vnitřkem, f : G R. Rovnice x(n) = f t, x, x , x , . . . , x(n-1) (1.5) se nazývá obyčejná diferenciální rovnice n-tého řádu rozřešená vzhledem k nejvyšší derivaci. Řešením této rovnice se rozumí n-krát diferencovatelná funkce x : J R, kde J R je interval, která splňuje podmínky t, x(t), x (t), x (t), . . . , x(n-1) (t) G, x(n) (t) = f t, x(t), x (t), x (t), . . . , x(n-1) (t) pro každé t J. Počáteční (Cauchyovu) podmínku pro rovnici (1.5) zadáváme ve tvaru x(t0) = x0, x (t0) = x1 0, x (t0) = x2 0, . . . , x(n-1) (t0) = xn-1 0 , (1.6) kde (t0, x1 0, x2 0, . . . , xn-1 0 ) G. Úplné řešení, obecné řešení, partikulární řešení rovnice (1.5) definujeme analogicky jako u rovnic prvního řádu. Obecné řešení závisí na n parametrech. 1.0.12 Poznámka Řešení počáteční úlohy (1.5), (1.6) je ekvivalentní s řešením počátečního problému pro systém n diferenciálních rovnic prvního řádu: x 1 = x2 x 2 = x3 ... x n-1 = xn x n = f(t, x1, x2, . . . , xn) (1.7) x1(t0) = x0, x2(t0) = x2 0, . . . , xn(t0) = xn-1 0 , (1.8) v tomto smyslu: Je-li x = x(t) řešením úlohy (1.5), (1.6), pak n-tice funkcí x1 = x, x2 = x , x3 = x , . . . , xn = x(n-1) je řešením úlohy (1.7), (1.8) a je-li n-tice funkcí x1 = x1(t), x2 = x2(t), . . . , xn = xn(t) řešením úlohy (1.7), (1.8), pak je funkce x = x(t) = x1(t) řešením úlohy (1.5), (1.6). 3 4 Kapitola 2 Elementární metody řešení 2.1 Rovnice typu x = f(t) Jedná se v podstatě o rovnost, jíž je definována primitivní funkce k dané funkci f. Obecné řešení této rovnice tedy je x(t) = f(t)dt a partikulární řešení splňující počáteční podmínku (1.2) je x(t) = x0 + t t0 f()d; samozřejmě za předpokladu, že příslušná primitivní funkce nebo určitý integrál existují. Příklady na užití rovnice tohoto typu jsou 6.1.1 a 6.1.2. 2.2 Exaktní rovnice x = f(t, x) g(t, x) , přičemž f(t, x) x = g(t, x) t Tuto rovnici lze přepsat na tvar dx dt = f(t, x) g(t, x) , neboli 0 = f(t, x)dt - g(t, x)dx . Za uvedených podmínek je f(t, x)dt - g(t, x)dx totálním diferenciálem nějaké funkce F dvou proměnných (sr. např. Z. Došlá, O. Došlý: Diferenciální počet funkcí více proměnných. MU, Brno 1999, kap. 4), přičemž platí dF(t, x) = 0. Obecné řešení dané rovnice je tedy implicitně zadáno rovností F(t, x) = C, kde C je reálná konstanta. 2.3 Rovnice se separovanými proměnnými x = f(t)g(x) Rovností g(x) = 0 je implicitně zadáno singulární (konstantní) řešení. Poznamenejme, že singulární řešení může, ale nemusí, být zahrnuto v řešení obecném pro nějakou volbu integrační konstanty. Rovností dx g(x) = f(t)dt 5 je implicitně zadáno obecné řešení dané rovnice. Rovností x x0 d g() = t t0 f()d je implicitně zadáno partikulární řešení této rovnice, které splňuje počáteční podmínku (1.2) takovou, že g(x0) = 0. Pokud g(x0) = 0, pak řešením dané rovnice s počáteční podmínkou (1.2) je konstantní funkce x(t) x0. Poznamenejme, že rovnici typu 2.1 lze považovat za zvláštní případ rovnice se separovanými proměnnými pro g(x) 1. 2.4 Homogenní rovnice x = f x t Zavedeme funkci u = u(t) = x(t) t . Pak x(t) = tu(t), x = u + tu . Dosazením do původní rovnice dostaneme u = f(u) - u t , což je rovnice se separovanými proměnnými pro neznámou funkci u. 2.5 Rovnice typu x = f at + bx + c t + x + 1. c = = 0. Pak f at + bx t + x = f a + bx t + x t a daná rovnice je homogenní. 2. c2 + 2 = 0, a = b = k. Zavedeme funkci u = u(t) = at + bx. Pak u = a + bx a tedy x = u - a b . Dosazením do původní rovnice dostaneme u = a + bf u + c ku + , což je rovnice se separovanými proměnnými pro neznámou funkci u. 3. c2 + 2 = 0, a = b. Nechť m a n jsou řešením soustavy lineárních algebraických rovnic am + bn = -c m + n = - . Zavedeme funkce u = u(t) = t - m v = v(t) = x - n. Pak dt = du, dx = dv, at + bx + c = a(u + m) + b(v + n) + c = au + bv + (am + bn) + c = au + bv, t + x + = (u + m) + (v + n) + = u + v + (m + bn) + = u + v Daná rovnice přejde na tvar dv du = f au + bv u + v , což je rovnice typu 1. pro neznámou funkci v = v(u). 6 2.6 Lineární rovnice x = a(t)x + b(t) 1. b(t) 0 (homogenní rovnice) Je to rovnice se separovanými proměnnými. Partikulární řešení počátečního problému (s podmínkou (1.2)) je: x x0 d = t t0 a()d ln x - ln x0 = t t0 a()d x = x0 exp t t0 a()d Obecné řešení homogenní lineární rovnice lze tedy zapsat x = C exp t t0 a()d . 2. b(t) 0 (nehomogenní rovnice) Řešení hledáme ve tvaru x(t) = C(t) exp t t0 a()d (metoda variace konstanty). Pak x = (C (t) + a(t)C(t)) exp t t0 a()d. Dosazením do dané rovnice dostaneme C (t) + a(t)C(t) exp t t0 a()d = a(t)C(t) exp t t0 a()d + b(t) C (t) = b(t) exp t0 t a()d C(t) - C(t0) = t t0 b() exp t0 a()d d Obecné řešení nehomogenní rovnice tedy je x(t) = const + t t0 b() exp t0 a()d d exp t t0 a()d . Partikulární řešení splňující počáteční podmínku (1.2) je x(t) = x0 + t t0 b() exp t0 a()d d exp t t0 a()d . 7 Jsou-li koeficienty konstantní, a(t) A, b(t) B, pak x(t) = x0 + B A eA(t-t0) - B A . Jiný postup při řešení nehomogenní rovnice: x - a(t)x = b(t) / e- R a(t)dt x e- R a(t)dt - a(t)xe- R a(t)dt = b(t)e- R a(t)dt d dt xe- R a(t)dt = b(t)e- R a(t)dt xe- R a(t)dt = b(t)e- R a(t)dt dt x = e R a(t)dt b(t)e- R a(t)dt dt 2.7 Bernoulliova rovnice x = a(t)x + b(t)xr , r R Zavedeme funkci u = u(t) = x(t)1-r . Pak x = u 1 1-r , x = 1 1 - r u 1 1-r -1 u . Dosadíme do dané rovnice: 1 1 - r u r 1-r u = a(t)u 1 1-r + b(t)u r 1-r / (1 - r)u r r-1 u = (1 - r)a(t)u + (1 - r)b(t) . To je lineární rovnice pro neznámou funkci u. Jsou-li koeficienty konstantní, a(t) A, b(t) B, lze použít substituci x = y - B A - 1 r-1 . pak x = - 1 r - 1 y - B A - 1 r-1 -1 y a tedy - 1 r - 1 y - B A - 1 r-1 -1 y = A + B y - B A -1 y - B A - 1 r-1 - 1 r - 1 y = A + B A Ay - B Ay - B A 1 1 - r y = A2 y - AB + AB Ay - B Ay - B A y = (1 - r)Ay, což je lineární homogenní rovnice. Příklad na užití Bernoulliovy rovnice je 6.1.5. 2.8 Rovnice nerozřešená vzhledem k derivaci F(t, x, x ) = 0 (implicitní rovnice) Zavedeme funkci p = p(t) = x (t). 8 1. Rovnice autonomní F(x, x ) = 0 Rovnicí F(x, p) = 0 může být implicitně zadána funkce p = p(x). Rovnici F(x, p(x)) = 0 derivujeme podle proměnné x: F x (x, p) + F p (x, p) dp dx = 0 dp dx = - F x (x, p) F p (x, p) , což je rovnice pro neznámou funkci p nezávisle proměnné x rozřešená vzhledem k derivaci. Je-li p = p(x) řešením poslední rovnice, pak rovnice se separovanými proměnnými x = p(x) je řešením původní rovnice; její obecné řešení je tedy implicitně zadáno rovnicí dx p(x) = t + const . 2. Rovnice nezávisející na neznámé funkci F(t, x ) = 0. Rovnici F(t, p) = 0 derivujeme podle proměnné t: F t (t, p) + F p (t, p) dp dt = 0 dp dt = - F t (t, p) F p (t, p) , což je rovnice pro neznámou funkci p = p(t) rozřešená vzhledem k derivaci. Je-li p = p(t) řešením poslední rovnice, je funkce x = x(t) = p(t)dt obecným řešením dané rovnice. 3. Clairautova rovnice x = tx + g(x ). Rovnici x = tp + g(p) derivujeme podle proměnné t: p = p + t dp dt + g (p) dp dt 0 = t + g (p) dp dt Musí tedy být dp dt = 0 nebo t = -g (p). Z první rovnosti a dané rovnice dostaneme obecné řešení x(t) = ct + g(c), kde c R je libovolná konstanta; z druhé rovnice dostaneme parametrické vyjádření singulárního řešení t = -g (p) x = -pg (p) + g(p) , kde p je parametr. 4. Lagrangeova rovnice x = tf(x ) + g(x ). Rovnici x = tf(p) + g(p) derivujeme podle proměnné t: p = f(p) + tf (p) dp dt + g (p) dp dt p - f(p) = tf (p) + g (p) dp dt 9 Má-li rovnice p - f(p) = 0 řešení p c, pak x(t) = ct + c1 je singulárním řešením dané rovnice. Konstantu c1 určíme dosazením do dané rovnice: ct + c1 = tf(c) + g(c) c1 = t f(c) - c + g(c) a poněvadž f(c) = c, je c1 = g(c). Singulární řešení Lagrangeovy rovnice tedy je x(t) = ct + g(c) , kde c je řešením rovnice c = f(c) (je pevným bodem funkce f). Pro p = f(p) dostaneme dt dp = tf (p) + g (p) p - f(p) , což je lineární rovnice pro neznámou funkci t nezávisle proměnné p. Označíme-li její řešení t = t(p) = (p), pak t = (p) x = f(p)(p) + g(p) je parametrickým vyjádřením obecného řešení Lagrangeovy rovnice. 2.9 Rovnice druhého řádu typu x = f(x) Rovnici vynásobíme 2x : 2x x = 2x f(x) d dt x2 = 2 dx dt f(x) Položíme-li p = x , máme d dt p2 = 2 dx dt f(x) dp2 dx dx dt = 2f(x) dx dt dp2 dx = 2f(x) p2 = 2 f(x)dx . Položíme dále F(x) = 2 f(x)dx a dostaneme p = F(x) dx dt = F(x) , což je rovnice prvního řádu se separovanými proměnnými. 2.10 Rovnice typu F(t, x(k) , x(k+1) , . . . , x(n) ) = 0, k {1, . . ., n - 1} Položíme y = y(t) = x(k) (t) a dostaneme rovnici F t, y, y , y , . . . , y(n-k) = 0 , což je rovnice řádu o k nižšího, než daná rovnice. 10 2.11 Autonomní rovnice typu F x, x , x , . . ., x(n) = 0 Položíme p = p(t) = x (t). Pak x = dp dt = dp dx dx dt = p dp dx x = d dt p dp dx = d dx p dp dx dx dt = dp dx 2 + p d2 p dx2 p . Postupujeme-li tak dále, vidíme, že x(k) = fk p, dp dx , . . . , dk-1 p dxk-1 pro každé k N. (fk je nějaká funkce k proměnných.) Dosazením do původní rovnice tedy dosta- neme F x, p, p dp dx , f3 p, dp dx , d2 p dx2 , . . . , fn p, dp dx , . . . , dn-1 p dxn-1 = 0 , neboli G x, p, dp dx , . . . , dn-1 p dxn-1 = 0 , což je rovnice řádu o jedna nižšího, než daná rovnice. 2.12 Rovnice homogenní v x, x , x , . . . , x(n) Nechť F je funkce n + 1 proměnných splňující podmínku F(z0, cz1, cz2, . . . , czn) = cF(z0, z1, z2, . . . , zn) (*) pro každé c R a každé (z0, z1, z2, . . . , zn) Dom F. Řešení rovnice F t, x, x , x . . . , x(n) = 0 lze hledat ve tvaru x(t) = e R y(t)dt , kde y = y(t) je nová neznámá funkce. Je totiž x = ye R y(t)dt x = y e R y(t)dt + y2 e R y(t)dt = (y + y2 )e R y(t)dt x = (y + 2yy )e R y(t)dt + (y + y2 )y e R y(t)dt = (y + 3yy + y3 )e R y(t)dt ... Dosadíme-li z těchto rovnic do dané rovnice, vypadne vzhledem k podmínce (*) faktor e R y(t)dt a dostaneme rovnici řádu o 1 nižšího, než byla daná rovnice. 2.13 Ekvidimensionální rovnice Řekneme, že rovnice F t, x, x , . . . , x(n) = 0 je ekvidimensionální v nezávisle proměnné, jestliže změna měřítka nezávisle proměnné t at pro každé a R \ {0} nezmění její tvar. Transformace t = e převede danou rovnici na rovnici autonomní (typ 2.11). 11 2.14 Cvičení Řešte rovnice (Cauchyovy úlohy) 1) 2t(2x - 3)dt + (t2 + 1)dx = 0 2) dx dt = et-x 3) tex dx + t2 + 1 x dt = 0 4) 1 + t2 dx + x2 - 1 dt = 0 5) t2 dx + (x2 - tx)dt = 0 6) dt dx = t + x x - t 7) t sin x t - x cos x t dt + t cos x t dx = 0 8) 2 dx dt - x = et/2 9) tdx + xdt = sin tdt 10) (t - 1)3 x + 4(t - 1)2 x = t + 1 11) e2x dt + 2(te2x - x)dx = 0 12) (x2 + 1)dt + (2tx + 1)dx = 0 13) (t + x)dt + (t + x2 )dx = 0 14) tdx - xdt + t3 dt = 0 15) (t2 + t - x)dt + tdx = 0 16) (cos t + x cos t)dt + dx = 0; x(/2) = 0 17) x + 2x = t; x(0) = 2 18) (t + 2x)dt + (x + 2t)dx = 0; x(1) = 1 19) Určete konstanty a, b, c tak, aby rovnice (at2 + bx2 )dt + ctx dx = 0 byla exaktní a vyřešte ji. 20) Řešte počáteční úlohu pro implicitní rovnici druhého řádu x x = t (x ) 2 , x(1) = 1, x (1) = -1. Výsledky: 1)x = 3 2 + C (t2 + 1)2 2)ex = et + C 3)ex (x - 1) + t2 2 + ln |t| = C 4)(x + x2 - 1 )(t + t2 + 1 ) = C 5)x = t ln |t| + C 6) 1 2 ln(t2 +x2 )+arctg x t = C 7)x = t arcsin C t 8)x = t + C 2 et/2 9)x = C - cos t t 10)x = t3 - 3t + C 3(t - 1)4 11)t = x2 + C 2 e-2x 12)t = C - x x2 + 1 13) t2 2 + tx + x3 3 = C 14)x = Ct - t3 2 15)x = Ct - t2 - t ln |t| 16)x = e1-sin t - 1 17)x = t 2 + 9 4 e-2t - 1 4 18)t2 + 4tx + x2 = 6 19)c = 2b; at3 3 + btx2 = C 20)x = exp 2 arctg 1 - t 2 12 Kapitola 3 Obecné vlastnosti obyčejných diferenciálních rovnic 3.1 Vektorové a maticové funkce 3.1.1 Normy vektorů a matic Normu vektoru x = x1 x2 ... xn definujeme předpisem |x| = n i=1 |xi|. Normu matice A = a11 a12 . . . a1n a21 a22 . . . a2n ... ... ... ... an1 an2 . . . ann definujeme předpisem |A| = n i=1 n j=1 |aij|. Na množině vektorů zavádíme metriku (x, y) = |x - y|, na množině matic zavádíme metriku (A, B) = |A - B|. Jedná se o součtovou neboli taxíkářskou metrikou, sr. Z. Došlá, O. Došlý: Metrické prostory. MU, Brno 1996, I.1.1.2.iii. * Platí: |A x| |A| |x|. Této vlastnosti se říká, že maticová norma je souhlasná s vektorovou normou. D.: Pro libovolné i, k {1, 2, . . ., n} je |aik| n j=1 |aij|. Odtud plyne |Ax| = n k=1 aikxk = n i=1 n k=1 aikxk n i=1 n k=1 |aik| |xk| = = n k=1 |xk| n i=1 |aik| n k=1 |xk| n i=1 n j=1 |aij| = n k=1 |xk| n i=1 n j=1 |aij| . 13 3.1.2 Spojitost, derivace a integrál vektorových a maticových funkcí x(t) = x1(t) x2(t) ... xn(t) = xi(t) , A(t) = a11(t) a12(t) . . . a1n(t) a21(t) a22(t) . . . a2n(t) ... ... ... ... an1(t) an2(t) . . . ann(t) = aij(t) Vektorová funkce x (resp. maticová funkce A) je spojitá v bodě t0 svého definičního oboru, jestliže ke každému kladnému číslu existuje kladné číslo tak, že pro všechna t z definičního oboru funkce x z nerovnosti |t - t0| < plyne nerovnost |x(t) - x(t0)| < (resp. |A(t) - A(t0)| < ). Vektorová (resp. maticová) funkce je spojitá právě tehdy, když všechny její složky jsou spojité. lim tt0 x(t) = lim tt0 xi(t) , dx dt (t) = x (t) = x i(t) , t t0 x(s)ds = t t0 xi(s)ds lim tt0 A(t) = lim tt0 aij(t) , dA dt (t) = A (t) = a ij(t) , t t0 A(s)ds = t t0 aij(s)ds 3.2 Existence a jednoznačnost řešení systému ODR Budeme se zabývat úlohou (1.3), (1.4). 3.2.1 Lemma Buď f spojitá na G Rn+1 . Funkce x = x(t) je řešením úlohy (1.3), (1.4) na intervalu J právě tehdy, když pro každé t J je t, x(t) G a x(t) = x0 + t t0 f s, x(s) ds . (3.1) D.: ,, Nechť x = x(t) je řešením úlohy (1.3), (1.4) na J. Pak d ds x(s) = f s, x(s) na J. Integrací této rovnosti podle s v mezích [t0, t] dostaneme: [x(s)]t t0 = t t0 f s, x(s) ds x(t) - x(t0) = t t0 f s, x(s) ds a vzhledem k (1.4) funkce x = x(t) splňuje (3.1). ,, Nechť funkce x = x(t) splňuje (3.1). Pak x(t0) = x0 + t0 t0 f s, x(s) ds = x0 , tedy je splněna podmínka (1.4). Derivováním (3.1) podle t dostaneme (1.3). 14 Nechť C1 (J) je množina (vektorových) funkcí x = x(t) diferencovatelných na uzavřeném intervalu J takových, že x(t0) = x0 . Na této množině zavedeme metriku (x, y) = max{|x(t) - y(t)| : t J} (metrika stejnoměrné konvergence). Prostor C1 (J), je úplný (sr. Z. Došlá, O. Došlý: Metrické prostory. MU, Brno 1996, I.1.1.2.iv a III.1.3.6.i). Dále definujme zobrazení F : C1 (J) C1 (J) předpisem: F(x)(t) = x0 + t t0 f s, x(s) ds . Řešení úlohy (1.3), (1.4), tedy funkce, která splňuje (3.1), je zřejmě pevným bodem zobrazení F. Podaří-li se tedy ukázat, že F je kontrakce úplného metrického prostoru C1 (J), , z Banachovy věty vyplyne, že existuje jediný pevný bod tohoto zobrazení, tedy že existuje jediné diferencovatelné řešení úlohy (1.3), (1.4) (sr. Z. Došlá, O. Došlý: Metrické prostory. MU, Brno 1996, IV.2. a V.1.). 3.2.2 Věta (Picard [1856­1941]­Lindelöf [1870­1946]) Buďte a, b R, a, b > 0, t0 R, x0 Rn . Označme ~J = [t0, t0 +a], D = {x Rn : |x-x0 | b}, m = max{|f(t, x)| : (t, x) ~J × D}, = min a, b m . Nechť funkce f : ~J × D Rn je spojitá a vzhledem k x Lipschitzovská (tj. existuje L R tak, že platí |f(t, x) - f(t, y)| L|x - y| pro všechna t ~J, x, y D). Pak existuje právě jedno řešení počátečního problému (1.3), (1.4) definované na intervalu J = [t0, t0 + ]. Toto řešení je (stejnoměrnou) limitou posloupnosti funkcí {xn(t)} n=1; tato posloupnost je definována rekurentně vztahy x1 = x0 xk+1 = x0 + t t0 f s, xk(s) ds, k = 1, 2, 3, . . . D.: Funkce y(t) = F(x)(t) je diferencovatelná, y (t) = f t, x(t) (sr. V Novák: Integrální počet v R., MU, Brno 2001, 2.4, věta 4.2). To znamená, že zobrazení F definované před větou zobrazuje C1 (J) do sebe. Buď K > L. Na C1 (J) zavedeme metriku vztahem (x, y) = max e-K(t-t0) |x(t) - y(t)| : t J . Tato metrika je na C1 (J) ekvivalentní s metrikou stejnoměrné konvergence , neboť e-K (x, y) (x, y) (x, y) . Prostor C1 (J), je tedy úplný. Položme P = {x C1 (J) : |x(t) - x0 | b pro každé t J}. Poněvadž P je uzavřená podmnožina C1 (J), je prostor (P, ) úplný (sr. Z. Došlá, O. Došlý: Metrické prostory. MU, Brno 1996, III.1.3.3). Zobrazení F zobrazuje množinu P do sebe, neboť pro každou funkci x P platí |F(x)(t) - x0 | = t t0 f(s, x(s))ds t t0 f s, x(s) ds (t - t0)m m b . 15 Ukážeme, že F je kontrakcí prostoru (P, ): Pro každé t J platí F(x), F(y) e-K(t-t0) |F(x)(t) - F(y)(t)| = = e-K(t-t0) t t0 f s, x(s) ds - t t0 f s, y(s) ds e-K(t-t0) t t0 f s, x(s) - f s, y(s) ds e-K(t-t0) t t0 L|x(s) - y(s)|ds = = L t t0 e-K(t-s) e-K(s-t0) |x(s) - y(s)|ds L t t0 e-K(t-s) (x, y)ds = = L (x, y) 1 K e-K(t-s) t s=t0 = = L K (x, y) 1 - e-K(t-t0) L K (x, y) . Poněvadž L < K, je L K < 1, což znamená, že F je kontrakce. 3.2.3 Poznámky 1. Posloupnost funkcí zavedená v 3.2.2 se nazývá Picardova posloupnost postupných aproximací. 2. Analogické tvrzení platí, nahradíme-li v 3.2.2 interval ~J = [t0, t0 + a] intervalem [t0 - a, t0] nebo intervalem [t0 - a, t0 + a]. 3. Má-li funkce f(t, x) = f1(t, x1, x2, . . . , xn) f2(t, x1, x2, . . . , xn) ... fn(t, x1, x2, . . . , xn) ohraničené parciální derivace všech složek podle každé z proměnných x1, x2, . . . , xn na množině ~J × D (zavedené v 3.2.2), pak jsou předpoklady Picardovy-Lindelöfovy věty splněny. D.: Množina ~J ×D jakožto uzavřená a ohraničená podmnožina prostoru Rn+1 je kompaktní (sr. Z. Došlá, O. Došlý: Metrické prostory. MU, Brno 1996, III.3.3.16). Z ohraničenosti parciálních derivací funkce f plyne existence čísla M = max fi(t, x) xj : i = 1, 2, . . . , n, j = 1, 2, . . ., n, (t, x) ~J × D . Podle věty o střední hodnotě pro funkce více proměnných (sr. Z. Došlá, O. Došlý: Diferenciální počet funkcí více proměnných. MU, Brno 1999, 3.4) pro všechna t ~J 16 a x = (x1, x2, . . . , xn), y = (y1, y2, . . . , yn) D existují čísla k ležící mezi xk a yk, k = 1, 2, . . ., n taková, že |f(t, x) - f(t, y)| = n i=1 |fi(t, x) - fi(t, y)| = = n i=1 n k=1 fi xk (t, x1, x2, . . . , xk-1, k, yk+1, . . . , yn)(xk - yk) n i=1 n k=1 fi xk (t, x1, x2, . . . , xk-1, k, yk+1, . . . , yn) |xk - yk| n i=1 n k=1 M|xk - yk| = n i=1 M|x - y| = nM|x - y| , takže funkce f je vzhledem k x Lipschitzovská s konstantou nM. 3.2.4 Důsledky 1. Má-li (vektorová) funkce f(t, x) = f1(t, x1, x2, . . . , xn) f2(t, x1, x2, . . . , xn) ... fn(t, x1, x2, . . . , xn) ohraničené parciální derivace všech složek podle každé z proměnných x1, x2, . . . , xn v jistém okolí bodu (t0, x0 ), pak počáteční problém (1.3), (1.4) má v okolí t0 jediné řešení. 2. Má-li (skalární) funkce f(t, x1, x2, . . . , xn) v jistém okolí bodu (t0, x0, x1 0, x2 0, . . . , xn 0 ) ohraničené parciální derivace podle každé z proměnných x1, x2, . . . , xn, pak počáteční problém (1.5), (1.6) má v okolí t0 jediné řešení. 3.2.5 Věta (Peano [1890]) Buďte a, b R, a, b > 0, t0 R, x0 Rn . Označme ~J = [t0, t0 +a], D = {x Rn : |x-x0 | b}, m = max{|f(t, x)| : (t, x) ~J × D}, = min a, b m . Nechť funkce f : ~J × D Rn je spojitá. Pak existuje alespoň jedno řešení počátečního problému (1.3), (1.4) definované na intervalu J = [t0, t0 + ]. D.: Viz Kalas, Ráb: Obyčejné diferenciální rovnice, str. 67­70. 3.3 Globální vlastnosti řešení systému ODR 3.3.1 Věta (o existenci úplného řešení) Nechť funkce f : G Rn je spojitá na otevřené množině G Rn+1 . Je-li x = x(t) řešení rovnice (1.3), pak je toto řešení buď úplné, nebo existuje úplné řešení y = y(t), které je prodloužením řešení x. D.: Viz J. Kalas, M. Ráb: Obyčejné diferenciální rovnice. MU, Brno 2001, str. 73­76. Důkaz využívá věty 3.2.5. 3.3.2 Definice Buď G Rn+1 . Řekneme, že funkce f : G Rn je lokálně lipschitzovská v G vzhledem k x, jestliže ke každému (, a) G existuje okolí O,a G a číslo L,a R tak, že pro všechna (t, x), (t, y) O,a platí |f(t, x) - f(t, y)| L,a|x - y|. 17 3.3.3 Věta (o globální jednoznačnosti) Nechť funkce f : G Rn je spojitá a lokálně lipschitzovská v G vzhledem k x a nechť funkce x = x(t), y = y(t) jsou dvě řešení rovnice (1.3). Jestliže existuje t0 takové, že x(t0) = y(t0), pak x(t) = y(t) pro všechna t, v nichž jsou řešení x, y definována. D.: Připusťme, že existuje b > t0 takové, že x(b) = y(b). Označme c = inf{t : x(t) = y(t)}. Funkce x, y jsou spojité (poněvadž jsou diferencovatelné). Ukážeme, že x(c) = y(c): Připusťme, že x(c) = y(c). Položme = |x(c) - y(c)|. K 4 > 0 existuje > 0 tak, že pro každé t (c-, c) je |y(t)-y(c)| < 4 a |x(t)-x(c)| < 4 . Poněvadž pro t (c - , c) je x(t) = y(t), platí pro t (c - , c) nerovnost = |x(c) - y(c)| = |x(c) - x(t) + y(t) - y(c)| |x(c) - x(t)| + |y(t) - y(c)| 4 + 4 = 2 , což je spor, neboť > 0, a tedy x(c) = y(c). Podle 3.2.2 nyní existuje takové, že pro t [c, c + ] je x(t) = y(t), což je spor. Analogicky vyloučíme možnost existence b < t0 takového, že x(b) = y(b). 3.3.4 Definice Buď x = x(t) úplné řešení rovnice (1.3) definované na intervalu (S, T ), kde - S < T . Řekneme, že Rn je -limitní bod řešení x, jestliže existuje posloupnost {tk} k=1 taková, že tk < T pro všechna k N, lim k tk = T a lim k x(tk) = . Řekneme, že Rn je -limitní bod řešení x, jestliže existuje posloupnost {tk} k=1 taková, že tk > S pro všechna k N, lim k tk = S a lim k x(tk) = . Řekneme, že Rn je limitní bod řešení x, jestliže je -limitním bodem nebo -limitním bodem. Množina všech -limitních bodů řešení x se nazývá -limitní množina řešení x, množina všech -limitních bodů řešení x se nazývá -limitní množina řešení x, množina všech limitních bodů řešení x se nazývá limitní množina řešení x. 3.3.5 Příklady 1. x = eat je úplné řešení rovnice x = ax definované na intervalu (-, ). Je-li a > 0, je 0 limitním bodem tohoto řešení a -limitní body toto řešení nemá. Je-li a < 0, je 0 -limitním bodem tohoto řešení a -limitní body toto řešení nemá. Je-li a = 0, je 1 - i -limitním bodem tohoto řešení. 2. x = sin 1 t je úplné řešení rovnice x = cos 1 t t2 definované na intervalu (-, 0). Interval [-1, 1] je -limitní množinou tohoto řešení. 3. x y = cos tg t sin tg t je úplné řešení soustavy rovnic x y = - y cos2 t x cos2 t definované na intervalu - 2 , 2 . Limitní množina tohoto řešení je {(x, y) R2 : x2 + y2 = 1}. 3.3.6 Věta Nechť funkce f : G Rn je spojitá na otevřené množině G Rn+1 a x = x(t) je úplné řešení rovnice (1.3) definované na intervalu (S, T ). Pak platí: T = nebo každý -limitní bod řešení x leží na hranici G. S = - nebo každý -limitní bod řešení x leží na hranici G. 18 D.: Buď x = x(t) úplné řešení rovnice (1.3) definované na intervalu (S, T ), T < a buď jeho -limitní bod. Kdyby (T, ) G, pak by existovalo okolí OT, bodu (T, ) takové, že OT, G, neboť G je otevřená. Podle 3.2.5 by existovalo řešení y = y(t) rovnice (1.3) s počáteční podmínkou y(T ) = definované na [T, T + ], kde je vhodné (malé) číslo. Funkce z = z(t) = x(t), S < t < T y(t), T + by byla řešením rovnice (1.3), které by bylo prodloužením řešení x, což by byl spor s úplností řešení x. Pro -limitní bod se důkaz provede analogicky s využitím ,,levostranné varianty věty 3.2.5. 3.3.7 Důsledek Nechť J = [t0, ), D = {x Rn : |x| < a}, kde t0 R a 0 < a a nechť funkce f : J × D Rn je spojitá. Jestliže existuje spojitá funkce g : J R taková, že úplné řešení x = x(t) rovnice (1.3) definované na (S, T ) splňuje pro každé t [t0, T ) podmínku |x(t)| g(t) < a, pak T = . 3.3.8 Důsledek Nechť J = [t0, ) a funkce f : J × Rn Rn je spojitá. Jestliže existuje m R takové, že pro každé (t, x) J × Rn Rn platí |f(t, x)| m, pak každé úplné řešení rovnice (1.3) je definováno pro všechna t t0. D.: Buď x = x(t) úplné řešení rovnice (1.3). Podle 3.2.1 je x(t) = x(t0) + t t0 f s, x(s) ds . Pro každé t, pro něž je x(t) definováno, platí |x(t)| = x(t0) + t t0 f s, x(s) ds |x(t0)| + t t0 |f(s, x(s))| ds |x(t0)| + m(t - t0) . Tvrzení tedy plyne z 3.3.7, položíme-li g(t) = |x(t0)| + m(t - t0). Toto tvrzení umožňuje rozhodnout, zda lze každé řešení rovnice (1.3) prodloužit do nekonečna, aniž bychom toto řešení znali. 3.3.9 Věta (o spojité závislosti řešení na počátečních podmínkách a pa- rametrech) Buď otevřená množina v R1+n+m a nechť funkce g : Rn je taková, že pro všechna R, Rn , Rm splňující podmínku (, , ) , má počáteční problém x = g(t, x, ), x() = právě jedno úplné řešení x = x(t) = x(t; , , ). Pak toto řešení, chápané jako zobrazení R1+1+n+m Rn , je spojité. 19 D.: Viz J. Kalas, M. Ráb: Obyčejné diferenciální rovnice. MU, Brno 2001, str. 82­83. Tato věta říká, že změní-li se málo funkce f v rovnici (1.3) a málo se změní počáteční podmínka (1.4), pak se řešení nového -- změněného -- problému liší na konečném intervalu málo od řešení původního problému. 3.4 Odhady řešení 3.4.1 Definice Řešení x = x (t) úlohy (1.1), (1.2) se nazývá maximální řešení, jestliže pro každé řešení x = x(t) této úlohy platí x(t) x (t) pro všechna t, v nichž jsou obě řešení definována. Řešení x = x(t) úlohy (1.1), (1.2) se nazývá minimální řešení, jestliže pro každé řešení x = x(t) této úlohy platí x(t) x(t) pro všechna t, v nichž jsou obě řešení definována. 3.4.2 Věta (srovnávací) Nechť t0 R, 0 < a , J = [t0, ), G = {(t, x) Rn+1 : t J, |x| < a}. Nechť dále f : G Rn je spojitá funkce a g : J×[0, a) [0, ) je spojitá funkce taková, že |f(t, x)| g(t, |x|) pro (t, x) G. Buď u0 |x0 | a u = u (t) maximální řešení úlohy u = g(t, u), u(t0) = u0 na intervalu J. Pak každé úplné řešení úlohy (1.3), (1.4) je definováno pro všechna t J a platí |x(t)| u (t) pro t J . D.: Viz J. Kalas, M. Ráb: Obyčejné diferenciální rovnice. MU, Brno 2001, str. 91. 3.4.3 Důsledek Nechť symboly t0, a, J, G mají stejný význam jako v 3.4.2. Nechť funkce f : G Rn je spojitá a nechť existuje spojitá funkce : J [0, ) taková, že |f(t, x)| (t)|x| pro (t, x) G. Pak pro každé x0 Rn takové, že |x0 | exp t t0 ()d < a pro všechna t J, jsou úplná řešení úlohy (1.3), (1.4) definována na celém intervalu J a platí |x(t)| |x0 | exp t t0 ()d pro t J . D.: Plyne z 3.4.2 volbou g(t, u) = (t)u, u0 = |x0 |. Jediné úplné (tedy maximální) řešení úlohy u = (t)u, u(t0) = u0 je podle 2.6 u(t) = u0 exp t t0 ()d . 20 Kapitola 4 Lineární rovnice 4.1 Systémy lineárních ODR Systém lineárních obyčejných diferenciálních rovnic je tvaru x 1(t) = a11(t)x1 + a12(t)x2 + . . . + a1n(t)xn + b1(t) x 2(t) = a21(t)x1 + a22(t)x2 + . . . + a2n(t)xn + b2(t) ... ... ... ... ... ... x n(t) = an1(t)x1 + an2(t)x2 + . . . + ann(t)xn + bn(t) Při označení A(t) = a11(t) a12(t) . . . a1n(t) a21(t) a22(t) . . . a2n(t) ... ... ... ... an1(t) an2(t) . . . ann(t) , b(t) = b1(t) b2(t) ... bn(t) , x(t) = x1(t) x2(t) ... xn(t) lze tento systém zapsat v maticovém tvaru x = A(t) x + b(t) . (4.1) Tuto rovnici nazýváme nehomogenní lineární rovnice. Spolu s rovnicí (4.1) budeme uvažovat počáteční podmínku x(t0) = x0 . (4.2) Rovnici x = A(t) x (4.3) nazýváme lineární homogenní rovnicí přidruženou k rovnici (4.1). 4.1.1 Věta o existenci a jednoznačnosti řešení Nechť maticová funkce A = A(t) a vektorová funkce b = b(t) jsou spojité na intervalu J R. Pak má počáteční problém (4.1), (4.2) právě jedno řešení, které existuje na celém intervalu J. D.: Vzhledem k 3.2.2 a 3.3.3 stačí ukázat, že ke každému J existuje okolí O J takové, že funkce A(t)x + b(t) je Lipschitzovská vzhledem k x na Rn × O . Je-li vnitřní bod intervalu J, existuje a R, a > 0 takové, že [ - a, + a] J. 21 Položme L = max{|A(t)| : t [ - a, + a]} (toto maximum existuje podle druhé Weierstrassovy věty), O = ( - a, + a). Pak pro každé t O a každé dva vektory x, y Rn platí podle tvrzení v 3.1.1 A(t)x + b(t) - A(t)y + b(t) = |A(t)x-A(t)y| = |A(t)(x-y)| |A(t)| |x-y| L |x-y| . Je-li pravý krajní bod intervalu J, položíme L = max{|A(t)| : - a t }, O = ( - a, ] a provedeme analogickou úvahu. Pro levý krajní bod intervalu J provedeme důkaz podobně. 4.1.2 Poznámka Řešení problému (4.1), (4.2) lze hledat jako limitu Picardovy posloupnosti postupných aproximací. Zejména pro A(t) = A = (aij) (konstantní matice) a b(t) 0 je (označíme-li E jednotkovou matici): x0 (t) = x0 = (x0 i ) x1 (t) = x0 + t t0 A x0 ds = x0 i + t t0 n j=1 aijx0 j ds = x0 i + n j=1 aijx0 j (t - t0) x2 (t) = x0 + t t0 n k=1 aik x0 k + n j=1 akjx0 j (s - t0) ds = = x0 + n k=1 aikx0 k(t - t0) + n k=1 n j=1 aikakjx0 j (t - t0)2 2 = = x0 + A x0 (t - t0) + (A A) x0 (t - t0)2 2 = E + A(t - t0) + A2 2 (t - t0)2 x0 ... xk (t) = E + A(t - t0) + A2 2! (t - t0)2 + + Ak k! (t - t0)k x0 ... Pokud bychom A považovali za konstantu a E za jedničku (tj. pro n = 1), je výraz v závorce k-tým částečným součtem Taylorovy řady funkce eA(t-t0) . Řešení úlohy x = A x, x(t0) = x0 lze proto zapsat ve tvaru x(t) = eA(t-t0) x0 . Matice eA(t-t0) je dána nekonečnou řadou eA(t-t0) = k=0 (t - t0)k k! Ak . Lze ukázat, že tato řada konverguje pro každé t R a tato konvergence je stejnoměrná. 4.1.3 Věta (princip superpozice) Jsou-li x = x(t), y = y(t) řešení rovnice (4.3), pak také c1x(t)+c2y(t), kde c1, c2 jsou konstanty, je řešením rovnice (4.3). D.: d dt (c1x+c2y) = c1x +c2y = c1A(t)x+c2A(t)y = A(t)(c1x)+A(t)(c2y) = A(t)(c1x+c2y) 22 4.1.4 Věta Je-li maticová funkce A = A(t) = aij(t) n i,j=1 spojitá na intervalu J, pak množina všech řešení rovnice (4.3) tvoří n-rozměrný vektorový prostor. D.: x(t) o je řešením rovnice (4.3). Podle 4.1.3 je libovolná lineární kombinace řešení rovnice (4.3) také řešením této rovnice. Buďte x1 , x2 , . . . , xn báze n-rozměrného vektorového prostoru a yk = yk (t) řešení rovnice (4.3) s počátečními podmínkami yk (t0) = xk , k = 1, 2, . . ., n. Vektory y1 (t), y2 (t), . . . , yn (t) jsou lineárně nezávislé pro každé t J: Kdyby existovalo t1 J a konstanty c1, c2, . . . , cn, ne všechny rovny nule, takové, že c1y1 (t1) + c2y2 (t1) + + cnyn (t1) = o, pak by podle 4.1.3 také funkce z = z(t) = c1y1 (t) + c2y2 (t) + + cnyn (t) byla řešením rovnice (4.3) s počáteční podmínkou z(t1) = o. Poněvadž ale x(t) o je řešením rovnice (4.3), vzhledem k jednoznačnosti řešení (sr. 4.1.1) by z(t) o. Zejména tedy c1y1 (t0) + c2y2 (t0) + + cnyn (t0) = c1x1 + c2x2 + + cnxn = z(t0) = o, což by byl spor s lineární nezávislostí vektorů x1 , x2 , . . . , xn . Dimense prostoru všech řešení je tedy alespoň n. Buď ~y = ~y(t) libovolné řešení rovnice (4.3). Pak ~y(t0) = 1x1 + 2x2 + + nxn pro vhodné konstanty 1, 2, . . . , n, neboť x1 , x2 , . . . , xn je báze. To znamená, že ~y(t0) = 1y1 (t0) + 2y2 (t0) + + nyn (t0). Podle principu superpozice je w(t) = ~y(t) - 1y1 (t) - 2y2 (t) - - nyn (t) také řešením rovnice (4.3); přitom w(t0) = o. Vzhledem k jednoznačnosti řešení rovnice (4.3) je w(t) = o pro všechna t J a tedy ~y(t) = 1y1 (t) + 2y2 (t) + + nyn (t) pro všechna t J. To znamená, že funkce y1 = y1 (t), y2 = y2 (t), . . . , yn = yn (t) tvoří bázi prostoru všech řešení rovnice (4.3). 4.1.5 Definice Libovolná báze prostoru všech řešení rovnice (4.3) se nazývá fundamentální systém řešení rovnice (4.3). Nechť y1 = y1 (t), y2 = y2 (t), . . . , yn = yn (t) je fundamentální systém řešení rovnice (4.3). Obecné řešení této rovnice je y = y(t) = c1y1 (t) + c2y2 (t) + + cnyn (t) . 23 Označme Y = Y(t) = y1 1(t) y2 1(t) . . . yn 1 (t) y1 2(t) y2 2(t) . . . yn 2 (t) ... ... ... ... y1 n(t) y2 n(t) . . . yn n(t) , c = c1 c2 ... cn . Matice Y = Y(t) se nazývá fundamentální matice řešení systému (4.3). Tato matice je v důsledku lineární nezávislosti sloupců regulární, det Y(t) = 0 pro každé t J. Obecné řešení rovnice (4.3) lze zapsat y = y(t) = Y(t) c . Partikulární řešení počátečního problému (4.1), (4.2) je y = y(t) = Y(t)c0 , kde c0 = Y(t0)-1 x0 . Tedy y = y(t) = Y(t) Y(t0)-1 x0 . Pro fundamentální matici řešení systému (4.3) Y = Y(t) zřejmě platí Y (t) = A(t) Y(t). 4.1.6 Věta Obecné řešení rovnice (systému) (4.1) je součtem obecného řešení přidružené homogenní rovnice (4.3) a nějakého partikulárního řešení rovnice (4.1): x(t) = Y(t) c + ~x(t) , kde Y(t) je fundamentální matice řešení rovnice (4.3) a ~x(t) je libovolné řešení rovnice (4.1). D.: x(t) = Y(t) c + ~x(t) je řešením rovnice (4.1), neboť x = Y c + ~x = A Y c + A ~x + b = A (Y c + ~x) + b. Přímým výpočtem ověříme, že x(t) = Y(t) c0 + ~x(t), kde c0 = Y(t0)-1 x0 - ~x(t0) je řešením problému (4.1), (4.2). 4.1.7 Nalezení partikulárního řešení rovnice (4.1) -- metoda variace konstant Řešení rovnice (4.1) hledáme ve tvaru ~x = ~x(t) = Y(t) c(t), kde c = c(t) je nějaká vektorová funkce. Pak ~x = Y c + Y c = A Y c + Y c ~x = A ~x + b = A Y c + b Odtud Y(t) c (t) = b(t) c (t) = Y(t)-1 b(t) c(t) = + t t0 Y(s)-1 b(s)ds , kde je konstantní vektor. Obecným řešením rovnice (4.1) je tedy x(t) = Y(t) + t t0 Y(t) Y(s)-1 b(s)ds . 24 Aby toto řešení splňovalo počáteční podmínku (4.2), musí platit x0 = Y(t0), tj. = Y(t0)-1 x0 . Řešení počátečního problému (4.1), (4.2) je tedy x(t) = Y(t) Y(t0)-1 x0 + t t0 Y(t) Y(s)-1 b(s)ds . 4.1.8 Řešení lineární homogenní rovnice s konstantní maticí Řešení rovnice x = A x (4.4) budeme hledat ve tvaru x(t) = et , kde Rn je konstantní vektor. Hodnota musí splňovat rovnosti x (t) = et , A x = A et , takže vzhledem k tomu, že et = 0, musí platit A = . To znamená, že je vlastní hodnotou matice A a je příslušný vlastní vektor. * Jsou-li 1, 2 různé vlastní hodnoty matice A a 1, 2 jsou příslušné vlastní vektory, pak 1e1t , 2e2t jsou lineárně nezávislá řešení rovnice (4.4). D.: Tvrzení plyne z předchozího výpočtu a faktu, že vlastní vektory příslušné k různým vlastním hodnotám jsou lineárně nezávislé. * Má-li rovnice (4.4) komplexní řešení x(t) = (t) + i(t), kde a jsou reálné vektorové funkce, a řešení x je lineárně nezávislé na libovolném nenulovém reálném řešení této rovnice, pak a jsou lineárně nezávislými reálnými řešeními rovnice (4.4). D.: Platí (t) + i (t) = (t) + i(t) = A (t) + i(t) = A (t) + iA (t). Porovnáním reálných a imaginárních částí této rovnosti dostaneme, že funkce a jsou řešeními rovnice (4.4). Kdyby vektorové funkce a byly lineárně závislé, tj. = c, pak by + i = (c+i) a řešení +i by bylo násobkem reálného nenulového řešení . To by byl spor s předpokladem tvrzení. * Je-li vlastní hodnota matice A, příslušný vlastní vektor, přičemž je k-násobným kořenem charakteristického polynomu det(A - E), pak funkce et , 1,0et + 1,1tet , . . . , k-1,0et + k-1,1tet + + k-1,k-1tk-1 et jsou pro vhodné vektory 1,0, 1,1, . . . , k-1,0, k-1,1, . . . , k-1,k-1 řešením rovnice (4.4). D.: Důkaz ukážeme pro k = 2. V případě vyšší násobnosti kořene charakteristické rovnice lze postupovat analogicky. Nechť je dvojnásobný kořen charakteristického polynomu. Buď B = B(s) diferencovatelná (a tedy spojitá) maticová funkce definovaná na okolí nuly taková, že B(0) = A, je pro každé s z definičního oboru funkce B jednoduchým kořenem charakteristického polynomu matice B(s) a existuje jednoduchý kořen (s) charakteristického polynomu, pro nějž platí lim s0 (s) = . (Matici A nepatrně porušíme tak, aby se dvojnásobný kořen rozdělil na dva různé jednoduché.) Označme (s) (resp. (s)) vlastní vektor matice B(s) příslušný k vlastní hodnotě (s) 25 (resp. ). Z diferencovatelnosti funkce B plyne podle věty o diferencovatelnosti implicitně zadané funkce (sr. Z. Došlá, O. Došlý: Diferenciální počet funkcí více proměnných. MU, Brno 1999, 8.1) diferencovatelnost funkcí a , zejména tedy existence limit lim s0 (s) = ,0, lim s0 (s) = ,0. Rovnice x = B(s) x má podle předchozí úvahy řešení (s)et a (s)e(s)t a podle principu superpozice 4.1.3 také řešení (s)et - (s)e(s)t - (s) . Poněvadž lim s0 B(s) = A, plyne z věty o spojité závislosti řešení na parametrech 3.3.9, že rovnice (4.4) má řešení lim s0 (s)et - (s)e(s)t - (s) = lim s0 (s)et - (s)e(s)t - (s)t (s)e(s)t -(s) = = ,0 - ,0 -(0) + (0)t et ; při výpočtu limity bylo využito de l'Hospitalovo pravidlo. Označíme-li nyní 1,0 = ,0 - ,0 (0) , 1,1 = (0), dostaneme tvrzení. 4.1.9 Převedení systému lineárních diferenciálních rovnic na rovnici vyššího řádu Uvažujme systém dvou rovnic x = a(t)x+b(t)y+c(t), y = (t)x+(t)y+(t); (4.5) o funkcích a, b, c, , , předpokládáme, že jsou definovány na nějakém intervalu J a jsou na něm diferencovatelné. V případě, že existuje t0 J takové, že b(t0) = 0, je funkce b na nějakém podintervalu I intervalu J nenulová. Pro t I z první rovnice vyjádříme druhou složku y = 1 b (x - ax - c) (4.6) a dosadíme ji do druhé rovnice: y = x + b (x - ax - c) + . (4.7) První z rovnic (4.5) zderivujeme x = a x + ax + b y + by + c , za y dosadíme (4.6) a za y dosadíme (4.7). Dostaneme x = a x + ax + b b (x - ax - c) + bx + (x - ax - c) + b + c = = a + + b b x - a - b + ab - a b b x + b + c - c - cb b . 26 První složka řešení systému (4.5) je tedy na intervalu I řešením rovnice druhého řádu x - a + + b b x - a - b + ab - a b b x = b + c - c - cb b , jeho druhá složka je pak dána rovností (4.6). V případě konstantních funkcí a, b, , a funkcí c, identicky rovných nule (homogenního systému s konstantní maticí), dostaneme rovnici x - (a + )x + (a - b)x = 0. (4.8) Povšimněme si, že charakteristický polynom konstantní matice a b je det a - b - = 2 - (a + ) + a - b. jeho koeficienty jsou tedy shodné s koeficienty na levé straně rovnice (4.8). Analogicky lze postupovat u systémů n rovnic. 4.2 Lineární diferenciální rovnice vyššího řádu Jedná se o rovnici tvaru x(n) + an-1(t)x(n-1) + an-2(t)x(n-2) + + a1(t)x + a0(t)x = f(t) , (4.9) kde funkce an-1, an-2, . . . , a1, a0, f jsou definované na nějakém intervalu J R. Je-li f(t) 0, rovnice se nazývá homogenní, v opačném případě nehomogenní. Rovnice x(n) + an-1(t)xn-1 + an-2(t)xn-2 + + a1(t)x + a0(t)x = 0 (4.10) se nazývá přidružená homogenní rovnice k rovnici (4.9). Spolu s rovnicí (4.9) uvažujeme počáteční podmínky x(t0) = x0, x (t0) = x1 0, . . . , x(n-1) (t0) = xn-1 0 . (4.11) Podle 1.0.12 je rovnice (4.9) ekvivalentní s vektorovou rovnicí (se systémem rovnic) x 1 = x2 x 2 = x3 ... x n-1 = xn x n = -an-1(t)xn - an-2(t)xn-1 - - a1(t)x2 - a0(t)x1 + f(t) Odtud a z 4.1.1, 4.1.3 a 4.1.4 plynou následující tři tvrzení: 4.2.1 Věta (o existenci a jednoznačnosti řešení) Jsou-li všechny funkce a0, a1, . . . , an-1, f spojité na intervalu J R a t0 J, pak má počáteční problém (4.9), (4.11) právě jedno řešení, které existuje na celém intervalu J. 4.2.2 Věta (princip superpozice) Jsou-li x = x(t), y = y(t) řešením homogenní lineární rovnice (4.10) a c1, c2 jsou libovolné konstanty, pak také z = z(t) = c1x(t) + c2y(t) je řešením této rovnice. 27 4.2.3 Věta Jsou-li všechny funkce a0, a1, . . . , an-1 spojité na intervalu J R, pak množina všech řešení rovnice (4.10) tvoří n-rozměrný vektorový prostor. 4.2.4 Definice Báze x1 = x1(t), x2 = x2(t), . . . , xn = xn(t) vektorového prostoru všech řešení rovnice (4.10) se nazývá fundamentální systém řešení rovnice (4.10). Řešení x1 = x1(t), x2 = x2(t), . . . , xk = xk(t) rovnice (4.9) jsou lineárně nezávislá na intervalu J, jsou-li vektory x1(t) x 1(t) ... x (n-1) 1 (t) , x2(t) x 2(t) ... x (n-1) 2 (t) , . . . , xk(t) x k(t) ... x (n-1) k (t) , lineárně nezávislé pro každé t J. Tvoří-li funkce x1 = x1(t), x2 = x2(t), . . . , xn = xn(t) fundamentální systém řešení rovnice (4.10), pak obecné řešení této rovnice je x = x(t) = c1x1(t)+c2x2(t)+ +cnxn(t), kde c1, c2, . . . , cn jsou konstanty. 4.2.5 Definice Buďte 1, 2, . . . , m funkce. Funkce W = W(t) = W(t; 1, 2, . . . , m) = det 1(t) 2(t) . . . m(t) 1(t) 2(t) . . . m(t) ... ... ... ... (m-1) 1 (t) (m-1) 2 (t) . . . (m-1) m (t) se nazývá wronskián funkcí 1, 2, . . . , m. 4.2.6 Věta Funkce x1 = x1(t), x2 = x2(t), . . . , xn = xn(t) tvoří fundamentální systém řešení rovnice (4.10) právě tehdy, když každé z nich je řešením této rovnice a jejich wronskián je v nějakém bodě intervalu J nenulový. D.: plyne z poznámek za 4.1.5. Z jednoznačnosti řešení rovnice (4.10) také plyne, že je-li wronskián řešení rovnice 4.10 nenulový v jednom bodě intervalu J, pak je nenulový ve všech bodech intervalu J. 4.2.7 Věta Buď x1, x2, . . . , xn fundamentální systém řešení rovnice (4.10). Obecné řešení rovnice (4.9) je x = x(t) = c1x1(t) + c2x2(t) + + cnxn(t) + ~x(t) , kde c1, c2, . . . , cn jsou konstanty a ~xn je libovolné partikulární řešení rovnice (4.9). D.: plyne z 4.1.6. 28 4.2.8 Metoda variace konstant Řešení rovnice (4.9) hledáme ve tvaru ~x = ~x(t) = c1(t)x1(t) + c2(t)x2(t) + + cn(t)xn(t) , kde c1(t), c2(t), . . . , cn(t) jsou funkce a x1(t), x2(t) . . . , xn(t) je fundamentální systém řešení rovnice (4.10). O funkcích c1, c2, . . . , cn budeme předpokládat, že splňují systém rovnic c 1x1 + c 2x2 + + c nxn = 0 c 1x 1 + c 2x 2 + + c nx n = 0 ... c 1x (n-2) 1 + c 2x (n-2) 2 + + c nx (n-2) n = 0 (4.12) Pak ~x = c 1x1 + c1x 1 + + c nxn + cnx n = c1x 1 + + cnx n ~x = c 1x 1 + c1x 1 + + c nx n + cnx n = c1x 1 + + cnx n ... ~x(n-1) = c1x (n-1) 1 + + cnx(n-1) n ~x(n) = c 1x (n-1) 1 + c1x (n) 1 + c nx(n-1) n + cnx(n) Současně ale platí ~x(n) = -an-1(t)~x(n-1) - an-2(t)~x(n-2) - - a1(t)~x - a0(t)~x + f(t) = = -an-1(c1x (n-1) 1 + + cnx(n-1) n ) - - a1(c1x 1 + + cnx n) - a0(c1x+ 1 + cnxn) + f = = c1(-an-1x (n-1) 1 - - a1x 1 - a0x1) + + cn(-an-1x(n-1) n - - a1x n - a0xn) + f = = c1x (n) 1 + + cnx(n) n + f . (Poslední rovnost plyne z toho, že x1, x2, . . . , xn jsou řešeními přidružené homogenní rovnice.) Celkem tedy dostaneme c 1x (n-1) 1 + + c nx(n-1) n = f (4.13) Systém rovnic (4.12) a (4.13) je soustava n lineárních rovnic pro n neznámých c 1, c 2, . . . , c n. Determinant této soustavy je wronskián fundamentálního systému řešení rovnice (4.10), je tedy podle 4.2.6 různý od 0 a systém rovnic (4.12), (4.13) má jediné řešení c 1(t) = 1(t), c 2(t) = 2(t), . . . , c n(t) = n(t). Integrací těchto rovnic určíme c1(t), c2(t), . . . , cn(t). 4.2.9 Homogenní lineární rovnice s konstatntními koeficienty x(n) + an-1x(n-1) + an-2x(n-2) + + a1x + a0x = 0 . (4.14) Řešení předpokládáme ve tvaru x(t) = et . Pak x (t) = et , x (t) = 2 et , . . . , x(n) (t) = n et a tedy n et + an-1n-1 et + an-2n-2 et + + a1et + a0et = 0 n + an-1n-1 + an-2n-2 + + a1 + a0 = 0 . 29 Poslední rovnice se nazývá charakteristická rovnice lineární diferenciální rovnice (4.14) (i) Jsou-li 1, 2 dva různé kořeny charakteristické rovnice, pak x1(t) = e1t , x2(t) = e2t jsou lineárně nezávislá řešení rovnice (4.14). D.: Že to jsou řešení, je vidět z předchozí úvahy, lineárně nezávislé jsou proto, že vektory 1 1 2 1 ... n-1 1 a 1 2 2 2 ... n-1 2 jsou pro 1 = 2 lineárně nezávislé. (ii) Je-li 0 k-násobný kořen charakteristické rovnice, pak funkce x1(t) = e0t , x2(t) = te0t , x3(t) = t2 e0t , . . . , xk(t) = tk-1 e0t jsou lineárně nezávislá řešení rovnice (4.14). D.: Z ekvivalence rovnice (4.14) s vektorovou lineární rovnicí s konstantní maticí a z úvah provedených v 4.1.8 plyne, že uvedené funkce jsou vskutku řešeními rovnice (4.14). Dále je W (0; e0t , te0t , t2 e0t , . . . , tk-1 e0t ) = = det 0 B B B B B B B @ 1 0 0 . . . 0 . . . 0 0 1 0 . . . 0 . . . 0 2 0 20 2 . . . 0 . . . 0 ... ... ... ... ... ... ... k-1 0 (k - 1)k-2 0 2 ,, k - 1 2 . . . ! ,, k - 1 k--1 0 . . . (k - 1)! 1 C C C C C C C A = = 2 6 24 (k - 1)! = 0 , takže funkce x1(t) = e0t , x2(t) = te0t , x3(t) = t2 e0t , . . . , xk(t) = tk-1 e0t jsou lineárně nezávislé. (iii) Jsou-li 1 = + i, 2 = - i k-násobné komplexně sdružené kořeny charakteristické rovnice, pak x1(t) = et cos t, x2(t) = tet cos t, x3(t) = t2 et cos t, . . . , xk(t) = tk-1 et cos t, xk+1(t) = et sin t, xk+2(t) = tet sin t, xk+3(t) = t2 et sin t, . . . , x2k(t) = tk-1 et sin t jsou lineárně nezávislá řešení rovnice (4.14). D.: t e(+i)t , t e(-i)t , 0 < k jsou lineárně nezávislá řešení rovnice (4.14) podle (ii). Tedy také x+1(t) = 1 2 t e(+i)t + t e(-i)t = t et cos t xk++1(t) = 1 2i t e(+i)t - t e(-i)t = t et sin t jsou řešení rovnice (4.14) a tato řešní jsou lineárně nezávislá, neboť funkce cos a sin jsou nezávislé. 4.2.10 Partikulární řešení nehomogenní lineární rovnice s konstantními koeficienty a se speciální pravou stranou x(n) + an-1x(n-1) + an-2x(n-2) + + a1x + a0x = f(t) . (4.15) 30 ˇ f(t) = Pm(t), kde Pm je polynom stupně m. Je-li nula k-násobným kořenem charakteristické rovnice (samozřejmě připouštíme i k = 0), lze partikulární řešení hledat ve tvaru ~x(t) = tk Qm(t), kde Qm je polynom stejného stupně jako Pm. * f(t) = et Pm(t). Substituce x(t) = et y(t) převede rovnici na lineární rovnici n-tého řádu s pravou stranou Pm (předchozí případ). * f(t) = cos(t)Pm(t) nebo f(t) = sin(t)Pm(t). Najdeme partikulární řešení rovnice x(n) + an-1x(n-1) + an-2x(n-2) + + a1x + a0x = eit Pm(t). Jeho reálná část je partikulárním řešením uvažované rovnice v prvním případě, imaginární část ve druhém. * f(t) = g(t) + h(t). Partikulární řešení je součtem partikulárních řešení rovnic x(n) + an-1x(n-1) + + a1x + a0x = g(t) a x(n) + an-1x(n-1) + + a1x + a0x = h(t). 4.2.11 Nalezení fundamentálního systému řešení rovnice x + P(t)x + Q(t)x = 0 v případě, že jedno řešení je známé Předpokládejme, že známe jedno nekonstantní řešení x1 = x1(t) dané rovnice. Zavedeme substituci x(t) = x1(t)y(t). Pak x = x 1y + x1y , x = x 1 y + 2x 1y + x1y , tedy x 1 y + 2x 1y + x1y + Px 1y + Px1y + Qx1y = 0 x1y + (2x 1 + Px1)y + (x 1 + Px 1 + Qx1) y = 0 y + P + 2 x 1 x1 y = 0 , což je rovnice typu 2.10 pro neznámou funkci y = y(t). Položíme z(t) = y (t). Pak z = - P(t) + 2 x 1(t) x1(t) z ln z = dz z = - P(t) + 2 x 1(t) x1(t) dt = - P(t)dt - ln(x1(t))2 z(t) = 1 (x1(t))2 e- R P (t)dt . Odtud plyne, že y(t) = 1 (x1(t))2 e- R P (t)dt dt a tedy druhé řešení dané rovnice je x2(t) = x1(t) 1 (x1(t))2 e- R P (t)dt dt . 31 Poněvadž platí W(t, x1, x2) = x1 x1 1 x2 1 e- R P x 1 x 1 1 x2 1 e- R P + 1 x1 e- R P = e- R P (t)dt > 0 , tak x1 = x1(t), x2 = x2(t) tvoří fundamentální systém řešení dané rovnice. 4.3 Eulerova a Riccatiho diferenciální rovnice 4.3.1 Eulerova rovnice tn x(n) + an-1tn-1 xn-1 + an-2tn-2 xn-2 + + a1tx + a0x = f(t) Zavedeme substituci t = e , tj. = ln t. Pak x = dx dt = dx d d dt = 1 t dx d x = d dt 1 t dx d = - 1 t2 dx d + 1 t d2 x d2 d dt = 1 t2 d2 x d2 - dx d x = d dt 1 t2 d2 x d2 - dx d = - 2 t3 d2 x d2 - dx d + 1 t2 d3 x d3 - d2 x d2 d dt = = 1 t3 d3 x d3 - 3 d2 x d2 + 2 dx d ... Dosadíme-li do dané rovnice, vypadnou faktory t, t2 , . . . , tn , takže dostaneme rovnici s konstantními koeficienty. 4.3.2 Riccatiho rovnice x = P(t)x2 + Q(t)x + R(t) Zavedeme substituci x(t) = - y (t) P(t)y(t) . Pak x = - Pyy - (P y + Py )y P2y2 = - y Py + P y P2y + y2 Py2 a tedy - y Py + P y P2y + y2 Py2 = Py2 P2y2 - Qy Py + R - y Py + 1 Py P P + Q y - R = 0 y - P P + Q y + PRy = 0 , což je lineární homogenní rovnice druhého řádu. 4.4 Cvičení Řešte rovnice 1) d2 x dt2 + dx dt = 0 2) x + tx = 0 3) tx - 2x = 0 32 Ukažte, že x = u(t) je řešením dané rovnice a rovnici vyřešte. 4) u = t2 ; t2 x - 2x = 0 5) u = t; x + x 4t2 = 0 (t > 0) Řešte rovnice (Cauchyovy úlohy) 6) x + 2x = 0 7) x + 6x + 5x = 0 8) x + 6x + 9x = 0 9) x - 2x + 4x = 0 10) x - x = 0; x(0) = 1, x (0) = -2 11) x + 4x = 0; x(0) = 0, x (0) = 2 12) x + x = t 13) x + x = sin t 14) x - x = et 15) x - 3x - 10x = -3 16) x - x = sin t 17) x - 3x = e3t - 12t 18) x + x = cotg t 19) x - 8x = e8t 20) x + 2x = t2 - et 21) t2 x - tx + x = t 22) t2 x - tx + 2x = (ln t)2 23) t3 x - t4 x2 - t2 x = 2 Řešte systémy rovnic 24) x = -2x + y 25) x = -x - 2 3 y + 1 3 et y = 3x - 4y y = 4 3 x + y - t 26) x + 3x + 2y = 5 sin t 27) 4x + 9y + 2x + 31y = et y - 2x + 7y = 8 cost 3x + 7y + x + 24y = 3 Výsledky: 1) x = C1e-t + C2 2) x = C1 e-t2 /2 dt + C2 3) x = C1t4 + C2t + C3 4) x = C1 t + C2t2 5) x = t (C1 ln |t| + C2) 6) x = C1 + C2e-2t 7) x = C1e-t + C2e-5t 8) x = (C1 + C2t)e-3t 9) x = et (C1 cos 3 t + C2 sin 3 t) 10) x = 3e-t - et 2 11 x = sin 2t 12) x = C1 + C2e-t + t2 2 - t 13) x = C1 cos t + C2 sin t t cos t 2 14) x = C1et + C2e-t + tet 2 15) x = C1e5t + C2e-2t + 3 10 16) x = C1 + C2et + cos t - sin t 2 17) x = C1 + C2e3t + 2t2 + te3t + 4t 3 18) x = C1 cos t + C2 sin t - sin t ln 1 + cos t sin t 19) x = C1 + C2 + t 8 e8t 20) x = C1 + C2e-2t + t3 6 - t2 4 + t 4 - et 3 21) x = t A ln t + B + 1 2 (ln t)2 22) At sin(B + ln t) + (1 + ln t)2 23) x = 2Ct - 1 t2(1 - Ct) 24) x = Ae-t + Be-5t , y = Ae-t - 3Be-5t 25) x = Aet/3 + Be-t/3 - 6t, y = -2Aet/3 - Be-t/3 + 9t + 1 2 et + 9 26) x = Ae-5t + Bte-5t + 365 338 sin t - 307 338 cos t, y = 2A-B 2 e-5t + Bte-5t + 144 338 sin t + 278 338 cos t 27) x = e-4t (A cos t + B sin t) + 31 26 et - 93 17 , y = e-4t ((B - A) cos t - (B + A) sin t) - 2 13 et + 6 17 33 34 Kapitola 5 Autonomní systémy Budeme uvažovat systém rovnic x = f(x) , (5.1) kde f = (f1, f2, . . . , fn) : Rn , Rn je otevřená množina. Systém (5.1) se nazývá autonomní systém obyčejných diferenciálních rovnic, definiční obor pravých stran se nazývá fázový prostor. V celé kapitole budeme předpokládat, že f je spojitá funkce taková, že počáteční problém: (5.1) s podmínkou x(t0) = x0 (5.2) má jediné řešení pro každé t0 R, x0 . 5.1 Fázový prostor, trajektorie, stacionární body 5.1.1 Věta Je-li x = x(t) řešením úlohy (5.1), (5.2), pak pro každé R je y = y(t) = x(t + ) řešením rovnice (5.1) s počáteční podmínkou y(t0 - ) = x0 . Je-li x definováno na intervalu (t1, t2), je y definováno na intervalu (t1 - , t2 - ). D.: y (t) = d dt x(t + ) = f(x(t + )) = f(y(t)). Tato věta ukazuje, že autonomní systémy popisují děje invariantní vzhledem k posunutí v čase. Bez újmy na obecnosti se tedy u autonomních systémů můžeme omezit na počátečními problémy s počáteční podmínkou x(0) = x0 . (5.3) Řešení x = x(t) problému (5.1), (5.3) lze interpretovat buďto jako graf zobrazení x : R , nebo jako křivku C = {x(t) : t R} ve fázovém prostoru zadanou parametricky. Tuto křivku nazveme trajektorií řešení x. Křivku C+ = {x(t) : t 0}, resp. C= {x(t) : t 0}, nazveme pozitivní, resp. negativní, polotrajektorií systému 5.1. Příklad: Systém x = y y = -x má řešení x(t) = r cos(t - ) y(t) = -r sin(t - ) , kde r = x2 0 + y2 0 , tg = y0 x0 nebo cotg = x0 y0 . Poněvadž x(t)2 + y(t)2 = r2 , jsou trajektorie řešení kružnice se středem v počátku. 35 5.1.2 Věta Jsou-li x = x(t), y = y(t) řešení systému (5.1), pak jejich trajektorie buďto splývají, nebo nemají žádný společný bod. D.: Nechť x(t1) = y(t2) pro nějaká t1, t2 0. Podle 5.1.1 je z = z(t) = x(t + (t1 - t2)) řešením rovnice (5.1) s počáteční podmínkou z(t2) = x(t1). Trajektorie řešení z a x zřejmě splývají. Současně ale z je řešením rovnice (5.1) s počáteční podmínkou z(t2) = y(t2) a z předpokládané jednoznačné řešitelnosti problému (5.1) s libovolnou počáteční podmínkou plyne z y. 5.1.3 Definice Bod x se nazývá stacionární bod (rovnovážný bod, ekvilibrium, kritický bod, singulární bod, degenerovaná trajektorie) rovnice (5.1), jestliže f(x ) = 0. Trajektorie rovnice (5.1) se nazývá cyklus, je-li uzavřenou křivkou. Trajektorie {x(t) : t R} rovnice (5.1) se nazývá homoklinická, jestliže existuje stacionární bod x takový, že lim t x(t) = lim tx(t) = x . Trajektorie {x(t) : t R} rovnice (5.1) se nazývá heteroklinická, jestliže existují stacionární body x , x takové, že lim t x(t) = x a lim tx(t) = x . 5.1.4 Věta Trajektorie řešení autonomního systému (5.1) jsou jednoho z typů: ­ Stacionární body (odpovídají konstantním řešením); ­ Cykly (odpovídají nekonstantním periodickým řešením); ­ Trajektorie, které samy sebe neprotínají. D.: Plyne z 5.1.2. 5.1.5 Definice Neprázdná podmnožina A fázového prostoru systému (5.1) se nazývá ­ pozitivně invariantní (invariantní vpřed, forward invariant), jestliže pro každé řešení x( ) systému (5.1) s počáteční hodnotou x(0) A platí, že x(t) A pro všechna t 0; ­ negativně invariantní (invariantní vzad, backward invariant), jestliže pro každé řešení x( ) systému (5.1) s počáteční hodnotou x(0) A platí, že x(t) A pro všechna t 0; ­ invariantní, je-li současně pozitivně i negativně invariantní. 5.1.6 Poznámky 1. Jsou-li množiny A, B pozitivně (resp. negativně) invariantní, pak také množiny A B a A B jsou pozitivně (resp. negativně) invariantní. 2. Libovolná trajektorie C systému (5.1) je invariantní množinou tohoto systému. 5.1.7 Definice Nechť A, B , B = a je nějaká metrika na ekvivalentní s euklidovskou. Řekneme, že ­ množina A atrahuje (přitahuje) množinu B (množina A je atraktorem množiny B), jestliže pro každé řešení systému (5.1) s počáteční hodnotou x(0) B platí, že lim t x(t), A = 0; ­ množina A je (globální) atraktor, jestliže A přitahuje ; 36 ­ množina A absorbuje množinu B, jestliže A je pozitivně invariantní a ke každému řešení x systému (5.1) s počáteční hodnotou x(0) B existuje T 0 takové, že x(T ) A; ­ množina A je globálně absorbující, jestliže absorbuje množinu . 5.1.8 Poznámky 1. Nechť x( ) je řešením systému (5.1) s počáteční podmínkou x(0) = x0 . Pokud množina A je -limitní množinou řešení x( ), pak je pozitivně invariantním atraktorem množiny {x0}. 2. Trajektorii C systému (5.1) nazveme limitní trajektorií, jestliže existuje množina B taková, že B ( \ C) = a C atrahuje množinu B. Je-li C navíc cyklem, nazveme ho limitním cyklem. 3. Buď i {1, 2, . . ., n}. Jestliže existují kladné konstanty K, takové, že pro každý bod x = (x1, x2, . . . , xn) splňující podmínku |xi| K platí sgn xifi(x) |xi|, pak množina Ai = {x = (x1, x2, . . . , xn) : |xi| K} je globálně absorbující množinou systému (5.1). D.: Nejprve ukážeme, že každá trajektorie protíná množinu Ai. Připusťme, že existuje řešení x( ) = x1( ), x2( ), . . . , xn( ) systému (5.1) takové, že pro všechna t 0 je |xi(t)| > K. Položme u(t) = |xi(t)|. Pak pro všechna t 0 je u (t) = d dt |xi(t)| = sgn xi(t)fi x(t) |xi(t)| = -u(t), neboli u (t) u(t) -. Integrací této nerovnosti v mezích od 0 po t dostaneme ln u(t) - ln u(0) -t, tj. |xi(t)| = u(t) = u(0)e-t = |xi(0)|e-t pro libovolné t 0. Odtud plyne, že lim t |xi(t)| = 0, což je ve sporu s předpokladem |xi(t)| > K > 0. Množina Ai má neprázdný průnik s libovolnou trajektorií, je tedy neprázdná. Ukážeme, že je navíc pozitivně invariantní. Připusťme, že existuje řešení x( ) = x1( ), x2( ), . . . , xn( ) systému (5.1) s počáteční hodnotou x(0) Ai takové, že pro jisté t1 > 0 je x(t1) Ai, tj. |xi(t1)| > K. Položme M = {t R : 0 t < t1, |xi(t)| K} . Pak 0 M, takže M je neprázdná shora ohraničená množina reálných čísel. Existuje tedy T = sup M. Ze spojitosti funkce xi( ) a z vlastností suprema plyne, že T < t1, xi(T ) = K a funkce xi( ) je v bodě T rostoucí. Avšak d dt |xi(t)| t=T = sgn xi(T )fi x(T ) -|xi(t)| = -K < 0, což je spor s faktem, že funkce xi( ) je v T rostoucí. 37 5.1.9 Definice Systém (5.1) se nazývá dissipativní, jestliže existuje ohraničená globálně atrahující množina. 5.1.10 Poznámky 1. Je-li systém (5.1) dissipativní, pak každé jeho řešení je ohraničené. 2. Jestliže existují kladné konstanty K, takové, že pro všechna i {1, 2, . . ., n} a všechna x taková, že |x| K platí sgn xifi(x) -, pak je systém (5.1) dissipativní a globálně absorbující je množina A = {x : |x| K}. D.: Nejprve ukážeme, že každá trajektorie protíná množinu A. Připusťme, že existuje řešení x( ) systému (5.1) takové, že x(t) A pro všechna t 0. Pak |x(t)| > K pro všechna t > 0, a tedy pro libovolné t > 0 platí d dt |x(t)| = n i=1 d dt |xi(t)| = n i=1 sgn xi(t)x i(t) = n i=1 sgn xi(t)fi x(t) n i=1 (-) = -n. Integrací této nerovnosti dostaneme |x(t)| |x(0)| - nt, takže pro t |x(0)| - K n je |x(t)| K, což je spor. Každá trajektorie má tedy s množinou A neprázdný průnik, což také znamená, že množina A je neprázdná. Ukážeme, že množina A je navíc pozitivně invariantní. Nechť x( ) je řešením systému (5.1) s počáteční podmínkou x(0) A. Připusťme, že existuje t1 > 0, pro něž x(t1) A. Pak |x(t1)| > K. Položme M = {t R : 0 t < t1, |x(t)| K} . Pak 0 M, takže M je neprázdná shora ohraničená množina reálných čísel. Existuje tedy T = sup M. Ze spojitosti funkce |x( )| a z vlastností suprema plyne, že |x(T )| = K a funkce |x( )| je v bodě T rostoucí. Avšak d dt |x(t)| t=T = n i=1 sgn xi(T )x i(T ) = n i=1 sgn xi(T )fi x(T ) -n < 0, což je spor s tím, že funkce |x( )| je v bodě T rostoucí. Pro všechna t > 0 je tedy x A a množina A je invariantní. 3. Nechť ke každému i {1, 2, . . ., n} existují kladné konstanty Ki, i takové, že pro všechna x = (x1, x2, . . . , xn) z nerovnosti |xi| K plyne nerovnost sgn xifi(x) -i|xi|. Pak je systém (5.1) dissipativní s globálně absorbující množinou A = {x = (x1, x2, . . . , xn) : |x1| K1, |x2| K2, . . . , |xn| Kn} . D.: Položme Ai = {x : |xi| Ki}. Pak každá z množin Ai je podle třetího z tvrzení 5.1.8 globálně absorbující množinou a A = A1 A2 An. Podle 5.1.6 je množina A pozitivně invariantní. Ukážeme, že je také globálně absorbující. Buď x(t) libovolné řešení systému (5.1). Podle třetího z tvrzení 5.1.8 existuje t1 0 takové, že |x1(t1)| K1 pro všechna t t1. Dále existuje t2 t1 takové, že pro všechna t t2 je |x2(t2)| K2 atd. Nakonec existuje tn tn-1 takové, že |xn(t)| Kn pro všechna t tn. Takže pro všechna t tn je |x1(t)| K1, |x2(t)| K2, . . . , |xn(t)| Kn, tj. x(t) A. 38 5.1.11 Věta (Poincaré [1854­1912]-Bendixson [1861-1935]) Jestliže rovnice (5.1) má trajektorii C+ = {x(t) : t 0}, která je ohraničená a její uzávěr neobsahuje stacionární body rovnice (5.1), pak existuje cyklus rovnice (5.1), který leží v C+. D.: P. Hartman: Ordinary Differential Equations. John Willey&Sons, New York-London-Sydney, 1964, kap. VII. 5.2 Autonomní systémy v rovině V tomto oddílu se budeme zabývat systémem (5.1) pro n = 2, tedy systémem x = f(x, y) y = g(x, y). (5.4) 5.2.1 Definice Křivka zadaná implicitně rovnicí f(x, y) = 0 (resp. g(x, y) = 0) se nazývá x-nulklina (resp. ynulklina) rovnice (5.4). Průsečík nulklin je stacionární bod, tečna k trajektorii v jejím průsečíku s x-nulklinou (resp. y-nulklinou) je rovnoběžná s osou y (resp. x). 5.2.2 Definice (typy stacionárních bodů v rovině) Stacionární bod (x , y ) systému (5.4) se nazývá bod rotace, jestliže v jeho libovolném okolí leží cyklus, obsahující (x , y ) ve svém vnitřku; střed, jestliže existuje jeho ryzí okolí U takové, že každá trajektorie s x(0), y(0) U je cyklem obsahujícím (x , y ) ve svém vnitřku (střed je speciálním případem bodu rotace); ohnisko, jestliže existuje jeho ryzí okolí U takové, že pro každou trajektorii s x(0), y(0) U platí lim t x(t), y(t) = (x , y ) nebo lim tx(t), y(t) = (x , y ) a pro orientovaný úhel (t), který svírá vektor x(t), y(t) - (x , y ) s nějakým pevným vektorem platí lim t (t) = nebo lim t (t) = - nebo lim t(t) = nebo lim t(t) = (trajektorie se přibližuje ke stacionárnímu bodu (nebo se od něho vzdaluje) po spirále); uzel, jestliže existuje jeho ryzí okolí U takové, že pro každou trajektorii s x(0), y(0) U platí lim t x(t), y(t) = (x , y ) nebo lim tx(t), y(t) = (x , y ) a pro orientovaný úhel (t), který svírá vektor x(t), y(t) - (x , y ) s nějakým pevným vektorem existuje vlastní lim t (t) nebo lim t- (t); sedlo, jestliže existuje jen konečný počet trajektorií (x, y) = x(t), y(t) takových, že lim t x(t), y(t) = (x , y ) nebo lim tx(t), y(t) = (x , y ) . 39 5.2.3 Stacionární body lineárního homogenního autonomního systému x = ax + by y = cx + dy. (5.5) Označme A = a b c d . Pokud det A = ad - bc = 0, má systém (5.5) jediný stacionární bod (0, 0). Vlastní čísla matice A jsou kořeny charakteristické rovnice a - b c d - = 2 - (a + d) + ad - cb = 2 - (tr A) + det A = 0, (5.6) tedy při označení D = (tr A)2 - 4 det A je 1,2 = 1 2 tr A D . (i) det A < 0 V tomto případě je D > (tr A) 2 > 0, což znamená, že rovnice (5.6) má dva reálné různé kořeny 1, 2. Poněvadž D > |tr A|, mají tyto kořeny opačná znaménka. Nechť pro určitost 1 < 0 < 2. Označme v1, resp. v2, vlastní vektor matice A příslušný k vlastní hodnotě 1, resp. 2. Obecné řešení systému (5.5) je podle 4.1.8 x(t) y(t) = v1e1t + v2e2t , kde , jsou nějaké konstanty. Pro = 0 = je lim t- x(t) y(t) = |v2| lim t- e2t = 0, pro = 0 = je lim t x(t) y(t) = |v1| lim t e1t = 0 a pro = 0 = je lim t- x(t) y(t) = |v1| lim t- e1t + |v2| lim t- e2t = + 0 = , lim t x(t) y(t) = |v1| lim t e1t + |v2| lim t e2t = 0 + = . To znamená, že stacionární bod (0, 0) je sedlo. (ii) det A > 0 (ii.1) tr A = 0 Nechť x(t), y(t) je řešením systému (5.5) a označme (t) úhel, který svírá přímka procházející body (0, 0), x(t), y(t) s vodorovnou osou. Platí tg (t) = y(t) x(t) , pokud x(t) = 0, cotg (t) = x(t) y(t) , pokud y(t) = 0. (ii.1.a) det A < 1 4 (tr A) 2 . V tomto případě je D > 0 a D < |tr A|. Charakteristická rovnice (5.6) má dva reálné různé kořeny 1, 2 takové, že oba mají stejné znaménko jako tr A. Nechť pro určitost 1 < 2 a u = u1 u2 , resp. v = v1 v2 , 40 je vlastní vektor příslušný k vlastní hodnotě 1, resp. 2. Alespoň jedna ze souřadnic každého z vlastních vektorů je nenulová. Obecné řešení systému (5.5) s počáteční podmínkou x(0), y(0) = (0, 0) je x(t) y(t) = u1 u2 e1t + v1 v2 e2t ; přitom alespoň jedna z konstant , je nenulová. Je-li tr A > 0, pak 0 < 1 < 2 a lim t- x(t) y(t) = lim t- u1 u2 e1t + v1 v2 e2t = u1 u2 0 + v1 v2 0 = 0 0 ; je-li tr A < 0, pak 1 < 2 < 0 a lim t x(t) y(t) = lim t u1 u2 e1t + v1 v2 e2t = u1 u2 0 + v1 v2 0 = 0 0 . Pokud = 0 = u1, pak lim t- y(t) x(t) = lim t- u2e1t + v2e2t u1e1t + v1e2t = lim tu2 + v2e(2-1)t u1 + v1e(2-1)t = u2 u1 a pokud = 0 = v1, pak lim t y(t) x(t) = lim t u2e1t + v2e2t u1e1t + v1e2t = lim t- u2e(1-2)t + v2 u1e(1-2)t + v1 = v2 v1 . Analogicky, pokud = 0, pak lim t- x(t) y(t) = u1 u2 když u2 = 0, lim t x(t) y(t) = v1 v2 když v2 = 0. Pokud = 0, pak lim t y(t) x(t) = v2 v1 když v1 = 0, lim t x(t) y(t) = v1 v2 když v2 = 0. Stacionární bod (0, 0) je v tomto případě uzel. Směrový vektor polotečny k trajektorii ve stacionárním bodě je vlastním vektorem matice A. (ii.1.b) det A = 1 4 (tr A) 2 V tomto případě je tr A = 0, neboť det A = 0, charakteristická rovnice (5.6) má dvojnásobný kořen 1 2 tr A a systém (5.5) s počáteční podmínkou x(0), y(0) = (0, 0) má řešení x(t) y(t) = + t + t e(1 2 tr A)t , kde , , , jsou nějaké konstnty, z nichž aspoň dvě jsou nenulové. Proto platí lim t y(t) x(t) = pro = 0, lim t x(t) y(t) = pro = 0, lim t y(t) x(t) = pro = 0 = . Je-li tr A < 0, pak lim t e(1 2 tr A)t = 0, lim t te(1 2 tr A)t = lim t t e-(1 2 tr A)t = lim t 1 -1 2 tr A e-(1 2 tr A)t = 0 a tedy lim t x(t) y(t) = 0 0 . Je-li tr A > 0 pak analogicky lim t- x(t) y(t) = 0 0 . 41 Stacionární bod (0, 0) je v tomto případě uzel. Nyní však již obecně neplatí, že směrový vektor polotečny k trajektorii ve stacionárním bodě je vlastním vektorem matice A; v případě = 0 = (tj. pokud vlastní hodnotě matice A přísluší dva lineárně nezávislé vlastní vektory) je každá přímka procházející bodem (0, 0) polotečnou nějaké trajektorie. (ii.1.c) det A > 1 4 (tr A) 2 V tomto případě je D < 0, charakteristická rovnice (5.6) má dva různé komplexně sdružené kořeny 1 2 tr A i -D a systém (5.5) s počáteční podmínkou x(0), y(0) = (0, 0) má řešení x(t) y(t) = cos -D 2 t - - - c b sin -D 2 t - - ~ e(1 2 tr A)t , kde , , ~ jsou vhodné konstanty, přičemž = 0, tg ~ = d - a -D , ~ 2k + 1 2 : k Z . Jedná se o parametrické vyjádření spirály, která se ,,navíjí na stacionární bod (0, 0) nebo se z něho ,,odvíjí . Pokud tr A > 0, pak lim t- x(t) y(t) = 0 0 , pokud tr A < 0, pak lim t x(t) y(t) = 0 0 . (ii.2) tr A = 0 V tomto případě je D = -4 det A < 0 a kořeny charakteristické rovnice (5.6) jsou ryze imaginární, 1 = i det A , 2 = -i det A , takže řešení systému (5.5) s počáteční podmínkou x(0), y(0) = (0, 0) je x(t) y(t) = cos det A t - - - c b sin det A t - + ~ pro vhodné konstanty , , přičemž = 0, tg ~ = a det A , ~ 2k + 1 2 : k Z . Jedná se o parametrické vyjádření elips se středem (0, 0), každá trajektorie je tedy cyklem se stacionárním bodem (0, 0) ve svém vnitřku. To znamená, že stacionární bod (0, 0) je střed. 5.2.4 Věta Uvažujme autonomní systém x = ax + by + P(x, y) y = cx + dy + Q(x, y) , (5.7) kde P, Q jsou funkce dvou proměnných spojité v okolí počátku. Nechť ad-bc = 0 a nechť existuje > 0 takové, že lim (x,y)(0,0) |P(x, y)| + |Q(x, y)| (|x| + |y|)1+ = 0 . Je-li bod (0, 0) uzlem nebo ohniskem pro systém (5.5), pak je stejného typu i pro systém (5.7). Je-li bod (0, 0) středem pro systém (5.5), pak je bodem rotace nebo ohniskem pro systém (5.7). Je-li bod (0, 0) sedlem pro systém (5.5) a funkce P, Q mají spojité parciální derivace podle obou proměnných v okolí počátku, pak je (0, 0) sedlem i pro systém (5.7). D.: P. Hartman: Ordinary Differential Equations. John Willey&Sons, New York-London-Sydney, 1964, kap. VIII. 42 5.2.5 Důsledek Nechť (x , y ) je stacionárním bodem systému (5.4) (tj. f(x , y ) = g(x , y ) = 0) a funkce f, g mají spojité druhé parciální derivace podle obou proměnných v okolí bodu (x , y ). Označme f1 = f x (x , y ), f2 = f y (x , y ), g1 = f x (x , y ), g2 = f y (x , y ) a nechť f1g2 - f2g1 = 0. Pak je bod (x , y ) uzlem, ohniskem nebo sedlem pro systém (5.4), je-li počátek stacionárním bodem stejného typu pro lineární homogenní systém x = f1x + f2y y = g1x + g2y . (5.8) Je-li počátek středem pro systém (5.8), je bod (x , y ) buďto ohniskem nebo bodem rotace pro systém (5.4). D.: Plyne z 5.2.4 a z Taylorovy věty pro funkce více proměnných (sr. Z. Došlá, O. Došlý: Diferenciální počet funkcí více proměnných. MU, Brno 1999, 5.2) 5.2.6 Věta (Dulacovo kritérium) Nechť funkce f, g jsou spojitě diferencovatelné a existuje funkce q : R, která je rovněž spojitě diferencovatelná a taková, že výraz q(x, y)f(x, y) x + q(x, y)g(x, y) y je v nějaké jednoduše souvislé oblasti B stále nezáporný nebo nekladný, přičemž není identicky roven nule v žádné otevřené podmnožině množiny B, pak v množině B neexistuje cyklus systému (5.4). D.: J. Kalas, Z. Pospíšil: Spojité modely v biologii. MU, Brno 2001, str. 9­10. 5.3 Stabilita 5.3.1 Definice (Persidskij [1903­1970]) Nechť x0 = x0(t) je řešení systému (5.1) definované na intervalu [0, ). Řešení x0 se nazývá stejnoměrně stabilní, jestliže ke každému > 0 existuje > 0 tak, že pro každé t0 0 všechna řešení x = x(t) systému (5.1) splňující podmínku |x(t0) - x0(t0)| < existují pro všechna t t0 a splňují pro ně nerovnost |x(t) - x0(t)| < . Není-li řešení x0 = x0(t) systému (5.1) stejnoměrně stabilní, nazývá se nestabilní. 5.3.2 Definice (Ljapunov [1857­1918]) Nechť x0 = x0(t) je řešení systému (5.1) definované na intervalu [0, ). Řešení x0 se nazývá stejnoměrně asymptoticky stabilní, je-li stejnoměrně stabilní a existuje > 0 tak, že pro každé t1 0 a všechna řešení x = x(t) systému (5.1) splňující podmínku |x(t1) - x0(t1)| < platí lim t |x(t) - x0(t)| = 0. Ze struktury prostoru řešení lineárního homogenního systému s konstantními koeficienty (sr. 4.1.8) plynou následující tři věty. 43 5.3.3 Věta Buď A konstantní matice. Jestliže všechny kořeny její charakteristické rovnice det(A - E) = 0 (vlastní čísla matice A) mají nekladnou reálnou část a ty s nulovou reálnou částí jsou jednoduché, pak řešení x0 = x0(t) 0 lineárního autonomního systému x = A x (5.9) je stejnoměrně stabilní. 5.3.4 Věta Jestliže alespoň jedno vlastní číslo matice A má kladnou reálnou část, pak řešení x0 = x0(t) 0 lineárního autonomního systému (5.9) je nestabilní. 5.3.5 Věta Řešení x0 = x0(t) 0 lineárního autonomního systému (5.9) je stejnoměrně asymptoticky stabilní právě tehdy, když každé vlastní číslo matice A má zápornou reálnou část. Uvažujme nyní perturbovaný lineární systém s konstantními koeficienty x = A x + g(x) . (5.10) 5.3.6 Věta Buď Y = Y(t) fundamentální matice řešení systému (5.9). Jestliže existují konstanty K > 0 a < 1 K takové, že t 0 |Y(t)Y(s) -1 |ds K pro t 0 (5.11) a |g(x)| |x| pro x , pak řešení x0 = x0(t) 0 systému (5.10) je stejnoměrně asymptoticky stabilní. D.: J. Kalas, M. Ráb: Obyčejné diferenciální rovnice. MU, Brno 2001, str. 130­131. Poznámka: Podmínka (5.11) zaručí stejnoměrnou asymptotickou stabilitu nulového řešení systému (5.9). Věta říká, že je-li perturbace g v jistém smyslu dostatečně malá, zůstává zachována stejnoměrná asymptotická stabilita nulového řešení rovnice (5.10). Z hlediska aplikací je důležité vyšetřování stejnoměrné asymptotické stability konstantních řešení (stacionárních bodů) rovnice (5.1). Je-li funkce f dvakrát spojitě diferencovatelná a f(x ) = o, pak podle Taylorovy věty pro funkce více proměnných platí f(x) = A (x - x ) + r1(x, x ) , kde A = (aij) = fi xj (x ) a r1 je příslušný Taylorův zbytek. Vyšetřování stejnoměrné asymptotické stability konstantních řešení rovnice (5.1) lze transformací y = x-x převést na vyšetřování stejnoměrné asymptotické stability nulového řešení rovnice y = A y + g(y), kde g(y) = r1(y + x , x ) a tu vyšetřit podle věty 5.3.6. 44 5.3.7 Definice Buď x stacionární bod systému (5.1). Jacobiho matice zobrazení f = (f1, f2, . . . , fn) v bodě x J(x ) = f1 x1 (x ) f1 xn (x ) ... ... ... fn x1 (x ) fn xn (x ) se nazývá variační matice systému (5.1) ve stacionárním bodě x a homogenní lineární autonomní systém z = J(x ) z se nazývá variační rovnice systému (5.1). 5.3.8 Věta Buď x stacionární bod systému (5.1) a nechť zobrazení f je spojitě diferencovatelné. Mají-li všechna vlastní čísla variační matice J(x ) záporné reálné části, pak konstantní řešení x(t) x systému (5.1) je stejnoměrně asymptoticky stabilní. Pokud existuje vlastní číslo variační matice J(x ) s kladnou reálnou částí, pak je konstantní řešení x(t) x systému (5.1) nestabilní. D.: J. Kalas, M. Ráb: Obyčejné diferenciální rovnice. MU, Brno 2001, str. 137­138. 5.3.9 Kvalitativní vlastnosti řešení dvourozměrného autonomního systému (5.4) Nechť (x , y ) je stacionární bod systému (5.4), funkce f, g jsou dvakrát spojitě diferencovatelné a J(x , y ) je variační matice tohoto systému v bodě (x , y ). Spojením výsledků z 5.2.2, 5.2.5 a 5.3.8 dostaneme dostatečné podmínky pro to, aby stacionární bod (x , y ) byl sedlem, stabilním nebo nestabilním uzlem a ohniskem; tyto podmínky jsou shrnuty v tabulce 5.1. det J(x , y ) < 0 sedlo tr J(x , y ) 2 4 det J(x , y ) nestabilní uzel tr J(x , y ) > 0 tr J(x , y ) 2 < 4 det J(x , y ) nestabilní ohnisko tr J(x , y ) 2 4 det J(x , y ) stabilní uzel det J(x , y ) > 0 tr J(x , y ) < 0 tr J(x , y ) 2 < 4 det J(x , y ) stabilní ohnisko Tabulka 5.1: Klasifikace stacionárních bodů systému (5.4). Uvedené podmínky jsou dostatečné pro to, aby stacionární bod (x , y ) byl typu uvedeného v posledním sloupci tabulky; J(x , y ) označuje variační matici systému (5.4) ve stacionárním bodě (x , y ). Do konce tohoto odstavce budeme symbolem (t; x0 ) = 1(t; x0 ), . . . , n(t; x0 ) označovat řešení počáteční úlohy (5.1), (5.2). 5.3.10 Definice Buď x a G okolí bodu x ve fázovém prostoru . Spojitá funkce V : G R se nazývá ljapunovská funkce systému (5.1) v bodě x , jestliže 45 (i) V (x ) = 0 a V (x) > 0 pro každé x G \ {x }. (ii) Pro každé G je složená funkce V (t; ) (chápaná jako funkce jedné reálné proměnné t) nerostoucí pro všechna t 0. 5.3.11 Věta (Přímá Ljapunova metoda) Existuje-li ljapunovská funkce systému (5.1) v bodě x , pak x je stacionárním bodem systému (5.1) a konstantní řešení x(t) x tohoto systému je stejnoměrně stabilní. Pokud navíc podmínku (ii) z definice 5.3.10 lze nahradit silnější podmínkou (ii ) Pro každé G je složená funkce V (t; ) (chápaná jako funkce jedné reálné proměnné t) klesající pro všechna t 0, pak je konstantní řešení x(t) x systému (5.1) stejnoměrně asymptoticky stabilní. D.: Pokud by existovalo > 0 takové, že (; x ) = x , pak by V (; x ) > 0 = V (0; x ) a funkce V (; x ) by nebyla nerostoucí. Bod x je tedy stacionárním bodem systému (5.1). Bez újmy na obecnosti lze předpokládat, že x = o. V opačném případě bychom totiž mohli systém (5.1) substitucí y = x - x transformovat na systém y = f(y + x ), pro který je o stacionárním bodem. Buď > 0 libovolné číslo takové, že {x Rn : |x| } G a označme = min {V (x) : |x| = } . Pak > 0 a V (x) pro každé x takové, že |x| = . Ze spojitosti funkce V plyne, že existuje > 0 takové, že |V (x) - 0| = V (x) < pro všechna x taková, že |x| < . Zřejmě je < . Buď dále takové, že || < . Pak V (0; ) < a poněvadž funkce V ( ; ) je nerostoucí, platí V (t; ) < pro všechna t > 0 z definičního oboru funkce ( ; ). (5.12) Kdyby nyní existovalo t1 > 0 takové, že |(t1, )| , pak by ze spojitosti funkce |( ; )| a z Bolzanovy věty plynula existence t0 (0, t1) takového, že |(t0; )| = a platilo by V (t0; ) , což by byl spor s (5.12). Pro všechna t > 0 z definičního oboru funkce ( ; ) tedy platí |(t; )| < . Odtud navíc podle 3.3.7 plyne, že ( ; ) je definována pro všechna t > 0. Tvrzení o stejnoměrné stabilitě je tedy dokázáno. V případě, že = x , platí: (t; ) = (t; x ) = x pro každé t 0, takže lim t (t; ) = x . Nechť = x = o a funkce V ( ; ) je klesající. Poněvadž funkce V ( ; ) je monotonní, existuje lim t V (t; ) = R (5.13) a z nezápornosti funkce V plyne 0. Připusťme > 0. Z toho, že lim xx V (x) = V (x ) = 0, plyne existence t2 > 0 a > 0 takových, že |(t; )| pro všechna t t2. Pro všechna t t2 je tedy |(t; )| . 46 Položme v(z) = V (1; z) - V (0; z) . Funkce v je podle 3.3.9 spojitá na kompaktní množině {z Rn : |z| } a je zde záporná. Podle Weierstrassových vět existuje = max {v(z) : |z| } ; je < 0 a pro k N platí V (t2 + k; ) = V (t2 + k; ) - V (t2 + k - 1; ) + V (t2 + k - 1; ) = = v (t2 + k - 1; ) + V (t2 + k - 1; ) = . . . = k i=1 v (t2 + i - 1; ) + V (t2; ) k + V (t2; ) . Poněvadž lim k k = -, je také lim k V (t2 + k; ) = -, což je spor s (5.13). Tento spor dokazuje, že = 0. Ze spojitosti funkce V , z faktu (t; ) = x pro t > 0 a = x , z podmínky (i) v definici 5.3.10 a ze vztahu (5.13) nyní plyne lim t (t; ) = x . Tím je dokázáno i tvrzení o stejnoměrné asymptotické stabilitě. 5.3.12 Důsledek Buď V : G R diferencovatelná funkce, která splňuje podmínku (i) z definice 5.3.10. Jestliže pro každé x G platí ˙V (x) = n i=1 V (x) xi fi(x) 0, (5.14) pak funkce V je ljapunovskou funkcí systému (5.1) v bodě x a tedy konstantní řešení x(t) x tohoto systému je stejnoměrně stabilní. Jestliže pro každé x G \ {x } platí ˙V (x ) < 0, pak funkce V splňuje podmínku (ii ) z věty 5.3.11 a tedy konstantní řešení x(t) x tohoto systému je stejnoměrně asymptoticky stabilní. D.: Pro každé G a každé t 0 je x = (t, ) G a podle věty o derivaci složené funkce platí d dt V (t; ) = n i=1 V xi (t; ) di(t; ) dt = n i=1 V xi (t; ) fi (t; ) = = n i=1 V (x) xi fi(x) = ˙V (x). Odtud a ze známých vět o vyšetřování průběhu funkce jedné proměnné pomocí derivace plynou obě tvrzení. Je-li U diferencovatelná funkce definovaná na G , pak výraz ˙U(x) definovaný vztahem ˙U(x) = n i=1 U(x) xi fi(x) se nazývá derivace funkce U vzhledem k systému (5.1). 47 48 Kapitola 6 Aplikace 6.1 Některé klasické elementární úlohy V tomto oddílu je uvedeno několik úloh vedoucích na obyčejné diferenciální rovnice, které lze vyřešit elementárními metodami z kapitoly 2. Úlohy vychází z různých oborů -- kinematiky (6.1.1, 6.1.4), geometrické optiky (6.1.3), dynamiky (6.1.2), makroekonomie (6.1.5). 6.1.1 Traktrisa Po stole táhneme hodinky na napjatém řetízku délky tak, že koncem řetízku sledujeme hranu stolu. Na počátku svírá řetízek a hrana stolu úhel (0, 1 2 ]. Úkolem je určit dráhu hodinek. Zvolíme orthonormální souřadnou soustavu tak, že svislá osa splývá s hranou stolu a je souhlasně orientovaná se směrem pohybu konce řetízku, viz obr. 6.1. Při této volbě budou hodinky na počátku v bodě (- sin , 0). Dráhu hodinek vyjádříme jako graf funkce y = y(x). Hodinky se pohybují ve směru působící síly, síla působí ve směru řetízku. To znamená, že přímka incidentní s řetízkem je tečnou ke grafu funkce y v každém bodě. Směrnice této tečny je tedy rovna y (x) = 2 - x2 -x . (6.1) Hledaná funkce je řešením této obyčejné diferenciální rovnice s počáteční podmínkou y(- sin) = 0. (6.2) Na pravé straně rovnice (6.1) se nevyskytuje hledaná funkce y, proto můžeme řešení úlohy (6.1), (6.2) bezprostředně psát ve tvaru určitého integrálu y(x) = x - sin 2 - 2 d = ln + 2 - 2 || - 2 - 2 x =- sin = = cos - 1 - x2 2 + ln + 2 - x2 -x sin 1 + cos . Úlohu o dráze hodinek tažených na řetízku po stole jako první řešil Gottfried Wilhelm von Leibniz (1646­1716). Křivku podrobně studoval v roce 1692 Christiaan Huygens, který jí také dal jméno tractrix (z latinského trahere, táhnout). 49 y x 2 - x2 Obrázek 6.1: Traktrisa 6.1.2 Ciolkovského rovnice Pohyb rakety budeme popisovat v souřadné soustavě takové, aby na raketu nepůsobily žádné vnější síly (tedy ve stavu beztíže). Nechť v čase t0 = 0 se raketa pohybuje rychlostí v0. V čase t0 se zažehne palivo, které rovnoměrně shoří za čas T a v podobě plynů proudí z trysky na zádi rakety rychlostí u vzhledem k raketě. Úlohou je určit rychlost rakety po provedení popsaného manévru, tedy její rychlost v čase T . Označme M . . . hmotnost rakety na počátku (v čase t0 = 0), . . . hmotnost paliva vyhořelého za čas T , m = m(t) . . . hmotnost rakety (s dosud nevyhořelým palivem) v čase t, v = v(t) . . . rychlost rakety v čase t. Předpoklad o rovnoměrném hoření paliva zapíšeme rovností m(t) = M - T t = MT - t T . (6.3) Rychlost v neznáme. Budeme však o ní předpokládat, že je spojitě diferencovatelnou funkcí svého argumentu (času). Hybnost rakety se zbývajícím palivem v čase t je p(t) = m(t)v(t). (6.4) Uvažujme krátký časový interval [t, t + t] [0, T ]. Během něho shoří palivo o hmotnosti = T t. (6.5) Rychlost vytékajících plynů v souřadné soustavě, v níž pohyb popisujeme, je v čase t rovna v(t)-u a v průběhu intervalu se mění v rozmezí od této hodnoty po hodnotu v(t + t) - u. Hybnost vyhořelého paliva vytrysklého v uvažovaném časovém intervalu proto vyjádříme jako pP (t, t) = w(t, t), (6.6) 50 kde w(t, t) je integrální průměr vytékajících plynů v časovém intervalu délky t, tj. w(t, t) = 1 t t+t t v() - u d = 1 t t+t t v()d - u. Podle první věty o střední hodnotě integrálního počtu existuje číslo (0, 1) takové, že t+t t v()d = v(t + t)t, takže w(t, t) = v(t + t) - u. S využitím této rovnosti a rovnosti (6.5) vyjádříme hybnost (6.6) vytékajícího plynu výrazem pP (t, t) = v(t + t) - u T t. (6.7) Hybnost rakety v čase t + t je vzhledem k (6.3) rovna pR(t + t) = m(t + t)v(t + t) = M - T (t + t) v(t + t) = m(t) - T t v(t + t). Podle Lagrangeovy věty o střední hodnotě platí v(t + t) = v(t) + v (t + t)t, kde (0, 1). Dosazením této rovnosti do předchozí dostaneme pR(t + t) = m(t) - T t v(t) + v (t + t)t = = m(t)v(t) - T v(t) - m(t)v (t + t) (t) - T v (t + t)(t)2 . (6.8) Souhrnná hybnost rakety a vyhořelého paliva je v čase t + t rovna p(t + t) = pR(t + t) + pP (t, t). Odtud a z (6.7), (6.8) dostaneme p(t + t) - p(t) t = v(t + t) - u - v(t) T + m(t)v (t + t) - T v (t + t)t. Limitním přechodem t 0 a jednoduchou úpravou vyjádříme derivaci hybnosti soustavy rakety s palivem ve tvaru p (t) = m(t)v (t) - u T . Podle zákona o zachování hybnosti je p (t) = 0, takže s využitím rovnosti (6.3) dostaneme diferenciální rovnici pro neznámou funkci v ve tvaru v (t) = u MT - t . Na její pravé straně se nevyskytuje hledaná funkce v, stačí tedy integrovat obě strany rovnice v mezích od 0 po t. S využitím počáteční podmínky v(0) = v0 dostaneme v(t) = v0 + t 0 u MT - d = v0 - u [ln |MT - |]t =0 = v0 + u ln MT MT - t = = v0 + u ln 1 + t MT - t . 51 r p q . P x0 y0 O y = y(x) Obrázek 6.2: K Archimédově úloze: y = y(x) ­ zrcadlo, p ­ přicházející paprsek, q ­ odražený paprsek, P ­ bod dopadu a odrazu paprsku, O ­ ohnisko, ­ tečna k zrcadlu v bodě dopadu přicházejícího paprsku, ­ normála k zrcadlu, ­úhel odrazu. Zejména pro t = T máme v(T ) = v0 + u ln 1 + M - . (6.9) Tato formule se nazývá Ciolkovského rovnice. Rovnici (6.9) odvodil William Moore ve výzkumné zprávě A Treatise on the Motion of Rockets pro Royal Military Academy, Woolwich, England, v roce 1813. Tato práce byla zapomenuta a nezávisle na ní rovnici objevil roku 1898 Konstantin Eduardovič Ciolkovskij. S její pomocí v článku Ó ÓÚ× ¸ Á× ÓÚ ÖÓÚÝ ÔÖÓ×Ö ×Ú Ö ÚÝ ÔÖ ÓÖ Ó Ó ÓÖ 1903, Ó X, No. 5 zdůvodnil, že rakety mohou létat naprosto nezávisle na okolním prostředí, a proto mohou být vhodným prostředkem pro lety do vesmíru. 6.1.3 Archimédova úloha Určete tvar zrcadla, které odráží rovnoběžné světelné paprsky do jediného bodu (ohniska). Zvolíme souřadnou soustavu tak, aby ohnisko bylo v jejím počátku O, přicházející paprsky byly rovnoběžné se svislou osou a směřovaly proti její orientaci (kreslete si obrázek 6.2). Uvažujme přicházející paprsek p, který se od zrcadla odráží v libovolném, ale pevně zvoleném bodě P = (x0, y0), x0 > 0, y0 < 0. Nechť tvar zrcadla je v okolí tohoto bodu popsán funkcí y = y(x); přitom samozřejmě y(x0) = y0. Označme , resp. , tečnu, resp. normálu, k zrcadlu v bodě P, q přímku incidentní s odraženým paprskem PO, r vodorovnou přímku procházející bodem P. Nechť dále = q je úhel, který svírá odražený paprsek s normálou . Úhel odrazu se rovná úhlu dopadu a tedy p = . Odtud plyne, že p = 1 2 -. Dále platí r = 1 2 -p = . Poněvadž je tečnou ke křivce o rovnici y = y(x), platí dy dx (x0) = tg(r) = tg . (6.10) Poněvadž přímky p a r jsou kolmé, je qr = 1 2 - p - q = 1 2 - 2 a tedy tg(qr) = tg 2 - 2 = cotg(2) = 1 - (tg )2 2 tg . (6.11) Současně tg(qr) = |y0| x0 = - y0 x0 . (6.12) 52 Spojením (6.10), (6.11) a (6.12) dostaneme rovnost - y0 x0 = 1 - dy dx (x0) 2 2 dy dx (x0) . Poněvadž bod P = (x0, y0) byl libovolný, dostáváme pro tvar zrcadla diferenciální rovnici dy dx 2 - 1 = 2 y x dy dx . To je rovnice nerozřešená vzhledem k derivaci. Jedná se však o jednoduchou kvadratickou rovnici pro neznámou derivaci, takže ji můžeme vyjádřit ve tvaru dy dx = y x y x 2 + 1 . Pro y < 0 a x > 0 je dy dx > 0, viz obrázek 6.2. Znaménko před odmocninou tedy musí být +. Dostáváme tak diferenciální rovnici pro tvar požadovaného zrcadla dy dx = y x + y x 2 + 1 . To je rovnice homogenní. Substitucí u = u(x) = y(x) x , tedy y(x) = xu(x), dy dx = u + x du dx dostaneme rovnici se separovanými proměnnými. Její řešení v implicitním tvaru je du u2 + 1 = dx x , tedy ln u + u2 + 1 = ln |x| + const. Odtud u(x) = 1 2 x C - C x , kde C je integrační konstanta. V původních proměnných dostaneme rovnost y(x) = 1 2 x2 C - C , neboli x2 = C(C + 2y). To je rovnice paraboly s ohniskem (0, 0) a řídící přímkou x = -C. Název ,,Archimédova úloha vychází z tradované historky, podle níž Archimédes při obléhání Syrakus armádou římského vojevůdce Marcella v letech 214­212 př. n. l. z vyleštěných štítů obránců města sestavoval zrcadla, kterými soustředil sluneční paprsky a tak zapaloval prosmolené lodě obléhatelů. 6.1.4 ,,Psí křivka Pes pronásleduje zajíce. Zajíc se pohybuje rovnoměrně přímočaře rychlostí u, pes běží ve směru k zajíci rovnoměrnou rychlostí v, v > u. Určete tvar dráhy psa a čas T , za který pes zajíce dohoní. Zvolíme orthonormální souřadnou soustavu tak, aby se zajíc pohyboval po druhé ose souhlasně s její orientací a na počátku, tj. v čase t0 = 0, se zajíc nacházel v bodě (0, b) a pes v bodě (-a, 0). Nechť pro určitost je a > 0; případ a = 0 je triviální a v případě a < 0 bude tvar dráhy zřejmě obrazem tvaru pro a > 0 v osové symetrii kolem druhé souřadné osy. 53 -a x b y b + ut u v y x b -a a) b) y x -a y x -a b c) d) Obrázek 6.3: a) K odvození rovnice ,,psí křivky . Vektor rychlosti zajíce u má v každém okamžiku souřadnice (0, u), vektor rychlosti psa v má v každém okamžiku velikost v a v čase t směřuje k zajíci, tj. je rovnoběžný s vektorem o souřadnicích (|x|, b + ut - y). b) ,,Psí křivka pro a < 0, b > 0. c) ,,Psí křivka pro a > 0, b = 0. d) ,,Psí křivka pro a > 0, b < 0. 54 Situace je znázorněna na obr. 6.3 a). Dráhu psa vyjádříme jako funkci y = y(x). V čase t = 0 je x = -a a y = 0, tj. y(-a) = 0. (6.13) Pes k zajíci směřuje od začátku, tj. y (-a) = b a . (6.14) V jistém čase t, t < T , se pes nachází v bodě (x, y), x (-a, 0), a pes v bodě (0, b + ut). Poněvadž pes stále směřuje k zajíci, platí y (x) = b + ut - y(x) |x| = y(x) - b - ut x , neboli ut = y - xy (x) - b. (6.15) Za čas t urazí pes dráhu délky vt. Této hodnotě tedy musí být rovna délka křivky (grafu funkce) y = y(x) od bodu (-a, 0) po bod (x, y), tedy vt = x -a 1 + y() 2 d. Z této rovnosti vyjádříme t a dosadíme do (6.15), v u x -a 1 + y() 2 d = y(x) - xy (x) - b. (6.16) Označíme s = u v . (6.17) Podle předpokladu je s < 1. Obě strany rovnosti (6.16) zderivujeme podle x. Dostaneme s 1 + y(x) 2 = y (x) - y (x) - xy (x) a po úpravě xy (x) + s 1 + y(x) 2 = 0. (6.18) Dráha psa je tedy řešením neautonomní nelineární diferenciální rovnice druhého řádu (6.18) s počátečními podmínkami (6.13), (6.14). Rovnice (6.18) je typu 2.10. Proto zavedeme novou neznámou funkci p = p(x) = y (x). Dosadíme ji do rovnice (6.18) a počáteční podmínky (6.14). Po snadné úpravě dostaneme počáteční úlohu p = - s x 1 + p2 , p(-a) = b a . Jedná se o rovnici se separovanými proměnnými. Řešení úlohy v implicitním tvaru tedy podle 2.3 je p b a d 1 + 2 = -s x -a d . Integrací dostaneme ln a p + 1 + p2 b + a2 + b2 = ln - a x s 55 a odtud p = 1 2C C2 - a x s - - x a s , kde C = b a + 1 + b a 2 . (6.19) Poněvadž p = y a funkce y splňuje podmínku (6.13), dostaneme řešení úlohy integrací poslední rovnosti, tedy y(x) = 1 2C x -a C2 - a s - - a s d = = Ca 2(1 - s) 1 - - x a 1-s - a 2C(1 + s) 1 - - x a 1+s . Za konstanty s a C dosadíme z rovností (6.17) a (6.19). Po úpravách dostaneme ,,psí křivku ve tvaru y(x) = v vb + u a2 + b2 v2 - u2 - v 2 b - a2 + b2 v + u x a 1+ u v + b + a2 + b2 v - u x a 1- u v . Nalezená funkce y je sudá, vyjadřuje tedy tvar dráhy psa pro a > 0 i pro a < 0; v prvním případě bychom za definiční obor považovali interval [-a, 0], ve druhém interval [0, -a]. Pes dostihne zajíce v bodě 0, y(0) , To znamená, že zajíc rychlostí u urazí dráhu délky y(0)-b a čas, za který pes zajíce dohoní, je tedy roven T = y(0) - b u = 1 u v vb + u a2 + b2 v2 - u2 - b = ub + v a2 + b2 v2 - u2 . ,,Psí křivku (,,courbe chien ) jako první studoval v roce 1732 francouzský matematik Pierre Bouger (který je známější jako účastník expedice do Peru v roce 1735, která změřila délku jednoho stupně zeměpisné délky na rovníku). Křivka je nejjednodušším případem křivek sledování (pursuit curves, pojem poprvé použil George Boole ve svém spisu ,,Treatise on Differential equations v roce 1859), které jsou definovány takto: jestliže body A a P se pohybují rovnoměrně, bod A po dané křivce a směr pohybu bodu P stále míří k bodu A, pak bod P opisuje křivku sledování. Úloha bývá někdy formulována tak, že pes sleduje svého pána, nebo že liška honí králíka. 6.1.5 Ekonomický růst (Solowův-Swanův neoklasický model) Budeme se snažit popsat dynamiku (vývoj v čase) základních makroekonomických ukazatelů v uzavřené ekonomice, tj. v ekonomice, v níž nedochází k žádné výměně s okolními ekonomikami (k exportu nebo importu). Za základní ukazatele budeme považovat: Y = Y (t) . . . hrubý domácí produkt v čase t. K = K(t) . . . kapitál v čase t. Kapitálem budeme rozumět nejen kapitál finanční, tj. peníze, ale také kapitál hmotný, tj. budovy, stroje, zařízení ap. Celkové množství kapitálu lze však vyjádřit pomocí peněžní jednotky. L = L(t) . . . disponibilní pracovní síla v čase t. Lze si ji představit jako množství práceschopného (nebo práceochotného) obyvatelstva. I = I(t) . . . investice v čase t, tj. peněžní prostředky použité k tvorbě nebo obnově kapitálu. S = S(t) . . . spotřeba v čase t. Za spotřebu budeme považovat peněžní prostředky k tvorbě nebo obnově kapitálu nevyužité, tj. nejen realizovanou spotřebu ale také např. úspory obyvatel- stva. 56 Vyjdeme z několika jednoduchých a z ekonomického hlediska přijatelných předpokladů: P1) Jedinými produkčními faktory jsou kapitál a práce. P2) Kapitál se vytváří investicemi. P3) Kapitál se amortizuje (znehodnocuje) tak, že poměr znehodnoceného kapitálu za jednotku času ke všemu kapitálu je konstantní. P4) Relativní přírůstek pracovní síly v čase je konstantní; v podstatě odpovídá přirozenému přírůstku obyvatel. P5) Sklon ke spotřebě, tj. podíl spotřebovaného produktu, je v čase konstantní. P6) Veškerý produkt se rozdělí na investice a spotřebu. P2) a P6) nejsou ve vlastním smyslu předpoklady, je jimi pouze specifikováno, co se rozumí pojmy ,,investice a ,,spotřeba . Základní makroekonomické ukazatele budeme považovat za nezáporné diferencovatelné funkce definované na intervalu [0, ), tj. zajímá nás vývoj od jistého okamžiku do budoucnosti. Nyní můžeme předpoklady vyjádřit matematicky: P1) Produkce je funkcí kapitálu a práce, tj. Y (t) = f K(t), l(t) , (6.20) kde f : [0, )2 [0, ) je nějaká funkce rostoucí v každé ze svých proměnných. Nazýváme ji produkční funkce. P2) Množství kapitálu vytvořeného za časový interval délky t, tj. od času t po čas t + t, je rovno I(t)t. P3) Označme a množství kapitálu amortizovaného za jednotku času. Pak pro každý čas t platí a K(t) = , kde je nějaká nezáporná konstanta, kterou nazveme mírou amortizace. Ta bývá vyjádřena pomocí odpisů. Za časový interval délky t se amortizuje at kapitálu, takže množství kapitálu znehodnoceného za interval délky t je rovno K(t)t. P4) Relativní přírůstek pracovní síly za jednotku času je konstantní, takže relativní přírůstek pracovní síly za časový interval délky t je roven t, tj. L(t + t) - L(t) L(t) = t. P5) Existuje konstanta s nazývaná mezní sklon ke spotřebě taková, že pro každé t 0 je S(t) Y (t) = s. P6) V každém čase t platí Y (t) = I(t) + S(t). (6.21) Dále potřebujeme specifikovat produkční funkci f. Budeme tedy ještě předpokládat: P7) Ke zdvojnásobení produkce je potřeba dvojnásobného kapitálu i dvojnásobné pracovní síly. Obecněji, ke zvětšení produkce o q% je potřeba zvětšit kapitál i pracovní sílu také o q%. Matematicky, f(K, L) = f(K, L) (6.22) pro každou konstantu > 0; jinak řečeno, produkční funkce f je homogenní prvního řádu. 57 P8) Není-li v ekonomice žádný kapitál, tak jakékoliv jeho přidání způsobí veliký nárůst produkce; přesněji, mezní produkt kapitálu roste do nekonečna, pokud se množství kapitálu v ekonomice přibližuje k nule. Je-li v ekonomice nadbytek kapitálu, tak jeho zvětšení již nezpůsobí významný nárůst produkce; přesněji, mezní produkt kapitálu klesá k nule, pokud jeho množství v ekonomice neomezeně roste. Analogické vztahy má produkce a práce. Tyto předpoklady zapíšeme ve tvaru lim K0+ f K (K, L) = lim L0+ f L (K, L) = , lim K f K (K, L) = lim L f L (K, L) = 0; nazýváme je Inadovy podmínky. Nyní již můžeme sestavit rovnice popisující vývoj makroekonomických ukazatelů. Podle P2) a P3) je kapitál v čase t + t roven K(t + t) = K(t) + I(t)t - K(t)t. Odtud dostaneme K(t + t) - K(t) t = I(t) - K(t) a limitním přechodem t 0 získáme diferenciální rovnici K = I - K. (6.23) Tímtéž limitním přechodem dostaneme po snadné úpravě z předpokladu P4) diferenciální rovnici L = L. (6.24) Z předpokladu P6) máme 1 = I(t) Y (t) + S(t) Y (t) a dále s využitím předpokladu P5) dostaneme I(t) = (1 - s)Y (t). (6.25) Systém dvou diferenciálních rovnic (6.23), (6.24) spolu s omezujícími rovnostmi (6.20), (6.25) lze považovat za matematický model dynamiky produkce, kapitálu a práce. Model (6.23), (6.24), (6.20), (6.25) můžeme ještě dále upravit. Zavedeme veličinu k = k(t) vztahem k(t) = K(t) L(t) ; (6.26) nazýváme ji míra vybavenosti práce kapitálem. Z rovností (6.20) a (6.22) dostaneme Y (t) = f K(t), L(t) = L(t)f K(t) L(t) , 1 = L(t)f k(t), 1 . Výraz f(k, 1) se nazývá intenzivní tvar produkční funkce. Derivováním vztahu K = Lk definujícího vybavenost práce kapitálem dostaneme K = L k + Lk . Po dosazení z rovnic (6.24), (6.23) máme I - K = Lk + Lk . Za proměnnou I nejprve dosadíme z rovnosti (6.25) a pak za proměnnou Y z rovnosti (6.20). Dostaneme tak (1 - s)Lf(k, 1) - K = Lk + Lk . Po vydělení hodnotou L máme (1-s)f(k, 1)-k = k+k . Odtud vyjádříme derivaci k a dostaneme základní dynamickou rovnici neoklasického modelu k = -( + )k + (1 - s)f(k, 1). (6.27) V základní rovnici stále zůstává neurčená produkční funkce f. Jednoduchá funkce, která splňuje podmínky P7) a P8) (a tedy může představovat jistý popis ekonomické reality) je CobbovaDouglasova produkční funkce, která je tvaru f(K, L) = AK L1, (6.28) 58 kde kladný koeficient A představuje produkci při jednotkovém kapitálu i práci a je nějaká konstanta taková, že 0 < < 1. Cobbova-Douglasova produkční funkce v intenzivním tvaru je f(k, 1) = Ak . (6.29) Dosazením této funkce do základní rovnice (6.27) dostaneme obyčejnou diferenciální rovnici pro neznámou funkci k ve tvaru k = -( + )k + A(1 - s)k . (6.30) To je rovnice Bernoulliova, sr. 2.7. Budeme ji řešit zavedením substituce r = k1, (6.31) tedy r = (1 - )k- k = (1 - )k(-( + )k + A(1 - s)k ) = = ( - 1)( + )k1- A(1 - s)( - 1), neboli r = ( - 1)( + )r - A(1 - s)( - 1). (6.32) To je nehomogenní lineární rovnice pro neznámou funkci k a její řešení s počáteční podmínkou r(0) = r0 je podle 2.6 rovno r(t) = r0 - t 0 A(1 - s)( - 1)e-(-1)(+) d e(-1)(+)t = = r0 + A 1 - s + e-(-1)(+)t - 1 e(-1)(+)t = = r0 - A 1 - s + e(-1)(+)t + A 1 - s + . (6.33) Míra amortizace je podle předpokladu P3) nezáporná. O relativním přírůstku pracovní síly jsme dosud nic nepředpokládali, může být kladný (obyvatelstva přibývá) i záporný (obyvatelstvo vymírá). Pro další úvahy budeme předpokládat, že + > 0 (materiál chátrá rychleji než obyvatelstvo). Poněvadž < 1, je funkce r daná rovností (6.33) monotonní (v případě (+)r0 > A(1-s) klesající, v případě ( + )r0 < A(1 - s) rostoucí) a platí lim t r(t) = A 1 - s + . (6.34) S využitím rovností (6.26), (6.29) a (6.31) můžeme psát r = k1- = K L 1- = AK AKL1= A K Y . (6.35) Z rovnosti (6.34) nyní plyne lim t K(t) Y (t) = 1 - s + ; (6.36) tento výsledek můžeme interpretovat tak, že kapitálová náročnost jednotky produkce (množství kapitálu potřebné k vytvoření jednotkového produktu) se v uzavřené ekonomice ustálí na jisté hodnotě, zvané mezní poměr kapitálu a produkce, která závisí pouze na sklonu ke spotřebě, míře amortizace a přirozeném přírůstku (nebo úbytku) obyvatel. Ve stabilizované uzavřené ekonomice tedy platí K Y = 1 - s + , neboli Y = 1 - s + K; 59 produkce je přímo úměrná kapitálu. Návratem k proměnné k = r1/(1+) můžeme vyjádřit řešení základní rovnice (6.30) s CobbovouDouglasovou produkční funkcí a s počáteční podmínkou k(0) = k0 ve tvaru k(t) = k1- 0 - A 1 - s + e(-1)(+)t + A 1 - s + 1 1. (6.37) Pro funkci k platí, že lim t k(t) = A 1 - s + 1 1; (6.38) její chování v dlouhém časovém úseku nezávisí na počáteční hodnotě. Ekonomika tedy směřuje ke konstantní (rovnovážné) vybavenosti práce kapitálem. Nyní můžeme pomocí řešení (6.33) a (6.37) rovnic (6.32) a (6.30) vyjádřit řešení rovnic původního modelu (6.20), (6.23), (6.24), (6.25), (6.21) v případě, že produkční funkce je CobbovaDouglasova tvaru (6.28). Budeme uvažovat počáteční podmínky K(0) = K0, L(0) = L0. Pak podle (6.28), (6.21) a (6.25) je Y (0) = Y0 = AK 0 L1- 0 , S(0) = S0 = sY0, I(0) = I0 = (1 - s)Y0. Rovnice (6.24) je lineární homogenní a její řešení splňující uvedenou počáteční podmínku je L(t) = L0et . Podle (6.26) je K(t) = k(t)L(t) = L0et K0 L0 1- A 1 - s + e(-1)(+)t + A 1 - s + 1 1- , podle (6.28), (6.20) je Y (t) = AK(t) L(t)1= AL0et K0 L0 1- A 1 - s + e(-1)(+)t + A 1 - s + 1- . Spotřeba S(t) je s-násobkem produkce a investice I(t) je jejím (1 - s)-násobkem. Z rovností (6.38), (6.36) plyne, že produkce, kapitál a pracovní síla jsou asymptoticky ekvivalentní funkce; zhruba řečeno, tyto makroekonomické charakteristiky rostou stejně rychle. Zejména platí lim t Y (t) K(t) = + 1 - s = const < . Základní dogma ekonomie však říká, že ekonomika roste tak, že se produkuje stále více se stále menšími náklady, tj. lim t Y (t) K(t) = . (6.39) Tato disproporce může být způsobena tím, že jsme neuvažovali technologický pokrok. Ten se projevuje tak, že efektivita práce v čase roste. To znamená, že pracovní síla nebude vyjádřena pouze množstvím pracujících. Předpoklad P4) tedy nahradíme předpokladem modifikovaným: P4') Relativní přírůstek pracovní síly za jednotku času v průběhu času roste, tj. L(t + t) - L(t) L(t) = (t)t, kde je rostoucí funkce. 60 Nyní můžeme zopakovat všechny úvahy. Těmi dojdeme k modifikované rovnici neoklasického mo- delu k = - (t) + )k + (1 - s)f(k, 1), (6.40) která se od rovnice (6.27) liší pouze v tom, že relativní přírůstek pracovní síly závisí na čase. Abychom tuto závislost specifikovali, přijmeme další předpoklad: P9) Relativní přírůstek pracovní síly vyjadřuje dosaženou technologickou úroveň. Technologický růst je proces kumulativní, tj. jeho přírůstek je úměrný úrovni dosažené a době rozvoje, tj. (t + t) - (t) = p(t)t, kde p je kladná konstanta. Uvedený předpoklad lze limitním přechodem t 0 přepsat ve tvaru homogenní diferenciální rovnice = p, (6.41) jejíž řešení je podle 2.6 rovno (t) = 0ept , kde 0 vyjadřuje počáteční technologickou úroveň. Do rovnice (6.40) nyní můžeme dosadit Cobbovu-Douglasovu produkční funkci (6.29) a vyjádření (6.41). Opět dostaneme Bernoulliovu rovnici k = - 0ept + k + A(1 - s)k , kterou substituce (6.31) transformuje na rovnici lineární r = ( - 1) 0ept + r + A(1 - s)(1 - ). Její řešení s počáteční podmínkou r(0) = r0 je podle 2.6 rovno r(t) = r0 + A(1 - )(1 - s) t 0 exp (1 - ) + 0 p (ep - 1) d × × exp -(1 - ) t + 0 p ept - 1 . Vzhledem k (6.35) je Y (t) K(t) = A r(t) , takže s použitím de l'Hospitalova pravidla můžeme vypočítat lim t Y (t) K(t) = A lim t exp (1 - ) t + 0 p (ept - 1) r0 + A(1 - )(1 - s) t 0 exp (1 - ) + 0 p (ept - 1) d = = A lim t (1 - ) + 0 p ept exp (1 - ) t + 0 p (ept - 1) A(1 - )(1 - s) exp (1 - ) t + 0 p (ept - 1) = lim t p + 0ept p(1 - s) = . Podmínka (6.39) je nyní splněna. Analogicky vypočítáme lim t K(t) Y (t) /0ept = 0. 61 Poslední výsledek lze interpretovat tak, že v uzavřené ekonomice s plynulým technologickým pokrokem klesá kapitálová náročnost jednotky produkce řádově rychleji, než technologická úroveň roste. Základní dynamickou rovnici neoklasického modelu sestavili nezávisle na sobě Robert M. Solow a Trevor W. Swan jako rozšíření modelu produktivity kapitálu, který nezávisle vytvořili Sir Roy F. Harrod (v roce 1939) a Ewsey Domar (v roce 1946). Publikovali ji v článcích Solow, R. M. A Contribution to the Theory of Economic Growth. Quaterly Journal of Economic. 1956, vol. 70, No. 1 (February), p. 65­94. Swan, T. W. Economic Growth and Capital Accumulation. Economic Record. 1956, No. 32 (November), p. 334­361. Robert Solow za rozpracování neoklasického modelu obdržel v roce 1987 Cenu Švédské národní banky za rozvoj ekonomické vědy na památku Alfreda Nobela (lidově zvanou Nobelova cena za ekonomii). 6.2 Epidemiologický model Daniela Bernoulliho Uvažujme chorobu, která trvá krátce, někteří pacienti na ni zemřou, jiní se uzdraví a získají vůči nákaze imunitu; typickým představitelem takové infekce byly neštovice. Budeme modelovat epidemii této choroby, tj. její šíření v nějaké kohortě. Kohortou rozumíme skupinu osob narozených ve stejnou dobu. Zavedeme označení: N počet osob tvořících kohortu, a jejich věk (tj. čas od počátku), S = S(a), resp. R = R(a) počet osob věku a, které neprodělaly, resp. prodělaly, chorobu. Další symboly zavedeme na základě následujících předpokladů: * Počet osob věku a, které zemřou z jiných příčin, než je uvažovaná infekce, je úměrná délce (krátkého) časového intervalu sledování a a počtu nenakažených osob S(a). Konstantu úměrnosti -- přirozenou úmrtnost ve věku a -- označíme (a). * Počet osob věku a, které se nakazí uvažovanou chorobou je úměrná délce sledování a a počtu nenakažených osob S(a). Konstantu úměrnosti -- incidenci choroby ve věku a -- označíme (a) * Počet nemocných osob věku a, které se uzdraví za časový interval a je úměrný počtu infikovaných osob tohoto věku a délce intervalu a. Konstantu úměrnosti -- ukazatel přežití choroby osobami věku a -- označíme s(a). Úmrtnost (a) lze interpretovat jako pravděpodobnost, že ,,zdravá osoba (tj. ta, která nemá uvažovanou chorobu) věku a zemře během časového intervalu délky a; incidenci (a) jako pravděpodobnost, že se ,,zdravá osoba věku a nakazí během časového intervalu délky a; ukazatel přežití s(a) jako pravděpodobnost, že nakažená osoba věku a se během časového intervalu délky a uzdraví. Budeme předpokládat, že onemocnění a uzdravení jsou jevy nezávislé, tj. že pravděpodobnost, že ,,zdravá osoba se během časového intervalu délky a nakazí a uzdraví, je rovna s(a)(a). Dále zavedeme letalitu choroby ve věku a vztahem c(a) = 1 - s(a); lze ji interpretovat jako pravděpodobnost, že nemocná osoba věku a během časového intervalu délky a zemře. Proměnné u = u(a) = S(a) N , resp. w = w(a) = R(a) N vyjadřují (klasickou) pravděpodobnost, že osoba se dožila věku a a neprodělala, resp. prodělala, chorobu. Novorozenec určitě chorobu neprodělal, tedy platí u(0) = 1, w(0) = 0. (6.42) 62 ? - ? ? &% '$ nemocní S(a) náchylní k chorobě R(a) imunní (a) (a)c(a) (a) s(a) Obrázek 6.4: Schema vývoje kohorty ohrožené chorobou Vývoj kohorty, v níž probíhá choroba, lze schematicky znázornit obrázkem 6.2 a předpoklady vyjádřit ve tvaru rovností: S(a + a) = S(a) - (a)S(a)a - (a)S(a)a = S(a) - (a) + (a) S(a)a, R(a+a) = R(a)+s(a)(a)S(a)a-(a)R(a)a = R(a)+ 1-c(a) (a)S(a)a-(a)R(a)a. V první z uvedených rovností převedeme na levou stranu S(a) a ve druhé z nich R(a), rovnosti vydělíme výrazem Na a provedeme limitní přechod a 0. Pro zjednodušení modelu budeme předpokládat, že funkce u a w jsou diferencovatelné; takový předpoklad je v případě velké kohorty dostatečně realistický. Dostaneme tak systém neautonomních diferenciálních rovnic du da = - (a) + (a) u, dw da = - 1 - c(a) (a)u - (a)w; (6.43) jejich řešení splňuje počáteční podmínky (6.42). První rovnice systému (6.43) je lineární homogenní rovnicí pro neznámou funkci u. Její řešení s počáteční podmínkou (6.42) je při označení M(a) = a 0 ()d, (a) = a 0 ()d (6.44) podle 2.6 rovno u(a) = e-(a)-M(a) . (6.45) Toto vyjádření dosadíme do druhé rovnice systému (6.43) a dostaneme dw da = -(a)w - 1 - c(a) (a)e-(a)-M(a) , což je lineární nehomogenní rovnice pro neznámou funkci w. Její řešení s počáteční podmínkou (6.42) je opět podle 2.6 rovno w(a) = e-M(a) 1 - a 0 c()()e-() d - e-(a)-M(a) . (6.46) Dosud provedené úvahy a výpočty lze shrnout: Pravděpodobnosti u(a), resp. w(a), že se osoba dožije věku a a neprodělá, resp. prodělá, chorobu, jsou řešením soustavy rovnic (6.43) s počátečními podmínkami (6.42) a jsou dány výrazy (6.45), resp. (6.46), kde funkce M a jsou dány výrazy (6.44). Pravděpodobnost, že se osoba dožije věku a za předpokladu, že choroba se v kohortě neobjevuje (tj. 0 a v důsledku toho také 0), je rovna 0(a) = u(a) = e-M(a) . 63 Pravděpodobnost, že se osoba dožije věku a pokud se choroba vyskytuje, je rovna (a) = u(a) + w(a) = e-M(a) 1 - a 0 c()()e-() d = 0(a) 1 - a 0 c()()e-() d . Pravděpodobnost dožití věku a je tedy součinem pravděpodobnosti dožití věku a při přirozené úmrtnosti a faktoru, který závisí pouze na incidenci a letalitě choroby. Označme dále x(a) = u(a) (a) , z(a) = w(a) (a) = (a) - u(a) (a) = 1 - x(a); Veličina x(a), resp. z(a), vyjadřuje podmíněnou pravděpodobnost, že osoba věku a neprodělala, resp. prodělala, chorobu za podmínky, že se věku a dožila. Poněvadž x(a) = e-(a)-M(a) e-M(a) 1 - a 0 c()()e-()d = e-(a) 1 - a 0 c()()e-()d , platí rovnost x(0) = 1 a dále dx(a) da = -(a)e-(a) 1 - a 0 c()()e-() d + e-(a) c(a)(a)e-(a) 1 - a 0 c()()e-()d 2 = = -(a) e-(a) 1 - a 0 c()()e-()d - c(a)e-2(a) 1 - a 0 c()()e-()d 2 = = -(a)x(a) 1 - c(a)x(a) . Relativní zastoupení osob věku a, které v uvažované kohortě neprodělaly chorobu, je tedy veličina x(a), která je řešením počáteční úlohy pro Bernoulliovu rovnici dx da = -(a)x 1 - c(a)x , x(0) = 1. 6.3 Model populace produkující škodlivé odpady Označme N = N(t) velikost nějaké populace v čase t. Specifická míra růstu nebo růstový koeficient p této populace je definován jako relativní změna velikosti populace, tj. p = N N . Vývoj populace je tedy modelován diferenciální rovnicí N = pN. (6.47) V případě konstantního růstového koeficientu dostaneme klasický Malthusův1 model růstu populace N(t) = N0ept , kde N0 = N(0) je počáteční velikost populace. V něm je exponenciální růst (pro p > 0) nebo úbytek (pro p < 0) velikosti populace nerealistický. 1Správněji malthusovský, Thomas Robert Malthus (1766­1834) model v takovém tvaru nikdy nepublikoval. 64 Model (6.47) se přiblíží realitě, pokud specifickou míru růstu p nebudeme považovat za nezávislou konstantu populace, ale za veličinu závislou na její velikosti, tedy p = p(N), nebo obecněji na nějakých ,,projevech její velikosti, tj. p = p F(N) , kde F je nějaký funkcionál, tedy zobrazení z množiny funkcí do množiny reálných čísel. V tomto oddílu budeme uvažovat populaci, která produkuje odpady svého metabolismu, které jsou toxické, nebo přinejmenším zmenšují schopnost přežívání populace. Tyto odpady se v prostředí hromadí, ale také rozkládají, mizí nebo přeměňují v něco, co populaci neomezuje. Budeme tedy předpokládat: 1. V čistém prostředí (bez uvažovaných odpadů) je specifická míra růstu rovna nějaké konstantě r (vnitřnímu koeficientu růstu, intrinsic growth rate). 2. V každém okamžiku populace produkuje odpad, jehož množství je úměrné velikosti populace. Množství odpadu vyprodukovaného v čase t označíme Pp(t); platí pro něho Pp(t) = cN(t), kde c je nějaká kladná konstanta. 3. Odpad se rozkládá konstantní relativní rychlostí > 0, tj. označíme-li P(t) množství odpadu v čase t a neuvažujeme jeho produkci, platí P (t) = -P(t). (6.48) 4. Specifická míra růstu populace klesá s rostoucím množstvím odpadu. Budeme uvažovat nejjednodušší možnost, že tato závislost je lineární. 5. Existuje jistá velikost populace K > 0, při které je populace se svým prostředím v dynamické rovnováze, její velikost se v čase nemění. Konstanta K představuje kapacitu prostředí (úživnost) pro uvažovanou populaci. Uvažujme na chvíli idealizovanou situaci, že pouze v čase s vzniklo množství Pp(s) odpadu a žádný další odpad není do prostředí dodáván. Množství odpadu v čase t > s tedy bude podle předpokladů 2. a 3. řešením rovnice (6.48) s počáteční podmínkou P(s) = Pp(s) = cN(s), tj. P(t) = cN(s)e-(t-s) . V reálné situaci se však odpad v prostředí kumuluje, v čase t ho tedy bude množství, které zůstalo ze všech odpadů vzniklých až do okamžiku t, tj. množství odpadu závislé na celé předchozí historii velikosti populace bude F(N) = t - cN(s)e-(t-s) ds. Předpoklad 4. lze nyní přepsat ve tvaru p = p F(N) = - F(N), kde > 0. Z předpokladu 1. plyne, že p(0) = r, tj. = r. Pro funkci ~N = ~N(t) K podle předpokladu 5. nyní platí 0 = p F( ~N) = r - F( ~N) = r - c t - Ke-(t-s) ds = r - cK e-(t-s) t s== r - cK . Odtud dostaneme, že c = r K a specifická míra růstu populace je p = r 1 - K t - N(s)e-(t-s) ds . 65 Model (6.47) je tedy nyní ve tvaru integrodiferenciální2 rovnice N (t) = rN(t) 1 - K t - N(s)e-(t-s) ds . (6.49) Zavedeme nové neznámé funkce x a y nové nezávisle proměnné následujícími vztahy: x() = rK N r , y() = K r - N(s)e-( r -s)ds. Pak x () = dx() d = rK N r 1 r = r2K rN r 1 - K r - N(s)e-( r -s)ds = = rK rK x() 1 - K r - N(s)e-( r -s)ds = x() 1 - y() , y () = dy() d = K d d r - N(s)e-( r -s)ds = = K 1 r N r e-( r - r ) - r r - N(s)e-( r -s)ds = = rK N r - r K r - N(s)e-( r -s)ds = x() - r y(). Rovnice (6.49) se tedy transformuje na dvourozměrný autonomní systém x = x(1 - y), y = x - r y. (6.50) Nejprve si všimněme, že uzavřený první kvadrant R2 + = (x, y) R2 ; x 0, y 0 je pozitivně invariantní množinou tohoto systému. Uzavřená polopřímka {(0, y) : y 0} je totiž pozitivně invariantní množinou (x(t) 0, y(t) = y0e-(/r)t je řešením systému (6.50) pro každé y0 0), pro řešení s počáteční podmínkou x(0) = x0 > 0, y(0) = 0 platí x (0) > 0, y (0) > 0 a tedy příslušná trajektorie směřuje dovnitř prvního kvadrantu. Systém (6.50) má stacionární body(0, 0) a (x , y ) = r , 1 a jeho variační matice je J(x, y) = 1 - y -x 1 - r . 2V této rovnici vystupuje neznámá funkce N za znakem integrálu i jako derivovaná. 66 Tedy J(0, 0) = 1 0 1 - r , det J(0, 0) = - r > 0, takže stacionární bod (0, 0) je sedlo. Dále J(x , y ) = 0 - r 1 - r , det J(x , y ) = r > 0, tr J(x , y ) = - r < 0, tr J(x , y ) 2 - 4 det J(x , y ) = r2 ( - 4r), takže v případě 4r je vnitřní stacionární bod (x , y ) stabilní uzel, v opačném případě se jedná o stabilní ohnisko. S využitím Dulacova kriteria 5.2.6 vyloučíme existenci cyklu v prvním kvadrantu. Položíme q(x, y) = 1 x . Pak x 1 x x(1 - y) + y 1 x x - r y = - rx < 0 pro všechna x > 0. Uvnitř prvního kvadrantu tedy neexistuje cyklus systému (6.50). Uvažujme nyní situaci, že na počátku (v čase t = 0) se dostane malá populace do nového prostředí. K rovnici (6.49) přidáme tedy počáteční podmínky N(0) = N0, N(t) = 0 pro t < 0. (6.51) Počáteční podmínky pro systém (6.50) v tomto případě budou x(0) = rK N0, y(0) = K 0 - N(s)es ds = 0 a z analýzy systému (6.50) plyne, že pro řešení N počáteční úlohy (6.49), (6.51) platí lim t N(t) = rK x = K; funkce N konverguje k hodnotě K v případě 4r monotonně, v opačném případě s tlumenými oscilacemi. 6.4 Lotkovy-Volterrovy systémy x i = xi bi - n j=1 aijxj , i = 1, 2, . . ., n. (6.52) Tyto systémy modelují vývoj společenstva (časové změny velikostí jednotlivých populací, z nichž se společenstvo skládá). Neznámé funkce a parametry interpretujeme následovně: xi = xi(t) . . . velikost i-té populace bi . . . růstový koeficient izolované i-té populace (vnitřní koeficient růstu i-té populace) bi > 0 . . . i-tá populace je soběstačná (producent) bi 0 . . . i-tá populace závisí na jiných populacích (konzument) aii . . . koeficient vnitrodruhových vztahů i-té populace aii > 0 . . . v i-té populaci se projevuje vnitrodruhová konkurence 67 aii < 0 . . . v i-té populaci se projevuje vnitrodruhová kooperace aij . . . koeficient vlivu j-té populace na i-tou min {aij, aji} > 0 . . . i-tá a j-tá populace jsou ve vztahu konkurence max {aij, aji} < 0 . . . i-tá a j-tá populace jsou ve vztahu mutualismu (symbiózy) aij < 0 < aji . . . j-tá populace je kořistí (hostitelem) i-té populace; i-tá populace je predátorem (parazitem) j-té populace aij > 0 . . . j-tá populace je amenzalistou i-té populace aij < 0 . . . j-tá populace je komenzalistou i-té populace Fázový prostor systému (6.52) je n-rozměrný uzavřený kladný orthant Rn + = {(x1, x2, . . . , xn) Rn : x1 0, x2 0, . . . , xn 0} . 6.4.1 Příklad (Vztah Lotkových-Volterrových systémů a Verhulstovy logistické rovnice ) Logistickou rovnici x = rx 1 - x K , v níž jsou oba parametry r (vnitřní koeficient růstu) a K (kapacita prostředí pro modelovanou populaci) kladné, lze považovat za jednorozměrný Lotkův-Volterrův systém s b1 = r a a11 = r/K, tedy za model soběstačné populace s vnitrodruhovou konkurencí (tak byla Verhulstova rovnice sestavena). Také platí K = b1/a11; odtud lze usoudit, že pro soběstačnou populaci s vnitrodruhovou konkurencí představuje podíl vnitřního koeficientu růstu a koeficientu vnitrodruhové konkurence kapacitu prostředí neovlivněnou ostatními populacemi společenstva. Jinou interpretaci logistické rovnice lze získat následující úvahou: Označme y = 1 - x K = K - x K . Poněvadž y = -x /K, dostaneme x = rxy y = - r K xy. (6.53) Jedná se o dvojrozměrný Lotkův-Volterrův systém s parametry b1 = b2 = 0, a11 = a22 = 0, a12 = r, a21 = - r K . Proměnnou y lze interpretovat jako relativní dostupnost zdrojů pro modelovanou populaci vzhledem k celkové kapacitě prostředí K. Velikost populace a relativní dostupnost zdrojů jsou tedy ve vztahu predace, obě tyto ,,složky společenstva nejsou ani producenty ani konzumenty a neprojevuje se u nich žádný vnitrodruhový vztah. Poznamenejme, že systém (6.53) nemá izolované stacionární body. 6.4.2 Příklad (dissipativita konkurenčních systémů) Uvažujme společenstvo n soběstačných populací, z nichž každá projevuje vnitrodruhovou konkurenci a každá z populací je amenzalistou jiné nebo ji neovlivňuje (zejména tedy každé dvě populace mohou být ve vztahu konkurence). Vývoj takového společenstva lze modelovat systémem (6.52) s kladnými parametry bi, aii, i = 1, 2, . . ., n a snezápornými parametry aij pro i = j. S využitím tvrzení 5.1.10 ukážeme, že takový systém je dissipativní, tedy že všechny složky jeho řešení jsou ohraničené: 68 Nechť > 0 a i {1, 2, . . . , n} jsou libovolná. Položme Ki = bi aii + , i = aii. Pak Ki > 0, i > 0 a pro všechna xj Kj, j {1, 2, . . . , n} platí xi bi - n j=1 aijxj xi(bi - aiixi) xi(bi - aiiKi) = xiaii bi aii - Ki = = xiaii(Ki - - Ki) = -xiaii = -ixi, takže předpoklady třetího z tvrzení 5.1.10 jsou splněny. Poněvadž kladná konstanta je libovolně malá, pro každé řešení x( ) = x1(), x2(), . . . , xn() systému (6.52) s bi > 0, aii > 0, aij 0, i, j = 1, 2, . . . , n existuje T 0 takové, že pro všechna t T je x1(t) b1 a11 , x2(t) b2 a22 , . . . , xn(t) bn ann . V dlouhém časovém horizontu populace nepřekračují velikost danou kapacitou prostředí pro populace izolované. Zavedeme označení b = b1 ... bn , A = a11 a1n ... ... ... an1 ann , matice A se nazývá matice interakcí společenstva. Pro libovolný vektor v = (v1, v2, . . . , vn)T polo- žíme diag v = v1 0 0 0 0 v2 0 0 ... ... ... ... ... 0 0 vn-1 0 0 0 0 vn a vektory ze standardní orthonormální báze n-rozměrného vektorového prostoru označíme ej , ej = 1,j 2,j ... n,j , kde ij = 0, i = j 1, i = j je Kroneckerův symbol. Systém (6.52) lze zapsat jako vektorovou rovnici x = diag x (b - A x) . (6.54) Je-li matice A regulární, existuje nejvýše jeden stacionární bod x = A-1 b systému (6.52) takový, že všechny jeho složky jsou kladné. Takový stacionární bod budeme nazývat vnitřní. Pokud vnitřní stacionární bod existuje, lze tuto skutečnost interpretovat jako možnou koexistenci všech populací společenstva, přičemž koexistující populace mají dynamicky stálé velikosti dané složkami vektoru x . 69 6.4.3 Příklad (Společenstvo se dvěma trofickými úrovněmi) Uvažujme společenstvo tvořené dvěma skupinami druhů -- producenty (kořistí) a konzumenty (predátory). Mezi druhy uvnitř jednotlivých trofických úrovní nejsou žádné interakce a konzumenti nemohou bez producentů přežít. Je-li takové společenstvo tvořeno n druhy producentů a m druhy konzumentů, lze jeho vývoj popsat systémem Lotkových-Volterrových rovnic tvaru x i = xi ri - m k=1 aikyk , i = 1, 2, . . ., n, y j = yj -sj + n k=1 bjkxk , j = 1, 2, . . . , m; (6.55) xi označuje velikost i-tého druhu producentů, yj velikost j-tého druhu konzumentů, parametry ri, sj, aij, bji jsou kladné. Systém (6.55) můžeme při zavedení vektorů x = (x1, x2, . . . , xn)T , y = (y1, y2, . . . , ym)T , r = (r1, r2, . . . , rn)T , s = (s1, s2, . . . , sm)T , a matic A = aij 1in 1jm = a11 a12 . . . a1m a21 a22 . . . a2m ... ... ... ... an1 an2 . . . anm , B = bji 1jm 1in = b11 b12 . . . b1n b21 b22 . . . b2n ... ... ... ... bm1 bm2 . . . bmn zapsat vektorově x = diag x (r - A y), y = diag y (-s + B x). Stavové proměnné xi, yj transformujeme na nové, které označíme ui, vj a definujeme rovnostmi ui = ln xi, vj = ln yj, i = 1, 2, . . . , n, j = 1, 2, . . ., m. Pak u i = x i xi = ri - m k=1 aikyk = ri - m k=1 aikevk a podobně v j = -sj + n k=1 bjkeuk . Zavedeme-li zobrazení f = (f1, f2, . . . , fn)T : Rm Rn a g = (g1, g2, . . . , gm)T : Rn Rm vztahy fi(v) = fi(v1, v2, . . . , vm) = ri - m k=1 aikevk , i = 1, 2, . . . , n, gj(v) = gi(u1, u2, . . . , un) = -sj + n k=1 bjkeuk , j = 1, 2, . . . , m, můžeme transformovaný systém (6.55) zapsat ve tvaru u i = fi(vi, v2, . . . , vm), i = 1, 2, . . ., n, v j = gj(ui, u2, . . . , un), j = 1, 2, . . . , m nebo vektorově u = f(v), v = g(u); derivace první sady proměnných závisí pouze na druhé sadě, derivace druhé sady proměnných závisí pouze na první sadě. Systémy tohoto tvaru se nazývají bipartitní. Hodnota aij vyjadřuje množství i-tého druhu kořisti, kterou za jednotku času zničí predátoři j-tého druhu za předpokladu, že populace i-tého druhu kořisti i j-tého druhu predátora měly jednotkovou velikost. Stručněji, aij je specifická úmrtnost i-tého druhu kořisti způsobená populací 70 j-tého druhu predátora o jednotkové velikosti. Hodnota bji je specifická porodnost j-tého druhu predátora po konzumaci jednotkového množství populace i-tého druhu kořisti. Poměr bji/aij lze tedy chápat jako efektivitu, s jakou se úbytek i-tého druhu kořisti přeměňuje do růstu populace j-tého druhu predátora. Předpokládejme nyní, že každý druh predátora využívá všechny druhy kořisti stejně efektivně, tj. že ke každému j = 1, 2, . . . , m existuje konstanta cj > 0 taková, že bji aij = cj pro všechny indexy i = 1, 2, . . . , n. Jinak řečeno, nechť existuje vektor c = (c1, c2, . . . , cm)T pro nějž platí BT = A diag c, neboli B = diag c AT . Předpokládejme dále, že existuje vnitřní stacionární bod systému (6.55), tj. existují vektory p = (p1, p2, . . . , pn)T , q = (q1, q2, . . . , qm)T se všemi složkami kladnými, takové že A q = r, B p = s, tj. ri = m k=1 aikqk pro i = 1, 2, . . . , n, sj = n k=1 bjkpk pro j = 1, 2, . . . , m. (Poznamenejme, že v případě m = n nemusí být některý z vektorů p, q určen jednoznačně.) Definujme nyní funkci H : Rn+m + R předpisem H(x, y) = H(x1, x2, . . . , xn, y1, y2, . . . , yn) = n i=1 (pi ln xi - xi) + m j=1 1 cj (qj ln yj - yj) . Pokud x, y jsou řešením systému (6.55), která mají všechny složky v každém čase kladné, pak platí d dt H(x, y) = n i=1 pi x i xi - x i + m j=1 1 cj qj y j yj - y j = = n i=1 (pi - xi) ri - m k=1 aikyk + m j=1 1 cj -sj + n k=1 bjkxk (qj - yj) = = n i=1 (pi - xi) m k=1 aikqk - m k=1 aikyk + m j=1 1 cj - n k=1 bjkpk + n k=1 bjkxk (qj - yj) = = n i=1 m k=1 (pi - xi)aik(qk - yk) - n k=1 m j=1 (pk - xk) bjk cj (qj - yj) = 0, neboť bjk/cj = akj. Jinak řečeno, funkce H je na trajektoriích systému (6.55) konstantní, funkce H je invariantem tohoto systému. Parciální derivace pravé strany rovnice (6.54) podle j-té proměnné je xj diag x (b - A x) = xj diag x (b - A x) + diag x xj (b - A x) = = diag ej (b - A x) + diag x -A ej a pro vnitřní stacionární bod x platí b - A x = O. Proto variační matice systému (6.52) ve vnitřním stacionárním bodě x je J(x ) = - diag x A. Odtud a z 5.3.8 plyne: 71 6.4.4 Věta Buď x stacionární bod systému (6.52), jehož všechny složky jsou nenulové. Mají-li všechna vlastní čísla matice diag x A kladnou reálnou část, pak konstantní řešení x(t) x systému (6.52) je stejnoměrně asymptoticky stabilní. Pokud existuje vlastní číslo matice diag x A které má zápornou reálnou část, pak je konstantní řešení x(t) x systému (6.52) nestabilní. 6.4.5 Příklad (koloběh dusíku v planktonu) Uvažujme proces schématicky znázorněný na obrázku 6.5: Ve fytoplanktonu probíhá fotosyntéza a při ní se dusík z okolního prostředí váže v jeho buňkách; fytoplankton slouží jako potrava pro zooplankton, takže dusík ze zkonzumovaného fytoplanktonu se stává součástí zooplanktonu. Plankton v důsledku svého metabolismu dusík opět vylučuje do okolního prostředí a také při rozkladu mrtvého planktonu se dusík uvolňuje. Dusík z prostředí není odebírán ani není nějakým způsobem do něho přidáván. Dusíku vylučovaného planktonem je tím více, čím je více planktonu, dusíku vázaného ve fytoplanktonu přibývá tím více, čím je více volného dusíku a fytoplanktonu; dusíku vázaného v zooplanktonu přibývá tím více, čím více je fytoplanktonu pozřeného zooplanktonem a toho je tím více, čím více je fytoplanktonu i zooplanktonu. Označme po řadě N, P a Z množství dusíku v prostředí, vázaného ve fytoplanktonu a vázaného v zooplanktonu. Všechny tyto veličiny se mění s časem, tj. N = N(t), P = P(t) a Z = Z(t). Celkové množství dusíku v systému je rovno V = N + P + Z. Koloběh dusíku lze nejjednodušeji modelovat systémem rovnic N = aP + bZ - cNP, P = cNP - dPZ - aP, Z = dPZ - bZ; všechny parametry a, b, c, d jsou kladné. - - 6 ? ? a b c dN dusík v prostředí P fytoplankton Z zooplankton Obrázek 6.5: Schéma koloběhu dusíku Nejprve si všimněme, že V = N + P + Z = 0, což znamená, že celkové množství dusíku V je konstantní. Proto lze množství dusíku v prostředí vyjádřit jako N(t) = V - P(t) - Z(t) a dosadit do druhé a třetí rovnice systému. Dostaneme P = (V c - a)P - cP2 - (c + d)PZ = P V c - a - cP - (c + d)Z , Z = -bZ + dPZ = Z - b + dP . (6.56) Jedná se tedy o Lotkův-Volterrův systém s vektorem růstových koeficientů a maticí interakcí b = V c - a -b a A = c -(c + d) d 0 . Vnitřní stacionární bod systému je tedy P Z = A-1 b = 1 d(c + d) 0 -(c + d) d c V c - a -b = b d c c + d V - a c - b d . Povšimněme si, že stacionární hodnota P je kladná a nezávisí na celkovém množství dusíku V . Pokud se zvětší přísun živin, nemá z toho užitek fytoplankton, ale jeho predátor zooplankton. 72 Pokud V > a c + b d , (6.57) pak Z > 0; v takovém případě je tedy možná koexistence fyto i zooplanktonu. Dále platí J(P , Z ) = - b d c c + d V - a c - b d 0 -(c + d) d c = = 0 b(c + d) d - cd c + d V - a c - b d - c2 c + d V - a c - b d . Je-li splněna podmínka (6.57), pak tr J(P , Z ) = - c2 c + d V - a c - b d < 0, det J(P , Z ) = bc V - a c - b d > 0, což znamená, že reálná část vlastních čísel variační matice J(P , Z ) je záporná, a tedy vnitřní stacionární řešení (P , Z ) systému (6.56) je stejnoměrně asymptoticky stabilní. Systém (6.56) má další stacionární body s0 = 0 0 , s1 = V - a c 0 , jeho variační matice v obecném bodě je J(P, Z) = V c - a - 2cP - (c + d)Z -(c + d)P dZ -b + dP , tedy J(s0) = V c - a 0 0 -b , tr J(s0) = c V - a c - b, det J(s0) = -bc V - a c , J(s1) = -c V - a c -(c + d) V - a c 0 d V - a c - b , tr J(s1) = d V - a c - b d - c V - a c , det J(s1) = -cd V - a c V - a c - b d . Odtud je vidět, že pokud V < a/c, pak s0 je stejnoměrně asymptoticky stabilní, pokud V > a c (6.58) a současně V < a/c + b/d, pak s0 je sedlo a s1 je stejnoměrně asymptoticky stabilní. Pokud platí nerovnost (6.57), pak s0 i s1 jsou sedla. Z dosud provedených úvah a výpočtů lze učinit závěr: je-li v prostředí dusíku málo (méně než a/c), plankton nepřežívá; je-li dusíku více, ale ne příliš mnoho (něco mezi a/c a a/c+b/d), přežívá fytoplankton, ale nikoliv zooplankton; je-li dusíku v prostředí dostatek (více než a/c+b/d), přežívá fyto i zooplankton. Nerovnosti (6.58) a (6.57) lze přepsat na tvar c > a V a c > ad V d - b = a V 1 - b V d . Koeficient c vyjadřuje, s jakou intenzitou je dusík z prostředí vázán do biomasy fytoplanktonu, což je proces fotosyntézy. Její intenzita roste s množstvím slunečního světla a to se mění s ročním obdobím. Dosažené výsledky lze tedy také interpretovat: v zimě se plankton nevyskytuje, na jaře se nejprve objeví fytoplankton a s přibýváním světla i zooplankton. 73 6.4.6 Poznámka Pro čtvercovou matici M položme SM = 1 2 M + MT . Matice SM je zřejmě symetrická. Pro každý n-rozměrný vektor v a čtvercovou matici M řádu n platí vT M v = vT SM v. D.: Poněvadž vT M v je číslo, tj. čtvercová matice řádu 1, platí vT M v = vT M v T = vT MT v a odtud plyne vT Mv = 2vT 1 2 (M + MT ) - 1 2 MT v = 2vT SMv-vT MT v = 2vT SMv-vT Mv a tato rovnost je již ekvivalentní s dokazovaným vztahem. 6.4.7 Věta Buď x = (x 1, x 2, . . . , x n) = A-1 b vnitřní stacionární bod systému (6.52). Jestliže existuje konstantní vektor c = (c1, c2, . . . , cn)T se všemi složkami kladnými a existuje okolí U bodu x takové, že pro všechna x U je výraz (x - x ) T S(diag c A) (x - x ) (6.59) nezáporný, pak funkce V (x) = V (x1, x2, . . . , xn) = n i=1 ci xi x i - x i d je ljapunovskou funkcí systému (6.52), tj. konstantní řešení x(t) x systému (6.52) je stejnoměrně stabilní. Pokud je výraz (6.59) pro všechna x U \{x } kladný, pak je toto řešení stejnoměrně asymptoticky stabilní. D.: Funkce V je definována pro všechna x1 > 0, x2 > 0, . . . , xn > 0. Platí V (x ) = n i=1 ci x i x i - x i d = 0. Pro každé xi > 0 je xi x i - x i d 0, neboť integrovaná funkce je kladná pro xi > x i (tj. v případě, že horní mez integrálu je větší, než dolní mez) a záporná pro xi < x i (horní mez integrálu menší než dolní mez). Rovnost nastane právě tehdy, když xi = x i . Odtud plyne, že pro x = x a takové, že všechny jeho složky jsou kladné, platí V (x) > 0. Dále podle věty o derivaci integrálu jako funkce horní meze platí V (x) xi = ci xi - x i xi , 74 a poněvadž b = A x , platí dále bi = n j=1 aijx j , takže derivace funkce V vzhledem k systému (6.52) je ˙V (x) = n i=1 ci xi - x i xi xi bi - n j=1 aijxj = n i=1 ci (xi - x i ) n j=1 aijx j - n j=1 aijxj = = - n i=1 n j=1 (xi - x i ) ciaij xj - x j = - (x - x ) T (diag c A) (x - x ) = = - (x - x ) T S (diag c A) (x - x ) (poslední rovnost plyne z poznámky 6.4.6). Věta nyní plyne z 5.3.11 a 5.3.12. 6.4.8 Důsledek Nechť systém (6.52) má vnitřní stacionární bod x = (x 1, x 2, . . . , x n). Jestliže existuje konstantní vektor c = (c1, c2, . . . , cn)T se všemi složkami kladnými takový, že matice S (diag c A) (6.60) je pozitivně semidefinitní, pak konstantní řešení x(t) x systému (6.52) je stejnoměrně stabilní. Pokud je matice (6.60) pozitivně definitní, pak konstantní řešení x(t) x systému (6.52) je stejnoměrně asymptoticky stabilní. 6.4.9 Příklad (trofický řetězec) Trofický řetězec je takové společenstvo, v němž je první druh producentem a každý jiný druh je nesoběstačným specializovaným predátorem právě jednoho dalšího druhu. Označíme x1 velikost populace producenta, x2 velikost populace jeho predátora, x3 velikost populace, která je predátorem populace o velikosti x2, atd. Každá z populací na některé trofické úrovni nemusí být tvořena jedním biologickým druhem, může jít o společenstvo organismů majících stejný způsob obživy. Trofický řetěz o n úrovních lze tedy modelovat systémem x 1 = x1(r - ax1) - p1x1x2 x 2 = -d2x2 + q2x1x2 - p2x2x3 ... x k = -dkxk + qkxk-1xk - pkxkxk+1 ... x n-1 = -dn-1xn-1 + qn-1xn-2xn-1 - pn-1xn-1xn x n = -dnxn + qnxn-1xn, (6.61) parametry r, d2, d3, . . . , dn, q2, q3, . . . , qn, p1, p2, . . . , pn-1 jsou kladné, parametr a je nezáporný (producent může, ale nemusí projevovat vnitrodruhovou konkurenci). V tomto případě je A = a p1 0 0 0 0 -q2 0 p2 0 0 0 0 -q3 0 0 0 0 ... ... ... ... ... ... ... 0 0 0 -qn-1 0 pn-1 0 0 0 0 -qn 0 , b = r -d2 -d3 ... -dn-1 -dn . 75 Položme c1 = 1, c2 = p1 q2 , c3 = p1p2 q2q3 , . . . , cn-1 = p1p2 pn-2 q2q3 qn-1 , cn = p1p2 pn-1 q2q3 qn . Pak je diag c A = a p1 0 0 0 0 -p1 0 p1p2 q2 0 0 0 0 - p1p2 q2 0 0 0 0 ... ... ... ... ... ... ... 0 0 0 p1p2 pn-2 q2q3 qn-2 0 p1p2 pn-1 q2q3 qn-1 0 0 0 0 p1p2 pn-1 q2q3 qn-1 0 , takže S (diag c A) = a 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ... ... ... ... ... ... ... 0 0 0 0 0 0 0 0 0 0 0 0 z čehož plyne (x - x ) T S(diag c A) (x - x ) = a(x1 - x 1)2 0. Pokud existuje vnitřní stacionární bod x uvažovaného systému, pak je příslušné konstantní řešení stejnoměrně stabilní. Hledejme nyní podmínky, které zaručí existenci takového stacionárního bodu x . Jeho souřadnice splňují n-rozměrný systém algebraických rovnic ax 1 + p1x 2 = r, qkx k-1 - pkx k+1 = dk, k = 2, 3, . . ., n - 1. qnx n-1 = dn. (6.62) ,,Prostřední rovnice tohoto systému lze přepsat ve tvaru rekurentních formulí x k-1 = 1 qk pkx k+1 + dk nebo x k+1 = 1 pk qkx k-1 - dk , k = 2, 3, . . . , n - 1. (6.63) Poněvadž všechny koeficienty pk, qk, dk jsou kladné, plyne z tohoto vyjádření: (i) je-li x 0 > 0 pro nějaké l0 {2, 3, . . ., n} , pak je x > 0 pro všechna 0, 0 - 2, 0 - 4, . . . , 1 2 3 + (-1)0 ; (ii) je-li x 1 0 pro nějaké 1 {1, 3, . . . , n - 2} , pak je x < 0 pro všechna 1 + 2, 1 + 4, . . . , n - 1 2 1 - (-1)1+n . Podle poslední rovnice systému (6.62) je x n-1 = dn qn . 76 Z první rekurentní formule (6.63) postupně vyjádříme x n-3 = 1 qn-2 pn-2x n-1 + dn-2 = pn-2 qn-2 dn qn + dn-2 qn-2 , x n-5 = 1 qn-4 pn-4x n-3 + dn-4 = pn-4 qn-4 pn-2 qn-2 dn qn + pn-4 qn-4 dn-2 qn-2 + dn-4 qn-4 , atd. Celkem dostaneme x n-(2+1) = i=0 dn-2i qn-2i j=i+1 pn-2j qn-2j , pro = 0, 1, . . . n 2 - 1, (6.64) kde [] označuje celou část z čísla a klademe k-1 j=k j = 1 pro libovolné přirozené k a každou posloupnost {j} j=0.3 Přímým výpočtem se lze přesvědčit, že (6.64) je skutečně řešením druhé až n-té rovnice systému (6.62). Nechť nejprve je n sudé. V tomto případě lze rovnost (6.64) přepsat na tvar x 2k-1 = x n-(2(n 2 -k)+1) = n 2 -k i=0 dn-2i qn-2i n 2 -k j=i+1 pn-2j qn-2j = n 2 i=k d2i q2i i-1 j=k p2j q2j , k = 1, 2, . . ., n 2 . Z tohoto vyjádření je vidět, že všechny souřadnice stacionárního bodu x s lichými indexy jsou kladné. Pro jeho první souřadnici platí x 1 = n 2 i=1 d2i q2i i-1 j=1 p2j q2j . (6.65) Z první rovnice systému (6.62) nyní dostaneme x 2 = r - ax 1 p1 , a ze druhé rekurentní formule (6.63) x 4 = 1 p3 (q3x 2 - d3) = q3 p3 r - ax 1 p1 - d3 p3 , x6 = 1 p5 (q5x 4 - d5) = q5 p5 q3 p3 r - ax 1 p1 - q5 p5 d3 p3 - d5 p5 , atd. Obecně x 2k = r - ax 1 p1 k-1 i=1 q2i+1 p2i+1 - k-1 i=1 d2i+1 p2i+1 k-1 j=i+1 q2j+1 p2j+1 , k = 1, 2, . . . , n 2 . Souřadnice x 1 je vyjádřena formulí (6.65). Tedy platí x n = x 2 n 2 = r - ax 1 p1 n 2 -1 =1 q2+1 p2+1 - n 2 -1 i=1 d2i+1 p2i+1 n 2 -1 j=i+1 q2j+1 p2j+1 = = n 2 -1 =1 q2+1 p2+1 r - ax 1 p1 - n 2 -1 i=1 d2i+1 p2i+1 i j=1 p2j+1 q2j+1 = = 1 p1 n 2 -1 =1 q2+1 p2+1 r - a n 2 i=1 d2i q2i i-1 j=1 p2j q2j - p1 n 2 -1 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 . 3Uvedená konvence je přirozeným rozšířením rovnosti kQ j=m j = k k-1Q j=m j , která platí pro libovolné k > m, také pro k = m. 77 Nutnou a dostatečnou podmínkou pro to, aby všechny souřadnice stacionárního bodu x byly kladné, je tedy podle tvrzení (i) a (ii) nerovnost r > p1 n 2 -1 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 + a n 2 i=1 d2i q2i i-1 j=1 p2j q2j . (6.66) Nechť nyní je n liché. V tomto případě lze rovnost (6.64) přepsat na tvar x 2k = xn-(2(n-1 2 -k)+1) = n-1 2 -k i=0 dn-2i qn-2i n-1 2 -k j=i+1 pn-2j qn-2j = n-1 2 i=k d2i+1 q2i+1 i-1 j=k p2j+1 q2j+1 , k = 1, 2, . . . , n - 1 2 . Z něho je vidět, že všechny souřadnice stacionárního bodu x se sudými indexy jsou kladné. Zejména jeho druhá souřadnice je x 2 = n-1 2 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 . (6.67) Je-li a = 0, dostaneme z první rovnice systému (6.62) x 1 = r - p1x 2 a , ze druhé rekurentní formule (6.63) nyní můžeme postupně vyjádřit x 3 = 1 p2 (q2x 1 - d2) = q2 p2 r - p1x 2 a - d2 p2 , x 5 = 1 p4 (q4x 3 - d4) = q4 p4 q2 p2 r - p1x 2 a - q4 p4 d2 p2 - d4 p4 , atd. Obecně dostaneme x 2k-1 = r - p1x 2 a k-1 i=1 q2i p2i - k-1 i=1 d2i p2i k-1 j=i+1 q2j p2j , k = 1, 2, . . ., n + 1 2 . Odtud s využitím (6.67) vyjádříme x n = x 2 n+1 2 -1 = r - p1x 2 a n-1 2 =1 q2 p2 - n-1 2 i=1 d2i p2i n-1 2 j=i+1 q2j p2j = = 1 a n-1 2 =1 q2 p2 r - p1 n-1 2 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 - a n-1 2 i=1 d2i q2i i-1 j=1 p2j q2j . Pro liché n a a = 0 tedy dostáváme jako nutnou a dostatečnou podmínku pro to, aby všechny souřadnice stacionárního bodu x byly kladné, nerovnost r > p1 n-1 2 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 + a n-1 2 i=1 d2i q2i i-1 j=1 p2j q2j . (6.68) Pokud je n liché a a = 0, dostaneme z první rovnice systému (6.62) rovnost x 2 = r p1 . 78 Současně však musí platit rovnost (6.67), takže soustava rovnic (6.62) má řešení (a to nekonečně mnoho řešení; stacionární bod není v takovém případě izolovaný) pouze tehdy, když r = p1 n-1 2 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 . Pravděpodobnost, že tato rovnost bude splněna pro systém (6.61) modelující reálné společenstvo, je však nulová. Povšimněme si ještě, že nerovnosti (6.66) a (6.68) lze zapsat jednotně ve tvaru r > p1 [n-1 2 ] i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 + a [n 2 ] i=1 d2i q2i i-1 j=1 p2j q2j . (6.69) Závěr: Je-li a > 0 (základní zdroj je omezený, v populaci producenta je vnitropopulační konkurence), pak vnitřní stacionární bod systému (6.61) existuje a je stejnoměrně asymptoticky stabilní (je možná koexistence všech populací tvořících trofický řetězec) právě tehdy, když je splněna podmínka (6.69) (vnitřní koeficient růstu producenta je dostatečně velký). Je-li a = 0 (základní zdroj je neomezený), pak vnitřní stacionární bod systému (6.61) existuje pouze pro sudé n (je možná koexistence pouze sudého počtu trofických úrovní); vnitřní stacionární bod v takovém případě je stejnoměrně stabilní a existuje právě tehdy, když je splněna podmínka r > p1 n 2 -1 i=1 d2i+1 q2i+1 i-1 j=1 p2j+1 q2j+1 . 6.4.10 Zobecněné Lotkovy-Volterrovy systémy (Grossberg 1978) Vlivy populací tvořících společenstvo na růst jednotlivých populací nemusí být tvaru přímé úměrnosti. Proto může být realističtější místo systému (6.52) uvažovat systém x i = gi(xi) bi - n j=1 aijfj(xj) , i = 1, 2, . . ., n. (6.70) Funkce fi, gi, i = 1, 2, . . ., n jsou definovány a spojité na intervalu [0, ) a splňují podmínky: * (i)gi(0) = 0 . . . je-li velikost i-té populace nulová (tj. i-tá populace ve společenstvu není), pak nulovou zůstane; uvažujeme tedy izolovaná společenstva, kde nedochází k imigraci nových druhů. * (i)( > 0)gi() > 0 . . . skutečnost, zda je i-tá populace soběstačná nebo ne, nezávisí na její velikosti; neuvažujeme tedy např. Alleeho efekt. * (j)fj(0) = 0 . . . není-li j-tá populace ve společenstvu přítomná, nijak neovlivňuje růst ostatních populací. * (j)fj je rostoucí . . . s rostoucí velikostí populace roste i její vliv na růst populací ostatních. Systém (6.70) lze zapsat vektorově: x = G(x) b - A f(x) , kde G(x) = G(x1, x2, . . . , xn) = diag g1(x1), g2(x2), . . . , gn(xn) , 79 f(x) = f1(x1), f2(x2), . . . , fn(xn) T . Poněvadž všechny složky zobrazení f : Rn Rn jsou rostoucí (tedy prosté) funkce, je toto zobrazení prosté a existuje k němu zobrazení inverzní f-1 = (f-1 1 , f-1 2 , . . . , f-1 n ). Je-li matice interakcí společenstva A regulární, existuje nejvýše jeden vnitřní stacionární bod x = (x 1, x 2, . . . , x n) = f-1 A-1 b systému (6.70), tj. takový bod, že x 1 > 0, x 2 > 0, . . . , x n > 0, který lze opět interpretovat jako dynamicky stálé velikosti všech populací koexistujících ve společenstvu. Analogicky jako v důkazu věty 6.4.7 ověříme, že pokud existuje okolí U vnitřního stacionárního bodu x a existuje konstantní vektor c = (c1, c2, . . . , cn)T se všemi složkami kladnými, pro něž je výraz f(x) - f(x ) T S(diag c A) f(x) - f(x ) nezáporný pro každé x U, pak je funkce V (x) = V (x1, x2, . . . , xn) = n i=1 ci xi x i fi() - fi(x i ) gi() d ljapunovskou funkcí systému (6.70) ve stacionárním bodě x . Odtud je vidět, že tvrzení důsledku 6.4.8 platí také pro systém (6.70). 80