Zagadnienia

8.1 Układ konkurujących gatunków
8.2 Modelowanie symbiozy

8. Modele oddziaływań między dwiema populacjami III

8.1. Układ konkurujących gatunków

Korzystając z doświadczeń związanych z omówionymi już modelami oddziaływań typu drapieżnik – ofiara przejdziemy teraz do opisu oddziaływań innego typu. Jak poprzednio zakładamy, że w środowisku występują dwa gatunki $\mathcal{P}_{1}$ i $\mathcal{P}_{2}$ , ale tym razem osobniki tych gatunków konkurują ze sobą o zasoby środowiska. Zwykle budując taki model przyjmujemy, że jest to konkurencja o pożywienie, ale może ona dotyczyć także dostępnej przestrzeni życiowej, a często mamy do czynienia z obydwoma rodzajami konkurencji. Ponieważ gatunki konkurują ze sobą, to nie widać powodu do wyróżniania któregokolwiek z nich, zatem opis powinien być symetryczny — jeśli zamienimy te gatunki i nazwiemy gatunek $\mathcal{P}_{i}$ ( $i=1$ , $2$ ) gatunkiem $\mathcal{P}_{j}$ ( $j=2$ , $1$ ) i na odwrót, to układ powinien pozostać taki sam.

Skoro mamy skorzystać z nabytych doświadczeń, zbudujemy model w oparciu o wewnętrzną dynamikę logistyczną dla każdego z gatunków $\mathcal{P}_{i}$ . Wobec tego, jeśli w środowisku występuje tylko jeden z tych gatunków, to jego zagęszczenie $N_{i}(t)$ w chwili $t$ opisuje równanie

$\dot{N}_{i}=r_{i}N_{i}\left(1-\frac{N_{i}}{K_{i}}\right),$

gdzie $r_{i}$ oznacza współczynnik rozrodczości netto dla gatunku $\mathcal{P}_{i}$ , zaś $K_{i}$ — pojemność środowiska dla tego gatunku. Oczywiście założenie ograniczonej pojemności środowiska implikuje też występowanie konkurencji między osobnikami tego samego gatunku. Tę konkurencję nazywamy w tym kontekście konkurencją wewnątrzgatunkową lub krócej — konkurencją wewnętrzną. Oprócz tego mamy też do czynienia z konkurencją zewnątrzgatunkową, albo konkurencją zewnętrzną, którą opisujemy podobnie do konkurencji wewnętrznej. Zakładamy więc, że osobniki konkurują ze sobą w trakcie spotkań między dwoma osobnikami różnego gatunku i liczba tych spotkań jest proporcjonalna do liczebności każdego z gatunków $\mathcal{P}_{i}$ , $i=1$ , $2$ . Wobec tego składnik konkurencji zewnętrznej w obu równaniach układu jest proporcjonalny do $N_{1}(t)N_{2}(t)$ .

Ostatecznie otrzymujemy dwóch układ równań

$\begin{array}[]{ccl}\dot{N}_{1}&=&r_{1}N_{1}\left(1-\frac{N_{1}}{K_{1}}-a_{{12}}\frac{N_{2}}{K_{2}}\right),\\ \dot{N}_{2}&=&r_{2}N_{2}\left(1-\frac{N_{2}}{K_{2}}-a_{{21}}\frac{N_{1}}{K_{1}}\right),\end{array}$

(8.1)

gdzie $a_{{ij}}$ są współczynnikami konkurencji zewnętrznej, przy czym dla uproszczenia obliczeń związanych z układem (8.1) współczynnik ten odnosi się do liczebności populacji $\mathcal{P}_{j}$ w stosunku do jej pojemności środowiska. Zauważmy, że dzięki takiemu zapisowi możemy łatwo przeprowadzić ubezwymiarowienie układu (8.1) wprowadzając zmienne $x_{i}=\tfrac{N_{i}}{K_{i}}$ , $i=1,2$ . Przypomnijmy, że taką zamianę zmiennych wprowadziliśmy w przypadku dyskretnego równania logistycznego, co umożliwiło ograniczenie rozważań do odcinka $[0,1]$ . W nowych zmiennych układ (8.1) przyjmuje postać

