• Nie Znaleziono Wyników

Układy równań różniczkowych zwyczajnych Wstęp

N/A
N/A
Protected

Academic year: 2021

Share "Układy równań różniczkowych zwyczajnych Wstęp"

Copied!
7
0
0

Pełen tekst

(1)

Układy równań różniczkowych zwyczajnych

Wstęp

W zastosowaniach bardzo często zamiast pojedynczego równania występują układy równao równao różniczkowych zwyczajnych (ang. System of Ordinary Differential Equations, ODEs). Oznacza to, że szukamy kilku funkcji, których pochodne spełniają pewne związki. W zagadnieniach kinetyki reakcji chemicznych funkcje y t1( ), y t2( ) mogą opisywad stężenia reagujących substancji. Na przykład pewien układ reakcji zwany Brusselatorem może prowadzid do następującego układu równao

2

1 1 2 1

2

2 1 1 2

1 4 ,

3 .

y y y y

y y y y

   



  

 (1)

Bardzo rzadko się zdarza, aby układy takie jak (1) miały rozwiązania, które można wyrazid „wzorami”

analitycznymi. W tej sytuacji równania takie możemy badad numerycznie, poprzez uzyskiwanie przybliżeo dla konkretnych warunków początkowych. Można też analizowad własności rozwiązao bez wyrażania ich wzorami, ale w oparciu o pewne twierdzenia i matematyczne teorie. Zajmuje się tym dział matematyki zwany teorią układów dynamicznych. Obie metody – numeryczna i jakościową – można łączyd.

Innym obszarem, który dostarcza przykładów układów równao różniczkowych zwyczajnych są różne modele biologiczne. Dobrze znanym jest tzw. układ LotkiVolterry, albo inaczej model drapieżnik- ofiara. Zakładamy, że na danym obszarze występują dwa gatunki: drapieżniki (np. lisy) i ofiary (np.

zające). Drapieżniki żywią się ofiarami. Jeżeli wprowadzimy oznaczenia:

1( )

y tliczebnośd populacji ofiar (można też mówid o gęstości populacji ofiar) w chwili t ,

2( )

y tliczebnośd populacji drapieżników (lub gęstości populacji drapieżników) w chwili t ,

to dośd prosty model opisujący jak zmienia się w czasie populacja ofiar i drapieżników, które na siebie wzajemnie oddziałują, zawarty jest w następującym układzie równao

1 2 1

2 1 2

( ) ,

( ) .

y b ay y y cy d y

  

   

(2)

W układzie tym występują dodatnie stałe a b c d, , ,  charakteryzujące jakośd środowiska oraz możliwości obu gatunków. Jeżeli przepiszemy pierwsze równanie układu (2) następująco

1

2 1

y , y b ay

   (3)

to widzimy, że względna zmiana w czasie liczebności ofiar (wszak 1 1

1 1

1 1/ 1 dydt y yyt) yy   jest proporcjonalna do b, a zatem parametr ten określa możliwości rozmnażania się ofiar: im więcej jest osobników tym więcej będzie ich przybywad. Gdyby w równości tej pominąd składnik a y2, to otrzymalibyśmy

(2)

1

1 1

1

( ) (0) bt,

y b y t y e

y

   

czyli wzrost wykładniczy zależny od b. Można powiedzied, że parametr ten określa zdolności reprodukcyjne ofiar oraz jakośd środowiska, które są stałe  nie zależą od liczebnośd (lub gęstości) żadnej z populacji. Z drugiej strony mamy jednak w równaniu (3) czynnik hamujący wzrost liczby ofiar – są nim drapieżnicy. Mianowicie im więcej drapieżników, tym więcej potrzebują pożywienia czyli ofiar. Czynnik ay2 zawiera parametr a, który charakteryzuje skutecznośd drapieżników w polowaniu na ofiary.

Podobnie można przeanalizowad drugie równanie układu (2) zapisując je w formie

2 1 2

y .

cy d y

   (4)

