{"cells": [{"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["***\n", "# **EDO -- L3 MEU 2022/2023 -- Universit\u00e9 Paris-Saclay**\n", "***\n", "\n", "\n", "\n", "\n", "\n", "\n", "#  TP 1 : Repr\u00e9sentation des solutions d'une EDO, r\u00e9solution approch\u00e9e par le sch\u00e9ma d'Euler explicite"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["Soit $f:\\mathbb{R}\\times\\mathbb{R}\\longrightarrow\\mathbb{R}$ une fonction de classe $C^1$ et consid\u00e9rons l'\u00e9quation diff\u00e9rentielle\n", "$$\n", "(E)\\ \\ \\ \\ \\ \\ \\ y'(t)=f(t,y(t)).\n", "$$\n", "On remarque que le signe de la fonction $f$ dans le plan $(t,y)$ donne le sens de variation des solutions de l'\u00e9quation $(E).$ \n", "\n", "On appelle **champ de vecteurs** associ\u00e9 \u00e0 $(E)$ l'application \n", "$$\n", "V:\\mathbb{R}\\times\\mathbb{R}\\longrightarrow\\mathbb{R}\\times\\mathbb{R},\\ \\ \\ (t,y)\\mapsto\\big(1,f(t,y)\\big).\n", "$$\n", "\n", "On repr\u00e9sente le champ de vecteurs en un point $P=(\\bar{t},\\bar{y})$ du plan par une fl\u00e8che correspondant au segment $[P,P+\\varepsilon V(P)],$ avec $\\varepsilon$ petit, partant de $P$ et pointant vers $P+\\varepsilon V(P)$.\n", "\n", "On remarque alors que si $y(t)$ est une solution de $(E)$ qui passe par le point $P$, le vecteur $V(P)$ est tangent au graphe de $y$ au point $P.$\n", "\n", "Dans la premi\u00e8re partie de ce TP nous allons observer ces propri\u00e9t\u00e9s.\n", "\n", "## 1. Repr\u00e9sentation des solutions d'une EDO\n", "\n", "Nous nous int\u00e9ressons aux trois \u00e9quations suivantes :\n", "\n", "\\begin{equation}\n", "(1)\\ \\ \\ \\ \\ \\ \\ y'=-2ty,\n", "\\end{equation}\n", "\n", "\\begin{equation}\n", "(2)\\ \\ \\ \\ \\ \\ \\ y'=y-t^2,\n", "\\end{equation}\n", "\n", "\\begin{equation}\n", "(3)\\ \\ \\ \\ \\ \\ \\ y'=ty^2.\n", "\\end{equation}\n", "\n", "**Q1)** D\u00e9terminer les solutions de chacune des \u00e9quations 1, 2 et 3.\n", "\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q2)** Consid\u00e9rons l'\u00e9quation 1. Tracer sur le m\u00eame graphique, dans l'intrevalle $[-3,3],$ la solution de 1 qui v\u00e9rifie $y(0)=1$ et la solution de 1 qui v\u00e9rifie la condition initiale $y(0)=2.$ \n"]}, {"cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": ["import numpy as np\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline\n", "\n", "# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["### Trac\u00e9 du champ de vecteurs.\n", "\n", "**Q3)** Dans la figure de l'exercice 1, tracer dans le plan $(t,y)$ le champ de vecteurs associ\u00e9 \u00e0 l'\u00e9quation 1, en utilisant les commandes :\n", "`t=np.linspace(-3,3,31)\n", "y=np.linspace(0,2.2,23)\n", "T,Y=np.meshgrid(t,y)\n", "U=np.ones(T.shape)/np.sqrt(1+f1(T,Y)**2)\n", "V=f1(T,Y)/np.sqrt(1+f1(T,Y)**2)\n", "plt.quiver(T,Y,U,V,angles='xy',scale=20,color='blue')`\n", "\n", "Ci-dessus, `f` d\u00e9signe la fonction d\u00e9finissant le second membre de 1.\n", "Limiter le graphique \u00e0 la fen\u00eatre $[-3,3]\\times[0,2.2]$ \u00e0 l'aide de la commande `plt.axis`\n", "Rajouter des l\u00e9gendes, des labels pour les axes et un titre \u00e0 votre graphique."]}, {"cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q4)** Refaire le m\u00eame exercice pour l'\u00e9quation 3, en repr\u00e9sentant les solutions qui v\u00e9rifient $y(0)=1$ et  $y(0)=-1$"]}, {"cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["### Calcul d'une solution approch\u00e9e d'une EDO : votre premier sch\u00e9ma num\u00e9rique.\n", "\n", "On consid\u00e8re le probl\u00e8me de Cauchy pour l'\u00e9quation diff\u00e9rentielle $(E),$ \n", "\\begin{equation}\n", "(4)\\ \\ \\ \\ \\left\\{\\begin{aligned}\n", "  y'(t) & =  f(t,y(t)), \\\\\n", "  y(t_0) & =  y_0,\n", "\\end{aligned} \\right.\n", "\\end{equation}\n", "o\u00f9 l'instant initial $t_0\\in\\mathbb{R}$ et la valeur initiale $y_0\\in\\mathbb{R}$ sont donn\u00e9s. \n", "\n", "On suppose que ce probl\u00e8me a une unique solution $y:J\\longrightarrow\\mathbb{R}$ d\u00e9finie dans un intervalle $J\\subseteq\\mathbb{R}.$ On s'int\u00e9resse ici \u00e0 trouver num\u00e9riquement une solution approch\u00e9e de $y$, dans un intervalle de la forme $[t_0, t_f]\\subseteq J.$ \n", "\n", "Pour ce faire, on se donne $\\Delta t>0$ *petit* et on consid\u00e8re une subdivision uniforme de pas $\\Delta t$ de l'intervalle $[t_0, t_f].$ On peut le faire en se donnant $N\\in\\mathbb{N}$ *grand*, en posant $\\Delta t=\\frac{t_f-t_0}{N}$ et en consid\u00e9rant les $(N+1)$ points  \n", "$$\n", "t_0,\\ t_1:=t_0+\\Delta t,\\ \\cdots,\\ t_N:=t_0+N\\Delta t=t_f. \n", "$$\n", "\n", "On cherche alors des valeurs $y^0,\\ y^1,\\cdots,\\ y^N$ qui approchent la solution exacte de (PC) respectivement aux points $t_0,\\ t_1,\\cdots,\\ t_N.$ On va d\u00e9finir ces valeurs par un *sch\u00e9ma num\u00e9rique*.\n", "\n", "On pose $y^0=y_0.$ Soit $n\\in\\{1,\\dots,N\\}.$ Une id\u00e9e pour calculer une valeur $y^n$ approchant $y(t_n)$ (la solution exacte \u00e0 l'instant $t_n$) est d'approcher la d\u00e9riv\u00e9e de $y$ par un taux d'accroissement au point $t_n$ et d'utiliser l'\u00e9quation $(E).$ La solution exacte de $(E)$ v\u00e9rifie\n", "$$\n", "y'(t)=f(t,y(t)),\n", "$$\n", "en chaque point $t$ et \n", "$$\n", "y'(t)\\approx \\frac{y(t+\\Delta t)-y(t)}{\\Delta t},\n", "$$\n", "on remplace alors l'\u00e9quation $y'(t)=f(t,y(t))$ par l'\u00e9quation\n", " $$\\frac{y(t+\\Delta t)-y(t)}{\\Delta t}=f(t,y(t)),$$\n", "en chaque point $t_n$ de la subdivision consid\u00e9r\u00e9e, puisque l'on a alors\n", " $$\n", "\\frac{y(t_{n+1})-y(t_n)}{\\Delta t}\\approx f(t_n,y(t_n)).\n", " $$\n", "\n", "\n", "On obtient ainsi tr\u00e8s simplement un premier sch\u00e9ma num\u00e9rique : le sch\u00e9ma d'Euler explicite, pour lequel $y^0,\\dots,y^N$ sont d\u00e9finis comme suit\n", "\n", "** Sch\u00e9ma d'Euler explicite** :   \n", "$$\n", "\\begin{cases}\n", "y^0=y_0,\\\\\n", "y^{n+1}=y^{n}+\\Delta t f(t_n,y^{n}),\\ \\ \\ \\textrm{pour }\\ n\\in\\{0,\\dots,N-1\\}.\n", "\\end{cases}\n", "$$  \n", "\n", "\n", "**Q5)** \u00c9crire sur papier le sch\u00e9ma d'Euler associ\u00e9 \u00e0 l'\u00e9quation (1).\n", "\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q6)** \u00c9crire une fonction python de la forme\n", "\n", "`euler_exp(t0, tf, y0, deltaT)`\n", "\n", "prenant en argument  les extr\u00e9mit\u00e9s **t0** et **tf** de l'intervalle de temps $[t_0,t_f],$ la donn\u00e9e initiale $y_0$, le pas de temps $\\Delta t$. Cette fonction `euler_exp` devra retourner deux vecteurs : \n", "\n", "- $[t_0,\\, t_1,\\, ...,\\, t_N],$  tableau `numpy` unidimensionnel de taille $(N+1)\\times 1,$ repr\u00e9sentant la subdivision de l'intervalle $[t_0,t_f]$ de pas $\\Delta t$ consid\u00e9r\u00e9e, \n", "\n", "- $[y_0,\\, y_1,\\, ...,\\, y_N],$ tableau `numpy` de taille $(N+1)\\times 1,$ repr\u00e9sentant la solution approch\u00e9e aux instants $t_n,\\ n=0,\\dots,N.$\n"]}, {"cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q7)** Tracer sur une m\u00eame fen\u00eatre graphique pour la donn\u00e9e initiale $y(0)=1$ :\n", "\n", "- La solution exacte sur l'intervalle $[0,1]$ discr\u00e9tis\u00e9 avec un pas de temps de $10^{-4}$ ;\n", "\n", "- Les 3 solutions approch\u00e9es sur le m\u00eame intervalle obtenues avec la m\u00e9thode d'Euler avec des pas de temps $\\Delta t$ \u00e9gaux \u00e0 $1/4,\\ 1/10,\\ 1/100.$\n", "Rajouter des l\u00e9gendes \u00e0 votre graphique."]}, {"cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q8)** Modifier la fonction `euler_exp` (ou d\u00e9finissez une nouvelle fonction) de fa\u00e7on \u00e0 ce qu'elle prenne la fonction $f$ d\u00e9finissant l'\u00e9quation diff\u00e9rentielle $y'=f(t,y)$ en argument.  \n", "`euler\\_exp(t0, tf, y0, deltaT,f)`"]}, {"cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": ["# <completer>\n", "\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q9)** Reprenez la question **Q7)** avec votre nouvelle fonction.\n", "\n", "**Remarque** Pour programmer le sch\u00e9ma d'Euler, on peut aussi construire une fonction python de la forme `euler_exp(t0, tf, y0, N,f)`, en donnant comme argument le nombre de points $N$ que l'on consid\u00e8re dans la subdivision de l'intervalle $[t_0,t_f].$ Dans ce cas, il faudra d\u00e9finir le pas $\\Delta t$ \u00e0 l'int\u00e9rieur de la fonction."]}, {"cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["## 2. Mod\u00e8les de croissance d'une population\n", "\n", "**Mod\u00e8le de Malthus** \n", "\n", "Consid\u00e9rons une population d'individus isol\u00e9e dont on note $P(t)$ sa taille \u00e0 l'instant $t.$ \n", "\n", "Le mod\u00e8le le plus simple pour d\u00e9crire l'\u00e9volution de cette population est donn\u00e9 par l'\u00e9quation diff\u00e9rentielle\n", "\\begin{equation}\n", "P'(t)=(n-d)P(t),\n", "\\tag{M}\n", "\\end{equation}\n", "o\u00f9 $n$ et $d$ sont des constantes positives qui repr\u00e9sentent respectivement les taux de naissances et de d\u00e9c\u00e8s de la population.\n", "\n", "**Q1)** Soit $P_0\\in\\mathbb{R}^+$ la population \u00e0 l'instant initial $t=0.$ Donner la solution $P(t)$ de $(M)$ qui v\u00e9rifie $P(0)=P_0.$ Comment \u00e9volue la population lorsque $n>d$ ? Et lorsque $d>n$ ? Interpr\u00e9ter. \n", "\n", "\n", "**Mod\u00e8le de Verhulst** Un deuxi\u00e8me mod\u00e8le suppose que la croissance de la population est inhib\u00e9e par la limitation des ressources et que la population ne d\u00e9passe pas une taille critique donn\u00e9e par une constante $K>0.$ Si on note $P_0\\ge0$ la population \u00e0 l'instant initiale, l'\u00e9volution de la population est mod\u00e9lis\u00e9e par les \u00e9quations \n", "\n", "\\begin{equation}\n", "\\begin{cases}\n", "P'(t)=(n-d)P(t)\\big(1-\\frac{P(t)}{K}\\big),\\\\\n", "P(0)=P_0.\n", "\\end{cases}\n", "\\tag{V}\n", "\\end{equation}\n", "\n", " **Q2)** Montrer que le probl\u00e8me de Cauchy $V$ admet une unique solution maximale d\u00e9finie dans un intervalle $J\\subseteq\\mathbb{R}$ tel que $0\\in J.$\n", "\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["\n", " **Q3)** Quelle est la solution de $V$ lorsque  $P_0=0$ et lorsque $P_0=K$ ?\n", "\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q4)** Soit $P_0\\in (0,K)$, montrer que la solution maximale est globale \u00e0 droite et calculer $\\lim_{t\\to\\infty}P(t)$\n", "\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Q5)** On consid\u00e8re dans cette question $n=1,\\ d=0.75$ et $K=200.$ Utiliser le sch\u00e9ma d'Euler explicite pour obtenir une solution approch\u00e9e du probl\u00e8me $V$, de donn\u00e9e initiale $P_0=10,$ dans l'intervalle de temps $[0,50],$ avec un pas $\\Delta t=1.$ Repr\u00e9senter graphiquement la solution obtenue. \n", "\n", "Illustrer, dans un nouveau graphique, l'\u00e9volution de la population lorsque la donn\u00e9e initiale est telle que $P_0>K.$ Interpr\u00e9ter le r\u00e9sultat obtenu.\n"]}, {"cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "markdown", "metadata": {"deletable": false}, "source": ["**Un mod\u00e8le d'\u00e9volution d'une population de saumons**\n", "Le mod\u00e8le suivant d\u00e9crit l'\u00e9volution du nombre d'individus d'une population de saumons :\n", "\n", "\\begin{equation}\n", "\\begin{cases}\n", "P'(t)=(2-\\cos(t))P(t)-\\frac 12 P^2(t)-1,\\\\\n", "P(0)=P_0,\n", "\\end{cases}\n", "\\tag{S}\n", "\\end{equation}\n", "\n", "o\u00f9 $P_0$ est la population \u00e0 l'instant initial. Le terme positif p\u00e9riodique $2-\\cos(t)$ correspond \u00e0 un taux de naissances saisonnier. Le mod\u00e8le suppose que les d\u00e9c\u00e8s sont proportionnels au carr\u00e9 du nombre d'individus dans la population et le terme ind\u00e9pendante avec signe n\u00e9gatif $-1$ mod\u00e9lise une proportion d'individus p\u00e9ch\u00e9s.\n", "\n", "**Q6)** Utiliser le sch\u00e9ma d'Euler explicite pour obtenir une solution approch\u00e9e du probl\u00e8me $S$, de donn\u00e9e initiale $P_0=5,$ dans l'intervalle de temps $[0,10],$ avec un pas $\\Delta t=0.01.$ Repr\u00e9senter graphiquement la solution obtenue. \n"]}, {"cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": ["# <completer>\n"]}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": []}, {"cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": []}], "metadata": {"kernelspec": {"display_name": "Python 3", "language": "python", "name": "python3"}, "language_info": {"codemirror_mode": {"name": "ipython", "version": 3}, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.0"}, "celltoolbar": "None"}, "nbformat": 4, "nbformat_minor": 2}