logo

Differential Stages in Numerical Analysis 📂Numerical Analysis

Differential Stages in Numerical Analysis

Definition

A Divided Difference of a function f:RRf : \mathbb{R} \to \mathbb{R} for distinct x1,,xnx_{1} , \cdots , x_{n} is defined as follows:

f[x0]:=f(x0)f[x0,x1]:=f(x1)f(x0)x1x0f[x0,x1,x2]:=f[x1,x2]f[x0,x1]x2x0f[x0,,xn]:=f[x1,,xn]f[x0,,xn1]xnx0 \begin{align*} f[x_{0}] :=& f( x_{0} ) \\ f [ x_{0} , x_{1} ] :=& {{ f ( x_{1} ) - f ( x_{0} ) } \over { x_{1} - x_{0} }} \\ f [ x_{0} , x_{1} , x_{2} ] :=& {{ f [ x_{1} , x_{2} ] - f [ x_{0} , x_{1} ] } \over { x_{2} - x_{0} }} \\ f [ x_{0} , \cdots , x_{n} ] :=& {{ f [ x_{1} , \cdots , x_{n} ] - f [ x_{0} , \cdots , x_{n-1} ] } \over { x_{n} - x_{0} }} \end{align*}

Theorem

  • [1]: f[x0,,xn]f [ x_{0} , \cdots , x_{n} ] is always the same regardless of the order of x0,,xnx_{0} , \cdots , x_{n}.
  • [2]: f[x0,,xn]=i=0nf(xi)j{0,,n}{i}(xixj)f [ x_{0} , \cdots , x_{n} ] = \sum_{i=0}^{n} {{ f( x_{i} ) } \over { \prod_{j \in \left\{ 0 , \cdots , n \right\} \setminus \left\{ i \right\} } ( x_{i} - x_{j} ) }}
  • [3]: There exists ξ\xi that satisfies the following.f(ξ)=f[x0,x1]f’ ( \xi ) = f [ x_0 , x_{1} ]
  • [4]: There exists ξ\xi that satisfies the following.1n!f(n)(ξ)=f[x0,,xn]{{1} \over {n!}} f^{(n)} ( \xi ) =f [ x_{0} , \cdots , x_{n} ]
  • [3]’: f[xi,xi]=f(xi)f [ x_{i} , x_{i} ] = f ' ( x_{i} )
  • [4]’: f[xi,,xin+1]=1n!f(n)(xi)f [ \underbrace{ x_{i} , \cdots , x_{i} }_{ n+1 } ] = {{1} \over {n!}} f^{(n)} ( x_{i} )

  • H{a,b,c,}\mathscr{H} \left\{ a,b,c, \cdots \right\} denotes the smallest interval that includes a,b,c,a,b,c, \cdots.

Explanation

Divided differences greatly save space when expressing various theories of numerical analysis.

Theorem [2] is especially useful in actual calculations, since it is expressed as f[x0,x1,x2]=f(x0)(x0x1)(x0x2)+f(x1)(x1x2)(x1x0)+f(x2)(x2x0)(x2x1) f [ x_{0} , x_{1} , x_{2} ] = {{ f( x_{0} ) } \over { ( x_{0} - x_{1} ) ( x_{0} - x_{2} ) }} + {{ f( x_{1} ) } \over { ( x_{1} - x_{2} ) ( x_{1} - x_{0} ) }} + {{ f( x_{2} ) } \over { ( x_{2} - x_{0} ) ( x_{2} - x_{1} ) }} , allowing for iterative calculations to be simplified compared to the original definition.

[3] is essentially Mean Value Theorem, and if you have ever considered the concept of ‘instantaneous rate of change’ when differentiating, this will help you understand why divided differences can replace differentiation and why it’s called Divided Difference. Theorem [4] can be generalized in this manner. The proof can be thought of as the variables of divided differences becoming infinitely close, or taking a limit, which can be seen in [3]’, [4]’ as derivatives of degree nn. Such a conceptual approach overcomes the limitations arising from the definition of divided differences and enriches the theory of numerical analysis.

