diff --git a/2.1 Interpolation.ipynb b/2.1 Interpolation.ipynb
new file mode 100644
index 0000000..bc8f1e1
--- /dev/null
+++ b/2.1 Interpolation.ipynb	
@@ -0,0 +1,347 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# 1 Interpolation et approximation de données\n",
+    "\n",
+    "## 1.1 Position du problème"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Interpolation de données\n",
+    "\n",
+    "Soit $n \\geq 0$ un nombre entier. Etant donnés <font color=blue>$(n+1)$ noeuds</font> \n",
+    "distincts $x_0$, $x_1$,$\\dots$ $x_n$ et $(n+1)$ valeurs $y_0$,\n",
+    "$y_1$,$\\dots$ $y_n$, on cherche un polynôme $p$\n",
+    " <font color=blue>de degré $n$</font>, tel que\n",
+    "\n",
+    "$$p(x_j)=y_j \\qquad \\text{ pour } \\, 0\\leq j\\leq n.$$\n",
+    "\n",
+    "**Exemple** On cherche le polynôme $\\Pi_n$ de degré $n=4$ tel que $\\Pi_n(x_j) = y_j, j =1,...,5$ avec \n",
+    "les données suivantes \n",
+    "\n",
+    "| $x_k$ | $y_k$ |\n",
+    "| --- | --- |\n",
+    "| 1 | 3 |\n",
+    "| 1.5 | 4 |\n",
+    "| 2 | 2 |\n",
+    "| 2.5 | 5 |\n",
+    "| 3 | 1 |\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {
+    "slideshow": {
+     "slide_type": "slide"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "# importing libraries used in this book\n",
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "\n",
+    "# Some data given: x=1, 1.5, 2, 2.5, 3 and y = 3,4,2,5,1 \n",
+    "x = np.linspace(1, 3, 5) # equivalent to np.array([ 1, 1.5, 2, 2.5, 3 ])\n",
+    "y = np.array([3, 4, 2, 5, 1])\n",
+    "\n",
+    "# Plot the points using matplotlib\n",
+    "plt.plot(x, y, 'ro')\n",
+    "\n",
+    "plt.xlabel('x'); plt.ylabel('y'); #plt.title('data')\n",
+    "plt.legend(['data'])\n",
+    "plt.show()\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Si ce polynôme existe, on le note $p=\\Pi_n$. On appelle $\\Pi_n$ le\n",
+    "<font color=blue>polynôme d'interpolation des valeurs $y_j$ aux noeuds $x_j,\\, j=0,\\dots,n$</font>."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Plot the interpolating function\n",
+    "\n",
+    "# Defining the polynomial function \n",
+    "def p(x):\n",
+    "    # coefficients of the interpolating polynomial\n",
+    "    a = np.array([-140.0, 343.0, -872./3.,  104.0,  -40./3.])\n",
+    "    \n",
+    "    # value of the polynomial in all the points t\n",
+    "    return  a[0] + a[1]*x + a[2]*(x**2) + a[3]*(x**3) + a[4]*(x**4)\n",
+    "\n",
+    "\n",
+    "# points used to plot the graph \n",
+    "z = np.linspace(1, 3, 100)\n",
+    "\n",
+    "plt.plot(x, y, 'ro', z, p(z))\n",
+    "plt.xlabel('x'); plt.ylabel('y'); #plt.title('data')\n",
+    "plt.legend(['data','$\\Pi_2(x)$'])\n",
+    "plt.show()\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Interpolation de fonctions\n",
+    "\n",
+    "\n",
+    "Soit $f\\in C^0(I)$ et $x_0,\\ldots,x_n\\in I$. \n",
+    "Si on prend $$y_j=f(x_j),\\quad 0\\leq j\\leq n,$$ \n",
+    "alors le polynôme d'interpolation\n",
+    "$\\Pi_n(x)$ est noté $\\Pi_n f(x)$ et est appelé <font color=blue>l'interpolant de $f$ aux\n",
+    "noeuds $x_0$,$\\dots$ $x_n$</font>.\n",
+    "\n",
+    "**Exemple** Soient $$x_1=1, x_2=1.75, x_3=2.5, x_4=3.25, x_5=4$$ les points d'interpolation et $$f(x) = x \\sin(2\\pi x).$$ On cherche l'interpolant $\\Pi_n f$ de degré $n=4$\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# defining the fonction that we want to interpolate\n",
+    "def f(x):\n",
+    "    return x*np.sin(x*2.*np.pi) \n",
+    "\n",
+    "# The interpolation must occour at points x=1, 1.75, 2.5, 3.25, 4\n",
+    "x = np.linspace(1, 4, 5)\n",
+    "\n",
+    "# points used to plot the graph \n",
+    "z = np.linspace(1, 4, 100)\n",
+    "\n",
+    "plt.plot(x, f(x), 'ro', z, f(z),':')\n",
+    "\n",
+    "# labels, title, legend\n",
+    "plt.xlabel('x'); plt.ylabel('$f(x)$'); #plt.title('data')\n",
+    "plt.legend(['$f(x_k)$','$f(x)$'])\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Plot the interpolating function\n",
+    "\n",
+    "# Defining the polynomial function \n",
+    "def p(x):\n",
+    "    # coefficients of the interpolating polynomial\n",
+    "    a = np.array([0,       7.9012,      -13.037,       5.9259,     -0.79012])\n",
+    "    \n",
+    "    # value of the polynomial in all the points x\n",
+    "    return  a[0] + a[1]*x + a[2]*(x**2) + a[3]*(x**3) + a[4]*(x**4)\n",
+    "\n",
+    "\n",
+    "# points where to evaluate the polynomial\n",
+    "z = np.linspace(1, 4, 100)\n",
+    "\n",
+    "plt.plot(x, f(x), 'ro', z, f(z),':', z, p(z))\n",
+    "plt.xlabel('$x$'); plt.ylabel('$f(x)$'); #plt.title('data')\n",
+    "plt.legend(['$f(x_k)$','$f(x)$','$\\Pi_n(x)$'])\n",
+    "\n",
+    "plt.show()\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Matrice de Vandermonde\n",
+    "\n",
+    "Il est possible d'écrire un système d'équations et de trouver les coefficients de manière directe.\n",
+    "Ce n'est pas toujours la meilleure solution.\n",
+    "\n",
+    "Nous cherchons les coefficients du polynôme $p(x) = a_0 + a_1 x + ... + a_n x^n$ qui satisfont les $(n+1)$ équations $p(x_k) = y_k, k=0,...,n$, c'est-à-dire\n",
+    "\n",
+    "$$a_0 + a_1 x_k + ... + a_n x_k^n = y_k, k=0,...,n$$\n",
+    "\n",
+    "Ce système s'écrit sous forme matricielle\n",
+    "\n",
+    "$$\\begin{pmatrix}\n",
+    "1 & x_0 & x_0^2 & \\cdots & x_0^n \\\\\n",
+    "\\vdots &   &   &  & \\vdots \\\\\n",
+    "1 & x_n & x_n^2 & \\cdots & x_n^n\n",
+    "\\end{pmatrix}\n",
+    "\\begin{pmatrix} a_0 \\\\ \\vdots \\\\ a_n \\end{pmatrix}\n",
+    "=\\begin{pmatrix} y_0 \\\\ \\vdots \\\\ y_n \\end{pmatrix}$$\n",
+    "\n",
+    "Pour construire cette matrice, vous pouvez utiliser la fonction\n",
+    "```python\n",
+    "# Defining the mxn Vandermonde matrix \n",
+    "def VandermondeMatrix(x):\n",
+    "    # Input\n",
+    "    # x : +1 array with interpolation nodes\n",
+    "    # Output\n",
+    "    # Matrix of Vandermonde of size n x n\n",
+    "```\n",
+    "que vous pouvez importer avec la commande\n",
+    "```python\n",
+    "from InterpolationLib import VandermondeMatrix \n",
+    "```"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "**Exemple** On cherche les coefficients du polynôme d'interpolation de degré $n=4$ \n",
+    "des valeurs suivantes \n",
+    "\n",
+    "| $x_k$ | $y_k$ |\n",
+    "| --- | --- |\n",
+    "| 1 | 3 |\n",
+    "| 1.5 | 4 |\n",
+    "| 2 | 2 |\n",
+    "| 2.5 | 5 |\n",
+    "| 3 | 1 |\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from InterpolationLib import VandermondeMatrix \n",
+    "\n",
+    "# Some data given: x=1, 1.5, 2, 2.5, 3 and y = 3,4,2,5,1 \n",
+    "x = np.linspace(1, 3, 5)\n",
+    "y = np.array([3, 4, 2, 5, 1])\n",
+    "n = x.size - 1\n",
+    "\n",
+    "A = VandermondeMatrix(x)\n",
+    "# print(A)\n",
+    "\n",
+    "# compute coefficients\n",
+    "# Resouds Ax = b avec b=y et rends x\n",
+    "# >>> COMPLETE HERE <<<\n",
+    "\n",
+    "# print the coefficients on screen\n",
+    "print('The coefficients a_0, ..., a_n are')\n",
+    "print(a)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Now we can define the polynomial\n",
+    "p = lambda x : a[0] + a[1]*x + a[2]*(x**2) + a[3]*(x**3) + a[4]*(x**4)\n",
+    "\n",
+    "# points used to plot the graph \n",
+    "z = np.linspace(1, 3, 100)\n",
+    "\n",
+    "# Plot points and function\n",
+    "# >>> COMPLETE HERE <<<\n",
+    "plt.xlabel('x'); plt.ylabel('y'); #plt.title('data')\n",
+    "plt.legend(['data','p(x)'])\n",
+    "plt.show()\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Alternatives : polyfit et polyval\n",
+    "\n",
+    "Les fonctions `polyfit` et `polyval` de `numpy` font essentiellement la même chose que les paragraphes ci-dessous. Plus tard nous verrons des méthodes plus performantes.\n",
+    "\n",
+    "`a = numpy.polyfit(x, y, n, ... )` :\n",
+    "\n",
+    " * input : $x,y$ les données à interpoler, $n$ le degré du polynôme recherché\n",
+    " * output : les coefficients du polynôme, _dans l'ordre inverse_ de ce que nous avons vu !\n",
+    " "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Some data given: x=1, 1.5, 2, 2.5, 3 and y = 3,4,2,5,1 \n",
+    "x = np.linspace(1, 3, 5)\n",
+    "y = np.array([3, 4, 2, 5, 1])\n",
+    "n = x.size - 1\n",
+    "\n",
+    "a = np.polyfit(x,y,n)\n",
+    "\n",
+    "# Now we can define the polynomial, with coeffs in the reverse order !\n",
+    "p = lambda x : a[4] + a[3]*x + a[2]*(x**2) + a[1]*(x**3) + a[0]*(x**4)\n",
+    "\n",
+    "# We can also use polyval instead !\n",
+    "# np.polyval(a,x)\n",
+    "\n",
+    "# points used to plot the graph \n",
+    "z = np.linspace(1, 3, 100)\n",
+    "\n",
+    "plt.plot(x, y, 'ro', z, p(z), '.', z, np.polyval(a,z))\n",
+    "plt.xlabel('x'); plt.ylabel('y'); #plt.title('data')\n",
+    "plt.legend(['data','p(x)','polyval'])\n",
+    "plt.show()\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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/2.2 Interpolation de Lagrange.ipynb b/2.2 Interpolation de Lagrange.ipynb
new file mode 100644
index 0000000..9559b95
--- /dev/null
+++ b/2.2 Interpolation de Lagrange.ipynb	
@@ -0,0 +1,191 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 1.2 Interpolation de Lagrange\n",
+    "\n",
+    "### Base de Lagrange\n",
+    "\n",
+    "On considère les polynômes $\\varphi_k$, $k=0,\\dots, n$ de\n",
+    "degré $n$ tels que\n",
+    "\n",
+    "$$\\varphi_k(x_j)=\\delta_{jk}, \\qquad k,j=0,\\dots, n,$$\n",
+    "\n",
+    "où $\\delta_{jk}=1$ si $j=k$ et $\\delta_{jk}=0$ si $j\\neq k$.\n",
+    "Explicitement, on a\n",
+    "\n",
+    "$$\\varphi_k(x)=\\prod_{j=0,j\\ne k}^n\\dfrac{(x-x_j)}{(x_k-x_j)}.$$"
+   ]
+  },
+  {
+   "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"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "**Exercice (théorique)** Vérifiez que\n",
+    "\n",
+    "1. $B = \\{\\varphi_k, k=0,\\dots, n\\}$ est une base de $P_n(\\mathbb R)$\n",
+    "2. Chaque polynôme $\\varphi_k$ est de degré $n$\n",
+    "3. $\\varphi_k(x_j)=\\delta_{jk}, \\qquad k,j=0,\\dots, n$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Defining the Lagrange basis functions\n",
+    "def phi(x,k,z):\n",
+    "    # the input variables are:\n",
+    "    # x : the interpolatory points\n",
+    "    # k : which basis function\n",
+    "    # z : where to evaluate the function\n",
+    "    \n",
+    "\n",
+    "    # careful, there are n+1 interpolation points!\n",
+    "    n = x.size - 1 \n",
+    "\n",
+    "    # init result to one, of same type and size as z\n",
+    "    result = np.zeros_like(z) + 1\n",
+    "\n",
+    "    # first few checks on k:\n",
+    "    if (type(k) != int) or (x.size < 1) or (k > n) or (k < 0):\n",
+    "        raise ValueError('Lagrange basis needs a positive integer k, smaller than the size of x')\n",
+    "        \n",
+    "    # loop on n to compute the product\n",
+    "    for j in range(0,n+1) :\n",
+    "        if (j == k) :\n",
+    "            continue\n",
+    "        if (x[k] == x[j]) :\n",
+    "            raise ValueError('Lagrange basis: all the interpolation points need to be distinct')\n",
+    "\n",
+    "        # >>> COMPLETE HERE <<<\n",
+    "        result = \n",
+    "\n",
+    "    return result\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "#### Exemple\n",
+    "\n",
+    "Pour $n=2$, $x_0=-1$, $x_1=0$, $x_2=1$, les polynômes de la base de Lagrange\n",
+    "sont\n",
+    "\n",
+    "\\begin{equation*}\n",
+    "\\begin{array}{lcl}\n",
+    "\\varphi_0(x)&=&\\displaystyle\\dfrac{(x-x_1)(x-x_2)}{(x_0-x_1)(x_0-x_2)}=\n",
+    "\\dfrac{1}{2}x(x-1),\\\\[2mm]\n",
+    "\\varphi_1(x)&=&\\displaystyle\\dfrac{(x-x_0)(x-x_2)}{(x_1-x_0)(x_1-x_2)}=\n",
+    "-(x+1)(x-1),\\\\[2mm]\n",
+    "\\varphi_2(x)&=&\\displaystyle\\dfrac{(x-x_0)(x-x_1)}{(x_2-x_0)(x_2-x_1)}=\n",
+    "\\dfrac{1}{2}x(x+1).\n",
+    "\\end{array}\n",
+    "\\end{equation*}\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# plot the Lagrange Basis functions \n",
+    "x = np.linspace(-1., 1, 3)\n",
+    "z = np.linspace(-1.1, 1.1, 100)\n",
+    "\n",
+    "plt.plot(z, phi(x,0,z), 'g', z, phi(x,1,z), 'r', z, phi(x,2,z),':')\n",
+    "\n",
+    "plt.xlabel('x'); plt.ylabel('$\\\\varphi_{k}(x)$'); plt.title('Lagrange basis functions')\n",
+    "plt.legend(['$\\\\varphi_{0}$','$\\\\varphi_{1}$','$\\\\varphi_{2}$'])\n",
+    "plt.grid(True)\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "#### Exercice\n",
+    "\n",
+    "Visualisez la base de Lagrange associée aux points $x_k = 1, 1.5, 2, 2.5, 3$.\n",
+    "Evaluez sur le graphique les valeurs de $\\varphi_k(x_j)$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# plot the Lagrange Basis functions \n",
+    "x = np.linspace(1., 3, 5)\n",
+    "z = np.linspace(0.9, 3.1, 100)\n",
+    "\n",
+    "plt.plot(z, phi(x,0,z), 'g', z, phi(x,1,z), 'r', z, phi(x,2,z),':', z, phi(x,3,z),':', z, phi(x,4,z),':')\n",
+    "\n",
+    "plt.xlabel('x'); plt.ylabel('phi(x)'); plt.title('Lagrange Basis functions')\n",
+    "plt.legend(['$\\\\varphi_0$','$\\\\varphi_1$','$\\\\varphi_2$',\n",
+    "            '$\\\\varphi_3$','$\\\\varphi_4$'])\n",
+    "plt.grid(True)\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Polynôme d'interpolation\n",
+    "\n",
+    "Le polynôme d'interpolation $\\Pi_n$ des\n",
+    "valeurs $y_j$ aux noeuds $x_j$, $j=0,\\dots,n$, s'écrit\n",
+    "\\begin{equation}\\label{e:int_lagr}\n",
+    "  \\Pi_n(x)=\\sum_{k=0}^n y_k \\varphi_k(x),\n",
+    "\\end{equation}\n",
+    "car il vérifie $\\Pi_n(x_j)=\\sum_{k=0}^n y_k \\varphi_k(x_j)=y_j$.\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/2.3 Interpolation de fonction.ipynb b/2.3 Interpolation de fonction.ipynb
new file mode 100644
index 0000000..b2ce518
--- /dev/null
+++ b/2.3 Interpolation de fonction.ipynb	
@@ -0,0 +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 $x_0,\\ldots,x_n\\in [a,b]$ des noeuds distincts. Le polynôme d'interpolation\n",
+    "$\\Pi_n(x)$ est noté $\\Pi_n f(x)$ et est appelé <font color=blue>l'interpolant de $f$ aux\n",
+    "noeuds $x_0,\\dots,x_n$</font>.\n",
+    "\n",
+    "\n",
+    "Si on prend $$y_k=f(x_k), k= 0,..., n,$$ \n",
+    "alors on aura\n",
+    "\\begin{equation*}\n",
+    "\\Pi_nf(x)=\\sum_{k=0}^n\n",
+    "f(x_k)\\varphi_k(x).\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 (x,y), evaluated at ordinate(s) z\n",
+    "def LagrangePi(x,y,z):\n",
+    "    # the input variables are:\n",
+    "    # x : the interpolatory points\n",
+    "    # y : the corresponding data at the points x\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ées 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 (x,y), evaluated at ordinate(s) z\n",
+    "def LagrangePi(x,y,z):\n",
+    "    # the input variables are:\n",
+    "    # x : the interpolatory points\n",
+    "    # y : the corresponding data at the points x\n",
+    "    # z : where to evaluate the function\n",
+    "    \n",
+    "    # {phi(x,k,.), k=0,...,n} is a basis of the polynomials of degree n\n",
+    "    # y represents the coordinates of the interpolating polynomial with respect to this basis.\n",
+    "    # Therefore LagrangePi(x,y,.) = y[0] phi(x,0,.) + ... + y[n] phi(x,n,.)\n",
+    "\n",
+    "    # >>> COMPLETE HERE <<<\n",
+    "\n",
+    "    return result"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# vecteur des points d'interpolation\n",
+    "x = 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(x, y, 'ro', z, LagrangePi(x,y,z))\n",
+    "\n",
+    "plt.xlabel('x'); plt.ylabel('y');\n",
+    "plt.legend(['data','p(x)'])\n",
+    "plt.show()\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Erreur d'interpolation\n",
+    "\n",
+    "Soient $x_0$, $x_1$, $\\ldots$, $x_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(x) - \\Pi_n f(x)|  \\leq \\dfrac{1}{4(n+1)}\n",
+    "\\left(\\frac{b-a}{n}\\right)^{n+1}\\max_{x\\in I}|f^{(n+1)}(x)|.\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 $$x_0=1, x_1=1.75, x_2=2.5, x_3=3.25, x_4=4$$ \n",
+    "et la fonction $$f(x) = x \\sin(2\\pi x)$$.\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",
+    "n'implique pas forcément que $\\max_{x \\in I} |E_nf(x)|$ tende vers zéro quand $n\\rightarrow\\infty$.\n",
+    "\n",
+    "Soit $f(x)=\\dfrac{1}{1+x^2}$, $x\\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 x : 1./(1+x**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",
+    "for n in Nrange :\n",
+    "\n",
+    "    # define x and y=f(x)\n",
+    "    # compute the Lagrange interpolant\n",
+    "    # plot it\n",
+    "\n",
+    "    # >>> COMPLETE HERE <<<\n",
+    "\n",
+    "\n",
+    "plt.plot(z,f(z), 'b')\n",
+    "plt.xlabel('x'); plt.ylabel('y'); plt.title('Runge function')\n",
+    "plt.legend(Nrange)\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "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{x}_i=-\\cos(\\pi i/n)\\in[-1,1]$$\n",
+    "**les points de Chebyshev** et on définit\n",
+    "\n",
+    "$$x_i=\\dfrac{a+b}{2}+\\dfrac{b-a}{2}\\hat{x}_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",
+    "$\\{x_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és $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 x : 1./(1+x**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",
+    "    hx = -np.cos(np.pi*np.linspace(0,n,n+1)/n)\n",
+    "    # mapped to [a,b]\n",
+    "    x =(a+b)/2 + (b-a)/2*hx\n",
+    "\n",
+    "    y = f(x);\n",
+    "\n",
+    "    plt.plot(z, LagrangePi(x,y,z), ':')\n",
+    "\n",
+    "plt.plot(z,f(z), 'b')\n",
+    "plt.xlabel('x'); 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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/2.4 Interpolation par intervalles.ipynb b/2.4 Interpolation par intervalles.ipynb
new file mode 100644
index 0000000..831760f
--- /dev/null
+++ b/2.4 Interpolation par intervalles.ipynb	
@@ -0,0 +1,287 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 1.4 Interpolation par intervalles\n",
+    "\n",
+    "### Interpolation linéaire par morceaux\n",
+    "\n",
+    "Soit $f: \\mathbf [a,b] \\rightarrow R$ continue et $a=x_0<\\ldots<x_n=b$. \n",
+    "\n",
+    "On choisit une partition de $[a,b]$ en $N$ sous-intervalles de la forme $[x_i, x_{i+1}], i=0,...,N-1]$.\n",
+    "Sur chaque sous-intervalle, on fait une interpolation de degré 1 avec les 2 noeuds $x_i, x_{i+1}$.\n",
+    "Sur chaque sous-intervalle $I_i=[x_i, x_{i+1}]$, on interpole $f_{\\mid I_i}$ par\n",
+    "un polynôme de degré 1. Le polynôme par morceaux (polynôme\n",
+    "  composite) qu'on obtient est noté $\\Pi_1^H f(x)$ et on a:\n",
+    "$$\\Pi_1^H f(x) = f(x_i) + \\dfrac{f(x_{i+1}) - f(x_i)}{x_{i+1} - x_i}\n",
+    "(x - x_i) \\quad \\text{pour } x\\in [x_i, x_{i+1}]$$\n",
+    "\n",
+    "Dans le cas de données $y$, cela s'écrit\n",
+    "$$\\Pi_1^H (x) = y_i + \\dfrac{y_{i+1} - y_i}{x_{i+1} - x_i}\n",
+    "(x - x_i) \\quad \\text{pour } x\\in [x_i, x_{i+1}]$$\n",
+    "\n",
+    "Le choix le plus simple est le suivant :\n",
+    "\n",
+    "* on pose $H = \\frac{b-a}{N}$\n",
+    "* ensuite $x_i = a + iH$ pour $i=0,...,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 PiecewiseLinearInterpolation as  PiH1\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "La fonction `PiH1` implémente l'interpolation linéaire par morceaux pour des points équidistribués\n",
+    "\n",
+    "```python\n",
+    "def PiecewiseLinearInterpolation(a,b,N,f,z):\n",
+    "    # the input variables are:\n",
+    "    # a,b : x[0] = a, x[n] = b\n",
+    "    # f : the corresponding data at the points x\n",
+    "    # z : where to evaluate the function\n",
+    "```\n",
+    "    \n",
+    "    "
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Exercice\n",
+    "\n",
+    "Soient $x_k, k=0,...,4$ des points équidistribués sur l'intervalle $[1,5]$ et \n",
+    "$y_k = (3.38, 3.86, 3.85, 3.59, 3.49)$\n",
+    "les valeurs d'une fonction en ces points. \n",
+    "\n",
+    "* Dessinez le graphe de l'interpolateur par morceaux de cette fonction\n",
+    "* Calculez numériquement (sur papier) la valeur de $\\Pi^H_1(4.5)$ et vérifiez le résultat sur le graphique\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# intervalle d'interpolation\n",
+    "a = 1; b = 5\n",
+    "\n",
+    "# define the vector of values at equidistributed points x\n",
+    "# N the number of intervals\n",
+    "# y the array of values\n",
+    "# then use PiH1(a,b,N,y,z) to compute and plot the piecewise linear interpolant\n",
+    "\n",
+    "# >>> COMPLETE HERE <<<\n",
+    "\n",
+    "\n",
+    "\n",
+    "plt.xlabel('x'); plt.ylabel('y');\n",
+    "plt.legend(['data','p(x)'])\n",
+    "plt.show()\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Exercice\n",
+    "\n",
+    "Dessinez le graphe de la fonction de Runge et l'interpolateur linéaire par morceaux \n",
+    "sur l'intervalle $[-5,5]$ pour $N=3,5,10$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# interval and function\n",
+    "a = -5; b = 5\n",
+    "f = lambda x : 1/(1+x**2)\n",
+    "\n",
+    "# >>> COMPLETE HERE <<<\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Erreur d'interpolation linéaire par morceaux \n",
+    "\n",
+    "#### Théorème \n",
+    "\n",
+    "Soient $f\\in C^{2}([a,b])$, $H = \\frac{b-a}{N}$, $x_i = a + iH$ pour $i=0,...,N$.\n",
+    "  \n",
+    "Soit $E^H_1f(x) = f(x) - \\Pi_1^{H} f(x)$, alors\n",
+    "$$\\max_{x\\in I}\\mid E^H_1f(x) \\mid\\leq \\dfrac{H^2}{8}\\max_{x \\in I} | f''(x)|.$$\n",
+    "\n",
+    "**Preuve**\n",
+    "D'après la formule, sur chaque intervalle $I_i=[x_i,x_{i+1}]$, on a\n",
+    "  \n",
+    "$$\\max_{x\\in[x_i,x_{i+1}]}\\mid E^H_1f(x)\\mid \\leq\n",
+    "      \\dfrac{H^2}{4(1+1)}\\max_{x\\in I_i}\\mid f''(x)\\mid.$$\n",
+    "\n",
+    "**Remarque**\n",
+    "On peut aussi montrer que, si l'on utilise un polynôme\n",
+    "de degré $n$ ($\\geq 1$) et si l'on dénote $E^H_nf(x) = f(x) - \\Pi_n^{H} f(x)$,\n",
+    " dans chaque sous-intervalle $I_i$, on\n",
+    "trouve\n",
+    "\n",
+    "$$\\max_{x\\in I}\\mid E^H_nf(x) \\mid \\leq \\dfrac{H^{n+1}}{4(n+1)} \\max_{x \\in I} |f^{(n+1)} (x)| \\, .$$\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "### Interpolation quadratique par morceaux\n",
+    "\n",
+    "Soit $f: \\mathbf [a,b] \\rightarrow R$ continue et $a=x_0<\\ldots<x_n=b$. \n",
+    "Sur chaque sous-intervalle $[x_i,x_{i+1}]$, on fait une interpolation de **degré 2** avec \n",
+    "les **3 noeuds** $x_i, x_{i+\\frac12}, x_{i+1}$, où $x_{i+\\frac12}$ est le milieu de $[x_i,x_{i+1}]$.\n",
+    "Sur chaque sous-intervalle $I_i=[x_i, x_{i+1}]$, on interpole $f_{\\mid I_i}$ par\n",
+    "un polynôme de degré 2. Le polynôme par morceaux (polynôme\n",
+    "  composite) qu'on obtient est noté $\\Pi_2^H f(x)$ et on a:\n",
+    "\n",
+    "$$\\Pi_2^H f(x) = f(x_i)\\varphi^{(i)}_0(x) + f(x_{i+\\frac12})\\varphi^{(i)}_1(x) \n",
+    "    + f(x_{i+1})\\varphi^{(i)}_2(x) \\quad \\text{pour } x\\in [x_i, x_{i+1}]$$\n",
+    "\n",
+    "où $\\varphi^{(i)}_0(x), \\varphi^{(i)}_1(x), \\varphi^{(i)}_0(x)$ sont les polynômes de la base de Lagrange\n",
+    "associés aux noeuds $(x_i,x_{i+\\frac12}, x_{i+1})$.\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from InterpolationLib import PiecewiseQuadraticInterpolation as  PiH2\n",
+    "\n",
+    "# intervalle et fonction\n",
+    "a = -5; b = 5\n",
+    "f = lambda x : 1/(1+x**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",
+    "for N in Nrange :\n",
+    "    plt.plot(z, PiH2(a,b,N,f,z), ':')\n",
+    "\n",
+    "plt.plot(z, f(z), 'b')\n",
+    "plt.xlabel('x'); plt.ylabel('y'); #plt.title('data')\n",
+    "plt.legend(Nrange)\n",
+    "plt.show()\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "#### Exercice\n",
+    "\n",
+    "Calculez l'erreur d'interpolation de la fonction de Runge dans l'intervalle $[-5,5]$ quand on utilise l'interpolation linéaire par morceaux\n",
+    "\n",
+    "1. Théoriquement\n",
+    "2. Faites un graphique qui montre la convergence en fonction de $H=\\frac{b-a}N=\\frac{10}{N}$\n",
+    "3. Refaites le même exercice pour l'interpolation quadratique par morceaux"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "#### Remarque\n",
+    "\n",
+    "On s'attend à ce que l'erreur soit quadratique en H. Pour le voir il faut utiliser \n",
+    "des axes logarithmiques dans les deux directions.\n",
+    "\n",
+    "Admettons que l'erreur converge quadratiquement vers 0 par rapport à $H$, e.g\n",
+    "```Python\n",
+    "err = lambda H : 3*H**2 + 2*H**3;\n",
+    "\n",
+    "```\n",
+    "Comment est le graphique de cette fonction dans l'intervalle $[10^{-6},1]$ ?\n",
+    "Que se passe-t-il si on change les axes avec `plt.xscale('log')` et  `plt.yscale('log')` ?\n",
+    "Quelle est la pente de `err` dans ce système ?\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "## Assume the error is quadratic\n",
+    "err = lambda H : 10*H**2 + 20*H**3;\n",
+    "\n",
+    "# take h has powers of 2: 2^(-20),2^(-19), ..., 2^(-1) \n",
+    "h = np.power(2,np.linspace(-20, -1, 20, endpoint=True) )\n",
+    "\n",
+    "plt.plot(h, err(h), 'b:.')\n",
+    "\n",
+    "plt.xlabel('h'); plt.ylabel('err');\n",
+    "\n",
+    "# plt.xscale('log') \n",
+    "# plt.yscale('log')\n",
+    "\n",
+    "plt.grid(True)\n",
+    "# try with and wothout equal axis\n",
+    "plt.axis('equal')\n",
+    "plt.show()\n",
+    "\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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/2.5 Approximations de donnee.ipynb b/2.5 Approximations de donnee.ipynb
new file mode 100644
index 0000000..8f31e33
--- /dev/null
+++ b/2.5 Approximations de donnee.ipynb	
@@ -0,0 +1,279 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 1.5 Approximation au sens des moindres carrés\n",
+    "\n",
+    "Soit $m \\geq 0$ un nombre entier. Etant donnés <font color=blue>$(n+1)$ points</font> \n",
+    "~~distincts~~ $x_0$, $x_1$,$\\dots$ $x_n$ et $(n+1)$ valeurs $y_0$,\n",
+    "$y_1$,$\\dots$ $y_n$, on cherche un polynôme $p$\n",
+    " <font color=blue>de degré $m<n$</font> , tel que\n",
+    "\n",
+    "$$\\sum_{j=0}^n |y_j - p(x_j)|^2 \\qquad \\text{ soit le plus petit possible}.$$\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "On appelle **polynôme aux moindres carrés de degré $m$**\n",
+    "le polynôme $\\hat{p}_n(x)$ de degré $m$ tel que\n",
+    "\\begin{equation}\n",
+    "  \\sum_{j=0}^n |y_j - \\hat{p}_n(x_j)|^2 \\leq \\sum_{j=0}^n |y_j - p_n(x_j)|^2\n",
+    "  \\qquad \\forall p_m(x) \\in \\mathbb{P}_m\n",
+    "\\end{equation}"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Nous cherchons les coefficients du polynôme $p(x) = a_0 + a_1 x + ... + a_m x^m$ qui satisfait *au mieux* les $(n+1)$ équations $p(x_k) = y_k, k=0,...,n$, c'est-à-dire\n",
+    "\n",
+    "$$a_0 + a_1 x_k + ... + a_m x_k^m = y_k, k=0,...,n$$\n",
+    "\n",
+    "Ce système s'écrit sous forme matricielle \n",
+    "\n",
+    "$$\\begin{pmatrix}\n",
+    "1 & x_0 & \\cdots & x_0^m \\\\\n",
+    "\\vdots   &   &  & \\vdots \\\\\n",
+    "1 & x_n & \\cdots & x_n^m\n",
+    "\\end{pmatrix}\n",
+    "\\begin{pmatrix} a_0 \\\\ \\vdots \\\\ a_m \\end{pmatrix}\n",
+    "=\\begin{pmatrix} y_0 \\\\ \\vdots \\\\ y_n \\end{pmatrix}$$\n",
+    "\n",
+    "Puisque $m<n$, on ne peut pas résoudre ce système de façon\n",
+    "classique.\n",
+    "\n",
+    "Il faut le résoudre au sens des moindres carrés, en considérant:\n",
+    "$$B^T B \\mathbf{a} = B^T \\mathbf{y}$$\n",
+    "\n",
+    "Où $B$ est la matrice $(m+1)\\times(n+1)$ du système, $\\mathbf{a}\\in\\mathbf R^{n+1}$ le vecteur des inconnues et $\\mathbf{y}\\in\\mathbf R^{m+1}$\n",
+    "le vecteur des données.\n",
+    "\n",
+    "Ce système linéaire est dit\n",
+    "**système d'équations normales**. On peut montrer que les\n",
+    "équations normales sont équivalentes au problème de minimisation.\n",
+    "\n",
+    "**Remarque:** *La résolution de ce système demande parfois des méthodes plus avancées que l'élimination de Gauss, comme la factorisation QR. Pour l'instant on va se contenter d'utiliser* `np.linalg.solve`"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "Pour construire cette matrice, vous pouvez utiliser la fonction\n",
+    "```python\n",
+    "def VandermondeMatrix(x, m=0):\n",
+    "    # Input\n",
+    "    # x : +1 array with interpolation nodes\n",
+    "    # m : degree of the polynomial. If empty, chooses m=size(x)-1 \n",
+    "    # Output\n",
+    "    # Matrix of Vandermonde of size m x n\n",
+    "```\n",
+    "que vous pouvez importer avec la commande\n",
+    "```python\n",
+    "from InterpolationLib import VandermondeMatrix \n",
+    "```"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "#### Exemple\n",
+    "\n",
+    "On considère un test mécanique pour\n",
+    "établir le lien entre contraintes ($MPa=100 N/cm^2$)\n",
+    "et déformations relatives (cm/cm) d'un échantillon de tissu biologique (disque intervertébral, selon P. Komarek,\n",
+    "Ch. 2 de *Biomechanics of Clinical Aspects of Biomedicine*, 1993, J. Valenta ed., Elsevier).\n",
+    "\n",
+    "![title](Figures/disque.png)\n",
+    "\n",
+    "\n",
+    "On cherche à approximer au sens des moindres carrés avec un polynôme $\\hat p_n$ de degré $n=1,2,3$.\n",
+    "Les mesures effectuées sont les suivantes \n",
+    "```Python\n",
+    "sigma = np.array([0.00, 0.06, 0.14, 0.25, 0.31, 0.47, 0.50, 0.70]);\n",
+    "epsilon = np.array([0.00, 0.08, 0.14, 0.20, 0.22, 0.26, 0.27, 0.29]);\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 VandermondeMatrix "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# data given:\n",
+    "sigma = np.array([0.00, 0.06, 0.14, 0.25, 0.31, 0.47, 0.50, 0.70]);\n",
+    "epsilon = np.array([0.00, 0.08, 0.14, 0.20, 0.22, 0.26, 0.27, 0.29]);\n",
+    "\n",
+    "# degree of the polynomial\n",
+    "m = 1;\n",
+    "\n",
+    "B = VandermondeMatrix(sigma,m)\n",
+    "\n",
+    "# print(B)\n",
+    "\n",
+    "# compute coefficients\n",
+    "# >>> COMPLETE HERE <<<\n",
+    "a = \n",
+    "\n",
+    "# print the coefficients on screen\n",
+    "print('The coefficients a_0, ..., a_m are')\n",
+    "print(a)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def polynomial(a,x):\n",
+    "    m = a.size-1\n",
+    "    # \\hat p = a_0 + a_1 x + ... + a_n x^m\n",
+    "    # is equal to the scalar product between the vectors a and (1, x, ..., x^m) :\n",
+    "    return np.power( np.tile(x, (m+1, 1)).T , np.linspace(0,m,m+1)).dot(a)\n",
+    "\n",
+    "\n",
+    "# points used to plot the graph, slightly larger than data\n",
+    "# >>> COMPLETE HERE <<<\n",
+    "z = \n",
+    "\n",
+    "plt.plot(sigma, epsilon, 'ro', z, polynomial(a,z),'b')\n",
+    "plt.xlabel('$\\sigma$'); plt.ylabel('$\\epsilon$');\n",
+    "plt.legend(['data','$\\hat p_n$'])\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "#### Exercice\n",
+    "\n",
+    "Depuis  https://hsso.ch/fr/2012/b/14 on a téléchargé les données de la population suisse dans\n",
+    "le fichier `Data/PopulationSuisse.csv`.\n",
+    "\n",
+    "Approximez l'évolution de la population avec un polynôme de degré $n=1,2,3,7$. \n",
+    "\n",
+    "Ensuite faites l'hypothèse de croissance exponentielle de la population, c'est-à-dire $p(x) = e^{a_1 x+a_0}$ où $a_0$ et $a_1$ sont des paramètres. Comment utiliser l'approximation polynômiale dans ce contexte ? "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Data of the population has been dowloaded from https://hsso.ch/fr/2012/b/14\n",
+    "# into the file Data/PopulationSuisse.csv\n",
+    "import pandas as pd \n",
+    "\n",
+    "# Read data from file 'PopulationSuisse.csv' \n",
+    "data = pd.read_csv(\"Data/PopulationSuisse.csv\") \n",
+    "# Preview the first 5 lines of the loaded data \n",
+    "data.head()\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "x = data['Année'].to_numpy()\n",
+    "y = data['Population'].to_numpy()\n",
+    "\n",
+    "m = 2\n",
+    "B = VandermondeMatrix(x,m)\n",
+    "\n",
+    "# >>> COMPLETE HERE <<<\n",
+    "# ... and solve the exercise\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def polynomial(a,x):\n",
+    "    m = a.size-1\n",
+    "    # \\hat p = a_0 + a_1 x + ... + a_n x^m\n",
+    "    # is equal to the scalar product between the vectors a and (1, x, ..., x^m) :\n",
+    "    return np.power( np.tile(x, (m+1, 1)).T , np.linspace(0,m,m+1)).dot(a)\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "On assume une croissance exponentielle : $population (x) = C * e^{a_1x} = e^{a_o + a_1 x}$\n",
+    "\n",
+    "En d'autres termes: $\\log (population) (x) = a_o + a_1 x$\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# \n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  },
+  {
+   "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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/InterpolationLib.py b/InterpolationLib.py
new file mode 100644
index 0000000..c1c26ad
--- /dev/null
+++ b/InterpolationLib.py
@@ -0,0 +1,165 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Thu Jan 30
+
+@author: Simone Deparis
+"""
+
+import numpy as np
+
+# Defining the mxn Vandermonde matrix 
+def VandermondeMatrix(x, m=0):
+    # Input
+    # x : +1 array with interpolation nodes
+    # m : degree of the polynomial. If empty, chooses m=size(x)-1 
+    # Output
+    # Matrix of Vandermonde of size m x n
+                
+    n = x.size - 1
+    if (m == 0):
+        m = n
+                      
+    # The following is a trick to put toghether the matrix A. Don't need to learn by heart ...
+    # np.tile(x, (n+1, 1)) : repete n+1 times of the array x 
+    #   .T : transpose
+    # np.linspace(0,n,n+1) : [0,...,n+1]
+    # np.tile : again: repete n+1 times the array  [0,...,n+1]
+    # np.power : element by element power funcion
+
+    A = np.power( np.tile(x, (m+1, 1)).T , np.tile(np.linspace(0,m,m+1), (n+1, 1)))
+
+    return A
+
+
+
+
+# Defining the Lagrange basis functions
+def LagrangeBasis(t,k,z):
+    # the input variables are:
+    # t : the interpolatory points
+    # k : which basis function
+    # z : where to evaluate the function
+    
+
+    # careful, there are n+1 interpolation points!
+    n = t.size - 1 
+
+    # init result to one, of same type and size as z
+    result = np.zeros_like(z) + 1
+
+    # first few checks on k:
+    if (type(k) != int) or (t.size < 1) or (k > n) or (k < 0):
+        raise ValueError('Lagrange basis needs a positive integer k, smaller than the size of t')
+        
+    # loop on n to compute the product
+    for j in range(0,n+1) :
+        if (j == k) :
+            continue
+        if (t[k] == t[j]) :
+            raise ValueError('Lagrange basis: all the interpolation points need to be distinct')
+
+        result = result * (z - t[j]) / (t[k] - t[j])
+
+    return result
+
+
+# Piece-wise linear interpolation, equidistributed points
+def PiecewiseLinearInterpolation(a,b,N,f,z):
+    # the input variables are:
+    # a,b : x[0] = a, x[n] = b
+    # f : the corresponding data at the points x
+    # z : where to evaluate the function
+
+    def localLinearInterpolation (a, H, x, y, zk):
+        # first find out in which interval we are. Easy here since equidistributed point
+        i = int( (zk-a)//H )
+
+        # if we are not in the interval, return constant function
+        if x[i] == zk :
+            return y[i]
+        
+        # use formula for local linear interpolation
+        return y[i] + (y[i+1] - y[i])/(x[i+1] - x[i])* (zk - x[i]) 
+    
+    # careful, there are N+1 points 
+    H = (b-a)/N
+    if np.isclose(H,0) :
+        raise ValueError('PiecewiseLinearInterpolation : a and b are too close')
+    
+    x = np.linspace(a,b,N+1)
+    
+    # if f is a function, we need to evaluate it at alle the locations
+    from types import FunctionType
+    if isinstance(f, FunctionType) :
+        y = f(x)
+    else :
+        y = f
+        
+    # if we want to evaluate the interpolating function in a single point
+    if np.ndim(z) == 0 :
+        return localLinearInterpolation (b, H, x, y, z)
+    
+    # if we are here, it means that we need to evaluate the interpolation on many points.
+    # init result to zero, of same type and size as z
+    result = np.zeros_like(z) 
+     
+    # The following part is not optimized
+    for k in range(z.size) :
+
+        result[k] = localLinearInterpolation (a, H, x, y, z[k])
+   
+    return result
+
+# Piece-wise quadratic interpolation, equidistributed points
+def PiecewiseQuadraticInterpolation(a,b,N,f,z):
+    # the input variables are:
+    # a,b : x[0] = a, x[n] = b
+    # f : the corresponding data at the points x
+    #     only works if fis a fonction
+    # z : where to evaluate the function
+
+    def localQuadraticInterpolation (a, H, x, f, zk):
+        # first find out in which interval we are. Easy here since equidistributed point
+        i = int( (zk-a)//H )
+
+        # if we are not in the interval, return constant function
+        if x[i] == zk :
+            return f(x[i])
+        
+        t = np.array([ x[i], (x[i]+x[i+1])/2, x[i+1]] )
+        
+        # use formula for local quadratic interpolation
+        return LagrangeInterpolation(t,f(t),zk)
+        
+    # careful, there are N+1 points 
+    H = (b-a)/N
+    if np.isclose(H,0) :
+        raise ValueError('PiecewiseQuadraticInterpolation : a and b are too close')
+    
+    x = np.linspace(a,b,N+1)
+    
+    # if f is a function, we need to evaluate it at alle the locations
+    # from types import FunctionType
+    # if ~isinstance(f, FunctionType) :
+    #    raise ValueError('PiecewiseQuadraticInterpolation : works only with a given function')
+        
+    # if we want to evaluate the interpolating function in a single point
+    if np.ndim(z) == 0 :
+        return localQuadraticInterpolation (b, H, x, f, z)
+    
+    # if we are here, it means that we need to evaluate the interpolation on many points.
+    # init result to zero, of same type and size as z
+    result = np.zeros_like(z) 
+     
+    # The following part is not optimized
+    for k in range(z.size) :
+
+        result[k] = localQuadraticInterpolation (a, H, x, f, z[k])
+   
+    return result
+
+
+   
+    
+