{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Barycentric Lagrange interpolation"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let us interpolate the Runge function on $[-1,+1]$\n",
"$$\n",
"f(x) = \\frac{1}{1 + 16 x^2}\n",
"$$\n",
"using the Barycentric formula\n",
"$$\n",
"p(x) = \\frac{ \\sum_{i=0}^N \\frac{w_i}{x - x_i} f_i }{ \\sum_{i=0}^N \\frac{w_i}{x - x_i} }\n",
"$$\n",
"where the weight $w_i$ is given by\n",
"$$\n",
"w_i = \\frac{1}{\\prod_{j=0,j\\ne i}^N (x_i - x_j)}\n",
"$$"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib inline\n",
"%config InlineBackend.figure_format = 'svg'\n",
"import numpy as np\n",
"from matplotlib import pyplot as plt"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This is the function we want to interpolate."
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"def fun(x):\n",
" f = 1.0/(1.0+16.0*x**2)\n",
" return f"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We next write a function that constructs and evaluates the Lagrange interpolation."
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [],
"source": [
"# X[nX], Y[nX] : Data points\n",
"# x[nx] : points where we want to evaluate\n",
"def LagrangeInterpolation(X,Y,x):\n",
" nx = np.size(x)\n",
" nX = np.size(X)\n",
" # compute the weights\n",
" w = np.ones(nX)\n",
" for i in range(nX):\n",
" for j in range(nX):\n",
" if i != j:\n",
" w[i] = w[i]/(X[i]-X[j])\n",
" # Evaluate the polynomial at x\n",
" num= np.zeros(nx)\n",
" den= np.zeros(nx)\n",
" eps=1.0e-14\n",
" for i in range(nX):\n",
" num = num + Y[i]*w[i]/((x-X[i])+eps)\n",
" den = den + w[i]/((x-X[i])+eps)\n",
" f = num/den\n",
" return f"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [],
"source": [
"xmin, xmax = -1.0, +1.0\n",
"N = 15 # Degree of polynomial"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We first interpolate on uniformly spaced points."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n",
"\n",
"\n"
],
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"X = np.linspace(xmin,xmax,N+1)\n",
"Y = fun(X)\n",
"x = np.linspace(xmin,xmax,100)\n",
"fi = LagrangeInterpolation(X,Y,x)\n",
"fe = fun(x)\n",
"plt.plot(x,fe,'b--',x,fi,'r-',X,Y,'o')\n",
"plt.title('Degree '+str(N)+' using uniform points')\n",
"plt.legend((\"True function\",\"Interpolation\",\"Data\"),loc='lower center');"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, we interpolate on Chebyshev points."
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n",
"\n",
"\n"
],
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"X = np.cos(np.linspace(0.0,np.pi,N+1))\n",
"Y = fun(X)\n",
"x = np.linspace(xmin,xmax,100)\n",
"fi = LagrangeInterpolation(X,Y,x)\n",
"fe = fun(x)\n",
"plt.plot(x,fe,'b--',x,fi,'r-',X,Y,'o')\n",
"plt.title('Degree '+str(N)+' using Chebyshev points')\n",
"plt.legend((\"True function\",\"Interpolation\",\"Data\"),loc='upper right');"
]
}
],
"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.6.6"
}
},
"nbformat": 4,
"nbformat_minor": 1
}