Divided differences

Given n + 1 data points

$${\displaystyle (x_{0},y_{0}),\ldots ,(x_{n},y_{n})}$$

where the ${\displaystyle x_{k}}$ are assumed to be pairwise distinct, the forward divided differences are defined as:

$${\displaystyle {\begin{aligned}{\mathopen {[}}y_{k}]&:=y_{k},&&k\in \{0,\ldots ,n\}\\{\mathopen {[}}y_{k},\ldots ,y_{k+j}]&:={\frac {[y_{k+1},\ldots ,y_{k+j}]-[y_{k},\ldots ,y_{k+j-1}]}{x_{k+j}-x_{k}}},&&k\in \{0,\ldots ,n-j\},\ j\in \{1,\ldots ,n\}.\end{aligned}}}$$

To make the recursive process of computation clearer, the divided differences can be put in tabular form, where the columns correspond to the value of j above, and each entry in the table is computed from the difference of the entries to its immediate lower left and to its immediate upper left, divided by a difference of corresponding x-values:

$${\displaystyle {\begin{matrix}x_{0}&y_{0}=[y_{0}]&&&\\&&[y_{0},y_{1}]&&\\x_{1}&y_{1}=[y_{1}]&&[y_{0},y_{1},y_{2}]&\\&&[y_{1},y_{2}]&&[y_{0},y_{1},y_{2},y_{3}]\\x_{2}&y_{2}=[y_{2}]&&[y_{1},y_{2},y_{3}]&\\&&[y_{2},y_{3}]&&\\x_{3}&y_{3}=[y_{3}]&&&\\\end{matrix}}}$$

Notation

Note that the divided difference ${\displaystyle [y_{k},\ldots ,y_{k+j}]}$ depends on the values ${\displaystyle x_{k},\ldots ,x_{k+j}}$ and ${\displaystyle y_{k},\ldots ,y_{k+j}}$, but the notation hides the dependency on the x-values. If the data points are given by a function f,

$${\displaystyle (x_{0},y_{0}),\ldots ,(x_{k},y_{n})=(x_{0},f(x_{0})),\ldots ,(x_{n},f(x_{n}))}$$

one sometimes writes the divided difference in the notation

$${\displaystyle f[x_{k},\ldots ,x_{k+j}]\ {\stackrel {\text{def}}{=}}\ [f(x_{k}),\ldots ,f(x_{k+j})]=[y_{k},\ldots ,y_{k+j}].}$$

Other notations for the divided difference of the function ƒ on the nodes x₀, ..., x_n are:

$${\displaystyle f[x_{k},\ldots ,x_{k+j}]={\mathopen {[}}x_{0},\ldots ,x_{n}]f={\mathopen {[}}x_{0},\ldots ,x_{n};f]=D[x_{0},\ldots ,x_{n}]f.}$$

Divided differences for ${\displaystyle k=0}$ and the first few values of ${\displaystyle j}$:

$${\displaystyle {\begin{aligned}{\mathopen {[}}y_{0}]&=y_{0}\\{\mathopen {[}}y_{0},y_{1}]&={\frac {y_{1}-y_{0}}{x_{1}-x_{0}}}\\{\mathopen {[}}y_{0},y_{1},y_{2}]&={\frac {{\mathopen {[}}y_{1},y_{2}]-{\mathopen {[}}y_{0},y_{1}]}{x_{2}-x_{0}}}={\frac {{\frac {y_{2}-y_{1}}{x_{2}-x_{1}}}-{\frac {y_{1}-y_{0}}{x_{1}-x_{0}}}}{x_{2}-x_{0}}}={\frac {y_{2}-y_{1}}{(x_{2}-x_{1})(x_{2}-x_{0})}}-{\frac {y_{1}-y_{0}}{(x_{1}-x_{0})(x_{2}-x_{0})}}\\{\mathopen {[}}y_{0},y_{1},y_{2},y_{3}]&={\frac {{\mathopen {[}}y_{1},y_{2},y_{3}]-{\mathopen {[}}y_{0},y_{1},y_{2}]}{x_{3}-x_{0}}}\end{aligned}}}$$

