Polynomial Interpolation

The Interpolation Problem

Given \(n+1\) distinct data points \((x_0,y_0),(x_1,y_1),\ldots,(x_n,y_n)\), the polynomial interpolation problem asks for a polynomial \(p(x)\) of degree at most \(n\) that passes exactly through every point:

The nodes \(x_i\) must be distinct (no repeated \(x\)-values). A foundational theorem guarantees that such a polynomial exists and is unique: it is the unique element of the \((n+1)\)-dimensional space \(\mathcal{P}_n\) of polynomials of degree \(\le n\) satisfying the \(n+1\) interpolation conditions. Different methods — Monomial, Lagrange, Newton, Hermite, splines — are simply different ways of representing or computing the same object (or, in the case of Hermite and splines, solving a related but distinct problem).

Monomial Basis (Vandermonde System)

Write \(p(x)=a_0+a_1 x+a_2 x^2+\cdots+a_n x^n\). Imposing the interpolation conditions \(p(x_i)=y_i\) gives the linear system \(Va=y\), where \(V\) is the Vandermonde matrix \(V_{ij}=x_i^j\). The system has a unique solution whenever all nodes are distinct (the Vandermonde determinant is \(\prod_{i>j}(x_i-x_j)\neq0\)). Although conceptually direct, the Vandermonde system is notoriously ill-conditioned for large \(n\) or widely spaced nodes: small perturbations in the data can produce large changes in the coefficients. The condition number grows exponentially in \(n\). Use this method for small \(n\) (say \(n\le6\)) where exactness and simplicity matter more than numerical stability.

Lagrange Form

The Lagrange representation avoids solving any linear system. Define the cardinal basis polynomials

Each \(\ell_i\) satisfies \(\ell_i(x_j)=\delta_{ij}\) (Kronecker delta), so \(p(x)=\sum_{i=0}^n y_i\,\ell_i(x)\) interpolates by construction. The Lagrange form is explicit and elegant, ideal for theoretical analysis and for adding a single new node (though adding a node requires recomputing all \(\ell_i\) from scratch). Each \(\ell_i\) is a degree-\(n\) polynomial involving a product of \(n\) linear factors.

Newton Divided Differences

The Newton form expresses \(p\) in the Newton basis \(\{1,\,(x-x_0),\,(x-x_0)(x-x_1),\,\ldots\}\):

The coefficients are divided differences, defined recursively: \(f[x_i]=y_i\) and \(f[x_i,\ldots,x_k]=(f[x_{i+1},\ldots,x_k]-f[x_i,\ldots,x_{k-1}])/(x_k-x_i)\). The divided difference table (shown highlighted in this tool) organises these computations in a triangular array; the first row of diagonal entries (marked in maroon) are precisely the Newton coefficients. A key advantage: adding a new node \(x_{n+1}\) requires only one new column appended to the existing table, making Newton’s form incremental and numerically efficient.

Hermite Interpolation

Hermite interpolation matches both function values \(f(x_i)=y_i\) and derivative values \(f'(x_i)=y'_i\) at each node, producing a polynomial of degree at most \(2n+1\) through \(n+1\) points. The construction uses a doubled-node divided difference table: each \(x_i\) appears twice in the sequence \(z_0\le z_1\le\cdots\le z_{2n+1}\), and coincident entries \(f[z_{2i},z_{2i+1}]\) are defined to equal \(y'_i\) rather than the usual divided difference (which would be \(0/0\)). The result is a smoother interpolant that respects the local slope at each data point — particularly valuable when fitting physical measurements where velocities or gradients are known alongside positions.

Cubic Splines — Natural and Clamped

For large data sets, a single high-degree polynomial oscillates wildly (Runge’s phenomenon, below). Cubic splines solve this by fitting a separate cubic polynomial \(S_i(x)\) on each subinterval \([x_i,x_{i+1}]\), subject to continuity conditions at the internal knots:

These \(C^2\) continuity conditions, combined with the interpolation requirements, leave exactly 2 free parameters — the boundary conditions. Two standard choices:

• Natural spline: \(S''(x_0)=S''(x_n)=0\). The spline is “free” at the ends — physically, a thin elastic beam supported but not clamped at the endpoints. This minimises \(\int(S'')^2\,dx\) among all interpolating functions, giving the minimum-curvature interpolant.
• Clamped spline: \(S'(x_0)=f'_0,\;S'(x_n)=f'_n\). The endpoint slopes are prescribed (from known derivatives or estimates). Clamped splines are generally more accurate than natural splines when derivative information is available.

Both variants reduce to solving a tridiagonal linear system for the second-derivative values \(M_i=S''(x_i)\), which can be solved in \(O(n)\) operations.

Runge’s Phenomenon and Chebyshev Nodes

The canonical cautionary example is \(f(x)=1/(1+25x^2)\) on \([-1,1]\) (Runge, 1901). With \(n+1\) equispaced nodes, the interpolation error \(\max|f(x)-p_n(x)|\to\infty\) as \(n\to\infty\), driven by wild oscillations near the endpoints despite \(f\) being perfectly smooth. The root cause is the Lebesgue function \(\Lambda_n(x)=\sum|\ell_i(x)|\): for equispaced nodes it grows exponentially near the boundary.

Chebyshev nodes are the minimisers of \(\max_{x\in[-1,1]}\prod_{i=0}^n|x-x_i|\):

They cluster near \(\pm1\), counteracting the endpoint oscillations. The Chebyshev interpolation error satisfies \(\max|f-p_n|\le\|f^{(n+1)}\|_\infty/(2^n(n+1)!)\), and for analytic \(f\) the error decreases geometrically in \(n\). Select Chebyshev in the node type to see the contrast with equispaced nodes on the Runge example.

Interpolation Error

For a function \(f\in C^{n+1}[a,b]\), the interpolation error at \(x\) is:

for some \(\xi\) between \(\min(x,x_i)\) and \(\max(x,x_i)\). The error depends on two factors: how large the \((n+1)\)th derivative of \(f\) is, and how small the nodal polynomial \(\omega_{n+1}(x)=\prod(x-x_i)\) can be made by the choice of nodes. Chebyshev nodes minimise the sup-norm of \(\omega_{n+1}\) over all choices of \(n+1\) nodes in \([a,b]\).

How to Use This Tool

• Choose a method using the six buttons at the top. The method determines how the polynomial is represented and displayed.
• Custom points: Enter \(x\) and \(y\) values as comma-separated lists. Expressions like pi/6, sin(pi/3), or exp(2) are accepted. For Hermite, also supply \(f'(x_i)\). For Clamped Spline, supply the two endpoint derivatives.
• From function: Enter a function expression, a domain \([a,b]\), and the number of interpolation nodes. Choose Equispaced or Chebyshev nodes. For Hermite and Clamped Spline in function mode, derivatives are computed automatically via symbolic differentiation (falling back to central differences if the symbolic attempt fails).
• Display options appear after computing: toggle the polynomial expression, the divided differences table, the Vandermonde matrix system, or the interpolation points on and off independently.
• Evaluate: Query the computed interpolant at any \(x\) values by entering them in the Evaluate field. When a function is known, absolute errors \(|f(x)-p(x)|\) are displayed alongside the interpolated values.
• Download: Save the graph as a PNG file using the Download button.