Interpretacja składnika c y1 jest taka, że im więcej ofiar tym szybciej rozmnażają się drapieżniki (jest dużo pożywienia). Dlatego c może oznaczad jakośd tego pożywienia (ofiary) oraz zdolności reprodukcyjne drapieżników.

W przypadku gdy mamy układ dwóch równao to zamiast numerowanie funkcji ( ( ),y t1 y t2( )) używa się często np. oznaczeo ( ),x t y t( ). Wtedy układ (2) zapiszemy następująco

( ) ,

( ) .

x b ay x y cx d y

  

   

(5)

Co ciekawe dokładnie taki sam układ uzyskuje się gdy rozważymy pewien teoretyczny mechanizm reakcji zaproponowany przez Alfreda Lotkę1 w 1920 roku, w którym występują dwa autokatalityczne etapy

1

2

3

(i) A X 2X

(ii) X Y 2Y

(iii) Y B

k

k

k

 

 

(6)

gdzie A i B to substrat i produkt, X i Y to formy pośrednie, a k1, k2, k3 to stałe szybkości odpowiednich reakcji. Jeżeli symbolem *X+ oznaczymy stężenie molowe (najczęściej wyrażane w jednostkach molm-

3 lub moldm-3), to równania chemiczne (i) oraz (ii) prowadzą do następującego układu równao różniczkowych zwyczajnych

1 2

2 2

[ ] [ ][ ] [ ][ ], [ ] [ ][ ] [ ],

d X k A X k X Y

dt

d Y k X Y k Y dt

 

 

(7)

1 Alfred Lotka (1880-1949)  urodzony we Lwowie polskoamerykaoski chemik, biofizyk i statystyk, twórca podstaw nauki od dynamice populacji.

(3)

gdzie stężenia *X+, *Y+ i *A+ są funkcjami czasu. Gdy stężenie substratu A jest stałe (np. składnik jest w nadmiarze i przez pewien czas ubytek można zaniedbad lub jest on uzupełniany w reaktorze chemicznym), to wtedy k1*A+ = const i układ równao (7) jest formalnie identyczny z układem (5). Inne jest tylko znaczenie symboli: ( ), ( )x t y t to ilości osobników dwóch gatunków w danym terenie, a *X+,

*Y+ to stężenia odpowiednich substancji w reaktorze chemicznym.

Rysunek 1. Przykładowe rozwiązanie układu LotkiVolterry. Jak widać rozwiązania są okresowe. W opisie chemicznym reakcji AB z mechanizmem (6) możemy powiedzieć, że stężenia form pośrednich X i Y wykazują oscylacje o stałej amplitudzie.

Układy równao różniczkowych zwyczajnych pierwszego rzędu (czyli takie, w których występują tylko pochodne pierwszego rzędu) pojawiają się także gdy przekształcamy równania różniczkowe wyższych rzędów. Rozważamy dla przykładu znane nam już równanie wahadła matematycznego w ośrodku stawiającym opór (Błąd! Nie można odnaleźd źródła odwołania.). Jeżeli ( ) t oznacza kąt wychylenia (liczony względem pionu), to z rysunku widad, że w kierunku stycznym do toru (który jest okręgiem o promieniu ) działają dwie siły: mg sin (składowa siły ciężkości mg w kierunku stycznym) oraz

k v (siła oporu ośrodka – zakładamy, że jest proporcjonalna do prędkości i oczywiście przeciwnie skierowana). Ponieważ przyspieszenie styczne as , a prędkośd styczna v , to z II zasady dynamiki Newtona otrzymujemy zależnośd

sin ,

m   k mg  (8)

czyli m k mgsin0, gdzie: mmasa kulki, długośd linki (pręta), kwspółczynnik oporu ośrodka, gprzyspieszenie grawitacyjne (g9,81 m/s ).2 Oznaczmy k m/ ,2g/ i przepisujemy równanie (8)

2sin .

    (9)

Jeżeli wprowadzimy teraz oznaczenia: y1, y2, to powyższe równanie można sprowadzid do następującego układu równao pierwszego rzędu

