chap10.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Programmation Python  pour les mathématiques"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Julien Guillod**, [Sorbonne Université](http://www.sorbonne-universite.fr/),\n",
    "Licence <a href=\"https://creativecommons.org/licenses/by-nc-nd/4.0/\">CC BY-NC-ND</a>\n",
    "\n",
    "L'entier des chapitres est disponible au format\n",
    "[HTML](https://python.guillod.org/) et [PDF](https://python.guillod.org/python.pdf).\n",
    "Ce notebook peut également être exécuté sur [GESIS](https://notebooks.gesis.org/binder/v2/gh/guillod/python/master?urlpath=lab/tree/chap10.ipynb)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Équations différentielles"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"ch:equadiff\"></div>\n",
    "\n",
    "Le but est d'introduire les méthodes de base permettant de résoudre des équations différentielles ordinaires du premier ordre du type:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{align*}\n",
    "\\dot{\\boldsymbol{x}}(t) &= f(t, \\boldsymbol{x}(t)) \\,, & \\boldsymbol{x}(0) &= \\boldsymbol{x}_0 \\,,\n",
    "\\end{align*}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où $f:\\mathbb{R}^+ \\times \\mathbb{R}^n \\to \\mathbb{R}^n$ est une fonction assez régulière et $\\boldsymbol{x}_0\\in\\mathbb{R}^n$ une donnée initiale.\n",
    "À noter que les équations différentielles ordinaires d'ordre supérieur peuvent être mises sous la forme précédente du premier ordre.\n",
    "\n",
    "**Concepts abordés:**\n",
    "\n",
    "* méthodes d'Euler\n",
    "\n",
    "* méthodes de Runge-Kutta\n",
    "\n",
    "* équation aux dérivées partielles non linéaire\n",
    "\n",
    "* différences finies\n",
    "\n",
    "* méthodes adaptatives"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.1: Méthodes d'Euler"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"exer:equdiff-euler\"></div>\n",
    "\n",
    "L'idée la plus simple pour résoudre de manière approchée une équation différentielle ordinaire est de discrétiser le temps avec un pas $h$ et d'approximer la dérivée temporelle sur chaque intervalle de longueur $h$. Il y a deux façons simples d'approximer la dérivée en temps. La première est l'approximation par différence finie avant:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\dot{\\boldsymbol{x}}(t) \\approx \\frac{\\boldsymbol{x}(t+h)-\\boldsymbol{x}(t)}{h} \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "la seconde par différence finie arrière:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\dot{\\boldsymbol{x}}(t) \\approx \\frac{\\boldsymbol{x}(t)-\\boldsymbol{x}(t-h)}{h} \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Les inconnues étant les évaluations de la solution $\\boldsymbol{x}$ aux temps $t_i = ih$ pour $i\\geq0$, c'est-à-dire $\\boldsymbol{x}_i = \\boldsymbol{x}(t_i)$. L'équation différentielle peut ainsi être approchée à l'aide des différences finies avant par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\frac{\\boldsymbol{x}_{i+1}-\\boldsymbol{x}_i}{t_{i+1}-t_i} = f(t_i, \\boldsymbol{x}_i) \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "ce qui donne la formule d'Euler explicite:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\boldsymbol{x}_{i+1} = \\boldsymbol{x}_i + (t_{i+1}-t_i) f(t_i, \\boldsymbol{x}_i) \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Avec l'approximation par différences finies arrière, on obtient la méthode d'Euler implicite:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\boldsymbol{x}_{i} = \\boldsymbol{x}_{i-1} + (t_i - t_{i-1}) f(t_i, \\boldsymbol{x}_i) \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La formule d'Euler explicite permet de calculer directement tous les $\\boldsymbol{x}_i$ par récurrence en connaissant $\\boldsymbol{x}_0$. En revanche la formule d'Euler implicite nécessite à chaque pas de temps la résolution d'une équation non linéaire pour $\\boldsymbol{x}_i$ par exemple avec la méthode de Newton.\n",
    "\n",
    "**a)**\n",
    "Écrire une fonction `euler_explicit(f,x0,t)` qui pour une donnée initiale `x0` retourne les valeurs $\\boldsymbol{x}_0, \\boldsymbol{x}_1, \\dots, \\boldsymbol{x}_m$ calculées avec la méthode d'Euler explicite aux temps $(t_i)_ {i=0}^m$ représentés par le vecteur `t`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**b)**\n",
    "Utiliser la méthode d'Euler explicite pour résoudre l'équation différentielle:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\dot{x}(t) + x(t) = \\sin(t) \\,, \\quad x(0)=1 \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "pour $t\\in[0,10]$. Comparer les résultats avec la solution exacte:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "x(t) = \\frac{1}{2}\\big( \\sin(t) - \\cos(t) + 3e^{-t} \\big) \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "pour différentes discrétisations du temps."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**c)**\n",
    "Résoudre le problème précédent avec la méthode d'Euler implicite.\n",
    "\n",
    "**Indication:**\n",
    "Vu que l'équation précédente est linéaire, on peut en fait rendre explicite la méthode d'Euler implicite en résolvant l'équation implicite à la main."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**d)**\n",
    "<span style=\"color:red\">!!</span> Définir une fonction `euler_implicit(f, Dxf, x0, t)` implémentant la méthode d'Euler implicite pour des équations non linéaires. À noter que pour résoudre le problème non linéaire avec la méthode de Newton, la dérivée de $f$ selon $\\boldsymbol{x}$ est nécessaire.\n",
    "\n",
    "**Indication:**\n",
    "Il est également possible d'utiliser l'algorithme de recherche de zéro de fonction `scipy.optimize.fsolve` qui ne nécessite pas de connaître la dérivée de $f$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**e)**\n",
    "<span style=\"color:red\">!</span> Utiliser les méthodes précédentes pour trouver une solution approchée du système:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{align*}\n",
    "\\dot{x}(t)+\\cos(y(t)) & =\\sin(t)\\,, & x(0) & =1\\,,\\\\ \n",
    "\\dot{y}(t)+\\cos(x(t)) & =0\\,, & y(0) & =0\\,.\n",
    "\\end{align*}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Réponse:**\n",
    "Voici le graphique de la solution de l'équation différentielle avec un pas de temps $h=0.1$:\n",
    "\n",
    "<center><img src=\"https://python.guillod.org/fig/euler-2d.png\" style=\"width:90%;max-width:800px;\"></center>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.2: Méthodes de Runge-Kutta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"exer:equadiff-RK\"></div>\n",
    "\n",
    "Le but de cet exercice est d'introduire une classe de méthodes plus précises que les méthodes d'Euler pour résoudre des équations différentielles ordinaires. Au lieu de faire une approximation au premier ordre en $h$ l'idée est de faire une approximation d'un ordre supérieur.\n",
    "\n",
    "L'idée de base est de construire une suite $\\boldsymbol{x}_i$ donnant une approximation de la solution de $\\dot{\\boldsymbol{x}}(t)=f(t,\\boldsymbol{x})$ au temps $t_i$ pour $i\\in\\mathbb{N}$. Cette suite est définie par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\boldsymbol{x}_{i+1} = \\boldsymbol{x}_i + M(t_i, \\boldsymbol{x}_i, t_{i+1}-t_i) \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "pour une certaine fonction $M$ appelée méthode.\n",
    "Par exemple pour la méthode d'Euler explicite, la fonction $M$ est donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "M(t, \\boldsymbol{x}, h) = h f(t,\\boldsymbol{x})\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Une méthode de Runge-Kutta d'ordre deux est donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "M(t,\\boldsymbol{x}, h) = h f\\bigg(t+\\frac{h}{2}, \\boldsymbol{x} + \\frac{h}{2} f(t,\\boldsymbol{x}) \\bigg) \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Une méthode de Runge-Kutta d'ordre quatre est donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "M(t,\\boldsymbol{x}, h) = \\frac{h}{6} \\big(\\boldsymbol{k}_1 + 2\\boldsymbol{k}_2 + 2\\boldsymbol{k}_3 + \\boldsymbol{k}_4\\big) \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{align*}\n",
    "\\boldsymbol{k}_{1} & =f(t,\\boldsymbol{x})\\,,\\\\ \n",
    "\\boldsymbol{k}_{2} & =f\\bigg(t+\\frac{h}{2},\\boldsymbol{x}+\\frac{h}{2}\\boldsymbol{k}_{1}\\bigg)\\,,\\\\ \n",
    "\\boldsymbol{k}_{3} & =f\\bigg(t+\\frac{h}{2},\\boldsymbol{x}+\\frac{h}{2}\\boldsymbol{k}_{2}\\bigg)\\,,\\\\ \n",
    "\\boldsymbol{k}_{4} & =f(t+h,\\boldsymbol{x}+h\\boldsymbol{k}_{3})\\,.\n",
    "\\end{align*}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "À noter que plus généralement, une méthode de Runge-Kutta d'ordre $s$ est donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "M(t,\\boldsymbol{x},h) = h \\sum_{i=1}^s b_i \\boldsymbol{k}_i \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{align*}\n",
    "\\boldsymbol{k}_{1} & =f(t,\\boldsymbol{x})\\,,\\\\ \n",
    "\\boldsymbol{k}_{2} & =f(t+c_{2}h,\\boldsymbol{x}+ha_{21}\\boldsymbol{k}_{1})\\,,\\\\ \n",
    "\\boldsymbol{k}_{3} & =f(t+c_{3}h,\\boldsymbol{x}+h(a_{31}\\boldsymbol{k}_{1}+a_{32}\\boldsymbol{k}_{2}))\\,,\\\\ \n",
    " & \\vdots\\\\ \n",
    "\\boldsymbol{k}_{s} & =f(t+c_{s}h,\\boldsymbol{x}+h(a_{s1}\\boldsymbol{k}_{1}+a_{s2}\\boldsymbol{k}_{2}+\\cdots+a_{s,s-1}\\boldsymbol{k}_{s-1}))\\,.\n",
    "\\end{align*}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Les coefficients $a_{ij}$ (pour $1 \\leq j < i \\leq s$), $c_i$ (pour $2 \\leq i \\leq s$), et $b_i$ (pour $1 \\leq i \\leq s$), sont souvent représentés dans un tableau de Butcher:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{c|ccccc}\n",
    "0\\\\ \n",
    "c_{2} & a_{21}\\\\ \n",
    "c_{3} & a_{31} & a_{32}\\\\ \n",
    "\\vdots & \\vdots &  & \\ddots\\\\ \n",
    "c_{s} & a_{s1} & a_{s2} & \\cdots & a_{s,s-1}\\\\ \n",
    "\\hline  & b_{1} & b_{2} & \\cdots & b_{s-1} & b_{s}\n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Par exemple, le tableau de Butcher de la méthode précédente d'ordre deux est:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{c|cc}\n",
    "0\\\\ \n",
    "\\frac{1}{2} & \\frac{1}{2}\\\\ \n",
    "\\hline  & 0 & 1\n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "et celui de la méthode d'ordre quatre:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{c|cccc}\n",
    "0\\\\ \n",
    "\\frac{1}{2} & \\frac{1}{2}\\\\ \n",
    "\\frac{1}{2} & 0 & \\frac{1}{2}\\\\ \n",
    "1 & 0 & 0 & 1\\\\ \n",
    "\\hline  & \\frac{1}{6} & \\frac{1}{3} & \\frac{1}{3} & \\frac{1}{6}\n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**a)**\n",
    "Définir une fonction `integrate(f, x0, t, M)` qui pour une liste de temps $(t_i)_ {i=0}^N$ donnée retourne les valeurs correspondantes $\\boldsymbol{x}_0,\\boldsymbol{x}_1,\\dots,\\boldsymbol{x}_N$ avec la méthode $M$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**b)**\n",
    "Implémenter les fonctions `M(f,t,x,h)` pour la méthode d'Euler explicite et la méthode de Runge-Kutta d'ordre deux. Comparer les deux méthodes."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**c)**\n",
    "Implémenter la fonction `M(f,t,x,h)` pour la méthode de Runge-Kutta d'ordre quatre. Comparer avec la méthode d'ordre deux."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.3: Mouvement d'une planète"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"exer:equadiff-planete\"></div>\n",
    "\n",
    "Le but est de simuler le mouvement bidimensionnel d'une planète gravitant autour d'une étoile fixe. L'étoile est supposée fixée à l'origine et la position de la planète dans le plan est décrite par le vecteur $\\boldsymbol{x}\\in\\mathbb{R}^2$. L'étoile est supposée interagir avec la planète avec le potentiel:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "V(\\boldsymbol{x}) = \\frac{1}{\\alpha} |\\boldsymbol{x}|^\\alpha \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "pour un certain $\\alpha\\in\\mathbb{R}$ où $|\\boldsymbol{x}|$ désigne la norme euclidienne du vecteur $\\boldsymbol{x}$. À noter que le potentiel gravitationnel correspond à $\\alpha=-1$. L'équation de la planète dans ce champ de force est donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\ddot{\\boldsymbol{x}} = -\\boldsymbol{\\nabla}V(\\boldsymbol{x}) = -\\boldsymbol{x} |\\boldsymbol{x}|^{\\alpha-2} \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**a)**\n",
    "Réécrire l'équation différentielle d'ordre deux comme une équation différentielle d'ordre un pour $\\boldsymbol{x}$ et $\\boldsymbol{p}=\\dot{\\boldsymbol{x}}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**b)**\n",
    "Implémenter la fonction `f(t,xp)` correspondant à l'équation trouvée au point précédent."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**c)**\n",
    "À l'aide de la méthode de Runge-Kutta d'ordre quatre, résoudre l'équation différentielle pour différentes données initiales et différentes valeurs de $\\alpha$ et tracer les trajectoires $\\boldsymbol{x}(t)$ dans le plan. Interpréter les résultats et expliquer en particulier pourquoi les cas $\\alpha=-1$ et $\\alpha=2$ sont différents des autres."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.4: <span style=\"color:red\">!!</span> Équation des ondes cubique"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"exer:equadiff-ondes\"></div>\n",
    "\n",
    "Le but est de résoudre numériquement l'équation des ondes non linéaire sur $\\mathbb{R}$:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{align*}-\\frac{\\partial^{2}u}{\\partial t^{2}}+\\frac{\\partial^{2}u}{\\partial x^{2}} & =u^{3}\\,, & u(0,\\cdot) & =u_{0}\\,, & \\frac{\\partial u}{\\partial t}(0,\\cdot) & =v_{0}\\,,\\end{align*}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "pour $u:\\mathbb{R}^+ \\times \\mathbb{R}$ avec $u_0, v_0:\\mathbb{R}\\to\\mathbb{R}$ deux fonctions données.\n",
    "\n",
    "**Remarque:**\n",
    "Les propriétés de cette équation d'apparence simple sont très mal comprises mathématiquement, voir l'article de recherche suivant pour plus de détails: [doi:10.2140/apde.2012.5.411](https://doi.org/10.2140/apde.2012.5.411).\n",
    "\n",
    "**a)**\n",
    "Réécrire l'équation précédente sous la forme de deux équations du premier ordre en temps pour $u$ et $v=\\frac{\\partial u}{\\partial t}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**b)**\n",
    "En approximant la seconde dérivée en espace par différences finies , montrer que l'équation peut s'approximer de la manière suivante:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{align*}\n",
    "\\frac{\\partial u_{n}}{\\partial t} & =v_{n} \\,, & u_{n}(0) & =u_{0}(x_{n})\\,,\\\\ \n",
    "\\frac{\\partial v_{n}}{\\partial t} & =\\frac{u_{n-1}-2u_{n}+u_{n+1}}{h^{2}}-u_{n}^{3} \\,, & v_{n}(0) & =v_{0}(x_{n})\\,,\n",
    "\\end{align*}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où $(x_n)_ {n=0}^{N}$ représente $N+1$ points équidistants de $h$ dans l'intervalle $[-L,L]$ et $u_n(t) = u(t,x_n)$ et $v_n(t) = v(t,x_n)$. Pour les conditions au bord du domaine, *i.e.* lorsque $n=0$ ou $n=N$, on prendra:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\frac{\\partial v_{0}}{\\partial t} = 0 \\,, \\qquad \\frac{\\partial v_{N}}{\\partial t} = 0 \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**c)**\n",
    "Déterminer la fonction $f:\\mathbb{R}\\times\\mathbb{R}^{2N+2} \\to \\mathbb{R}^{2N+2}$ permettant de mettre l'approximation précédente sous la forme $\\dot{\\boldsymbol{u}} = f(t,\\boldsymbol{u})$ pour $\\boldsymbol{u}=(u,v)$ et implémenter cette fonction."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**d)**\n",
    "Résoudre l'équation différentielle donnée par $\\dot{\\boldsymbol{u}} = f(t,\\boldsymbol{u})$ par exemple avec la méthode de Runge-Kutta d'ordre quatre. Un bon choix de paramètre est $L=100$, $N=1 000$, pour la donnée initiale $u_0(x)=e^{-x^2}$ et $v_0(x)=0$. La vitesse de propagation de l'onde est un et donc, après un temps plus grand que $L$, l'onde sort de la boite $[-L,L]$ et ne correspond plus à une bonne approximation de l'équation initiale."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**e)**\n",
    "À l'aide du module `animation` de Matplotlib réaliser une vidéo montrant l'évolution de l'onde en fonction du temps.\n",
    "\n",
    "**Indication:**\n",
    "Utiliser par exemple la fonction `FFMpegWriter`.\n",
    "\n",
    "**Réponse:**\n",
    "Voici la forme de la solution:\n",
    "\n",
    "<center><img src=\"https://python.guillod.org/fig/nonlinear-wave.png\" style=\"width:90%;max-width:800px;\"></center>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercice 10.5: <span style=\"color:red\">!!!</span> Méthodes de Bogacki-Shampine"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<div id=\"exer:equadiff-bogacki-shampine\"></div>\n",
    "\n",
    "En combinant deux méthodes de Runge-Kutta d'ordres différents (par exemple (2,3) ou (4,5)), on obtiendra une estimation empirique de l'erreur sur un pas de temps. En utilisant cette estimation d'erreur, il est possible d'adapter le pas de temps, soit en l'augmentant soit en le diminuant et ainsi de s'adapter à l'équation.\n",
    "\n",
    "Pour une méthode de Runge-Kutta d'ordre $s$, une méthode interne d'ordre moins élevé (généralement $s-1$) est donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "M^*(t,\\boldsymbol{x},h) = h \\sum_{i=1}^s b_i^* \\boldsymbol{k}_i \\,,\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "où les $\\boldsymbol{k}_i$ sont identiques à ceux de la méthode d'ordre $s$.\n",
    "Une estimation de l'erreur est alors donnée par:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "E(t,\\boldsymbol{x},h) = M(t,\\boldsymbol{x},h) - M^*(t,\\boldsymbol{x},h) =h \\sum_{i=1}^s (b_i - b_i^*) \\boldsymbol{k}_i \\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Une telle méthode est donnée par un tableau de Butcher étendu:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{array}{c|ccccc}\n",
    "0\\\\ \n",
    "c_{2} & a_{21}\\\\ \n",
    "c_{3} & a_{31} & a_{32}\\\\ \n",
    "\\vdots & \\vdots &  & \\ddots\\\\ \n",
    "c_{s} & a_{s1} & a_{s2} & \\cdots & a_{s,s-1}\\\\ \n",
    "\\hline  & b_{1} & b_{2} & \\cdots & b_{s-1} & b_{s}\\\\ \n",
    "& b_{1}^* & b_{2}^* & \\cdots & b_{s-1}^* & b_{s}^*\n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**a)**\n",
    "Implémenter la méthode de Bogacki-Shampine d'ordre (4,5). L'article original est disponible à l'adresse [doi:10.1016/0898-1221(96)00141-1](https://doi.org/10.1016/0898-1221(96)00141-1).\n",
    "\n",
    "**Indication:**\n",
    "Les coefficients des tables de Butcher sont implémentés dans un module `nodepy` donc la documentation est disponible [ici](https://nodepy.readthedocs.io/). Le nom de la méthode dans ce package est «BS5»."
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}