Example

Let’s calculate some divided differences for f(x):=x2+1f(x) := x^2 + 1 simply.

For one point, f[3]=f(3)=10 f[3] = f(3) = 10

For two points, f[0,5]=f(0)f(5)05=1265=5 f[0,5] = {{ f(0) - f(5) } \over { 0 - 5}} = {{1 - 26 } \over { -5 }} = 5

For three points, f[2,3,1]=f(2)f(3)23f(3)f(1)3(1)2(1)=510110243=523=1 \begin{align*} f[2,3,-1] =& {{ \displaystyle {{ f(2) - f( 3) } \over { 2 - 3 }} - {{ f(3) - f(-1) } \over { 3 - (-1) }} } \over { 2 - (-1) }} \\ =& {{ \displaystyle {{ 5 - 10 } \over { -1 }} - {{ 10 - 2 } \over { 4 }} } \over { 3}} \\ =& {{ 5 -2 } \over {3}} \\ =& 1 \end{align*} Divided differences can obviously be calculated for a finite vector, but usually, most applications go up to n=2n=2, and beyond that, it’s rarely used.

Implementation

Here is an implementation of divided differences in R code, which takes a function ff and a vector x=(x0,,xn)\mathbf{x} = ( x_{0} , \cdots , x_{n} ) to return f[x0,,xn]f [ x_{0} , \cdots , x_{n} ].

f<-function(x) {x^2 + 1}
dd<-function(f,X){
  if(length(X)==1) {return(f(X))}
  temp<-numeric(0)
  for(i in 1:length(X)) {temp[i]<-prod(X[i]-X[-i])}
  return(sum(f(X)/temp))
}
dd(f,c(3))
dd(f,c(0,5))
dd(f,c(2,3,-1))

The result of executing the above code is as follows, and it’s confirmed to exactly match the hand calculation for f(x)=x2+1f(x) = x^2 +1.

20190319\_161740.png

Proof

[2]

Strategy: Compare the degrees of Lagrange’s formula and Newton’s divided difference formula.


Ψn(x):=(xx0)(xxn) \Psi_{n} (x) := (x-x_{0}) \cdots (x - x_{n}) Differentiating Ψn\Psi_{n} with respect to xx and substituting x=xjx = x_{j} gives Ψn(xj):=(xx0)(xxj1)(xxj+1)(xxn) \Psi’_{n} (x_{j}) := (x-x_{0}) \cdots (x - x_{j-1})(x - x_{j+1}) \cdots (x - x_{n})

    Ψ(x)n(xxj)Ψn(xj)(xx0)(xxj1)(xxj+1)(xxn)(xjx0)(xjxj1)(xjxj+1)(xjxn)=lj(x) \implies {{ \Psi (x)_{n} } \over { (x - x_{j}) \Psi’_{n}(x_{j}) }} \equiv {{ (x - x_{0} ) \cdots (x - x_{j-1} ) (x - x_{j+1} ) \cdots (x - x_{n} ) } \over { (x_{j} - x_{0} ) \cdots (x_{j} - x_{j-1} ) (x_{j} - x_{j+1} ) \cdots (x_{j} - x_{n} ) }} = l_{j} (x) According to Lagrange’s formula, the polynomial interpolation of ff is pn(x)=j=0nΨn(x)(xxj)Ψn(xj)f(xj) p_{n} (x) = \sum_{j=0}^{n} {{ \Psi_{n} (x) } \over { (x - x_{j}) \Psi’_{n} (x_{j}) }} \cdot f( x_{j}) Therefore, since the coefficient of the highest order term of pnp_{n} is j=0nf(xj)Ψn(xj)\displaystyle \sum_{j=0}^{n} {{ f( x_{j}) } \over { \Psi’_{n} (x_{j}) }}, and in Newton’s divided difference formula, the coefficient of the highest order term of pnp_{n} is f[x0,,xn]f [x_{0} , \cdots , x_{n} ], f[x0,,xn]=i=0nf(xi)j{0,,n}{i}(xixj) f [ x_{0} , \cdots , x_{n} ] = \sum_{i=0}^{n} {{ f( x_{i} ) } \over { \displaystyle \prod_{j \in \left\{ 0 , \cdots , n \right\} \setminus \left\{ i \right\} } ( x_{i} - x_{j} ) }}

