Partie II : Systèmes bidimensionnels : écoulements dans le plan - Chapitre 9

Exemples avec des cycles limites

Toutes les versions de cet article : [English] [Español] [français]

Nous présentons divers exemples où des cycles limites attractifs apparaissent.

Cycles de proies et de prédateurs : le modèle Rosenzweig-McArthur

Nous revenons aux proies et aux prédateurs. Nous avons précédemment étudié un modèle de base pour l’interaction des requins et des sardines, à savoir .

$$ \begin{cases} \dot{x}= x(1-x) - xy\\ \dot{y}=\beta y ( x-\alpha) \end{cases} $$

$x$ représente la densité de sardines, $y$ la densité de requins, et $\alpha,\beta$ sont des paramètres positifs. Nous avons observé et prouvé qu’il ne peut y avoir de cycles limites dans ce modèle, ce qui signifie que les populations de requins et de sardines ne peuvent pas osciller périodiquement de manière robuste.

Une question naturelle se pose : quels facteurs biologiques créent des cycles limites dans un modèle proie-prédateur ? L’utilisation d’une réponse fonctionnelle plus réaliste est une possibilité. La réponse fonctionnelle est la vitesse à laquelle chaque prédateur capture ses proies. Dans le modèle ci-dessus, la réponse fonctionnelle est une fonction linéairement croissante de la densité des proies, ce qui signifie que les prédateurs ne sont jamais rassasiés ! Une façon simple de remédier à cela est de considérer une réponse fonctionnelle de la forme

$$ \frac{cx}{a+x} $$

$a,c$ sont des paramètres positifs.

Cette modification conduit au modèle suivant, que nous mettons sous une forme adimensionnée :

$$ \begin{cases} \dot{x}= x\left( 1-\frac{x}{\gamma}\right) - \frac{xy}{1+x}\\ \dot{y}=\beta y \left( \frac{x}{1+x}-\alpha\right) \end{cases} $$

$\alpha,\beta,\gamma$ sont des paramètres positifs. Ce modèle a été introduit par Rosenzweig et McArthur.

Invariance du quadrant positif. Nous ne nous intéressons qu’au quadrant positif pour la raison évidente que les densités de population sont des nombres non négatifs. Pour que le modèle soit bien défini, nous devons donc vérifier qu’en partant du quadrant positif, nous avons $(x(t),y(t))$ qui y demeure pour tout $t>0$. Nous avons

$$ \frac{\dot{x}}{x} = \frac{\text{d}}{\text{d} t}\ln x=1-\frac{x}{\gamma} -\frac{y}{1+x}. $$


En intégrant les deux côtés entre le temps $0$ et le temps $t$, on obtient

$$ x(t)=x_0 \, \exp\left( \int_0^t \left(1-\frac{x(s)}{\gamma} -\frac{y(s)}{1+x(s)}\right) \mathrm{d}s\right). $$

On en déduit que si $x_0=0$, alors $x(t)=0$ pour tout $t>0$, et que si $x_0>0$, alors $x(t)>0$ pour tout $t>0$. On obtient la même conclusion pour $y$. Ceci prouve plus que l’invariance du quadrant positif : elle prouve l’invariance de son intérieur et celle de sa frontière.

Avant de manipuler le modèle, déterminons les nulles-clines et les points fixes.

Nulles-clines.
En fixant $\dot{x}=0$, on trouve les $x$-nullclines :

$$ x=0\quad\text{and}\quad y=g(x)=(1+x)\left( 1-\frac{x}{\gamma}\right) $$

où nous avons introduit la fonction $g(x)$ pour plus de commodité. Cette fonction est une parabole dont le sommet est atteint pour

$$ x=\frac{\gamma-1}{2}. $$

En fixant $\dot{y}=0$, on trouve les nulles-clines pour $y$ :

$$ y=0\quad\text{and}\quad x=\frac{\alpha}{1-\alpha}. $$


Points fixes.
Pour toutes les valeurs des paramètres, $(0,0)$ et $(\gamma,0)$ sont des points fixes. Le premier correspond à l’absence des deux populations. Le second correspond à l’absence de prédateurs, alors que les proies ont une densité $\gamma$. Il existe éventuellement un troisième point fixe :

$$ (\bar{x}, g(\bar{x})) $$



$$ \bar{x}=\frac{\alpha}{1-\alpha}. $$

Comme les densités de population sont des nombres positifs, ce point fixe n’existe que si $\bar{x}\leq \gamma$. Cette condition signifie simplement que les nulles-clines $x=\alpha/(1-\alpha)$ et $y=g(x)$ se coupent à l’intérieur du quadrant positif.