• Linearity
  $${\displaystyle {\begin{aligned}(f+g)[x_{0},\dots ,x_{n}]&=f[x_{0},\dots ,x_{n}]+g[x_{0},\dots ,x_{n}]\\(\lambda \cdot f)[x_{0},\dots ,x_{n}]&=\lambda \cdot f[x_{0},\dots ,x_{n}]\end{aligned}}}$$

  
• Leibniz rule
  $${\displaystyle (f\cdot g)[x_{0},\dots ,x_{n}]=f[x_{0}]\cdot g[x_{0},\dots ,x_{n}]+f[x_{0},x_{1}]\cdot g[x_{1},\dots ,x_{n}]+\dots +f[x_{0},\dots ,x_{n}]\cdot g[x_{n}]=\sum _{r=0}^{n}f[x_{0},\ldots ,x_{r}]\cdot g[x_{r},\ldots ,x_{n}]}$$

  
• Divided differences are symmetric: If ${\displaystyle \sigma$ :\{0,\dots ,n\}\to \{0,\dots ,n\}} is a permutation then
  $${\displaystyle f[x_{0},\dots ,x_{n}]=f[x_{\sigma (0)},\dots ,x_{\sigma (n)}]}$$

  
• Polynomial interpolation in the Newton form: if ${\displaystyle P}$ is a polynomial function of degree ${\displaystyle \leq n}$, and ${\displaystyle p[x_{0},\dots ,x_{n}]}$ is the divided difference, then
  $${\displaystyle P_{n-1}(x)=p[x_{0}]+p[x_{0},x_{1}](x-x_{0})+p[x_{0},x_{1},x_{2}](x-x_{0})(x-x_{1})+\cdots +p[x_{0},\ldots ,x_{n}](x-x_{0})(x-x_{1})\cdots (x-x_{n-1})}$$

  
• If ${\displaystyle p}$ is a polynomial function of degree ${\displaystyle <n}$, then
  $${\displaystyle p[x_{0},\dots ,x_{n}]=0.}$$

  
• Mean value theorem for divided differences: if ${\displaystyle f}$ is n times differentiable, then
  $${\displaystyle f[x_{0},\dots ,x_{n}]={\frac {f^{(n)}(\xi )}{n!}}}$$

  
  for a number ${\displaystyle \xi }$ in the open interval determined by the smallest and largest of the ${\displaystyle x_{k}}$'s.

The divided difference scheme can be put into an upper triangular matrix:

$${\displaystyle T_{f}(x_{0},\dots ,x_{n})={\begin{pmatrix}f[x_{0}]&f[x_{0},x_{1}]&f[x_{0},x_{1},x_{2}]&\ldots &f[x_{0},\dots ,x_{n}]\\0&f[x_{1}]&f[x_{1},x_{2}]&\ldots &f[x_{1},\dots ,x_{n}]\\0&0&f[x_{2}]&\ldots &f[x_{2},\dots ,x_{n}]\\\vdots &\vdots &&\ddots &\vdots \\0&0&0&\ldots &f[x_{n}]\end{pmatrix}}.}$$

Then it holds

• ${\displaystyle T_{f+g}(x)=T_{f}(x)+T_{g}(x)}$
• ${\displaystyle T_{\lambda f}(x)=\lambda T_{f}(x)}$ if ${\displaystyle \lambda }$ is a scalar
• ${\displaystyle T_{f\cdot g}(x)=T_{f}(x)\cdot T_{g}(x)}$
  This follows from the Leibniz rule. It means that multiplication of such matrices is commutative. Summarised, the matrices of divided difference schemes with respect to the same set of nodes x form a commutative ring.
