Solution (source code)

= Solution

For distinct nodes, define the divided difference by
$$
f[x_0,\ldots,x_k]
=\sum_{j=0}^k
\frac{f(x_j)}{\prod_{\substack{0\leq i\leq k\\i\ne j}}(x_j-x_i)}.
$$
This is the coefficient of $x^k$ in the <Lagrange interpolation polynomial> through the first $k+1$ data points.

Let $p_k$ denote that interpolating <polynomial>. The difference $p_k-p_{k-1}$ vanishes at $x_0,\ldots,x_{k-1}$, so
$$
p_k(x)-p_{k-1}(x)=c_k\prod_{i=0}^{k-1}(x-x_i).
$$
The Lagrange formula shows that the leading coefficient of $p_k$ is $f[x_0,\ldots,x_k]$, whereas $p_{k-1}$ has degree at most $k-1$. Hence $c_k=f[x_0,\ldots,x_k]$. Starting with $p_0=f(x_0)$ and iterating gives the <Newton interpolation polynomial>
$$
\boxed{
p_n(x)=f(x_0)+\sum_{k=1}^n f[x_0,\ldots,x_k]
\prod_{i=0}^{k-1}(x-x_i)}.
$$

The divided-difference recurrence is
$$
\boxed{
f[x_0,\ldots,x_k]
=\frac{f[x_1,\ldots,x_k]-f[x_0,\ldots,x_{k-1}]}
{x_k-x_0}}.
$$
For three nodes, the triangular table is
$$
\begin{array}{ccc}
f[x_0]&\longrightarrow&f[x_0,x_1]\longrightarrow f[x_0,x_1,x_2]\\
f[x_1]&\longrightarrow&f[x_1,x_2]\\
f[x_2]
\end{array}
$$
where each entry in a new column uses the two adjacent entries to its left. There are $n$ first differences, $n-1$ second differences, and so on, each requiring one division. The exact total is
$$
\boxed{n+(n-1)+\cdots+1=\frac{n(n+1)}2}.
$$