$\begin{array}[]{ccl}\dot{x}_{1}&=&r_{1}x_{1}\left(1-{x_{1}}-a_{{12}}{x_{2}}\right),\\ \dot{x}_{2}&=&r_{2}x_{2}\left(1-{x_{2}}-a_{{21}}{x_{1}}\right),\end{array}$

(8.2)

przy czym $x_{i}$ odzwierciedla procentową eksploatację środowiska przez gatunek $i$ .

Przeprowadzimy teraz analizę dynamiki modelu konkurujących gatunków opisanego przez układ (8.2) podobnie jak dla układu (6.1). Zaczniemy od zbadania podstawowych własności takich jak istnienie, jednoznaczność, nieujemność rozwiązań dla nieujemnych warunków początkowych.

Stwierdzenie 8.1

Niech $x_{i}(0)=x_{i}^{0}\geq 0$ , $i=1$ , $2$ , będzie warunkiem początkowym dla układu (8.2). Dla danego warunku początkowego istnieje jednoznaczne, nieujemne rozwiązanie określone dla wszystkich $t\geq 0$ . Co więcej, rozwiązanie jest ograniczone, $x_{i}(t)\leq\max\{ x_{i}^{0},1\}$ , $i=1$ , $2$ .

Istnienie i jednoznaczność łatwo wykazać korzystając z tego, że prawa strona układu jest funkcją klasy $\mathbf{C}^{1}$ . Nieujemność rozwiązań wynika np. z postaci logarytmicznej

$\frac{\dot{x}_{i}}{x_{i}}=r_{i}\left(1-{x_{i}}-a_{{ij}}{x_{j}}\right)\ \Longrightarrow\ x_{i}(t)=x_{i}^{0}\exp\left(r_{i}\int _{0}^{t}(1-{x_{i}(s)}-a_{{ij}}{x_{j}(s)})ds\right)\geq 0,$

natomiast przedłużalność wnioskujemy bezpośrednio z nieujemności, gdyż

$\dot{x}_{i}\leq r_{i}x_{i},$

a jeśli prawa strona układu ma liniowe oszacowanie, to rozwiązania istnieją dla wszystkich $t\geq 0$ . Co więcej, ograniczoność rozwiązań także wynika z ich nieujemności, gdyż zamiast powyższego liniowego oszacowania możemy wziąć

$\dot{x}_{i}\leq r_{i}x_{i}\left(1-{x_{i}}\right)$

i korzystając z nierówności różniczkowych oszacować z góry rozwiązania układu (8.1) przez rozwiązania równania logistycznego z pojemnością środowiska $K=1$ , analogicznie jak w przypadku układu (7.1).

∎

Po zbadaniu podstawowych własności możemy przejść do przeanalizowania portretu fazowego. Izokliny zerowe układu (8.2) wyznaczamy jako

$\dot{x}_{1}=0\Longleftrightarrow x_{1}=0$ lub ${x_{2}}=\frac{1}{a_{{12}}}\left(1-{x_{1}}\right)$ ;
$\dot{x}_{2}=0\Longleftrightarrow x_{2}=0$ lub ${x_{2}}=1-a_{{21}}{x_{1}}$ .

Wzajemne położenie izoklin nietrywialnych zależy od wielkości współczynników $a_{{ij}}$ . Mamy trzy istotnie różne przypadki generyczne

jeden ze współczynników jest większy, a drugi mniejszy niż 1 — założymy, że $a_{{12}}>1$ i $a_{{21}}<1$ , drugi układ parametrów implikuje taką samą dynamikę z dokładnością do zamiany miejscami gatunków $\mathcal{P}_{1}$ i $\mathcal{P}_{2}$ ;
$a_{{12}}\,,\ a_{{21}}>1$ ;
$a_{{12}}\,,\ a_{{21}}<1$ .

