diff --git a/README.md b/README.md index bb8289d692c2f5bc0f9d3f77e75667cf8715ff3d..f2114897bf36b5d55265c3c1a46d684f6830cd45 100644 --- a/README.md +++ b/README.md @@ -1,4 +1,4 @@ -# Calcul différentiel et EDO +# Calcul différentiel et équations différentielles ordinaires Cours de calcul différentiel et équations différentielles ordinaires pour la formation ModIA. diff --git a/tp/runge-kutta.ipynb b/tp/runge-kutta.ipynb index c44e03a4138887f1316f4111c32fe5821f21fca1..19b5220ef40b2cf25e2e2801cd48348b899c697e 100644 --- a/tp/runge-kutta.ipynb +++ b/tp/runge-kutta.ipynb @@ -9,12 +9,49 @@ "\n", "# Méthodes de Runge-Kutta\n", "\n", - "- Date : 2023-2024\n", + "- Date : 2024-2025\n", "- Durée approximative : 1h15\n", "\n", - "## Introduction\n", + "<div style=\"width:95%;\n", + " margin:10px;\n", + " padding:8px;\n", + " background-color:#afa;\n", + " border:2px solid #bbffbb;\n", + " border-radius:10px;\n", + " font-weight:bold;\n", + " font-size:1.5em;\n", + " text-align:center;\">\n", + "Introduction\n", + "</div>\n", + "\n", + "**Méthode à un pas explicite.** On désire calculer une approximation de la solution, notée $x(\\cdot, t_0, x_0)$, sur l'intervalle \n", + "$I := [t_0, t_f]$ du problème de Cauchy\n", + "\\begin{equation*}\n", + " \\dot{x}(t) = f(t, x(t)), \\quad x(t_0) = x_0,\n", + "\\end{equation*}\n", + "où $f \\colon \\mathbb{R} \\times \\mathbb{R}^n \\to \\mathbb{R}^n$ est une application suffisamment régulière. \n", + "On considère pour cela une subdivision de l'intervalle $I$ : \n", + "$$\n", + " t_0 < t_1 < \\cdots < t_N := t_f\n", + "$$\n", + "\n", + "On note $h_i := t_{i+1} - t_i$ pour $i \\in \\llbracket 0, N-1 \\rrbracket$, les **pas** de la subdivision et \n", + "$h_{\\max} := \\max_i(h_i)$ le pas le plus long. L'idée est de calculer une approximation de \n", + "la solution en les points de discrétisation. On souhaite donc approcher $x(t_i, t_0, x_0)$ pour \n", + "$i \\in \\llbracket 1, N \\rrbracket$. \n", + "\n", + "On appelle *méthode à un pas explicite*, toute méthode calculant un N-uplet $(x_1, \\ldots, x_N) \\in {(\\mathbb{R}^n)}^N$ tel que la valeur $x_{i+1}$ est calculée en fonction de $t_i$, $h_i$ et $x_i$ de la forme :\n", + "\n", + "$$\n", + " x_{i+1} = x_i + h_i\\, \\Phi(t_i, x_i, h_i).\n", + "$$\n", + "\n", + "Pour approcher la solution, nous utiliserons une méthode à un pas. Si on note $x_h^* := (x_1^*, \\ldots, x_N^*)$ une solution au problème discrétisé (obtenue par la méthode à un pas considérée), alors on souhaite que \n", + "$$\n", + " x_i^* \\approx x(t_i, t_0, x_0) \\quad \\text{pour tout } i \\in \\llbracket 1, N \\rrbracket.\n", + "$$\n", "\n", - "Nous allons dans ce TP, implémenter quelques méthodes de Runge-Kutta et étudier leur convergence. On considère un pas de temps $h$ uniforme. On rappelle qu'ue méthode à un pas explicite est **convergente** si pour toute solution $x(\\cdot, t_0, x_0)$, la suite ${(x_i)}_i$ définie par $x_{i+1} = x_i + h\\, \\Phi(t_i, x_i, h)$ vérifie \n", + "**Ordre de convergence.** On considère à partir de maintenant un pas de temps $h$ **uniforme**. Une méthode à un pas explicite est **convergente** si pour toute solution $x(\\cdot, t_0, x_0)$, la suite ${(x_i)}_i$ définie par $x_{i+1} = x_i + h\\, \\Phi(t_i, x_i, h)$ vérifie \n", "$$\n", " \\max_{1 \\le i \\le N}\\, \\|{x(t_i, t_0, x_0) - x_i}\\| \\to 0 \n", " \\quad\\text{quand}\\quad h \\to 0.\n", @@ -32,7 +69,22 @@ " \\log(E) = \\log(M) + p\\, \\log(h).\n", "$$\n", "\n", - "Nous en déduisons que si on trace $\\log(E)$ en fonction de $\\log(h)$, on doit obtenir une droite de pente $p$. C'est ce que nous allons vérifier dans ce TP." + "Nous en déduisons que si on trace $\\log(E)$ en fonction de $\\log(h)$, on doit obtenir une droite de pente $p$. C'est ce que nous allons vérifier dans ce TP.\n", + "\n", + "**Objectif.** Nous allons dans ce TP, implémenter quelques méthodes (à un pas) de [Runge-Kutta](https://fr.wikipedia.org/wiki/Méthodes_de_Runge-Kutta) explicites et \n", + "étudier leur convergence. On appelle *méthode de Runge-Kutta explicite à $s$ étages*, la méthode définie par le schéma\n", + "\n", + "$$\n", + " \\begin{aligned}\n", + " k_1 &= f(t_0, x_0) \\\\\n", + " k_2 &= f(t_0 + c_2 h, x_0 + h a_{21} k_1) \\\\\n", + " \\vdots & \\\\[-1em]\n", + " k_s &= f(t_0 + c_s h, x_0 + h \\sum_{j=1}^{s-1} a_{sj} k_j) \\\\\n", + " x_1 &= x_0 + h \\sum_{i=1}^{s} b_i k_i,\n", + " \\end{aligned}\n", + "$$\n", + "où les coefficients $c_i$, $a_{ij}$ et $b_i$ sont des constantes qui définissent précisément le schéma. \n", + "On supposera toujours dans la suite que $c_1 = 0$ et $c_i = \\sum_{j=1}^{i-1} a_{ij}$, pour $i = 2, \\ldots, s$." ] }, { @@ -96,7 +148,7 @@ " # Ecriture des choix de paramètres\n", " println(\"Méthode d'intégration : \", method_name)\n", "\n", - " plts = [] # Liste des graphiques\n", + " plts = [] # Liste des graphiques\n", "\n", " for N ∈ Nspan\n", "\n", @@ -131,7 +183,7 @@ "function ordre(method, f, x0, tspan, sol, hspan, nfespan; \n", " xlims=xlims_, ylims=ylims_, xlims_nfe=xlims_nfe_, ylims_nfe=ylims_nfe_)\n", "\n", - " plts = [] # Liste des graphiques\n", + " plts = [] # Liste des graphiques\n", "\n", " #\n", " plt1 = plot(xaxis=:log, yaxis=:log); push!(plts, plt1)\n", @@ -224,31 +276,53 @@ "end;\n", "\n", "hspan_ = 10 .^ range(-4, stop=0, length=20)\n", - "Nfespan_ = 10 .^ range(0, stop=4, length=20);" + "Nfespan_ = 10 .^ range( 0, stop=4, length=20);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## La méthode d'Euler\n", - "\n", - "On rappelle que la méthode d'Euler est donnée par :\n", + "<div style=\"width:95%;\n", + " margin:10px;\n", + " padding:8px;\n", + " background-color:#afa;\n", + " border:2px solid #bbffbb;\n", + " border-radius:10px;\n", + " font-weight:bold;\n", + " font-size:1.5em;\n", + " text-align:center;\">\n", + "L'exemple d'étude\n", + "</div>\n", + "\n", + "On s'intéresse dans ce TP au problème de Cauchy\n", "\n", "$$\n", - "\\left\\{\n", - "\\begin{array}{l}\n", - "x_{n+1} = x_n + h f(t_n, x_n) \\\\\n", - "x_0 = x(t_0)\n", - "\\end{array}\n", - "\\right.\n", + " \\dot{x}(t) = (1-2t) x(t), \\quad x(0) = x_0 = 1\n", "$$\n", "\n", + "sur l'intervalle $[0, 3]$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Définition du problème de Cauchy\n", + "f(t, x) = (1-2t)x # Second membre f(t, x)\n", + "x0 = 1.0 # Condition initiale\n", + "tspan = (0.0, 3.0); # Intervalle de temps" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ "### Exercice 1\n", "\n", - "1. Calculer la solution analytique et la coder dans la fonction `sol` ci-dessous.\n", - "2. Implémenter la méthode d'Euler dans la fonction `euler` ci-dessous.\n", - "3. Vérifier la convergence de la méthode d'Euler sur l'exemple donné dans la cellule d'après." + "* Calculer la solution analytique de l'exemple ci-dessus et la coder dans la fonction `sol` ci-dessous." ] }, { @@ -263,6 +337,39 @@ "end;" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "<div style=\"width:95%;\n", + " margin:10px;\n", + " padding:8px;\n", + " background-color:#afa;\n", + " border:2px solid #bbffbb;\n", + " border-radius:10px;\n", + " font-weight:bold;\n", + " font-size:1.5em;\n", + " text-align:center;\">\n", + "La méthode d'Euler\n", + "</div>\n", + "\n", + "On rappelle que la méthode d'Euler est donnée par :\n", + "\n", + "$$\n", + "\\left\\{\n", + "\\begin{array}{l}\n", + "x_{n+1} = x_n + h f(t_n, x_n) \\\\\n", + "x_0 = x(t_0)\n", + "\\end{array}\n", + "\\right.\n", + "$$\n", + "\n", + "### Exercice 1\n", + "\n", + "1. Implémenter la méthode d'Euler dans la fonction `euler` ci-dessous.\n", + "2. Vérifier la convergence de la méthode d'Euler sur l'exemple donné dans la cellule d'après." + ] + }, { "cell_type": "code", "execution_count": null, @@ -287,18 +394,6 @@ "end;" ] }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Définition du problème de Cauchy\n", - "f(t, x) = (1-2t)x # Second membre f(t, x)\n", - "x0 = 1.0 # Condition initiale\n", - "tspan = (0.0, 3.0); # Intervalle de temps" - ] - }, { "cell_type": "code", "execution_count": null, @@ -313,16 +408,14 @@ "plts = convergence(method, f, x0, tspan, sol, Nspan)\n", "\n", "# Affichage des graphiques\n", - "plot(plts..., layout=(1, length(Nspan)), size=(900, 500))" + "plot(plts..., layout=(1, length(Nspan)), size=(800, 400))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "### Exercice 2 \n", - "\n", - "1. Afficher l'erreur globale de la méthode d'Euler en fonction de $h$ et vérifier que l'erreur est bien en $O(h)$." + "3. Afficher l'erreur globale de la méthode d'Euler en fonction de $h$ et vérifier que l'erreur est bien en $O(h)$." ] }, { @@ -340,14 +433,24 @@ "plt_order_euler = ordre(method, f, x0, tspan, sol, hspan, Nfespan)\n", "\n", "# Affichage du graphique\n", - "plot(plt_order_euler..., layout=(1, 2), size=(900, 500))" + "plot(plt_order_euler..., layout=(1, 2), size=(800, 400))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## La méthode de Runge \n", + "<div style=\"width:95%;\n", + " margin:10px;\n", + " padding:8px;\n", + " background-color:#afa;\n", + " border:2px solid #bbffbb;\n", + " border-radius:10px;\n", + " font-weight:bold;\n", + " font-size:1.5em;\n", + " text-align:center;\">\n", + "La méthode de Runge \n", + "</div>\n", "\n", "On rappelle que la méthode de Runge est donnée par :\n", "\n", @@ -406,7 +509,7 @@ "plts = convergence(method, f, x0, tspan, sol, Nspan)\n", "\n", "# Affichage des graphiques\n", - "plot(plts..., layout=(1, length(Nspan)), size=(900, 500))" + "plot(plts..., layout=(1, length(Nspan)), size=(800, 400))" ] }, { @@ -424,14 +527,24 @@ "plt_order_runge = ordre(plt_order_euler, method, f, x0, tspan, sol, hspan, Nfespan)\n", "\n", "# Affichage du graphique\n", - "plot(plt_order_runge..., layout=(1, 2), size=(900, 500))" + "plot(plt_order_runge..., layout=(1, 2), size=(800, 400))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## La méthode de Heun \n", + "<div style=\"width:95%;\n", + " margin:10px;\n", + " padding:8px;\n", + " background-color:#afa;\n", + " border:2px solid #bbffbb;\n", + " border-radius:10px;\n", + " font-weight:bold;\n", + " font-size:1.5em;\n", + " text-align:center;\">\n", + "La méthode de Heun \n", + "</div>\n", "\n", "On rappelle que la méthode de Heun est donnée par :\n", "\n", @@ -494,7 +607,7 @@ "plts = convergence(method, f, x0, tspan, sol, Nspan)\n", "\n", "# Affichage des graphiques\n", - "plot(plts..., layout=(1, length(Nspan)), size=(900, 500))" + "plot(plts..., layout=(1, length(Nspan)), size=(800, 400))" ] }, { @@ -512,14 +625,24 @@ "plt_order_heun = ordre(plt_order_runge, method, f, x0, tspan, sol, hspan, Nfespan)\n", "\n", "# Affichage du graphique\n", - "plot(plt_order_heun..., layout=(1, 2), size=(900, 500))" + "plot(plt_order_heun..., layout=(1, 2), size=(800, 400))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## La méthode de Runge-Kutta d'ordre 4\n", + "<div style=\"width:95%;\n", + " margin:10px;\n", + " padding:8px;\n", + " background-color:#afa;\n", + " border:2px solid #bbffbb;\n", + " border-radius:10px;\n", + " font-weight:bold;\n", + " font-size:1.5em;\n", + " text-align:center;\">\n", + "La méthode de Runge-Kutta d'ordre 4\n", + "</div>\n", "\n", "On rappelle que la méthode de Runge-Kutta d'ordre 4 est donnée par :\n", "\n", @@ -584,7 +707,7 @@ "plts = convergence(method, f, x0, tspan, sol, Nspan)\n", "\n", "# Affichage des graphiques\n", - "plot(plts..., layout=(1, length(Nspan)), size=(900, 500))" + "plot(plts..., layout=(1, length(Nspan)), size=(800, 400))" ] }, { @@ -602,7 +725,7 @@ "plt_order_rk4 = ordre(plt_order_heun, method, f, x0, tspan, sol, hspan, Nfespan)\n", "\n", "# Affichage du graphique\n", - "plot(plt_order_rk4..., layout=(1, 2), size=(900, 500))" + "plot(plt_order_rk4..., layout=(1, 2), size=(800, 400))" ] } ],