# Calculating coefficients of interpolating polynomial using Neville's algorithm

First of all, sorry for my bad math terminology as it's not my native language and I may misuse some terms in English.

I've been tasked with writing an application which calculates the general formula of interpolating polynomial (being given a set of $X_i$ and $Y_i$), however, I have to use Neville's algorithm for that. I know Lagrange's method can be used, but I've been told to use Neville's algorithm specifically.

So, my question is - is it even possible? If so, how should I calculate the coefficients of the interpolating polynomial using this algorithm $(a_0, a_1,\dots ,a_n)$?

I feel like I'm missing a point somewhere (maybe because I'm not that good at maths), but I can't find a solution to this.

I'll be honest with it - this is homework. But I'm not asking for a solved code, just need help with understanding this.

• If you need help with mathjax see here. Commented Oct 29, 2014 at 21:50
• Did you look at the Wikipedia page? en.wikipedia.org/wiki/Neville%27s_algorithm There's a nice explanation. Commented Oct 29, 2014 at 22:37
• @jflipp I did. But in this article it is only described how to calculate interpolated values. What I want to calculate are the values of polynomial coefficients only $(a_0, a_1,\dots ,a_n)$, to make it easier to code. Commented Oct 30, 2014 at 14:11

I think you can calculate the coefficients of the interpolation polynomial as follows.

The algorithm given on the Wikipedia page http://en.wikipedia.org/wiki/Neville%27s_algorithm works not only for the value of the interpolation polynomial at a given point x, but also for the complete interpolation polynomial (in a mathematical sense, regarding the interpolation polynomial as an element of the polynomial ring $\mathbb R[X]$). More specifically, the initial polynomials $p_{i,i}(X) = y_i$ are constant. From there, we can calculate the intermediate polynomials $$p_{i,j}(X) = ((x_j - X) p_{i,j-1}(X) + (X - x_i) p_{i+1,j}(X))/(x_ j - x_i)$$until we reach the complete interpolation polynomial.

Now, all you have to do is to translate this recursive scheme for polynomials into a form suitable for mainstream programming languages :-). You can do this by identifying a polynomial with its coefficient vector. We know that, given n + 1 interpolation points $(x_i, y_i)$, the interpolation polynomial has degree n, and more generally, the polynomial $p_{i,j}$ from the recursive scheme with $0 \leq i < j \leq n$ has degree $j - i \leq n$. So, all the polynomials we have to consider, and their corresponding coefficient vectors are of the form $$p(X) = a_0 + a_1 X + \cdots + a_n X^n <-> (a_0, \ldots, a_n).$$ And of course, in a programming language, a coefficient vector is just an array :-). With this correspondence, we can translate the recursive scheme from above into a recursive scheme for coefficient vectors.

First, the initial polynomial $p_{i,i}(X) = y_i$ has coefficient vector $(y_i, 0, \ldots, 0)$.

Next, we have to figure out the coefficient vector of $Xp_{i,j-1}(X)$ given the coefficient vector of $p_{i,j-1}(X)$. But that's easy. Multiplying a polynomial by $X$ shifts its coefficient vector to the right. So, if we denote the coefficient vector of the polynomial $p(X)$ by $V(p(X))$ and on coefficient vectors define the shift operator $$S((a_0, \ldots, a_n)) = (0, a_0, \ldots, a_{n - 1}),$$ then we have $$V(Xp(X)) = S(V(p(X)).$$Check that :-).

Now we can rewrite the recursion for polynomials from above as a recursion for their coefficient vectors and get $$V(p_{i,j}(X)) = (x_jV(p_{i,j-1}(X)) - S(V(p_{i,j-1}(X))) + S(V(p_{i+1,j}(X))) - x_iV(p_{i+1,j}(X))) / (x_j - x_i).$$

Now you have a recursive scheme for coefficient vectors of polynomials which takes you from the coefficient vectors of the initial, constant polynomials $p_{i,i}$ all the way to the coefficient vector of the complete interpolation polynomial.

• This seems to be very much what I've been looking for. Can you explain though how in $$V(p_{i,j}(X)) = (x_jV(p_{i,j-1}(X)) - S(V(p_{i,j-1}(X))) + S(V(p_{i+1,j}(X))) - x_iV(p_{i+1,j}(X))) / (x_j - x_i).$$ I should get $x_iV(p_{i+1,j}(X))$? For example, when calculating $V(p_{0,1}(X))$ it would be $x_iV(p_{1,1}(X))$ but it is like the next step I should do. Hope my mathjax is ok now ;) Commented Oct 31, 2014 at 18:49
• $x_i$ is a constant (as opposed to $X$, which is the variable in the polynomials). Multiplying a polynomial by a constant just multiplies its coefficient vector by that same constant. In other words, if $V(p_{i+1,j}(X)) = (a_0, \ldots, a_n)$, then $V(x_i p_{i+1,j}(X)) = x_i V(p_{i+1,j}(X)) = (x_i a_0, \ldots, x_i a_n).$ Commented Oct 31, 2014 at 23:20
• Either I'm really dumb or I miss the point. Or maybe both. I know it's much to ask, but could you calculate $V(p_{0,1}(X))$ with, for example, this set of points ${(1,0),(2,2),(4,12)}$? Commented Nov 1, 2014 at 8:31
• Ok, so we have $(x_0, y_0) = (1, 0)$, $(x_1, y_1) = (2, 2)$, $(x_2, y_2) = (4, 12)$. According to the recursion, in order to calculate $p_{0,1}(X)$, we need $p_{0,0}(X)$ and $p_{1,1}(X)$. We have $p_{0,0}(X) = y_0 = 0$ and $p_{1,1}(X) = y_1 = 2$. The coefficient vectors are $V(p_{0,0}(X)) = (0, 0, 0)$ and $V(p_{1,1}(X)) = (2, 0, 0)$. So, the recursion for the coefficient vectors gives $V(p_{0,1}(X)) = \left(2 (0,0,0) - S(0,0,0) + S(2,0,0) - 1 (2,0,0)\right) / (2 - 1) = \left( (0,0,0) - (0,0,0) + (0,2,0) - (2,0,0)\right) / 1 = (-2,2,0).$ Commented Nov 1, 2014 at 12:12
• Here, for example, the expression $-1(2,0,0)$ is the vector $(2,0,0)$ multiplied by the scalar $-1$ which is $-(2, 0, 0)$. All in all, we find $p_{0,1}(X) = V^{-1}((-2,2,0)) = -2 + 2X$. And this polynomial indeed interpolates the points (1,0) and (2,2), as expected. Commented Nov 1, 2014 at 12:17