[1]

By Theorem [2], changing the order within [x0,,xn][x_{0} , \cdots , x_{n}] is just changing the order of addition when calculating i=0nf(xi)j{0,,n}{i}(xixj)\displaystyle \sum_{i=0}^{n} {{ f( x_{i} ) } \over { \displaystyle \prod_{j \in \left\{ 0 , \cdots , n \right\} \setminus \left\{ i \right\} } ( x_{i} - x_{j} ) }}, so it’s always the same.

[4]

Strategy: Compare the degrees of Lagrange’s formula and Newton’s divided difference formula.


If we say pn(x)=pn1(x)+C(x)p_{n}(x) = p_{n-1}(x) + C(x), then degC=n\deg C = n and for i=0,,(n1)i= 0, \cdots, (n-1), pn(xi)=f(xi)=pn1(xi)p_{n}(x_{i}) = f(x_{i}) = p_{n-1}(x_{i}), so for some constant an0a_{n} \ne 0 we have C(x)=an(xx0)(xxn1)C (x) = a_{n} (x - x_{0}) \cdots (x - x_{n-1}). Substituting C(x)C(x) with x=xnx=x_{n} gives C(xn)=pn(xn)pn1(xn) C(x_{n} ) = p_{n} (x_{n}) - p_{n-1}(x_{n})

C(xn)=an(xnx0)(xnxn) C (x_{n}) = a_{n} (x_{n} - x_{0}) \cdots (x_{n} - x_{n}) Since pnp_{n} is the polynomial interpolation of ff, it becomes pn(xn)=f(xn)p_{n} (x_{n}) = f(x_{n}), an(xnx0)(xnxn1)=f(xn)pn1(xn) a_{n} (x_{n} - x_{0}) \cdots (x_{n} - x_{n-1}) = f (x_{n}) - p_{n-1}(x_{n})

an=f(xn)pn1(xn1)(xnx0)(xnxn1) a_{n} = {{ f (x_{n}) - p_{n-1}(x_{n-1}) } \over { (x_{n} - x_{0}) \cdots (x_{n} - x_{n-1}) }}

Error between polynomial interpolation and actual function: For some ξH{x0,,xn1}\xi \in \mathscr{H} \left\{ x_{0} , \cdots , x_{n-1} \right\}, the polynomial interpolation pn1p_{n-1} of ff for nn-times differentiable f:RRf : \mathbb{R} \to \mathbb{R} is for some xnRx_{n} \in \mathbb{R} f(xn)pn1(xn)=(xnx0)(xnxn1)n!f(n)(ξ) f(x_{n}) - p_{n-1} (x_{n}) = {{ (x_{n} - x_{0}) \cdots (x_{n} - x_{n-1}) } \over { n! }} f^{(n)} ( \xi )

According to the formula of Error between polynomial interpolation and actual function, an=f(n)(ξ)n! a_{n} = {{ f^{(n)} ( \xi ) } \over { n! }} Meanwhile, since the coefficient of the highest order term of pnp_{n} in Newton’s divided difference formula is an=f[x0,,xn]a_{n} = f [x_{0} , \cdots , x_{n} ], 1n!f(n)(ξ)=f[x0,,xn] {{1} \over {n!}} f^{(n)} ( \xi ) =f [ x_{0} , \cdots , x_{n} ]

[3]

It’s trivial by Theorem [4].