Oprócz wymienionych powyżej przypadków generycznych mamy też przypadki niegeneryczne, gdy co najmniej jeden ze współczynników konkurencji jest równy $1$ . Analizę tych przypadków pomijamy. Od strony biologicznej nieprawdopodobne wydaje się, żeby wartość jakiegokolwiek parametru utrzymywała się na stałym poziomie — zwykle występują drobne wahania, a co za tym idzie przypadki niegeneryczne są nieistotne z biologicznego punktu widzenia.

Zbadajmy najpierw stany stacjonarne i ich stabilność. Mamy zawsze trzy stany stacjonarne na brzegu przestrzeni fazowej $(\mathbb{R}^{+})^{2}$ : $(0,0)$ , $(1,0)$ i $(0,1)$ , które odzwierciedlają wymarcie co najmniej jednego z gatunków. Jeśli nietrywialne izokliny przecinają się, to mamy dodatni stan stacjonarny

$(\bar{x}_{1},\bar{x}_{2})=\left(\tfrac{1-a_{{12}}}{1-a_{{12}}a_{{21}}},\tfrac{1-a_{{21}}}{1-a_{{12}}a_{{21}}}\right)\,.$

Widzimy, że dodatni stan stacjonarny istnieje, jeśli albo oba współczynniki $a_{{12}}\,,\ a_{{21}}$ są poniżej $1$ , albo oba są powyżej. W celu zbadania stabilności wyznaczamy macierz Jacobiego

$MJ(x_{1},x_{2})=\left(\begin{array}[]{cc}r_{1}(1-2x_{1}-a_{{12}}x_{2})&-r_{1}a_{{12}}x_{1}\\ -r_{2}a_{{21}}x_{2}&r_{2}(1-2x_{2}-a_{{21}}x_{1})\end{array}\right)$

i obliczamy wartości własne

dla $(0,0)$ mamy $\lambda _{i}=r_{i}>0$ , $i=1,2$ ;
dla $(1,0)$ : $\lambda _{1}=-r_{1}<0$ , $\lambda _{2}=r_{2}(1-a_{{21}})$ i dla $(0,1)$ — symetrycznie;
dla $(\bar{x}_{1},\bar{x}_{2})$ wartości własne są rozwiązaniami równania charakterystycznego

$\lambda^{2}+\lambda(r_{1}\bar{x}_{1}+r_{2}\bar{x}_{2})+r_{1}r_{2}\bar{x}_{1}\bar{x}_{2}(1-a_{{12}}a_{{21}}).$

Wynika stąd, że $(0,0)$ jest zawsze węzłem niestabilnym, stabilność $(1,0)$ zależy od parametru $a_{{21}}$ — jeśli $a_{{21}}>1$ , to mamy węzeł stabilny, a jeśli $a_{{21}}<1$ , to siodło (analogicznie dla $(0,1)$ w zależności od $a_{{12}}$ ). Natomiast stabilność dodatniego stanu stacjonarnego zależy od iloczynu $a_{{12}}a_{{21}}$ — jeśli jest on większy niż 1, to $(\bar{x}_{1},\bar{x}_{2})$ jest siodłem, a dla nierówności przeciwnej — stabilnym węzłem. We wszystkich tych przypadkach następuje bifurkacja siodło – węzeł dla krytycznej wartości parametrów (czyli gdy są równe $1$ ).

Analizę zachowania rozwiązań w przestrzeni fazowej $(\mathbb{R}^{+})^{2}$ zaczniemy od przypadku 1, tj. $a_{{21}}<1<a_{{12}}$ . Przy takim układzie parametrów nietrywialne izokliny dla zmiennych $x_{1}$ i $x_{2}$ nie przecinają się w omawianej przestrzeni fazowej — nie istnieje dodatni stan stacjonarny. Wiemy już, że $(1,0)$ jest siodłem, a $(0,1)$ węzłem stabilnym. Izokliny dzielą przestrzeń fazową na 3 obszary: obszar $A$ pod izokliną dla $x_{1}$ , obszar $B$ pomiędzy izoklinami oraz $C$ — powyżej izokliny dla $x_{2}$ , por. rys. 8.1.