(4)

1 2

2

2 2 1

,

sin .

y y

yyy

  

    

 (10)

Rysunek 2. Wahadło matematyczne o długości .

W ogólnym przypadku może występowad n niewiadomych funkcji y1(t),…, yn(t), których ewolucje czasowe są powiązane. Oznacza to, że danych jest n równao różniczkowych zwyczajnych (czyli układ równań):

1 1 1

1

( , , , ),

( , , , ),

n

n n n

y f t y y

y f t y y

  



  

(11)

z warunkami początkowymi y t1( )0y10, ,y tn( )0yn0, gdzie liczby t0, y10, ,yn0 są dane.

Układy równań liniowych

W ogólnym przypadku układ liniowy pierwszego rzędu ma postad

1 11 1 1 1

1 1

( ) ( ) ( ),

( ) ( ) ( ).

n n

n n nn n n

y a t y a t y f t

y a t y a t y f t

    



     

(12)

W układzie tym funkcje: aija tij( ), fif ti( ) :  J są dane i określone na odcinku .J W praktyce najczęściej zakładamy, że są one przynajmniej ciągłe. Wtedy problem początkowy dla układu (12) ma jednoznaczne rozwiązanie dla każdych warunków początkowych.

Układ (12) można zapisad w sposób zwarty używając notacji macierzowej (algebra liniowa):

( ) ( ),

y A t yf t (13)

gdzie występują macierz ( )A t i wektory kolumnowe ( ),y t f t( ) :

(5)

11 1 1 1

1

( ) ( ) ( ) ( )

( ) , ( ) , ( ) .

( ) ( ) ( ) ( )

n

n nn n n

a t a t y t f t

A t y t f t

a t a t y t f t

     

     

     

     

     

(14)

Dokładniej zajmiemy się szczególnym przypadkiem: dwa równania (n=2), stałe współczynniki (tzn.

funkcje ai j( )t nie zależą od czasu). Zaczniemy od równania jednorodnego ( ( ) 0).f t  Dla wygody funkcje niewiadome będziemy oznaczad ( ), ( )x t y t zamiast y t1( ), y t2( ), a współczynniki , , ,a b c d zamiast a11,a12,a21,a22. Mamy zatem układ

, , x ax by y cx dy

  

   

 (15)

którego macierzą jest a b .

A c d

 

  

 

W przypadku układu (15) można zawsze wyrazid rozwiązanie prostymi wzorami. Sytuacja jest tutaj bardzo podobna do równania skalarnego liniowego drugiego rzędu ze stałymi współczynnikami (które było rozważane wcześniej), aybycy0. Musimy tym razem znaleźd wartości własne macierzy układu, czyli rozwiązad równanie

det a b 0,

c d

  

  

  (16)

czyli (a)(d)cb2 (a d)adbc0.

Oczywiście (16) jest równaniem kwadratowym na , które ma zawsze rozwiązanie (w dziedzinie zespolonej, ). Równanie (16) nazywamy równaniem charakterystycznym (macierzy, układu równao), a sam wielomian zmiennej  – wielomianem charakterystycznym. Pierwiastki wielomianu charakterystycznego nazywamy wartościami własnymi (macierzy, układu). Na przykład dla układ

2 3 , ,

x x y

y x y

  

    

 (17)

ma macierz I równanie charakterystyczne:

2 3 2 3

, det (2 )(1 ) 3 0.

1 1 1 1

A

 

    

         

Skąd (2)(1) 3 23 5 0, czyli pierwiastki charakterystyczne układu (17) wynoszą

1 2 1 2

3 11 3 11 3 11 3 11

, , .

2 2 2 i 2 2 i 2

       

      

Aby podad rozwiązanie układu (15) można obliczyd wektory własne macierzy układu dla każdej wartości własnej, lub bezpośrednio rozpatrzyd trzy zasadnicze przypadki: (a)  0 (dwa rzeczywiste

(6)

pierwiastki), (b)  0 (dwa zespolone pierwiastki, (c)  0 (jeden rzeczywisty pierwiastek podwójny (musi byd rzeczywisty, bo współczynniki równania charakterystycznego są liczbami rzeczywistymi)).

Przypadek a) Δ > 0  12 .