• Since ${\displaystyle T_{f}(x)}$ is a triangular matrix, its eigenvalues are obviously ${\displaystyle f(x_{0}),\dots ,f(x_{n})}$.
• Let ${\displaystyle \delta _{\xi }}$ be a Kronecker delta-like function, that is
  $${\displaystyle \delta _{\xi }(t)={\begin{cases}1&:t=\xi ,\\0&:{\mbox{else}}.\end{cases}}}$$

  
  Obviously ${\displaystyle f\cdot \delta _{\xi }=f(\xi )\cdot \delta _{\xi }}$, thus ${\displaystyle \delta _{\xi }}$ is an eigenfunction of the pointwise function multiplication. That is ${\displaystyle T_{\delta _{x_{i}}}(x)}$ is somehow an "eigenmatrix" of ${\displaystyle T_{f}(x)}$: ${\displaystyle T_{f}(x)\cdot T_{\delta _{x_{i}}}(x)=f(x_{i})\cdot T_{\delta _{x_{i}}}(x)}$. However, all columns of ${\displaystyle T_{\delta _{x_{i}}}(x)}$ are multiples of each other, the matrix rank of ${\displaystyle T_{\delta _{x_{i}}}(x)}$ is 1. So you can compose the matrix of all eigenvectors of ${\displaystyle T_{f}(x)}$ from the ${\displaystyle i}$-th column of each ${\displaystyle T_{\delta _{x_{i}}}(x)}$. Denote the matrix of eigenvectors with ${\displaystyle U(x)}$. Example
  $${\displaystyle U(x_{0},x_{1},x_{2},x_{3})={\begin{pmatrix}1&{\frac {1}{(x_{1}-x_{0})}}&{\frac {1}{(x_{2}-x_{0})(x_{2}-x_{1})}}&{\frac {1}{(x_{3}-x_{0})(x_{3}-x_{1})(x_{3}-x_{2})}}\\0&1&{\frac {1}{(x_{2}-x_{1})}}&{\frac {1}{(x_{3}-x_{1})(x_{3}-x_{2})}}\\0&0&1&{\frac {1}{(x_{3}-x_{2})}}\\0&0&0&1\end{pmatrix}}}$$

  
  The diagonalization of ${\displaystyle T_{f}(x)}$ can be written as
  $${\displaystyle U(x)\cdot \operatorname {diag} (f(x_{0}),\dots ,f(x_{n}))=T_{f}(x)\cdot U(x).}$$

  

Polynomials and power series

The matrix

$${\displaystyle J={\begin{pmatrix}x_{0}&1&0&0&\cdots &0\\0&x_{1}&1&0&\cdots &0\\0&0&x_{2}&1&&0\\\vdots &\vdots &&\ddots &\ddots &\\0&0&0&0&\;\ddots &1\\0&0&0&0&&x_{n}\end{pmatrix}}}$$

contains the divided difference scheme for the identity function with respect to the nodes ${\displaystyle x_{0},\dots ,x_{n}}$, thus ${\displaystyle J^{m}}$ contains the divided differences for the power function with exponent ${\displaystyle m}$. Consequently, you can obtain the divided differences for a polynomial function ${\displaystyle p}$ by applying ${\displaystyle p}$ to the matrix ${\displaystyle J}$: If

$${\displaystyle p(\xi )=a_{0}+a_{1}\cdot \xi +\dots +a_{m}\cdot \xi ^{m}}$$

and

$${\displaystyle p(J)=a_{0}+a_{1}\cdot J+\dots +a_{m}\cdot J^{m}}$$

then

$${\displaystyle T_{p}(x)=p(J).}$$

This is known as Opitz' formula.[2][3]

Now consider increasing the degree of ${\displaystyle p}$ to infinity, i.e. turn the Taylor polynomial into a Taylor series. Let ${\displaystyle f}$ be a function which corresponds to a power series. You can compute the divided difference scheme for ${\displaystyle f}$ by applying the corresponding matrix series to ${\displaystyle J}$: If

$${\displaystyle f(\xi )=\sum _{k=0}^{\infty }a_{k}\xi ^{k}}$$

and

$${\displaystyle f(J)=\sum _{k=0}^{\infty }a_{k}J^{k}}$$

then

$${\displaystyle T_{f}(x)=f(J).}$$

• Difference quotient
• Neville's algorithm
• Polynomial interpolation
• Mean value theorem for divided differences
• Nörlund–Rice integral
• Pascal's triangle

• Louis Melville Milne-Thomson (2000) [1933]. The Calculus of Finite Differences. American Mathematical Soc. Chapter 1: Divided Differences. ISBN 978-0-8218-2107-7.
• Myron B. Allen; Eli L. Isaacson (1998). Numerical Analysis for Applied Science. John Wiley & Sons. Appendix A. ISBN 978-1-118-03027-1.
• Ron Goldman (2002). Pyramid Algorithms: A Dynamic Programming Approach to Curves and Surfaces for Geometric Modeling. Morgan Kaufmann. Chapter 4:Newton Interpolation and Difference Triangles. ISBN 978-0-08-051547-2.

• An implementation in Haskell.