diff --git a/tp/rk-implicites.ipynb b/tp/rk-implicites.ipynb new file mode 100644 index 0000000000000000000000000000000000000000..83507d31a869ec600e3a875234215cb4a79daa48 --- /dev/null +++ b/tp/rk-implicites.ipynb @@ -0,0 +1,311 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "<img src=\"https://gitlab.irit.fr/toc/ens-n7/texCoursN7/-/raw/main/LOGO_INP_N7.png\" alt=\"N7\" height=\"80\"/>\n", + "\n", + "<img src=\"https://gitlab.irit.fr/toc/ens-n7/texCoursN7/-/raw/main/logo-insa.png\" alt=\"INSA\" height=\"80\"/>\n", + "\n", + "# Méthodes de Runge-Kutta implicites - TP Projet\n", + "\n", + "- Date : 2023-2024\n", + "- Durée approximative : inconnue" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**Nom** : \"Entrez votre nom ici\"\n", + "\n", + "**Prénom** : \"Entrez votre prénom ici\"" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Rendu et consignes\n", + "\n", + "Une fois le travail terminé, vous enverrez directement le fichier `.ipynb` par email à l'adresse : `olivier.cots@toulouse-inp.fr`.\n", + "\n", + "- **Date limite du rendu : mercredi 15/11/2023 à 19h.** Attention, à partir de 19h, 1 point est enlevé de la note finale toutes les 30 minutes.\n", + "- **Attention :** Le fichier doit être renommé de la façon suivante : `rk_implicites_NOM_Prenom.ipynb`\n", + "- **Documents autorisés :** vous pouvez utiliser votre polycopié et les TPs précédents." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Introduction\n", + "\n", + "Nous allons dans ce TP, implémenter quelques méthodes de Runge-Kutta **implicites** (voir **polycopié Section 8.2**) et étudier leur convergence. On considère un pas de temps $h$ uniforme. Une méthode à un pas implicite est convergente si pour toute solution $x(\\cdot, t_0, x_0)$ d'un problème de Cauchy, la suite approximante ${(x_i)}_i$ donnée par la méthode à un pas implicite 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", + "$$\n", + "\n", + "Si la convergence est d'ordre $p$ alors il existe une constante $C \\ge 0$ telle que l'**erreur globale** $E$ vérifie\n", + "\n", + "$$\n", + " E := \\max_{1 \\le i \\le N}\\, \\|{x(t_i, t_0, x_0) - x_i}\\| \\le C\\, h^p.\n", + "$$\n", + "\n", + "Faisons l'hypothèse que $E = M\\, h^p$ pour un certain $M \\ge 0$. En passant au logarithme, on obtient\n", + "\n", + "$$\n", + " \\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." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# activate local project\n", + "using Pkg\n", + "Pkg.activate(\".\")\n", + "\n", + "# load packages\n", + "using LinearAlgebra\n", + "using Plots\n", + "using Plots.PlotMeasures\n", + "using Polynomials\n", + "\n", + "#\n", + "px = PlotMeasures.px;" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Fonctons auxiliaires \n", + "\n", + "# VOUS POUVEZ METTRE VOS FONCTIONS AUXILIAIRES ICI" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## L'exemple d'étude\n", + "\n", + "On s'intéresse (pour les exercices 1, 2 et 3) au problème de Cauchy\n", + "\n", + "$$\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": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Solution analytique\n", + "function sol(t)\n", + " return exp(t-t^2)\n", + "end;" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Estimation de la constante de Lipschitz de f sur [0, 3]\n", + "# Voir Théorème 8.2.2 pour l'utilité de cette estimation\n", + "function dfx(t)\n", + " return 1-2t\n", + "end \n", + "L = maximum(abs.(dfx.(range(0, stop=3, length=1000))))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## La méthode d'Euler implicite\n", + "\n", + "La méthode d'Euler implicite est donnée par :\n", + "\n", + "$$\n", + "\\left\\{\n", + "\\begin{array}{l}\n", + "x_{n+1} = x_n + h f(t_n, x_{n+1}) \\\\\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 implicite avec le point fixe (penser à voir le polycopié Section 8.2).\n", + "2. Pourquoi si $h \\ge 0.2$, la méthode d'Euler implicite ne marche pas ?\n", + "3. Tracer la solution approchée et la solution exacte sur le même graphique pour différentes valeurs de $h$ que vous choisirez pour illustrer la convergence de la méthode.\n", + "4. Tracer l'erreur globale de la méthode d'Euler implicite en fonction de $h$ et vérifier que l'erreur est bien en $O(h)$.\n", + "\n", + "**Attention** : pour l'algorithme du point fixe, faites attention aux critères d'arrêts (il y en a 2) ! Voir votre polycopié Section 8.2. Vous fixerez la valeur de la tolérance à $10^{-6}$ et le nombre maximum d'itérations à $1000$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# AJOUTER VOTRE CODE ICI" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## La méthode des trapèzes\n", + "\n", + "La méthode des trapèzes est donnée par le tableau de Butcher :\n", + "\n", + "$$\n", + " \\begin{array}{c | c c}\n", + " 0 & 0 & 0 \\\\[0.2em]\n", + " 1 & 1/2 & 1/2 \\\\[0.2em]\n", + " \\hline\n", + " & 1/2 & 1/2 \\\\\n", + " \\end{array}\n", + "$$\n", + "\n", + "### Exercice 2\n", + "\n", + "1. Implémenter la méthode des trapèzes avec le point fixe.\n", + "2. Tracer la solution approchée et la solution exacte sur le même graphique pour différentes valeurs de $h$ que vous choisirez pour illustrer la convergence de la méthode.\n", + "3. Tracer l'erreur globale de la méthode des trapèzes. Quel est l'ordre de convergence de la méthode des trapèzes ?\n", + "4. Tracer l'erreur globale en fonction du nombre d'évaluation de la fonction $f$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# AJOUTER VOTRE CODE ICI" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## La méthode de Gauss à 2 étages\n", + "\n", + "La méthode de Gauss à 2 étages est donnée par le tableau de Butcher :\n", + "\n", + "$$\n", + " \\begin{array}{c | c c}\n", + " 1/2 - \\sqrt{3}/6 & 1/4 & 1/4 - \\sqrt{3}/6 \\\\[0.2em]\n", + " 1/2 + \\sqrt{3}/6 & 1/4 + \\sqrt{3}/6 & 1/4 \\\\[0.2em]\n", + " \\hline\n", + " & 1/2 & 1/2 \\\\\n", + " \\end{array}\n", + "$$\n", + "\n", + "### Exercice 3\n", + "\n", + "1. Implémenter la méthode de Gauss à 2 étages avec le point fixe.\n", + "2. Tracer la solution approchée et la solution exacte sur le même graphique pour différentes valeurs de $h$ que vous choisirez pour illustrer la convergence de la méthode.\n", + "3. Tracer l'erreur globale de la méthode de Gauss à 2 étages. Quel est l'ordre de convergence de la méthode de Gauss à 2 étages ?\n", + "4. Tracer l'erreur globale en fonction du nombre d'évaluation de la fonction $f$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# AJOUTER VOTRE CODE ICI" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Un autre exemple\n", + "\n", + "On considère à partir de maintenant l'équation différentielle en dimension 2 :\n", + "\n", + "$$\n", + " \\dot{x}_1(t) = x_2(t), \\quad \\dot{x}_2(t) = - x_1(t).\n", + "$$\n", + "\n", + "On peut montrer facilement que la norme de $x(t) = (x_1(t), x_2(t))$ est constante le long des solutions :\n", + "\n", + "$$\n", + " \\frac{\\mathrm{d}}{\\mathrm{d} t} \\|x(t)\\|^2 = 2\\, \\left( x(t) \\,|\\, \\dot{x}(t) \\right) = 2 \\left( x_1(t) x_2(t) - x_2(t) x_1(t) \\right) = 0.\n", + "$$\n", + "\n", + "### Exercice 4\n", + "\n", + "On considère le problème de Cauchy associé de condition initiale $x_0 = (1, 0)$.\n", + "\n", + "1. Afficher l'approximation de la solution sur $[0, 10]$ pour les méthodes :\n", + "- Euler explicite ;\n", + "- Euler implicite ;\n", + "- Trapèzes ;\n", + "- Gauss à 2 étages.\n", + "2. Commentaires.\n", + "\n", + "**Attention :** vous ferez un affichage dans le plan $(x_1, x_2)$. Vous fixerez le nombre de pas à $N=100$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# AJOUTER VOTRE CODE ICI" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.9.0", + "language": "julia", + "name": "julia-1.9" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.9.0" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +}