If the first derivatives of the function are known as well as the function value at each of the node points , i.e., we have available a set of values , then the function can be interpolated by a polynomial of degree :

(46) |

(47) |

In practice, the *Hermite interpolation* can be used in such a
case. First, we assume , and represent the Hermite polynomial
as a linear combination of
basis polynomials
(
) of degree :

(48) |

(49) |

(50) |

(51) |

and | (52) |

(53) |

(54) |

**Example**

The Hermite interpolation is carried out to the same function used in previous examples, with the result shown in the figure below, together with the basis polynomials . As the first order derivative is available as well as the function value at each node point , the interpolation matches the given function very well (almost identical on the plots), with an error , which is much reduced from of all methods previously discussed based only on .

The Matlab code that implements the Hermite interpolation method is listed below.

function [H a b]=HIL(u,x,y,dy) % Hermite interpolation (Lagrange) % u: discrete data points; % vector x: [x_1,...,x_n] % vector y: [y_1,...,y_n] % vector dy: [y'_1,...,y'_n] n=length(x); % number of interpolating points k=length(u); % number of discrete data points li=ones(n,k); % Lagrange basis polynomials a=zeros(n,k); % basis polynomials alpha(x) b=zeros(n,k); % basis polynomials beta(x) H=zeros(1,k); % Hermie interpolation polynomial H(x) for i=1:n dl=0; % derivative of Lagrange basis for j=1:n if j~=i dl=dl+1/(x(i)-x(j)); li(i,:)=li(i,:).*(u-x(j))/(x(i)-x(j)); end end l2=li(i,:).^2; b(i,:)=(u-x(i)).*l2; % basis polynomial alpha(x) a(i,:)=(1-2*(u-x(i))*dl).*l2; % basis polynomial beta(x) H=H+a(i,:)*y(i)+b(i,:)*dy(i); % Hermite polynomial H(x) end end

We next consider the case when . Now we have known values at each of the node points . Such a dataset can be represented by another variable that takes the values of each repeated times, as shown in the table:

(55) |

This function with
sample points can then be
interpolated by the Newton polynomial method. For example, if
, then the Newton's polynomial of degree
can be found to be:

(56) |

It can be verified that indeed for all and . Similar to the Newton polynomial method discussed previously, the divided difference coefficients can be obtained recursively, with the only difference that there exist repeated copies at each point , where the divided difference can be found by

(57) |

(58) |

(59) |

**Example:**

The Hermite interpolation based Newton's polynomials is again carried out to the same function used before. Now we assume both the first and second order derivatives and are available as well as at the points. The resulting Hermite interpolation is plotted together with in the figure below. We see that they are almost identical, with an error .

The Matlab code that implements this algorithm is listed below.

function [v]=HIN(u,x,dy) % Hermite interpolation (Newton) % u: discrete data points; % vector x: [x_1,...,x_n] % matrix dy contains m derivatives at each of the n points [n m]=size(dy); k=length(u); % number of discrete data points v=zeros(1,k); % interpolation results dd=DividedDifference2(x,dy); % get the divided difference array w=ones(1,k); for i=1:n p=u-x(i); for j=1:m l=(i-1)*m+j; % index of the coefficient v=v+dd(l,l).*w; % which is on the diagnal of array dd w=w.*p; end end end function dd=DividedDifference2(x,dy) % generate array of divided differences [n m]=size(dy); % n data points, m derivatives (0 to m-1) dd=zeros(n*m); % matrix of divided differences z=zeros(1,n*m); k=1; for i=1:n % n data points for j=1:m % m derivatives (0 to m-1) at each point k=(i-1)*m+j; % row index z(k)=x(i); dd(k,1)=dy(i,1); % 0th divided difference in first column fprintf('%6.3f\t%6.3f\t',z(k),dd(k,1)); for l=2:k % column index for the remaining columns %fprintf('(%f %f)\n',dd(k,l-1),dd(k-1,l-1)); if dd(k,l-1)==dd(k-1,l-1) % left and top-left neighbors are repeated dd(k,l)=dy(i,l)/factorial(l-1); fprintf('k=%d, l=%d\n',k,l); pause else dd(k,l)=(dd(k,l-1)-dd(k-1,l-1))/(z(k)-z(k-l+1)); end fprintf('%6.3f\t',dd(k,l)); end fprintf('\n'); end end end

The array of divided differences generated by the function `DividedDifference2`

is given below, the elements along the diagonal are the coefficients in the Hermite
polynomials.

(60) |