$Szkic portretu fazowego układu~\eqref{konk-x} w przypadku gdy $a_{{22}}<1<a_{{12}}$.$

Rys. 8.1. Szkic portretu fazowego układu (8.2) w przypadku gdy $a_{{22}}<1<a_{{12}}$ .

Zauważmy, że izokliny trywialne $x_{1}=0$ oraz $x_{2}=0$ stanowią rozwiązania układu. Rzeczywiście, jeśli $x_{1}^{0}=0$ , to $\dot{x}_{1}\equiv 0$ i stąd $x_{1}(t)=0$ dla $t\geq 0$ . Wtedy $\dot{x}_{2}=r_{2}x_{2}(1-x_{2})$ i $x_{2}$ spełnia równanie logistyczne, zatem na osi pionowej rozwiązanie zbiega do punktu $(0,1)$ . Dokładnie takie samo rozumowanie pokazuje, że na osi poziomej rozwiązanie zbiega do $(1,0)$ . Ponieważ ten ostatni punkt jest siodłem, więc wyznaczyliśmy w ten sposób rozmaitość stabilną. Z przebiegu pola wektorowego w otoczeniu punktu $(1,0)$ wnioskujemy, że interesujący nas fragment rozmaitości niestabilnej przebiega w obszarze $B$ . W obszarze $A$ obie zmienne rosną, nie mogą pozostać w tym obszarze, bo musiałyby być zbieżne do punktu o dodatnich współrzędnych na izoklinie dla zmiennej $x_{1}$ — ale taki punkt musiałby wtedy być stanem stacjonarnym, a nie ma dodatniego stanu stacjonarnego w tym przypadku. Wobec tego rozwiązanie przechodzi do obszaru $B$ . Z kolei jeśli punkt początkowy należy do $C$ , to obie zmienne maleją i analogiczne rozumowanie pokazuje, że albo rozwiązanie zostaje w $C$ i wtedy zbiega do stanu stacjonarnego $(0,1)$ , albo przechodzi do obszaru $B$ . Przebieg pola wektorowego pokazuje, że rozwiązanie nie może wyjść z obszaru $B$ , a ponieważ w tym obszarze jest także monotoniczne — $x_{1}$ maleje a $x_{2}$ rośnie — to zbiega do stanu stacjonarnego, czyli do $(0,1)$ . Rozmaitość niestabilna stanowi separatrysę, która oddziela rozwiązania zaczynające się w $A$ (dokładniej, jeśli przeprowadzimy analizę przebiegu takich rozwiązań dla $t<0$ , to zauważymy, że zbiegają one do $(0,0)$ przy $t\to-\infty$ ) i w związku z tym mające taką własność, że $x_{2}$ stale rośnie (pod izokliną dla zmiennej $x_{2}$ ), natomiast $x_{1}$ najpierw rośnie, osiąga wartość maksymalną przekraczając izoklinę, a następnie maleje — musi więc maleć do $0$ , od rozwiązań zaczynających się w $C$ i przechodzących do $B$ — takie rozwiązanie ma malejącą pierwszą współrzędną, natomiast druga współrzędna osiąga minimum na izoklinie, a następnie rośnie do $1$ . Ostatecznie wszystkie rozwiązania oprócz $x_{2}^{0}=0$ zbiegają asymptotycznie do stanu stacjonarnego $(0,1)$ i są od pewnego miejsca monotoniczne. Wnioskujemy, że w tym przypadku gatunek $\mathcal{P}_{1}$ ginie, natomiast $\mathcal{P}_{2}$ stabilizuje się na poziomie pojemności środowiska.

$Portrety fazowe dla układu~\eqref{konk} (po lewej) wraz z przykładowymi rozwiązaniami (po prawej).$