On suppose en outre que $\alpha<1$, car sinon toutes les solutions partant de densités initiales strictement positives tendent vers $(\gamma,0)$.

La raison en est simple : (…)

La raison en est simple : pour tout $x>0$, $x/(1+x)<1$, et si $\alpha>1$ alors $\dot{y}<0$, d’où $y(t)$ à $0$ lorsque $t$ à $+\infty$. Par conséquent, l’équation pour $x$ se réduit à l’équation logistique $\dot{x}=x\left( 1-\frac{x}{\gamma}\right)$, d’où $x(t)\to\gamma$ lorsque $t$ à $+\infty$.


Une bifurcation de Hopf supercritique. Bien qu’il ne fasse aucun doute qu’un cycle limite existe dans ce modèle, appliquons le le théorème de bifurcation d'Andronov-Hopf, avec $\alpha$ comme paramètre de bifurcation. (Nous pourrions également appliquer le théorème de Poincaré-Bendixson.) On se contentera de vérifier que la trace de la matrice jacobienne au point fixe $(\bar{x},g(\bar{x}))$ change de signe lorsque $\alpha$ passe par la valeur $(\gamma-1)/(\gamma+1)$, alors que le déterminant est toujours positif. Au passage, observons que pour voir la bifurcation se produire dans le quadrant positif, $\gamma$ doit être strictement supérieur à $1$, ce qui est en accord avec l’expérience numérique. Nous supposons donc que $\gamma>1$.

La matrice jacobienne à $(\bar{x},g(\bar{x}))$ est

$$ A= \begin{pmatrix} \alpha g’(\bar{x}) & -\alpha \\ \frac{\beta}{(1+\bar{x})^2} g(\bar{x}) & 0 \end{pmatrix}. $$


