diff --git a/1.3 Interpolation de fonction.ipynb b/1.3 Interpolation de fonction.ipynb index 077f2de..cd8a001 100644 --- a/1.3 Interpolation de fonction.ipynb +++ b/1.3 Interpolation de fonction.ipynb @@ -1,401 +1,401 @@ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## 1.3 Interpolation d'une fonction continue\n", "\n", "Soit $f: \\mathbf [a,b] \\rightarrow R$ continue et $t_0,\\ldots,t_n\\in [a,b]$ des noeuds distincts. Le polynôme d'interpolation\n", "$\\Pi_n(t)$ est noté $\\Pi_n f(t)$ et est appelé l'interpolant de $f$ aux\n", "noeuds $t_0,\\dots,t_n$.\n", "\n", "\n", "Si on prend $$y_k=f(t_k), k= 0,..., n,$$ \n", "alors on aura\n", "\\begin{equation*}\n", "\\Pi_nf(t)=\\sum_{k=0}^n\n", "f(t_k)\\varphi_k(t).\n", "\\end{equation*}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Exercice\n", "\n", "Ecrivez une fonction Python qui a la définition suivante, en utilisant la fonction `phi` définie plus haut.\n", "Ecrivez aussi un petit test sur la base de l'exercice précédent.\n", "\n", "```python\n", "# Lagrange Interpolation of data (t,y), evaluated at ordinate(s) z\n", "def LagrangePi(t,y,z):\n", " # the input variables are:\n", " # t : the interpolatory points\n", " # y : the corresponding data at the points t\n", " # z : where to evaluate the function\n", " \n", "```\n", "\n", "Utilisez le fait que $\\{\\varphi_k, k = 0,...,n\\}$ est une base des polynômes de degré $\\leq n$\n", "et que le vecteur $y$ représente les coordonnéee du polynôme d'interpolation recherché par rapport à cette base,\n", "c'est-à-dire\n", "$$\\Pi_n (z) = y_0 \\varphi_0 + ... + y_n \\varphi_n$$\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# importing libraries used in this book\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "from InterpolationLib import LagrangeBasis as phi" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Lagrange Interpolation of data (t,y), evaluated at ordinate(s) z\n", "def LagrangePi(t,y,z):\n", " # the input variables are:\n", " # t : the interpolatory points\n", " # y : the corresponding data at the points t\n", " # z : where to evaluate the function\n", " \n", " # (Complete the code below)\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", " return result " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# vecteur des points d'interpolation\n", "t = np.linspace(0, 1, 5)\n", "\n", "# vecteur des valeurs\n", "y = np.array([3.38, 3.86, 3.85, 3.59, 3.49])\n", "\n", "z = np.linspace(-0.1, 1.1, 100)\n", "plt.plot(t, y, 'ro', z, LagrangePi(t,y,z))\n", "\n", "plt.xlabel('t'); plt.ylabel('y');\n", "plt.legend(['data','p(t)'])\n", "plt.show()\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Erreur d'interpolation\n", "\n", "Soient $t_0$, $t_1$, $\\ldots$, $t_n$,\n", "$(n+1)$ nœuds *équirépartis* dans $I=[a,b]$ et soit $f\\in C^{n+1}(I)$.\n", "Alors\n", "\\begin{equation}\n", " \\max_{x\\in I} |f(t) - \\Pi_n f(t)| \\leq \\dfrac{1}{2(n+1)}\n", "\\left(\\frac{b-a}{n}\\right)^{n+1}\\max_{x\\in I}|f^{(n+1)}(t)|.\n", "\\end{equation}\n", "\n", "On remarque que l'erreur d'interpolation dépend de la dérivée\n", "$(n+1)$-ième de $f$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Exercice\n", "\n", "On considère les points d'interpolation $$t_0=1, t_1=1.75, t_2=2.5, t_3=3.25, t_4=4$$ \n", "et la fonction $$f(t) = t \\sin(2\\pi t)$$.\n", "\n", "1. Calculez la base de Lagrange associée à ces points. D'abord sur papier, ensuite utilisez Python pour en dessiner le graphique.\n", "\n", "2. Calculez le polynôme d'interpolation $\\Pi_n$ à l'aide de la base de Lagrange. D'abord sur papier, ensuite avec Python.\n", "\n", "4. Quelle est l'erreur théorique d'interpolation ? \n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Comportement pour $n$ grand: eg, la fonction de Runge.\n", "\n", "Le fait que\n", - "$$\\lim_{n\\rightarrow\\infty} \\dfrac{1}{4(n+1)}\\left(\\dfrac{b-a}{n}\\right)^{n+1}=0$$\n", + "$$\\lim_{n\\rightarrow\\infty} \\dfrac{1}{2(n+1)}\\left(\\dfrac{b-a}{n}\\right)^{n+1}=0$$\n", "n'implique pas forcément que $\\max_{t \\in I} |E_nf(t)|$ tende vers zéro quand $n\\rightarrow\\infty$.\n", "\n", "Soit $f(t)=\\dfrac{1}{1+t^2}$, $t\\in [-5,5]$. Si on l'interpole dans des\n", "noeuds équirépartis, l'interpolant présente des oscillations au voisinage des extrémités de l'intervalle.\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Runge fonction\n", "f = lambda t : 1./(1+t**2)\n", "\n", "# Values of N to use\n", "Nrange = [3,5,10]\n", "# plotting points\n", "z = np.linspace(-5, 5, 100)\n", "\n", "# (Complete the code below)\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`sympy` est une librairie pour le calcul symbolique en Python. Nous allons l'utiliser pour étudier le comportement de l'erreur d'interpolation de la fonction de Runge." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Using symbolic python to compute derivatives\n", "import sympy as sp\n", "\n", "# define x as symbol\n", "x = sp.symbols('x')\n", "\n", "# define Runge function\n", "f = 1/(1+x**2)\n", "\n", "# pretty print the 5th derivative of f\n", "f5 = sp.diff(f, x,5)\n", "# display f5\n", "sp.init_printing(use_unicode=True)\n", "display(f5)\n", "\n", "\n", "# evalf can be used to compute the value of a function at a given point\n", "print('5th derivative evaluated at 3 :')\n", "print( f5.evalf(subs={x: 3}) )\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Pour définir une fonction qui accepte un array de valeur, il faut utiliser les lignes suivantes.\n", "\n", "Ensuite on peut aussi dessiner le graphe ..." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "\n", "# to evaluate a function at many given points, we need the following trick\n", "diff_f_func = lambda t: float(sp.diff(f,x,k).evalf(subs={x: t}))\n", "diff_f = np.vectorize(diff_f_func)\n", "\n", "# the derivative can be set with k (not very elegant...)\n", "k = 4\n", "print(diff_f(4.5))\n", "\n", "# plotting points\n", "z = np.linspace(-5, 5, 100)\n", "# plot the derivative between -5 and 5\n", "\n", "plt.plot(z,diff_f(z), 'b')\n", "plt.xlabel('t'); plt.ylabel('y'); plt.title('Derivatives of Runge function')\n", "plt.legend(Nrange)\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "... ou évaluer le maximum pour plusieurs $n$ et voir le comportement de $\\max |f^{(n)}|$ en fonction de $n$." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Plot max(abs(fn)) in the range -5,5\n", "z = np.linspace(-5, 5, 100)\n", "\n", "Nmax = 10\n", "maxValFn = np.zeros(Nmax)\n", "for k in range(Nmax):\n", " maxValFn[k] = np.max(np.abs(diff_f(z)))\n", " \n", "plt.plot(range(10), maxValFn)\n", "\n", "plt.yscale('log')\n", "plt.xlabel('n'); plt.ylabel('$\\max|\\partial f|$');\n", "plt.title('Max of $|\\partial^n f|$');\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Interpolation de Chebyshev\n", "\n", "\n", "Pour chaque entier positif $n\\geq 1$, pour $i=0,\\dots n$, on note\n", "$$\\hat{t}_i=-\\cos(\\pi i/n)\\in[-1,1]$$\n", "**les points de Chebyshev** et on définit\n", "\n", "$$t_i=\\dfrac{a+b}{2}+\\dfrac{b-a}{2}\\hat{t}_i\\in[a,b],$$\n", "\n", "pour un intervalle arbitraire $[a,b]$. Pour une fonction continue $f\\in\n", "C^1([a,b])$, le polynôme d'interpolation $\\Pi_n f$ de degré $n$ aux noeuds\n", "$\\{t_i,i=0,\\ldots,n\\}$ converge uniformément vers $f$ quand\n", "$n\\rightarrow \\infty$.\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Chebichev points on the interval [-1,1]\n", "\n", "z = np.linspace(0,1, 100)\n", "plt.plot(np.cos(np.pi*z), np.sin(np.pi*z))\n", "\n", "n =5\n", "z = np.linspace(0,1, n+1)\n", "plt.plot(np.cos(np.pi*z), np.sin(np.pi*z), 'o')\n", "plt.plot(np.cos(np.pi*z), 0*z, 'x')\n", "\n", "for k in range(0,n+1) :\n", " plt.plot([np.cos(np.pi*z[k]),np.cos(np.pi*z[k])],[0,np.sin(np.pi*z[k])],':')\n", "\n", "plt.axis('equal')\n", "\n", "plt.xlabel('t'); \n", "plt.title('Chebyshev points')\n", "plt.show()\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Exemple\n", "\n", "On reprend le\n", " même exemple mais on interpole la fonction de Runge dans les points de\n", "Chebyshev. La figure montre les polynômes de\n", "Chebyshev de degré $n=5$ et $n=10$. On remarque que les oscillations\n", "diminuent lorsqu'on augmente le degré du polynôme.\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Runge fonction\n", "f = lambda t : 1./(1+t**2)\n", "\n", "# Values of N to use\n", "Nrange = [3,5,10]\n", "# plotting points\n", "[a,b] = [-5,5]\n", "z = np.linspace(a,b, 100)\n", "\n", "for n in Nrange :\n", "\n", " # Chebyshev points on [-1,1]\n", " ht = -np.cos(np.pi*np.linspace(0,n,n+1)/n)\n", " # mapped to [a,b]\n", " t =(a+b)/2 + (b-a)/2*ht\n", "\n", " y = f(t);\n", "\n", " plt.plot(z, LagrangePi(t,y,z), ':')\n", "\n", "plt.plot(z,f(z), 'b')\n", "plt.xlabel('t'); plt.ylabel('y'); plt.title('Chebyshev interpolation')\n", "plt.legend(Nrange)\n", "plt.show()" ] }, { "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 }