diff --git a/serie07.ipynb b/serie07.ipynb deleted file mode 100644 index 942c9a7..0000000 --- a/serie07.ipynb +++ /dev/null @@ -1,745 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Exercice: Formule du rectangle et du trapèze" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "import matplotlib.pyplot as plt\n", - "import numpy as np" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "On considère le calcul de l'intégrale\n", - "$$\n", - "I=\\int_{0}^1 f(x) \\, dx,\n", - "$$\n", - "où $f(x)$ est une fonction continue sur $[0,1]$.\n", - "\n", - "1. Ecrivez une fonction `midpoint` qui implémente la formule composite du rectangle (point\n", - "milieu) pour l'approximation de l'intégrale ci-dessus. Pour permettre le choix d'un intervalle d'intégration de la forme $[a,b]$, le nombre de sous-intervalles $N$\n", - "et la fonction $f(x)$, définie par la commande `f = lambda x : ... `\n", - " utilisez la structure suivante:\n", - "\n", - "```python\n", - "def midpoint( a, b, N, f ) :\n", - " # [a,b] : inteval\n", - " # N : number of subintervals\n", - " # f : fonction to integrate using the midpoint rule\n", - "``` \n", - " \n", - "2. Testez le code en considérant la fonction $f(x)=x^2$ et\n", - "$M=10$. Tracez ensuite le graphe de l'erreur $|I(f)-I_{\\mathrm{\\textit{pm}}}^c (f)|$\n", - "en fonction de $N$ (utilisez à cette fin les commandes `plt.xscale('log'); plt.yscale('log')`), en choisissant $a=0$, $b=1$,\n", - "et $M = 10^1, 10^2, 10^3, \\ldots 10^5$. Quel est l'ordre de\n", - "convergence de la formule composite du rectangle par rapport à $H=(b-a)/M$? Donnez une\n", - "interprétation des résultats d'après la théorie.\n", - "\n", - "3. Modifiez le code du point 1. pour permettre le calcul\n", - "de l'intégrale avec la formule composite du trapèze.\n", - "Tracez le graphe de l'erreur $|I(f)-I_{\\mathrm{\\textit{t}}}^c (f)|$\n", - "pour les mêmes valeurs de $M$.\n", - "Comparez les résultats avec ceux obtenus\n", - "au point 2.\n", - "\n", - "```python\n", - "def trapeziodal( a, b, N, f ) :\n", - " # [a,b] : inteval\n", - " # N : number of subintervals\n", - " # f : fonction to integrate using the trapezoidal rule\n", - "``` \n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 1" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "def midpoint( a, b, N, f ) :\n", - " # [a,b] : inteval\n", - " # M : number of subintervals\n", - " # f : fonction to integrate using the midpoint rule\n", - "\n", - " # (completez)\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 2" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "f = lambda x : x**2\n", - "\n", - "a = 0; b = 1; N = 10;\n", - "\n", - "intmp = midpoint( a, b, N, f )\n", - "\n", - "print(intmp)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "En sachant que la valeur exacte de l'intégrale est $I(f)=1/3$,\n", - "on peut écrire les commandes suivantes:" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "N = 10**np.linspace(1,5,5).astype(int)\n", - "\n", - "errmp = []\n", - "for i in range(5) :\n", - " # (completez)\n", - "\n", - "# (completez)\n", - " \n", - "\n", - " \n", - " \n", - "slope = ( np.log(errmp[4]) - np.log(errmp[0]) ) / ( np.log(N[4]) - np.log(N[0]) ) \n", - "\n", - "print(f'La pente est de {slope:.6f}')\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Donc, l'erreur décroît comme la puissance -2 du paramètre $N$.\n", - "Par conséquent, l'ordre de convergence par rapport à $H = 1/N$ est $2$.\n", - "\n", - "\n", - "#### Partie 3 \n", - "\n", - "On peut modifier le code du point milieu pour\n", - "implémenter la méthode du trapèze comme suit." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "def trapezoidal( a, b, N, f ) :\n", - " # [a,b] : inteval\n", - " # M : number of subintervals\n", - " # f : fonction to integrate using the trapezoidal rule\n", - "\n", - " # (completez)\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "errtrap = []\n", - "\n", - "# (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "On observe que la précision est encore $\\mathcal{O}(h^2)$, mais\n", - "que l'erreur est plus grande (le double, précisément)\n", - "que dans le cas du point milieu.\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Exercice: Formule de Simpson" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "import matplotlib.pyplot as plt\n", - "import numpy as np" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "On considère la fonction $f \\, : \\, [a,b] \\rightarrow \\mathbb R$ dans $C^0([a,b])$;\n", - "on est intéressé à approcher l'intégrale $I(f) = \\int_a^b f(x) \\, dx$. \n", - "\n", - "1. Ecrivez une fonction `Simpson` qui implémente la formule composite de Simpson pour l'approximation de l'intégrale ci-dessus. \n", - "Utilisez la structure suivante:\n", - "\n", - "```python\n", - "def Simpson( a, b, N, f ) :\n", - " # [a,b] : inteval\n", - " # N : number of subintervals\n", - " # f : fonction to integrate using the Simpson rule\n", - "```\n", - "\n", - "2. Testez ensuite pour quels monômes $f(x) = x^d$ cette formule intègre exactement la fonction, pour $d=0,1,...$ sur l'intervalle $[1,4]$, avec $N=1$ et ensuite $N=10$.\n", - " \n", - "3. Vérifiez numériquement pour quelques polynômes que la fonction ainsi écrite est linéaire en $f$ pour $N=10$\n", - "\n", - "*Suggestion*: En utilisant une fonction `lambda` il est possible de décider le paramètre $d$ d'un monôme à un moment ultérieur:\n", - "\n", - "```python\n", - "monomial = lambda x : x**d\n", - "```\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 1" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# This function just provides the Simpson quadrature rule.\n", - "def Simpson( a, b, N, f ) :\n", - " # [a,b] : inteval\n", - " # N : number of subintervals\n", - " # f : fonction to integrate using the trapezoidal rule\n", - "\n", - " # (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 2\n", - "\n", - "$$\\int_1^4 1 dx= 3$$\n", - "$$\\int_1^4 x^d dx = \\frac1{d+1} x^{d+1} \\vert_1^4 = \\frac{4^{d+1} - 1}{d+1}$$" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# Checking simpson fonction\n", - "a = 1; b = 4;\n", - "\n", - "# with lambda fonctions, it is possible to determine a parameter (here d)\n", - "# at a later moment\n", - "monomial = lambda x : x**d\n", - "\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "N = 1\n", - "for d in range(7) :\n", - " # (completez)\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "N = 10\n", - "\n", - "for d in range(7) :\n", - " # (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 3\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# (completez)\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Exercice: Formule du rectangle, du trapèze et de Simpson" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "On considère la fonction $f \\, : \\, [a,b] \\rightarrow \\mathbb R$ dans $C^0([a,b])$;\n", - "on est intéressé à approcher l'intégrale $I(f) = \\int_a^b f(x) \\, dx$. \n", - "\n", - "Plus précisément on prend\n", - "$f(x) = \\sin( \\frac72 \\, x ) + e^x- 1$ avec $a=0$ et $b=1$ ($f \\in C^{\\infty}([a,b])$) et on\n", - "peut calculer $I(f) = \\frac27 \\, (1-\\cos(\\frac72 )) + e - 2$.\n", - "\n", - "1. Calculez une approximation de l'intégrale en utilisant les formules du rectangle, du trapèze et\n", - "de Simpson *simples*, c-à-d avec un seul intervalle. \n", - "\n", - "2. Calculez une approximation de l'intégrale en utilisant les fonctions `Midpoint`, `Trapezoidal` et `Simpson` (déjà codées).\n", - " avec $N=10$ sous-intervalles de même taille.\n", - " On notera les valeurs approchées de l'intégrale\n", - " $I_{mp}^c(f)$, $I_{t}^c(f)$, and $I_{s}^c(f)$, respectivement.\n", - "\n", - "3. Répétez le point 2. avec $N=2^{k}$ pour $k=2,\\ldots,7$ et calculez les erreurs $E^c_{mp}(f) := | I(f) - I_{mp}^c(f)|$, $E^c_{t}(f) := | I(f) - I_{t}^c(f)|$, et $E^c_{s}(f) := | I(f) - I_{s}^c(f)|$.\n", - " Dessinez les erreurs en fonction de $H=(b-a)/N$ sur une échelle logarithmique sur les deux axes.\n", - " Quel est l'ordre de convergence de ces méthodes ? Est-ce en accord avec la théorie ?\n", - " Motivez votre réponse.\n", - "\n", - "4. On prend maintenant $f(x) = x^d$, $a=0$ et $b=1$, avec $d \\in \\mathbb N$. L'intégrale de $f$ vaut $I(f) = 1/(d+1)$.\n", - " Vérifiez numériquement les degrés d'exactitude de chacune des formules de quadrature du point 1.\n", - " Pour cela, il faut choisir plusieurs valeurs de $d = 0,1,2,\\ldots$. \n", - " Motivez votre réponse.\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "import matplotlib.pyplot as plt\n", - "import numpy as np\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 1" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "f = lambda x : np.sin(7/2*x) + np.exp(x) - 1\n", - "\n", - "a = 0; b = 1\n", - "Iexact = 2/7*(1-np.cos(7/2) ) + np.exp(1) - 2\n", - "\n", - "x = np.linspace(a,b,1000)\n", - "y=f(x)\n", - "\n", - "plt.plot(x, y, 'b')\n", - "plt.show()" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "intMD = \n", - "intTrap = \n", - "intSimp = \n", - "\n", - "print(f'exact \\t rect \\t\\t trap \\t\\t Simpson')\n", - "print(f'{Iexact:.4f} \\t {intMD:.4f} \\t {intTrap:.4f} \\t {intSimp:.4f}')\n", - " " - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 2" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "N = 10\n", - "intMD = Midpoint(a,b,N,f)\n", - "intTrap = Trapezoidal(a,b,N,f)\n", - "intSimp = Simpson(a,b,N,f)\n", - "\n", - "print(f'exact \\t rect \\t\\t trap \\t\\t Simpson')\n", - "print(f'{Iexact:.4f} \\t {intMD:.4f} \\t {intTrap:.4f} \\t {intSimp:.4f}')\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 3" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "errmp = []\n", - "errtrap = []\n", - "errSim = []\n", - "\n", - "# (completez)\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 4\n", - "\n", - " The function $f(x)=x^d$ is a polynomial of degree $d$ for $d \\in \\mathbb N$. The simple midpoint, trapezoidal, and Simpson quadrature formulas possesses degree of exactness equal to $1$, $1$, and $3$, respectively." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "N = 1\n", - "# with lambda fonctions, it is possible to determine a parameter (here d)\n", - "# at a later moment\n", - "monomial = lambda x : x**d\n", - "\n", - "# (completez)\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# (completez)\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "De ces simulations nous vérifions que les monômes $f(x)=x^d$ sont intégrées exactement pour les degrées $d=0$ et $d=1$ dans le cas des formules du rectangle et du trapèze, et pour $d=0,1,2,3$ dans le cas de la formule de Simpson." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "# Exercice: Formule du rectangle et du trapèze" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "import matplotlib.pyplot as plt\n", - "import numpy as np\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "On considère, sur l'intervalle $[-1,1]$, la fonction $f$ suivante:\n", - "\n", - "$$f(x) = \\left\\{ \\begin{array}{cc}\n", - "\\frac{}{}e^x & \\mbox{ si } x\\leq 0 \\\\\n", - "1 & \\mbox{ si } x > 0 \\\\\n", - "\\end{array} \\right.$$\n", - "\n", - "On peut définir une telle fonction en utilisant la commande\n", - "\n", - "```python\n", - "f = lambda x : np.exp(x)*(x<=0) + (1)*(x>0)\n", - "```\n", - "\n", - "\n", - "1. Utilisez $1000$ points équirépartis dans l'intervalle $[-1,1]$ pour afficher la fonction $f$ (utilisez la commande `axis` pour recadrer l'image).\n", - "\n", - "2. On s'intéresse à présent à l'intégrale $I=\\int_{-1}^1 f(x) \\ dx$. On peut calculer la valeur de $I$ analytiquement et on trouve $I=2-\\frac{1}{e} \\cong 1.6321$. Calculez des valeurs approchées de $I$ en considérant les formules du point milieu, du trapèze et de Simpson avec $N=1,9,99,999$ où $N$ est le nombre de sous-intervalles de formules composites. Utilisez les fonctions `midpoint`, `trapezoidal` et `simpson` (déjà codées).\n", - "\n", - "3. Reportez les erreurs calculées au point (b) dans un graphe montrant l'erreur en fonction de $H$ avec des échelles logarithmiques.\n", - "\n", - "4. Estimez, à partir des graphes obtenus au point précédent, l'ordre de chacune des méthodes. Comparez-les avec les ordres donnés au cours.\n", - "Y a-t-il des différences? Pourquoi? (Regardez les dérivées de $f$ ).\n", - "\n", - "5. Pourquoi obtient-on de bien meilleurs résultats pour la méthode de Simpson avec un nombre pair d'intervalles qu'avec un nombre impair (essayez avec $N=99$ et ensuite $N=100$ )?\n", - "\n", - "6. Refaites l'exercice avec $N=2,10,100,1000$.\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 1" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "f = lambda x : np.exp(x)*(x<=0) + (1)*(x>0)\n", - "\n", - "a = -1; b = 1\n", - "x = np.linspace(a,b,1000)\n", - "y=f(x)\n", - "\n", - "plt.plot(x, y, 'b')\n", - "plt.show()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 2\n", - "\n", - "On calcule tout d'abord la valeur exacte de l'intégrale, puis on calcule les erreurs, que l'on stocke dans des vecteurs:" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "Iexact = 2 - np.exp(-1)\n", - "\n", - "# (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 3" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Les graphes pour toutes les méthodes sont des droites. On remarque que lorsque $H$ est divisé par $10$, les \n", - "erreurs sont environ divisées par $100$. On peut donc supposer que les ordres sont $2$ par rapport à $H$. Pour le confirmer,\n", - "on peut ajouter sur le graphe la pente représentant l'ordre $2$ et vérifier qu'elle est parallèle aux\n", - "droites d'erreur." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 4" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# (completez)\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Les graphiques sont parallèles à la droite $H^2$ et donc l'ordre est $2$ pour toutes les méthodes. Pour la\n", - "méthode du point milieu et du trapèze, c'est l'ordre auquel on s'attend d'après la théorie. Par contre, pour Simpson, on\n", - "pouvait s'attendre à un ordre $4$. Le problème vient du fait que la fonction n'est pas très régulière, puisque \n", - "$f \\in C^0([-1,1])$ mais $f \\not\\in C^1([-1,1])$ puisque sa dérivée n'est pas continue en $x=0$. Dans la théorie, on demande\n", - "$f \\in C^4([-1,1])$ pour assurer l'ordre $4$. On obtient tout de même l'ordre $2$ car, mis à part en $x=0$, la fonction $f$\n", - "est très régulière." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 5\n", - "\n", - "Si on regarde l'erreur pour la méthode de Simpson avec $N=99$ et $N=100$ sous-intervalles" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "print (np.abs( Simpson( a, b, 99, f) - Iexact ) )\n", - "print (np.abs( Simpson( a, b, 100, f) - Iexact ) )\n" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "On obtient des valeurs d'environ $1.70 \\ 10^{-5}$ et $3.51 \\ 10^{-11}$. Les erreurs sont très différentes !\n", - "Voilà pourquoi:\n", - "\n", - "* si $N$ est pair, ...\n", - "\n", - "* si $N$ est impair, ..." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "#### Partie 6\n", - "\n", - "Il faut relancer les parties 2 et 3 avec `N= [2,10,100,1000]` à la place de `N = [1,9,99,999]`\n" - ] - }, - { - "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.6" - } - }, - "nbformat": 4, - "nbformat_minor": 4 -}