Puisque $g(x)>0$ chaque fois que $x\in [0,\gamma[$, et $\bar{x}\in\thinspace ]0,\gamma[$, nous avons

$$ \text{det}(A)=\frac{\alpha\beta}{(1+\bar{x})^2} g(\bar{x})>0. $$

On a

$$ \text{tr}(A)=\alpha g’(\bar{x})=\frac{\gamma-1-2\bar{x}}{\gamma}. $$

D’où

$$ \text{tr}(A) \begin{cases} <0 & \text{if}\quad \bar{x}>\frac{\gamma-1}{2}\\ =0 & \text{if}\quad \bar{x}=\frac{\gamma-1}{2}\\ >0 & \text{if}\quad \bar{x}<\frac{\gamma-1}{2}. \end{cases} $$

Le paramètre de bifurcation $\alpha^*$ se trouve en résolvant l’équation suivante

$$ \frac{\alpha}{1-\alpha}=\frac{\gamma-1}{2} $$

ce qui donne

$$ \alpha^*=\frac{\gamma-1}{\gamma+1}. $$

Un modèle de glycolyse

Nous examinons un modèle simple d’un processus biochimique fondamental appelé glycolyse dans lequel les cellules obtiennent de l’énergie en décomposant le sucre. Sous forme adimensionnelle, les équations sont

$$ \begin{cases} \dot{x}= -x+ay+x^2 y\\ \dot{y}=b-ay-x^2 y \end{cases} $$

$x$ et $y$ sont, respectivement, les concentrations d’ADP (adénosine diphosphate) et de F6P (fructose-6-phosphate), et $a,b>0$ sont des paramètres cinétiques. Il existe un seul point fixe

$$ \left(b\, ,\frac{b}{a+b^2}\right) $$


situé à l’intersection des nulles-clines

$$ y=\frac{x}{a+x^2}\quad\text{and}\quad y=\frac{b}{a+x^2}. $$

Nous pouvons observer qu’un cycle limite apparaît pour certaines valeurs de $(a,b)$. (Nous expliquons ci-dessous la forme de la région du plan $ab$ dans laquelle nous avons un cycle limite.)


Le théorème de Poincaré-Bendixson en action. Nous allons appliquer ce théorème pour prouver rigoureusement qu’il existe un cycle limite dans ce système. Pour atteindre ce but, nous construisons une région de piégeage. Nous affirmons que la « boîte » délimitée par la ligne rouge représentée sur la figure est une telle région.


Il faut montrer que tous les vecteurs de la frontière de la boîte pointent dans celle-ci. Sur les côtés horizontaux et verticaux, c’est laissé au lecteur. La partie délicate consiste à traiter la diagonale de pente $-1$ s’étendant du point $(b,b/a)$ à la nullité $y=x/(a+x^2)$. Considérons la quantité $\dot{x}-(-\dot{y})$. Nous trouvons

$$ \dot{x}-(-\dot{y})= -x+ay+x^2y+(b-ay-x^2y)=b-x. $$

D’où

$$ -\dot{y}>\dot{x}\quad\text{if}\quad x>b. $$

Cela implique que le champ de vecteurs pointe vers l’intérieur sur la diagonale de la figure, car

$$ \frac{\text{d}y}{\text{d}x}=\frac{\text{d}y/\text{d}t}{\text{d}x/\text{d}t} =\frac{\dot{y}}{\dot{x}} $$

est plus négatif que $-1$. Cela signifie que les vecteurs sont plus raides que la ligne diagonale. Nous concluons que la boîte est une région de piégeage.

Peut-on maintenant appliquer le théorème de Poincaré-Bendixson et conclure qu’il existe un cycle limite dans la boîte ? Pas tout à fait, car il y a un point fixe à l’intérieur ! Mais supposons que le point fixe soit un répulseur (source nodale ou spirale), alors nous sommes en bonne voie en considérant une version « perforée » de la boîte, comme le montre la figure.

Le répulseur entraîne toutes les solutions voisines dans la région hachurée, et cette région est exempte de points fixes : Le théorème de Poincaré-Bendixson s’applique !

Maintenant nous regardons sous quelles conditions le point fixe est un répulseur. La matrice jacobienne au point fixe est

$$ A= \begin{pmatrix} -1+\frac{2b^2}{a+b^2} & a+b^2 \\ -\frac{2b^2}{a+b^2} & -(a+b^2) \end{pmatrix}. $$

Nous avons $\text{det}(A)=a+b^2>0$, donc le point fixe n’est jamais une selle. Ensuite, nous trouvons

$$ \text{tr}(A)=-\frac{b^4+(2a-1)b^2+a+a^2}{a+b^2}. $$

Le point fixe est répulsif pour $\text{tr}(A)>0$ et attractif pour $\text{tr}(A)<0$. La ligne de partage $\text{tr}(A)=0$ dans le plan $ab$ est formée par les deux courbes

$$ b=\sqrt{\frac{1}{2}\left(1-2a\pm\sqrt{{\scriptsize 1-8a}}\right)}\, . $$

Un modèle pour une réaction chimique qui oscille

La réaction Belousov-Zhabotinsky a constitué un tournant majeur dans l’histoire de la chimie, car elle a montré que certaines réactions chimiques pouvaient osciller pendant de longues périodes avant de se stabiliser à l’équilibre. Elle a été découverte par Belousov dans les années 1950, mais ses travaux n’ont été diffusés qu’à la fin des années 1960. Avant cela, les chimistes pensaient que toutes les réactions chimiques tendaient de manière monotone vers l’équilibre.

Il existe désormais une grande catégorie de réactions chimiques présentant des oscillations. Une réaction particulièrement simple est celle du dioxyde de chlore et de l’acide iode-malonique (ClO$_2$-I$_2$-MA). Les équations différentielles exactes qui modélisent cette réaction sont extrêmement compliquées. Heureusement, il existe un moyen d’approximer la réaction car la concentration de certains des réactifs varie beaucoup plus lentement que d’autres. L’idée est de traiter la concentration des réactifs lents comme constante. Cette approximation permet de capturer les phénomènes d’oscillations. Après un adimensionnement approprié, on obtient le modèle bidimensionnel suivant :

$$ \begin{cases} \dot{x}= a-x-\frac{4xy}{1+x^2}\\ \dot{y}=bx\left( 1-\frac{y}{1+x^2}\right) \end{cases} $$

$x$ et $y$ sont les concentrations sans dimension de $\text{I}^-$ et de $\text{ClO}_2^-$, respectivement, et $a,b$ sont des paramètres positifs.

On a $\dot{x}=0$ sur la courbe

$$ y=\frac{(a-x)(1+x^2)}{4x} $$

et $\dot{y}=0$ sur l’axe $y$ et sur la parabole

$$ y=1+x^2. $$

Ces deux courbes se coupent en un seul point, nous n’avons donc qu’un seul point fixe :

$$ \left( \frac{a}{5}\, , \,1+\left(\frac{a}{5}\right)^2\right). $$

Vous pouvez observer que pour certains paramètres, il existe un cycle limite attractif. On peut également observer une bifurcation de Hopf supercritique.


Le lecteur peut chercher une région piégeuse, comme dans l’exemple précédent, pour appliquer le théorème de Poincaré-Bendixson. On peut montrer que la courbe

$$ b=\frac{3a}{4}-\frac{25}{a} $$

divise le plan $ab$ en deux régions : pour $(a,b)$ en dessous de la courbe, on a un cycle limite attractif entourant le point fixe qui est un repoussoir (source spirale), et au-dessus de la courbe, il n’y a pas de cycle limite et le point fixe est attractif (c’est un puits spirale).

Un modèle proie-prédateur avec effet Allee pour les proies

Nous considérons un modèle proie-prédateur qui est obtenu à nouveau à partir du modèle

$$ \begin{cases} \dot{x}= x(1-x) - xy\\ \dot{y}=\beta y ( x-\alpha). \end{cases} $$

Cette fois, nous ne modifions pas le terme d’interaction $xy$ (comme nous l’avons fait précédemment pour obtenir le modèle Rosenzweig-McArthur), mais le terme logistique $x(1-x)$. Nous supposons qu’il existe un seuil de croissance pour la population de proies (en l’absence de prédateurs). La manière la plus simple de le faire est de prendre

$$ \dot{x} = x\left(\frac{x}{\gamma}-1\right)(1-x) $$

$0<\gamma<1$. Ce modèle est facile à analyser graphiquement.

Nous arrivons donc au modèle suivant :

$$ \begin{cases} \dot{x} = x\left(\frac{x}{\gamma}-1\right)(1-x)-xy\\ \dot{y} = \beta y (x-\alpha) \end{cases} $$

$\alpha,\beta,\gamma$ sont des paramètres positifs avec $\gamma<1$. En procédant comme pour le modèle de Rosenzweig-McArthur, on peut vérifier que le quadrant positif est invariant (en fait, son intérieur et sa frontière sont tous deux invariants).

Points fixes Il y a trois points fixes qui existent pour toutes les valeurs des paramètres :

$$ (0,0),\ ; (\gamma,0),\ ; (1,0). $$

Un quatrième peut entrer en jeu si les lignes nulles $x=\alpha$ et $y=g(x)=\big(\frac{x}{\gamma}-1\big)(1-x)$ se coupent à l’intérieur du quadrant positif, à savoir

$$ (\alpha,g(\alpha)), $$

ce qui est le cas si et seulement si $\gamma\leq \alpha\leq 1$. Cela correspond à la coexistence des deux populations. Lorsque $\alpha$ varie, nous allons observer une bifurcation de Hopf (supercritique).


Le comportement observé des points fixes peut être confirmé par linéarisation. On peut également vérifier le mécanisme de base de la bifurcation de Hopf.


La matrice jacobienne en un point $(x,y)$ est .

$$ A= \begin{pmatrix} xg’(x)+g(x)-y & -x \\ \beta y & \beta(x-\alpha) \end{pmatrix}. $$

A $(x,y)=(0,0)$ on a $\text{det}(A)=\alpha\beta>0$ et $\text{tr}(A)=-1-\alpha\beta<0$, donc c’est un puits.
Le point fixe $(\gamma,0)$ est une source si $\alpha<\gamma$, et une selle si $\alpha>\gamma$.
Le point fixe $(1,0)$ est une selle si $\alpha<1$, et un puits si $\alpha>1$.
Le point fixe $(\alpha,g(\alpha))$, qui existe seulement si $\alpha\in [\gamma,1]$, est un puits si $\alpha>\frac{\gamma+1}{2}$, et une source si $\alpha>\frac{\gamma+1}{2}$. On peut vérifier que $\text{det}(A)=\alpha\beta g(\alpha)>0$ pour tout $\alpha\in ]\gamma,1[$, alors que $\text{tr}(A)$ change de signe de positif à négatif lorsque $\alpha$ augmente
et passe par la valeur $(\gamma+1)/2$.

Neurodynamique : le modèle de FitzHugh-Nagumo.

Le modèle de FitzHugh-Nagumo est un modèle du modèle de Hodgkin-Huxley. Bien que ce modèle ne soit pas aussi biologiquement précis que le modèle original, il capture néanmoins le comportement essentiel des impulsions nerveuses, y compris le phénomène appelé excitabilité, qui sera décrit ci-dessous. Le modèle est défini par les équations

$$ \begin{cases} \dot{x}= x-\frac{x^3}{3} -y +I \\ \dot{y}=\frac{1}{\tau}(x+\alpha-\beta y). \end{cases} $$

Dans ces équations, $x$ est similaire à la tension et représente l’excitabilité du système ; la variable $y$ représente une combinaison d’autres forces qui tendent à ramener le système au repos (rétroaction). Il existe quatre paramètres positifs : $\alpha,\beta,\tau,I$. Le paramètre $I$ est un paramètre de stimulation qui entraîne l’excitation du système ; $I$ est comme un courant de stimulation appliqué. Le paramètre $\tau$ est une échelle de temps qui permet de faire évoluer $x$ plus rapidement que $y$ ; en gros $\dot{x}/\dot{y}\approx \tau$.

La courbe sur laquelle $\dot{x}=0$ est la fonction cubique

$$ y=x-\frac{x^3}{3}+I $$

et $\dot{y}=0$ sur la courbe linéaire

$$ y=\frac{x+\alpha}{\beta}. $$

Il existe un seul point fixe $(\bar{x},\bar{y})$ qui est la solution d’une équation cubique.

Les valeurs typiques des paramètres $\alpha,\beta,\tau$ sont respectivement de 0,7$, 0,8,13$. Nous ne réglerons que le paramètre $I$. Nous observerons ce que l’on appelle le phénomène de bloc d’excitation.

Un exemple avec une infinité de cycles limites.

Nous considérons le système

$$ \begin{cases} \dot{x}= y + (x^2+y^2)^2 x\sin\left(\frac{1}{x^2+y^2}\right)\\\\ \dot{y}= -x + (x^2+y^2)^2 y\sin\left(\frac{1}{x^2+y^2}\right). \end{cases} $$

L’expérience numérique suggère clairement qu’il existe une infinité de cycles limites circulaires.


En passant en coordonnées polaires, on peut analyser le système assez facilement :

$$ \begin{cases} \dot{r}=r^5 \sin\left(\frac{1}{r^2}\right)=f(r)\\ \dot{\theta}=-1. \end{cases} $$

Nous obtenons des équations découplées. La seconde est simplement une rotation dans le sens inverse des aiguilles d’une montre à vitesse constante. L’équation pour $r(t)$ a des points fixes en

$$ \frac{1}{r^2}=n\pi,\, n\, \text{entier positif quelconque}. $$


Appliquons le test de la dérivée au système unidimensionnel $\dot{r}=f(r)$. Nous avons

$$ f’(r) = r^2\left( 5 r^2 \sin\left(\frac{1}{r^2}\right)-2 \cos\left(\frac{1}{r^2}\right)\right). $$


Aux points $\frac{1}{r^2}=n\pi$, on a $\sin\left(\frac{1}{r^2}\right)=0$ et

$$ f’(r)=-2r\cos(n\pi)\thinspace\thinspace\thinspace\thinspace\thinspace \begin{cases} \text{positif} & \text{si}\thinspace\thinspace \thinspace\text{impair}\\ \text{négatif} & \text{si}\thinspace\thinspace \thinspace \text{pair}. \end{cases} $$


Par conséquent, $\frac{1}{\sqrt{2\pi}}, \frac{1}{\sqrt{4\pi}}$,... sont des points fixes stables, et $\frac{1}{\sqrt{\pi}}, \frac{1}{\sqrt{3\pi}}$,... sont des points fixes instables. En revenant aux coordonnées cartésiennes, cela signifie que nous avons une infinité de cycles limites qui sont des cercles de rayon $\frac{1}{\sqrt{n\pi}}$, où ceux avec $n$ pair sont des cycles limites attractifs, ceux avec $n$ impair sont des cycles limites répulsifs.

Les cycles limites ne sont pas tous des courbes convexes.

Tous les exemples de cycles limites dans le plan que nous avons vus jusqu’à présent sont des courbes convexes. Une courbe est convexe si elle se trouve complètement d’un côté de chacune de ses tangentes. L’exemple suivant montre que ce n’est pas forcément le cas, bien que de tels exemples ne soient pas si faciles à découvrir !


Voulez-vous vraiment voir les équations de ce système ? Les voilà :

$$ \begin{cases} \dot{x}= y\\ \dot{y}=-x-y\big( a_2 x^2+a_4 x^4+a_6 x^6+ a_8 x^8 + a_{10} x^{10} + a_{12} x^{12} + a_{14} x^{14}\big) \end{cases} $$

$$ a_{2}=90, a_4=-882, a_6=2598.4, a_8=-3359.997, a_{10}=2133.34, a_{12}=-651.638, a_{14}=76.38. $$