Rys. 8.2. Pole wektorowe układu (8.1) (po lewej) wraz z przykładowymi rozwiązaniami (po prawej) w przypadku gdy $a_{{21}}<1<a_{{12}}$ (u góry), $a_{{12}},a_{{21}}<1$ (po środku) lub $a_{{12}},a_{{21}}>1$ (u dołu).

Wspomnieliśmy już, że jeśli układ parametrów $a_{{12}}$ i $a_{{21}}$ zmienia się, tj. gdy $a_{{12}}<1<a_{{21}}$ to rozwiązania zachowują się symetrycznie — dostajemy więc zbieżność wszystkich dodatnich rozwiązań do $(1,0)$ . Opisany przebieg rozwiązań wiąże się z tym, że konkurencja zewnętrzna okazuje się znacznie bardziej niebezpieczna ( $a_{{12}}>1$ lub $a_{{21}}>1$ ) dla jednego z gatunków, a mało wpływa na drugi. Wtedy ten, dla którego jest niebezpieczna — ginie.

Jeśli oba parametry $a_{{12}}\,,\ a_{{21}}>1$ lub $a_{{12}}\,,\ a_{{21}}<1$ , to nietrywialne izokliny przecinają się tworząc dodatni stan stacjonarny. Jeśli jednak oba parametry są powyżej $1$ , stan ten jest siodłem, natomiast w przeciwnym przypadku — węzłem stabilnym. W przestrzeni fazowej $(\mathbb{R}^{+})^{2}$ możemy teraz wydzielić 4 obszary i analogicznie jak w przypadku 1., analizując przebieg pola wektorowego wnioskujemy, że w przypadku 3. dodatni stan stacjonarny jest globalnie stabilny we wnętrzu przestrzeni fazowej, natomiast w przypadku 2. rozmaitość stabilna tworzy separatrysę rozdzielającą baseny przyciągania stanów $(1,0)$ i $(0,1)$ — w zależności od warunku początkowego jeden gatunek przeżywa, a drugi ginie.

Podsumowując przeprowadzoną analizę zauważmy, że jeśli tylko $a_{{12}}a_{{21}}>1$ , co oznacza, że co najmniej jeden ze współczynników $a_{{12}}$ , $a_{{21}}$ przewyższa $1$ , to dynamika modelu odzwierciedla znaną zasadę ekologiczną mówiącą o konkurencyjnym wykluczaniu gatunków — jeśli nisze ekologiczne dwóch konkurujących gatunków zbytnio się pokrywają, to jeden z gatunków wypiera drugi. Tylko przy bardzo szczególnych układach parametrów możliwe jest współistnienie takich gatunków w tym samym siedlisku.

8.2. Modelowanie symbiozy

Kolejnym typem oddziaływań występujących między dwoma gatunkami $\mathcal{P}_{1}$ , $\mathcal{P}_{2}$ jest symbioza. Analogicznie jak w przypadku modelu konkurujących gatunków budując model heurystyczny zakładamy, że wewnętrzną dynamikę każdego z gatunków można opisać za pomocą modelu logistycznego, ze względu na ograniczoną pojemność środowiska, a co za tym idzie — występowanie konkurencji wewnątrzgatunkowej. Oddziaływania międzygatunkowe pojawiają się oczywiście wtedy, gdy spotykają się osobniki różnych gatunków i — ponownie podobnie do poprzednich przypadków — zakłada się, że wpływ tych oddziaływań na dynamikę każdego z gatunków zależy od liczby spotkań między osobnikami gatunków $\mathcal{P}_{1}$ oraz $\mathcal{P}_{2}$ .

Niech $N_{i}(t)$ oznaczają liczbę osobników (zagęszczenia) gatunków $\mathcal{P}_{i}$ , $i=1,2$ . Układ symbiotyczny możemy opisać za pomocą następujących równań

$\begin{array}[]{lcl}\dot{N}_{1}&=&r_{1}N_{1}\left(1-\frac{N_{1}}{K_{1}}+b_{{12}}N_{2}\right),\\ \dot{N}_{2}&=&r_{2}N_{2}\left(1-\frac{N_{2}}{K_{2}}+b_{{21}}N_{1}\right),\end{array}$

