{ "metadata": { "name": "", "signature": "sha256:876f7be49420e51f271caf228705f59460c38237302a77c9d53b69349f4c3f79" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "heading", "level": 1, "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", "collapsed": true, "input": [ "import numpy as np\n", "from matplotlib import pyplot as plt" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 32 }, { "cell_type": "code", "collapsed": true, "input": [ "def fun(x):\n", " f = 1.0/(1.0+16.0*x**2)\n", " return f" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 33 }, { "cell_type": "code", "collapsed": true, "input": [ "def BaryInterp(X,Y,x):\n", " nx = np.size(x)\n", " nX = np.size(X)\n", " f = 0*x\n", " w = (-1.0)**arange(0,nX)\n", " w[0] = 0.5*w[0]\n", " w[nX-1] = 0.5*w[nX-1]\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" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 34 }, { "cell_type": "code", "collapsed": false, "input": [ "xmin, xmax = -1.0, +1.0\n", "N = 20" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 35 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let us interpolate on Chebyshev points." ] }, { "cell_type": "code", "collapsed": true, "input": [ "X = cos(np.linspace(0.0,pi,N))\n", "Y = fun(X)\n", "x = np.linspace(xmin,xmax,100)\n", "fi = BaryInterp(X,Y,x)\n", "fe = fun(x)\n", "plt.plot(x,fe,'b--',x,fi,'r-',X,Y,'o')\n", "plt.legend((\"True function\",\"Interpolation\",\"Data\"),'upper right')\n", "plt.axis([-1.0,+1.0,0.0,1.1])" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 36, "text": [ "[-1.0, 1.0, 0.0, 1.1]" ] }, { "metadata": {}, "output_type": "display_data", "svg": [ "\n", "\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n" ], "text": [ "" ] } ], "prompt_number": 36 } ], "metadata": {} } ] }