Rozwiązanie ogólne ma postad

1 2

1 2

1 2

1 1 2 2

( ) ,

( ) ( ) ( ) ,

t t

t t

x t C be C be

y t C a e C a e

 

 

    (18)

gdzie C C1, 2 są dowolnymi stałymi.

Przypadek b) Δ < 0  12 :    1 i , 2  i . Rozwiązanie ogólne ma postad

1 2

1 2 1 2

( ) ( sin cos ),

( ) ((( ) )sin ( ( ) ) cos ),

t t

x t be C t C t

y t e a C C t C a C t

 

     

 

      (19)

gdzie C C1, 2 są dowolnymi stałymi.

Przypadek c) Δ = 0   12  .

Rozwiązanie ogólne zależy dodatkowo od relacji pomiędzy a b c d, , , . Pojawiają się funkcje liniowe. Na przykład, gdy zachodzi ad, to

2

2

2

1 2

1 2 2

( ) 2 ,

( ) (( ) ( ) ) ,

a d

a d

t C

a d

t

x t b C C t e

y t d a C C d a C t e

 

     (20)

gdzie C C1, 2 są dowolnymi stałymi.

Zadania

Zad. 1) Zamieo na układ równao podane równania skalarne drugiego rzędu.

a) y y20. b) y 1 cosy3 .y b) y2yt y2 0. c) ya t y( ) b t y( ) c t( ).

Zad. 2) Podaj rozwiązanie ogólne układów ze stałymi współczynnikami.

a) x x y y x y

  

   

 b) 2

3

x x y

y x y

  

    

 c) 2

4 6

x x y

y x y

   

   

Zad. 3) Niech ( ),x t y t( ) będzie dowolnym rozwiązanie układu LotkiVolterry. Pokazad, że wtedy funkcja ( )f tcx t( )  dln ( )x ta y t( )  bln ( )y t jest stała względem czasu.

(7)

Zad. 4) Podaj rozwiązanie problemu wahadła matematycznego przyjmując założenie małych wychyleo. Oznacza to, że przybliżamy sinus dla małych kątów znanym wzorem (szereg Taylora!):

sin  . Tak więc zamiast (8) rozważamy równanie

0.

m k mg

Dla uproszczenia rozważyd następujący warunek początkowy: (0) 0, (0)0. Przedyskutowad charakter rozwiązao w zależności od relacji pomiędzy parametrami , ,m k g.

Cytaty

Powiązane dokumenty

[r]

Z powodów niezależnych od organizacji, które ukierunkowują jakiekol‑ wiek działanie człowieka poprzez złożone systemy oddziaływania i pobudzania zachowań

Wówczas, aby rozwiązać równanie wystarczy podać wszystkie jego rozwiązania integralne, gdyż każde inne rozwiązanie jest obcięciem pewnego rozwiązania integralnego do

Zanim przejdziemy do dalszej części wykładu przypomnijmy, że jedynymi zbiorami spój- nymi na prostej R są: zbiór pusty, zbiory jednoelementowe i dowolne przedziały.. Jest

Jak pokazaliśmy w przykładzie 1.3.1., każde rozwią- zanie tego równania określone jest na pewnym przedziale zawartym w dziedzinie jednego z powyższych rozwiązań, więc

Jest to równanie o zmiennych rozdzielonych.... Jest to równanie o

temperatury, natomiast, co już może dziwić, czasami widać, że zarejestrowane stężenie tlenu jest wyższe niż stężenie nasycenia, ale i to jest normalne i zdarza się,

W rozwiązaniu powinien znaleźć się skrypt rozwiązujący dane równanie w Matlabie oraz wyświetlający pole kierunkowe wraz z przykładowymi rozwiązaniami, jak również link do