(8.3)

gdzie tradycyjnie $r_{i}$ oznacza współczynnik rozrodczości gatunku $i$ , pojemności środowiska dla danego gatunku oznaczamy $K_{i}$ , natomiast $b_{{ij}}$ opisuje siłę wpływu oddziaływań symbiotycznych dla gatunku $i$ w symbiozie z gatunkiem $j$ .

Zauważmy, że lokalne istnienie, jednoznaczność i nieujemność rozwiązań układu (8.3) dla nieujemnego warunku początkowego wykazujemy dość łatwo — stosując podobne argumenty jak w przypadku układu (8.1). Natomiast określoność rozwiązań dla wszystkich $t>0$ staje się tu własnością pożądaną, ale nie zawsze występującą. Rozpatrzmy następującą sytuację. Niech oba gatunki charakteryzują się tymi samymi parametrami, $r_{1}=r_{2}$ , $K_{1}=K_{2}$ i $b_{{12}}=b_{{21}}$ . Załóżmy też, że na początku liczebność obu gatunków jest jednakowa. Wobec tego, ze względu na całkowitą symetrię — liczebność tych gatunków jest jednakowa dla dowolnej chwili $t>0$ , dla której istnieje rozwiązanie układu (8.3). Skoro $N_{1}(t)=N_{2}(t)$ , to

$\dot{N}_{1}=r_{1}N_{1}\left(1-\frac{N_{1}}{K_{1}}+b_{{12}}N_{1}\right)$

i jeśli $b_{{12}}>\tfrac{1}{K_{1}}$ , to $N_{1}$ stale rośnie. Co więcej, $\dot{N}_{1}>\gamma N_{1}^{2}$ , gdzie $\gamma=b_{{12}}-\tfrac{1}{K_{1}}>0$ . Wiadomo natomiast, że rozwiązania równania $\dot{x}=\gamma x^{2}$ dla dodatniego warunku początkowego $x_{0}>0$ są określone w skończonym przedziale czasu $[0,\bar{t})$ , ponieważ dla $\bar{t}=\tfrac{1}{\gamma x_{0}}$ zachodzi $\lim _{{t\to\bar{t}^{{-}}}}x(t)=+\infty$ . Rzeczywiście, rozwiązując to równanie metodą rozdzielenia zmiennych dostajemy $x(t)=\tfrac{x_{0}}{1-\gamma x_{0}t}$ i dla $t=\bar{t}$ mianownik rozwiązania zeruje się. Tego typu zachowania możemy się spodziewać za każdym razem, gdy w modelu występuje składnik kwadratowy z dodatnim współczynnikiem.

Przeprowadzimy teraz analizę portretu fazowego układu (8.3) w zależności od parametrów. Izokliny zerowe zadane są za pomocą następujących prostych

dla zmiennej $N_{1}$ : $N_{1}=0$ lub $N_{2}=\frac{1}{b_{{12}}}\left(\frac{N_{1}}{K_{1}}-1\right)$ ;
dla zmiennej $N_{2}$ : $N_{2}=0$ lub $N_{2}=K_{2}\left(1+b_{{21}}{N_{1}}\right)$ .

Widzimy więc, że niezerowe izokliny są prostymi o dodatnim współczynniku kierunkowym, zatem w zależności od parametrów mogą się przecinać lub nie w pierwszej ćwiartce $\left(\mathbb{R}^{+}\right)^{2}$ , czyli naszej przestrzeni fazowej. Wobec tego układ ma albo 3 stany stacjonarne $(0,0)$ , $(K_{1},0)$ i $(0,K_{2})$ , jeśli izokliny się nie przecinają, albo 4 takie stany, jeśli się przecinają i wtedy mamy dodatni stan stacjonarny

$\bar{N}_{1}=\frac{1+b_{{12}}K_{2}}{1-b_{{12}}b_{{21}}K_{1}K_{2}}K_{1},\quad\bar{N}_{2}=\frac{1+b_{{21}}K_{1}}{1-b_{{12}}b_{{21}}K_{1}K_{2}}K_{2}.$

