{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Barycentric Lagrange interpolation on Chebyshev points"
]
},
{
"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\\limits_{i=0}^N{}' \\frac{(-1)^i}{x - x_i} f_i }{ \\sum\\limits_{i=0}^N{}' \\frac{(-1)^i}{x - x_i} }\n",
"$$\n",
"where the prime on the summation means that the first and last terms must be multiplied by a factor of half."
]
},
{
"cell_type": "code",
"execution_count": 1,
"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": [
"Define the function to be interpolated."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"def fun(x):\n",
" f = 1.0/(1.0+16.0*x**2)\n",
" return f"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The next function evaluates the Lagrange interpolation using Chebyshev points."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"def BaryInterp(X,Y,x):\n",
" nx = np.size(x)\n",
" nX = np.size(X)\n",
" f = 0*x\n",
" # Compute weights\n",
" w = (-1.0)**np.arange(0,nX)\n",
" w[0] = 0.5*w[0]\n",
" w[nX-1] = 0.5*w[nX-1]\n",
" # Evaluate barycentric foruma at x values\n",
" for i in range(nx):\n",
" num, den = 0.0, 0.0\n",
" for j in range(nX):\n",
" if np.abs(x[i]-X[j]) < 1.0e-15:\n",
" num = Y[j]\n",
" den = 1.0\n",
" break\n",
" else:\n",
" num += Y[j]*w[j]/((x[i]-X[j]))\n",
" den += w[j]/(x[i]-X[j])\n",
" f[i] = num/den\n",
" return f"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"xmin, xmax = -1.0, +1.0\n",
"N = 19 # degree of polynomial"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let us interpolate on Chebyshev points."
]
},
{
"cell_type": "code",
"execution_count": 5,
"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 = BaryInterp(X,Y,x)\n",
"fe = fun(x)\n",
"plt.figure(figsize=(8,5))\n",
"plt.plot(x,fe,'b--',x,fi,'r-',X,Y,'o')\n",
"plt.legend((\"True function\",\"Interpolation\",\"Data\"),loc='upper right')\n",
"plt.axis([-1.0,+1.0,0.0,1.1]);"
]
}
],
"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
}