Zauważmy, że stan ten istnieje przy założeniu $b_{{12}}b_{{21}}K_{1}K_{2}<1$ , co w praktyce oznacza, że przy zadanych pojemnościach środowiska dla rozważanych gatunków zysk płynący z oddziaływań symbiotycznych nie może być zbyt duży, aby gatunki mogły pozostawać w pewnej równowadze, opisanej przez stan stacjonarny. W takim układzie mamy $\bar{N}_{i}>K_{i}$ , $i=1,2$ , czyli współrzędne dodatniego stanu stacjonarnego przewyższają pojemności środowiska.

Wyznaczmy teraz wartości własne odpowiadające poszczególnym stanom stacjonarnym. Macierz Jacobiego układu (8.3) jest równa

$MJ(N_{1},N_{2})=\left(\begin{array}[]{cc}r_{1}\left(1-2\frac{N_{1}}{K_{1}}+b_{{12}}N_{2}\right)&r_{1}b_{{12}}N_{1}\\ r_{2}b_{{21}}N_{2}&r_{2}\left(1-2\frac{N_{2}}{K_{2}}+b_{{21}}N_{1}\right)\end{array}\right)$

(8.4)

Podstawiając do wzoru (8.4) współrzędne poszczególnych stanów stacjonarnych otrzymujemy wartości własne

dla $(0,0)$ mamy $\lambda _{i}=r_{i}>0$ , $i=1,2$ ;
dla $(K_{1},0)$ : $\lambda _{1}=-r_{1}<0$ , $\lambda _{2}=r_{2}(1+b_{{21}}K_{1})>0$ i symetrycznie dla $(0,K_{2})$ : $\lambda _{1}=r_{1}(1+b_{{12}}K_{2})>0$ , $\lambda _{2}=-r_{2}<0$ ;
dla $(\bar{N}_{1},\bar{N}_{2})$ : $\lambda _{i}$ , $i=1,2$ są rozwiązaniami równania charakterystycznego

$\lambda^{2}+\left(r_{1}\frac{\bar{N}_{1}}{K_{1}}+r_{2}\frac{\bar{N}_{2}}{K_{2}}\right)\lambda+r_{1}r_{2}\frac{\bar{N}_{1}\bar{N}_{2}}{K_{1}K_{2}}(1-b_{{12}}b_{{21}}K_{1}K_{2})\Longrightarrow\lambda _{1},\ \lambda _{2}<0.$

Wobec tego $(0,0)$ jest węzłem niestabilnym, $(K_{1},0)$ i $(0,K_{2})$ są siodłami, natomiast dodatni stan stacjonarny, jeśli istnieje, to jest węzłem stabilnym.

$Szkic portretów fazowych oraz portrety fazowe układu~\eqref{symbioza}.$

Rys. 8.3. Szkic portretów fazowych układu (8.3) z zaznaczonymi izoklinami zerowymi i kierunkami przebiegu poszczególnych zmiennych (na górze) oraz portrety fazowe z naniesionymi przykładowymi krzywymi fazowymi (na dole) w przypadku gdy istnieje dodatni stan stacjonarny (po lewej) i gdy nie istnieje (po prawej).

Na rysunku 8.3 widzimy przebieg pola wektorowego $(\dot{N}_{1},\dot{N}_{2})$ w obu przypadkach. Jeśli izokliny nie przecinają się w analizowanej przestrzeni fazowej, to łatwo możemy wywnioskować, że obie współrzędne rozwiązania są nieograniczone i albo są stale rosnące, albo na początku maleją osiągając w pewnej chwili $t_{m}$ minimum, a potem rosną dla $t>t_{m}$ . Rzeczywiście, przestrzeń fazową możemy podzielić na trzy obszary, oznaczone na portrecie fazowym $A_{1}$ (ponad izokliną zerową dla $N_{2}$ ), $A_{2}$ (pod izokliną zerową dla $N_{1}$ ) i $B$ (pomiędzy izoklinami). W obszarze $A_{1}$ zmienna $N_{2}$ maleje, a $N_{1}$ rośnie. Wewnątrz tego obszaru nie ma stanów stacjonarnych, a stan stacjonarny na brzegu jest siodłem, przy czym jego rozmaitość stabilna pokrywa się z osią $N_{2}$ , zatem rozwiązanie nie może zbiegać do tego stanu. Wnioskujemy więc, że rozwiązanie opuszcza ten obszar i przechodzi do obszaru $B$ . Symetrycznie przebiegają orbity w obszarze $A_{2}$ , zatem także wchodzą do $B$ . Z kolei w $B$ obie zmienne są rosnące, przy czym rosną nieograniczenie — gdyby wzrost którejś z nich był ograniczony, to musiałaby być zbieżna, co oznacza, że pochodna też byłaby ograniczona, a nawet dążyłaby do $0$ , co jest niemożliwe ze względu na przebieg pola wektorowego.

Wróćmy teraz do problemu określoności rozwiązań dla wszystkich $t>0$ . Zauważmy, że gdyby np. zmienna $N_{1}(t)$ rosła do $+\infty$ w skończonym czasie, to oznaczałoby, że pewnej chwili musiałaby przeciąć izoklinę zerową zmiennej $N_{2}$ , co znów nie zgadza się z przebiegiem pola. Zatem obie zmienne są określone dla wszystkich $t>0$ i $\lim _{{t\to\infty}}N_{i}(t)=+\infty$ , $i=1,2$ , oczywiście poza przypadkiem omówionym powyżej, gdy poszczególne współczynniki dla obu gatunków pokrywają się.

Jeśli izokliny zerowe przecinają się w pierwszej ćwiartce układu współrzędnych, to przestrzeń fazowa dzieli się na 4 obszary, $A_{1}$ — ponad obydwoma izoklinami, $A_{2}$ — pod nimi, $B_{1}$ i $B_{2}$ — pomiędzy nimi. W obszarze $B_{1}$ obie zmienne rosną, natomiast w $B_{2}$ — obie maleją. Zauważmy, że każdy prostokąt o jednym wierzchołku w początku układu współrzędnych i dwóch bokach wzdłuż osi oraz naprzeciwległym wierzchołku w obszarze $B_{2}$ jest zbiorem niezmienniczym, gdyż pole wektorowe wchodzi do wewnątrz. Wobec tego każde rozwiązanie pozostaje ograniczone i analizując pole wektorowe dochodzimy do wniosku, że każde rozwiązanie układu (8.3) z dodatnim warunkiem początkowym zbiega do punktu krytycznego $(\bar{N}_{1},\bar{N}_{2})$ .

Zauważmy, że symbioza zawsze działa korzystnie na liczebność gatunków, które w taki sposób na siebie oddziałują, ale tylko w przypadku drugim, gdy efekt tych oddziaływań nie jest nadmierny, gatunki mogą współistnieć w środowisku. Natomiast jeśli oddziaływania symbiotyczne przekraczają pewien próg, to prowadzi do nieograniczonego wzrostu obu populacji, co oczywiście w warunkach rzeczywistych nie może mieć miejsca i musi zakończyć się katastrofą, por. rys. 8.4.

$Przebieg przykładowych rozwiązań układu~\eqref{symbioza}.$

Rys. 8.4. Przebieg przykładowych rozwiązań układu (8.3) w przypadku gdy istnieje dodatni stan stacjonarny (po lewej) i gdy ten stan nie istnieje (po prawej).

Treść automatycznie generowana z plików źródłowych LaTeXa za pomocą oprogramowania wykorzystującego LaTeXML.

Modelowanie matematyczne w biologii i medycynie

Zagadnienia

8. Modele oddziaływań między dwiema populacjami III

8.1. Układ konkurujących gatunków

Stwierdzenie 8.1

8.2. Modelowanie symbiozy