Module for the Newton Polynomial

Numerical Methods Using Mathematica 
Complementary Supplement 
for 
Numerical Analysis - Numerical Methods 
(c) 
John H. Mathews and Kurtis Fink, 2005

Background.

We have seen how to expand a function  f (x) in a Maclaurin polynomial about (Subscript[x, 0]) = 0 involving the powers (x^k) and a Taylor polynomial about (Subscript[x, 0]) ≠ 0 involving the powers (((x - (Subscript[x, 0])))^k).  These polynomials have a single "center" (Subscript[x, 0]).  Polynomial interpolation can be used to construct the polynomial of  degree  ≤n  that passes through the n+1 points  ((Subscript[x, k]), (Subscript[y, k]))  =  ((Subscript[x, k]), f ((Subscript[x, k]))),  for  k = 0, 1, ⋯, n.  If multiple "centers"    (Subscript[x, 0]), (Subscript[x, 1]),  ... , (Subscript[x, n])  are used, then the result is the so called Newton polynomial.  We attribute much of the founding theory to Sir Isaac Newton (1643-1727).

Theorem (Newton Polynomial).  Assume that  f∈(C^(n + 1))[a, b] and  (Subscript[x, k]) ∈ [a, b] for  k = 0, 1, ⋯, n  are distinct  values.  Then 

	f (x) = (Subscript[P, n]) (x) + (Subscript[R, n]) (x),
	
where (Subscript[P, n]) (x) is a polynomial that can be used to approximate  f (x), 

	(Subscript[P, n]) (x)  =  
        (Subscript[a, 0])  +  (Subscript[a, 1]) (x - (Subscript[x, 0]))  +  (Subscript[a, 2]) (x - 
                (Subscript[x, 0])) (x - (Subscript[x, 1]))  +  (Subscript[a, 3]) (x - (Subscript[x, 0])) (x - 
                (Subscript[x, 1])) (x - (Subscript[x, 2]))
                         +  ... 
\(\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \(\(+\ \ a\_n\) \((x - x\_0)\)\ \((x - x\_1)\)\ \((x - x\_2)\)\  ... \) \((x - 
            x\_\(n - 1\))\)\)  

and we write  

	f (x) ≈ (Subscript[P, n]) (x).

The Newton polynomial goes through the n + 1 points  (Subscript[{((Subscript[x, k]),(Subscript[y, k]))}, k = 0]^n),  i.e.

	(Subscript[P, n]) ((Subscript[x, k]))  =  f ((Subscript[x, k]))    for   k = 0,1,.... ,n.  

The remainder term  (Subscript[R, n]) (x) has the form

	(Subscript[R, n]) (x)=((Derivative[n + 1][f](c))/((n + 
                      1)!))(x - (Subscript[x, 0])) (x - (Subscript[x, 1])) (x - (Subscript[x, 2])) ...(x - (Subscript[x, n - 1]))(x - (Subscript[x, n])),

for some value c = c (x) that lies in the interval [a, b].  The coefficients  (Subscript[a, i])  are constructed using divided differences.

Proof  Newton Polynomials  Newton Polynomials

Definition. Divided Differences.

The divided differences for a function f[x] are defined as follows:

	\(f[x\_\(i - 1\), 
          x\_i] = \(f[x\_i]\  - \ f[x\_\(i - 1\)]\)\/\(\(x\_i\  - \ x\_\(i - 1\)\)\(
\)
          \)\)
\(f[x\_\(i - 2\), x\_\(i - 1\), 
          x\_i] = \(f[x\_\(i - 1\), x\_i]\  - \ f[x\_\(i - 2\), x\_\(i - 1\)]\)\/\(\(x\_i\  - \ x\_\(i - 2\)\)\(
\)
          \)\)
\(f[x\_\(i - 3\), x\_\(i - 2\), x\_\(i - 1\), 
          x\_i] = \(f[x\_\(i - 2\), x\_\(i - 1\), x\_i\ ] - \ f[x\_\(i - 3\), x\_\(i - 2\), x\_\(i - 1\)]\)\/\(\(x\_\(\(i\)\(\ \)\) - \ x\_\(i - 3\)\)\(
\)
          \)\)
\(f[x\_\(i - j\), x\_\(i - j + 1\),  ... , 
          x\_i] = \(f[x\_\(i - j + 1\),  ... , x\_i]\  - \ f[x\_\(i - j\),  ... , x\_\(i - 1\)]\)\/\(x\_i\  - \ x\_\(i - j\)\)\)

The divided difference formulae are used to construct the divided difference table:

({{(Subscript[x, i]), f[(Subscript[x, i])]}, {({{(Subscript[x, 0])}, {(Subscript[x, 1])}, {(Subscript[x, 2])}, {(Subscript[x, 3])}, {(Subscript[x, 4])}}), ({{f[(Subscript[x, 0])]}, {f[(Subscript[x, 1])]}, {f[(Subscript[x, 2])]}, {f[(Subscript[x, 3])]}, {f[(Subscript[x, 4])]}})}})({{f[(Subscript[x, i - 1]), (Subscript[x, i])], f[(Subscript[x, i - 2]), (Subscript[x, i - 1]), 
              (Subscript[x, i])]}, {({{  }, {f[(Subscript[x, 0]), (Subscript[x, 1])]}, {f[(Subscript[x, 1]), (Subscript[x, 2])]}, {f[(Subscript[x, 2]), (Subscript[x, 3])]}, {f[(Subscript[x, 3]), (Subscript[x, 4])]}}), ({{  }, {   }, {f[(Subscript[x, 0]), (Subscript[x, 1]), (Subscript[x, 2])]}, {f[(Subscript[x, 1]), (Subscript[x, 2]), (Subscript[x, 3])]}, {f[(Subscript[x, 2]), (Subscript[x, 3]), (Subscript[x, 4])]}})}})({{f[(Subscript[x, i - 3]), (Subscript[x, i - 2]), (Subscript[x, i - 1]), (Subscript[x, i])], f[
              (Subscript[x, i - 4]), (Subscript[x, i - 3]), (Subscript[x, i - 2]), (Subscript[x, i - 1]), 
              (Subscript[x, i])]}, {({{  }, {   }, {   }, {f[(Subscript[x, 0]), (Subscript[x, 1]), (Subscript[x, 2]), (Subscript[x, 3])]}, {f[(Subscript[x, 1]), (Subscript[x, 2]), (Subscript[x, 3]), (Subscript[x, 4])]}}), ({{  }, {   }, {   }, {   }, {f[(Subscript[x, 0]), (Subscript[x, 1]), (Subscript[x, 2]), (Subscript[x, 3]), (Subscript[x, 4])]}})}})

The coefficient (Subscript[a, i]) of the Newton polynomial  (Subscript[P, n]) (x) is  (Subscript[a, i]) = f[(Subscript[x, 0]), (Subscript[x, 1]),  ... , (Subscript[x, i])]  and it is the top element in the column of the i-th divided differences.

The Newton polynomial of  degree  ≤n  that passes through the  n + 1 points  ((Subscript[x, k]), (Subscript[y, k]))  =  ((Subscript[x, k]), f ((Subscript[x, k]))),  for  k = 0, 1, ⋯, n  is  

	(Subscript[P, n]) (x)  =  
        (Subscript[a, 0])  +  (Subscript[a, 1]) (x - (Subscript[x, 0]))  +  (Subscript[a, 2]) (x - 
                (Subscript[x, 0])) (x - (Subscript[x, 1]))  +  (Subscript[a, 3]) (x - (Subscript[x, 0])) (x - 
                (Subscript[x, 1])) (x - (Subscript[x, 2]))
                         +  ... 
\(\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \(\(+\ \ a\_n\) \((x - x\_0)\)\ \((x - x\_1)\)\ \((x - x\_2)\)\  ... \) \((x - 
            x\_\(n - 1\))\)\)

Proof  Newton Polynomials  Newton Polynomials

The cubic curve in the figure below illustrates a Newton polynomial of degree n = 3.

Clear[f, p, t, x, y, z]; 
f[x_, y___, z_] := 
      Module[{}, Return[((f[y, z] - f[x, y])/(z - x))]]
p[t_] = 
        f[(Subscript[x, 0])]  +  
          f[(Subscript[x, 0]), (Subscript[x, 1])] (t - (Subscript[x, 0]))  +  
f[(Subscript[x, 0]), (Subscript[x, 1]), 
              (Subscript[x, 2])] (t - (Subscript[x, 0])) (t - (Subscript[x, 1]))  +  
f[
              (Subscript[x, 0]), (Subscript[x, 1]), (Subscript[x, 2]), 
              (Subscript[x, 3])] (t - (Subscript[x, 0])) (t - (Subscript[x, 1])) (t - (Subscript[x, 2]));
Print["p[x] = ", p[x]]; 
\(\(Print[\*""p[\!\(x\_0\)] = "", 
        p[x\_0]];\)\)
\(\(Print[\*""p[\!\(x\_1\)] = "", p[x\_1]];\)\)
\(\(Print[\*""p[\!\(x\_2\)] = "", 
        p[x\_2]];\)\)
\(\(Print[\*""p[\!\(x\_2\)] = "", 
        Together[p[x\_2]]];\)\)
\(\(Print[\*""p[\!\(x\_3\)] = "", 
        p[x\_3]];\)\ \ \)
\(\(Print[\*""p[\!\(x\_3\)] = "", 
        Together[p[x\_3]]\ ];\)\ \ \)

Theorem.  (Error Bounds for Newton Interpolation, Equally Spaced Nodes)  Assume that  f (x)  defined on [a, b],  which contains the equally spaced nodes  (Subscript[x, k]) = (Subscript[x, 0]) +  kh.  Additionally, assume that    f (x)  and the derivatives of  f (x)  up to the order  n + 1  are continuous and bounded on the special subintervals  [(Subscript[x, 0]),(Subscript[x, 1])], [(Subscript[x, 0]),(Subscript[x, 2])], [(Subscript[x, 0]),(Subscript[x, 3])], [(Subscript[x, 0]),(Subscript[x, 4])], and [(Subscript[x, 0]),(Subscript[x, 5])], respectively;  that is, 

	|(f^((n + 1))) (x)| ≤  (Subscript[M, n + 1])   for   (Subscript[x, 0])<x<(Subscript[x, n]),  

for  n = 1,2,3,4,5.  The error terms corresponding to these three cases have the following useful bounds on their magnitude  

(i).	|(Subscript[R, 1]) (x)|≤  ((Subscript[M, 2])/8) (h^2)   is valid for  x∈[(Subscript[x, 0]),(Subscript[x, 1])],  

(ii).	\(\(\(|\)\(\(R\_2\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((Subscript[M, 3])/(9 (Sqrt[3]))) (h^3)   is valid for  x∈[(Subscript[x, 0]),(Subscript[x, 2])],  

(iii).	\(\(\(|\)\(\(R\_3\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((Subscript[M, 4])/24) (h^4)   is valid for  x∈[(Subscript[x, 0]),(Subscript[x, 3])],  

(iv).	\(\(\(|\)\(\(R\_4\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((Sqrt[4750 + 290 (Sqrt[145])])/3000) (Subscript[M, 5]) (h^5)   is valid for  x∈[(Subscript[x, 0]),(Subscript[x, 4])],  

(v).	\(\(\(|\)\(\(R\_5\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((10 + 7 (Sqrt[7]))/1215) (Subscript[M, 6]) (h^6)   is valid for  x∈[(Subscript[x, 0]),(Subscript[x, 5])].

Proof  Newton Polynomials  Newton Polynomials

Animations (Newton Polynomials  Newton Polynomials).  Internet hyperlinks to animations.

Computer Programs  Newton Polynomials  Newton Polynomials

Algorithm (Newton Interpolation Polynomial).   To construct and evaluate the Newton polynomial of  degree  ≤n  that passes through the n+1 points  ((Subscript[x, i]), (Subscript[y, i]))  =  ((Subscript[x, i]), f ((Subscript[x, i]))),  for  i = 0, 1, ⋯, n     

	\(\(P\_n\) \((x)\)\  = \ \(\(d\_\(0, 0\) + \(d\_\(1, 1\)\) \((x - 
                  x\_0)\) + \(d\_\(2, 2\)\) \((x - x\_0)\)\ \((x - x\_1)\)

          		\ \ \ \ \ \ \ \ \(\(\(+d\_\(3, 3\)\) \((x - x\_0)\)\ \((x - 
                    x\_1)\)\ \((x - x\_2)\)\)\(\ \ \)\(+\)\)\ \  ... \)

          \ \ \ \ \ \ \ \ \ \ \ \ \ \ \(+d\_\(n, n\)\) \((x - x\_0)\)\ \((x - 
                  x\_1)\)\ \((x - x\_2)\)\  ... \) \((x - 
              x\_\(n - 1\))\)\)  

where  
	(Subscript[d, i, 0]) = (Subscript[y, i])   for   i = 0, 1,  ... , n  
and  
	(Subscript[d, i, j])=(((Subscript[d, i, j - 1])  -  
                    (Subscript[d, i - 1, j - 1]))/((Subscript[x, i])  -  (Subscript[x, i - j])))   for  i=1,2,...,n
	                                   and  j=1,2,...,i  

Remark 1.  Newton polynomials are created "recursively." 

	(Subscript[P, n]) (x)  =  (Subscript[P, n - 1]) (x)  +  (Subscript[d, n, 
                    n]) (x - (Subscript[x, 0])) (x - (Subscript[x, 1])) (x - 
                    (Subscript[x, 2]))  ...  (x - (Subscript[x, n - 1])).

Remark 2.  Mathematica's arrays are lists and the subscripts must start with  1  instead of  0.

Mathematica Subroutine (Newton Polynomial).

\(\(\(NewtonPoly[
        XY_]\)\(:=\)\(
\)\(Module[{d, j, i, n, X, 
          Y}, 
X\  = \ \(Transpose[XY]\)\_\(\([[\)\(1\)\(]]\)\); \ 
Y\  = \ \(Transpose[XY]\)\_\(\([[\)\(2\)\(]]\)\); \ 
n\  = \ 
          Length[XY] - 1; \ 
d\  = \ 
          Table["", {n + 1}, {n + 
                1}\ ]; \ 
d\_\(\([[\)\(All, 1\)\(]]\)\)\  = \ 
          Y\_\(\([[\)\(All\)\(]]\)\); \ 
For[\ j = 1, \ 
          j ≤ 
            n, \ \(j++\), 
\(For[\ i = j, \ 
              i ≤ 
                n, \ \(i++\), 
d\_\(\([[\)\(i + 1, j + 1\)\(]]\)\)\  = \ \(d\_\(\([[\)\(i + 1, j\)\(]]\)\)\  - \ d\_\(\([[\)\(i, j\)\(]]\)\)\)\/\(X\_\(\([[\)\(i + 1\)\(]]\)\)\  - \ X\_\(\([[\)\(i + 1 - j\)\(]]\)\)\)\ ];\)\ ]; \ 
For[\ i = 0, \ 
          i ≤ 
            n, \ \(i++\), 
p[i + 1, 
              x_]\  = \ ∏\+\(j = 1\)\%i\((x - 
                  X\_\(\([[\)\(j\)\(]]\)\))\)\ ]; \ 
Dtable = 
          d; 
Return[\ \ ∑\+\(i = 0\)\%n d\_\(\([[\)\(i + 1, i + 1\)\(]]\)\)\ p[i + 1, 
                x]\ ];\ ]\)\(\ \)\)\)

Example 1.   Form the Newton polynomials of degree  n = 1,2, 3, 4, and 5  for the function  f[x] = Cos[x]  over the interval  [(Subscript[x, 0]), (Subscript[x, n])]  using equally spaced nodes selected from the following list  
{{(Subscript[x, k]), (Subscript[y, k])}_(k = 0)^5} = {{0, 1}, {(1/5), Cos[(1/5)]}, {(2/5), 
            Cos[(2/5)]}, {(3/5), Cos[(3/5)]}, {(4/5), Cos[(4/5)]}, {1, 
            Cos[1]}}

Solution 1 (a).

Use the nodes  {{(Subscript[x, k]), (Subscript[y, k])}_(k = 0)^1} to construct the Newton interpolation polynomial  (Subscript[P, 1])[x],  of degree n = 1, and compare it to the polynomial constructed with Mathematica's InterpolatingPolynomial procedure.

f[x_] = Cos[x]; 
n = 1; 
XY = Table[{(1.0/5) k, f[(1.0/5) k]}, {k, 0, n}]; 
(Subscript[p, 1])[x_] = NewtonPoly[XY]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.2].", "\n\
"];  
\(\(Print[\*""\!\(p\_1\)[x] = "", p\_1[x]\ ];\)\ \ \)

The polynomial obtained with Mathematica's InterpolatingPolynomial procedure is the nested form of the Newton polynomial.

\(\(\(Print[\*""\!\(p\_1\)[x] = "", 
        InterpolatingPolynomial[XY, x]\ ];\)\(\ \ \)\)\)

Now graph the function and polynomial, and interpolation nodes.

Needs["Graphics`Colors`"]; 
dots = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
gr1 = 
        Plot[{f[x], (Subscript[p, 1])[x]}, {x, -1.5, 2.5}, 
          PlotStyle -> {Magenta, Red}, 
          DisplayFunction -> Identity]; 
Show[gr1, dots, 
        PlotRange -> {{-0.5, 2.1}, {-0.5, 1.05}}, 
        DisplayFunction -> $DisplayFunction]; 
Print["f[x] = ", f[x]]; 
Print["Nodes = ", XY]; 
Print["The interval for interpolation is [0.0,0.2]."]; 
\(\(Print[\*""\!\(p\_1\)[x] = "", p\_1[x]];\)\ \)

Solution 1 (b).

Use the nodes  {{(Subscript[x, k]), (Subscript[y, k])}_(k = 0)^2}  to construct the Newton interpolation polynomial  (Subscript[P, 2])[x],  of degree n = 2, and compare it to the polynomial constructed with Mathematica's InterpolatingPolynomial procedure.

f[x_] = Cos[x]; 
n = 2; 
XY = Table[{(1.0/5) k, f[(1.0/5) k]}, {k, 0, n}]; 
(Subscript[p, 2])[x_] = NewtonPoly[XY]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.4].", "\n\
"]; 
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        p\_2[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        Expand[p\_2[x]]\ ];\)\ \ \)

The polynomial obtained with Mathematica's InterpolatingPolynomial procedure is the nested form of the Newton polynomial.

\(\(\(Print[\*""\!\(p\_2\)[x] = "", 
        InterpolatingPolynomial[XY, x]\ ];\)\(\ \ \)\)\)

Notice that  (Subscript[P, 3])[x]  is obtained from  (Subscript[P, 2])[x]  by adding one more term.

\(\(Print[\*""\!\(p\_1\)[x] = "", 
        p\_1[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_2\)[x] = "", p\_2[x]\ ];\)\ \ \)

Now graph the function and polynomial, and interpolation nodes.

Needs["Graphics`Colors`"]; 
dots = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
gr2 = 
        Plot[{f[x], (Subscript[p, 2])[x]}, {x, -1.5, 2.5}, 
          PlotStyle -> {Magenta, Brown}, 
          DisplayFunction -> Identity]; 
Show[gr2, dots, 
        PlotRange -> {{-0.5, 2.1}, {-0.5, 1.05}}, 
        DisplayFunction -> $DisplayFunction]; 
Print["f[x] = ", f[x]]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.4].", "\n\
"]; 
\(\(Print[\*""\!\(p\_2\)[x] = "", p\_2[x]\ ];\)\ \ \)

Solution 1 (c).

Use the nodes  {{(Subscript[x, k]), (Subscript[y, k])}_(k = 0)^3}  to construct the Newton interpolation polynomial  (Subscript[P, 3])[x],  of degree n = 3, and compare it to the polynomial constructed with Mathematica's InterpolatingPolynomial procedure.

f[x_] = Cos[x]; 
n = 3; 
XY = 
        Table[{(1.0/5) k, f[(1.0/5) k]}, {k, 0, 
            n}]; 
(Subscript[p, 3])[x_] = NewtonPoly[XY];  
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.6].", "\n\
"]; 
\(\(Print[\*""\!\(p\_3\)[x] = "", 
        p\_3[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_3\)[x] = "", 
        Expand[p\_3[x]]\ ];\)\ \)

The polynomial obtained with Mathematica's InterpolatingPolynomial procedure is the nested form of the Newton polynomial.

\(\(\(Print[\*""\!\(p\_3\)[x] = "", 
        InterpolatingPolynomial[XY, x]\ ];\)\(\ \ \)\)\)

Notice that  (Subscript[P, 3])[x]  is obtained from  (Subscript[P, 2])[x]  by adding one more term.

\(\(Print[\*""\!\(p\_1\)[x] = "", 
        p\_1[x]];\)\ \ \)
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        p\_2[x]];\)\ \ \)
\(\(Print[\*""\!\(p\_3\)[x] = "", p\_3[x]];\)\ \)

Now graph the function and polynomial, and interpolation nodes.

dots = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
gr3 = 
        Plot[{f[x], (Subscript[p, 3])[x]}, {x, -1.5, 2.5}, 
          PlotStyle -> {Magenta, Cyan}, 
          DisplayFunction -> Identity]; 
Show[gr3, dots, 
        PlotRange -> {{-0.5, 2.1}, {-0.5, 1.05}}, 
        DisplayFunction -> $DisplayFunction]; 
Print["f[x] = ", f[x]]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.6].", "\n\
"]; 
\(\(Print[\*""\!\(p\_3\)[x] = "", p\_3[x]\ ];\)\ \ \)

Solution 1 (d).

Use the nodes  {{(Subscript[x, k]), (Subscript[y, k])}_(k = 0)^4}  to construct the Newton interpolation polynomial  (Subscript[P, 4])[x],  of degree n = 4, and compare it to the polynomial constructed with Mathematica's InterpolatingPolynomial procedure.

f[x_] = Cos[x]; 
n = 4; 
XY = Table[{(1.0/5) k, f[(1.0/5) k]}, {k, 0, n}]; 
(Subscript[p, 4])[x_] = NewtonPoly[XY]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.8].", "\n\
"]; 
\(\(Print[\*""\!\(p\_4\)[x] = "", 
        p\_4[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_4\)[x] = "", 
        Expand[p\_4[x]]\ ];\)\ \)

The polynomial obtained with Mathematica's InterpolatingPolynomial procedure is the nested form of the Newton polynomial.

\(\(\(Print[\*""\!\(p\_4\)[x] = "", 
        InterpolatingPolynomial[XY, x]\ ];\)\(\ \ \)\)\)

Notice that  (Subscript[P, 4])[x]  is obtained from  (Subscript[P, 3])[x]  by adding one more term.

\(\(Print[\*""\!\(p\_1\)[x] = "", 
        p\_1[x]];\)\ \ \)
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        p\_2[x]];\)\ \ \)
\(\(Print[\*""\!\(p\_3\)[x] = "", 
        p\_3[x]];\)\ \)
\(\(Print[\*""\!\(p\_4\)[x] = "", p\_4[x]];\)\ \)

Now graph the function and polynomial, and interpolation nodes.

dots = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
gr4 = 
        Plot[{f[x], (Subscript[p, 4])[x]}, {x, -1.5, 2.5}, 
          PlotStyle -> {Magenta, DarkGreen}, 
          DisplayFunction -> Identity]; 
Show[gr4, dots, 
        PlotRange -> {{-0.5, 2.1}, {-0.5, 1.05}}, 
        DisplayFunction -> $DisplayFunction]; 
Print["f[x] = ", f[x]]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.8].", "\n\
"]; 
\(\(Print[\*""\!\(p\_4\)[x] = "", p\_4[x]\ ];\)\ \ \)

Solution 1 (e).

Use the nodes  {{(Subscript[x, k]), (Subscript[y, k])}_(k = 0)^5}  to construct the Newton interpolation polynomial  (Subscript[P, 5])[x],  of degree n = 5, and compare it to the polynomial constructed with Mathematica's InterpolatingPolynomial procedure.

Clear[f, P5, Q5, x]; 
f[x_] = Cos[x]; 
n = 5; 
XY = 
        Table[{(1.0/5) k, f[(1.0/5) k]}, {k, 0, 
            n}]; 
(Subscript[p, 5])[x_] = NewtonPoly[XY]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,0.8].", "\n\
"]; 
\(\(Print[\*""\!\(p\_5\)[x] = "", 
        p\_5[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_5\)[x] = "", 
        Expand[p\_5[x]]\ ];\)\ \)

The polynomial obtained with Mathematica's InterpolatingPolynomial procedure is the nested form of the Newton polynomial.

\(\(\(Print[\*""\!\(p\_5\)[x] = "", 
        InterpolatingPolynomial[XY, x]\ ];\)\(\ \ \)\)\)

Notice that  (Subscript[P, 5])[x]  is obtained from  (Subscript[P, 4])[x]  by adding one more term.

\(\(Print[\*""\!\(p\_1\)[x] = "", 
        p\_1[x]];\)\ \ \)
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        p\_2[x]];\)\ \ \)
\(\(Print[\*""\!\(p\_3\)[x] = "", 
        p\_3[x]];\)\ \)
\(\(Print[\*""\!\(p\_4\)[x] = "", 
        p\_4[x]];\)\ \)
\(\(Print[\*""\!\(p\_5\)[x] = "", p\_5[x]];\)\)

Now graph the function and polynomial, and interpolation nodes.

dots = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
gr5 = 
        Plot[{f[x], (Subscript[p, 5])[x]}, {x, -1.5, 2.5}, 
          PlotStyle -> {Magenta, Blue}, 
          DisplayFunction -> Identity]; 
Show[gr5, dots, 
        PlotRange -> {{-0.5, 2.1}, {-0.5, 1.05}}, 
        DisplayFunction -> $DisplayFunction]; 
Print["f[x] = ", f[x]]; 
Print["\n", "Nodes = ", 
        XY, "\n"];  
Print["The interval for interpolation is [0.0,1.0].", "\n\
"]; 
\(\(Print[\*""\!\(p\_5\)[x] = "", p\_5[x]\ ];\)\ \ \)

Solution 1 (f).

Show[gr1, gr2, gr3, gr4, gr5, dots, 
        PlotRange -> {{-1.5, 2.5}, {-1.2, 1.2}}, 
        DisplayFunction -> $DisplayFunction];  
Print["f[x] = ", f[x] ];  
\(\(Print[\*""\!\(p\_1\)[x] = "", 
        p\_1[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        p\_2[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_3\)[x] = "", 
        p\_3[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_4\)[x] = "", 
        p\_4[x]\ ];\)\ \ \)
\(\(Print[\*""\!\(p\_5\)[x] = "", p\_5[x]\ ];\)\ \ \)

Solution 1 (g).

Comparison with the Lagrange polynomial.

When we constructed the Lagrange polynomial of degree n = 5 through the six equally spaced nodes we obtained  

		(Subscript[P, 5])[
          x] =  - 
            26.0417 (-1 + x) (-0.8 + x) (-0.6 + 
              x) (-0.4 + x) (-0.2 + x)
                + 
          127.613 (-1 + x) (-0.8 + x) (-0.6 + 
            x) (-0.4 + x) x
                - 
          239.86 (-1 + x) (-0.8 + x) (-0.6 + 
            x) (-0.2 + x) x
                + 
          214.931 (-1 + x) (-0.8 + x) (-0.4 + 
            x) (-0.2 + x) x
                - 
          90.717 (-1 + x) (-0.6 + x) (-0.4 + 
            x) (-0.2 + x) x
                + 
          14.0704 (-0.8 + x) (-0.6 + x) (-0.4 + 
            x) (-0.2 + x) x 
 
 which involves 30 multiplications and  possibly 35 additions or subtractions (for this example 30 because the first node was 0).

The Newton polynomial is
  
		(Subscript[P, 5])[x] = 
        1  -  0.09966711079379187` x  -  
          0.48840202099497687` (-0.2` + x) x  +  
          0.04900763384891306` (-0.4` + x) (-0.2` + x) x
                    + 
          0.03812245731240869` (-0.6000000000000001` + 
            x) (-0.4` + x) (-0.2` + x) x
                    - 
          0.003962046706969001` (-0.8` + 
            x) (-0.6000000000000001` + x) (-0.4` + 
            x) (-0.2` + x) x 
and it involves 15 multiplications and possibly 20 additions or subtractions (for this example 15 because the first node was 0).

The Newton polynomial is
  
		(Subscript[P, 5])[x] = 
        1  -  0.09966711079379187` x  -  
          0.48840202099497687` (-0.2` + x) x  +  
          0.04900763384891306` (-0.4` + x) (-0.2` + x) x
                    + 
          0.03812245731240869` (-0.6000000000000001` + 
            x) (-0.4` + x) (-0.2` + x) x
                    - 
          0.003962046706969001` (-0.8` + 
            x) (-0.6000000000000001` + x) (-0.4` + 
            x) (-0.2` + x) x 
and it involves 15 multiplications and possibly 20 additions or subtractions (for this example 15 because the first node was 0).

The nested Newton polynomial obtained with Mathematica's InterpolatingPolynomial procedure is 

		\(p\_5[x]\  = 
        1 + \((\(-0.09966711079379187`\) + \((\(-0.48840202099497687`\) + \((\(\(0.04900763384891306`\)\(\)\) + \((\(\(0.038122457312408806`\)\(\)\) - 
                                  0.003962046706969931`\ \((\(-0.8`\) + 
                                        x)\))\)\ \((\(-0.6000000000000001`\) + x)\))\)\ \((\(-0.4`\) + x)\))\)\ \((\(-0.2`\) + x)\))\)\ x\)

and it involves 5 multiplications and possibly 10 additions or subtractions ( for this example 9 because the first node was 0).

In general the nth degree Lagrange polynomial requires  n (n  +  1)  multiplications and  n (n  +  2) additions and subtractions. 

The nth degree Newton polynomial requires  ((n (n  +  1))/2)  multiplications and  ((n (n  +  3))/2) additions and subtractions. 

And the nth degree Mathematica InterpolatingPolynomial requires n multiplications and  2  n additions and subtractions.

Observation.  Consider the Newton polynomial  (Subscript[P, 5])[x]  and Lagrange polynomial of degree 5 that use six equally spaced nodes in the interval  [0,1].  Are they the same ?  The answer is "yes"  even though there were two "routes" for constructing the polynomial.  It is a fact that the polynomial of degree n that passes through n+1 points is "unique."

Example 2.  Error Analysis.  Investigate the error for the Newton polynomial approximations in Example 1.

Solution 2 (a).

Investigate the error over the interval  [(Subscript[x, 0]), (Subscript[x, 1])]  =  [0.0,  0.2]  for the Newton interpolation polynomial  (Subscript[P, 1])[x],  of degree n = 1.

(Subscript[e, 1])[x_] = f[x] - (Subscript[p, 1])[x];
extrema1  =  
        FindMinimum[-(Subscript[e, 1])[x], {x, 0.1}];
max1 = -(Subscript[extrema1, [[1]]]);
extrema  =  {max1};
maxerror1 =  Max[Abs[{max1}]];
ge1 = 
        Plot[{max1, (Subscript[e, 1])[x]}, {x, 0.0, 0.2}, 
          PlotRange -> {{0, 1}, {0, 0.005}}, 
          PlotStyle -> {Cyan, Blue}];
Print[" f[x] = ", f[x]];
\(\(Print[\*""\!\(p\_1\)[x] = "", 
        p\_1[x]];\)\)
Print["The interval for interpolation is [0.0, 0.2]."];
\(\(Print[\*""Graph of the error (Subscript[e, 1])[x] = f[x]-(Subscript[p, 1])[x]""];\)\)
\(\(Print[\*""Extrema for \!\(e\_1\)[x] is "", 
        extrema\ ];\)\)
\(\(Print[\*""|\!\(e\_1\)[x]| ≤ "", 
        maxerror1\ ];\)\)

Use formula(i).	|(Subscript[R, 1]) (x)|≤  ((Subscript[M, 2])/8) (h^2)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 1])]  =  [0.0,  0.2],  and find the error bound for this example.

h = (1/5);
(Subscript[M, 2]) =Abs[Derivative[2][f][0]];
Plot[Derivative[2][f][x] ,{x, 0, 0.2},PlotRange -> {{0, 1}, {-1, 0}},PlotStyle -> Green];
Print[" 
(f^((2)))[x]  = ",Derivative[2][f][x]]; 
\(Print[\*""|(f^((2)))(x)| ≤  (Subscript[M, 2]) = "", 
        M\_2];\)
Print["h = ", 
        h];
\(Print[\*""The remainder term (Subscript[R, 2])(x) has the form""];\)
\(Print[\*""|(Subscript[R, 2])(x)| ≤ "", \*""((Subscript[M, 2])/8)(h^2) = "", \(M\_2\/8\) h\^2, " = ", 
        N[\(M\_2\/8\) h\^2]];\)

Thus,  \(\(\(|\)\(\(R\_2\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\(0.005`\)\)\)\)  is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 1])]  =  [0.0,  0.2],  which is a little bit larger than the maximum error  0.00497089.  After all, it is an error bound.

Solution 2 (b).

Investigate the error over the interval  [(Subscript[x, 0]), (Subscript[x, 2])]  =  [0.0,  0.4]  for the Newton interpolation polynomial  (Subscript[P, 2])[x],  of degree n = 2.

(Subscript[e, 2])[x_] = f[x] - (Subscript[p, 2])[x];
extrema1  =  
        FindMinimum[-(Subscript[e, 2])[x], {x, 0.1}];
max1 = -(Subscript[extrema1, [[1]]]);
extrema2  =  
        FindMinimum[(Subscript[e, 2])[x], {x, 0.3}];
min1 = 
        (Subscript[extrema2, [[1]]]);
extrema  =  {max1, min1};
maxerror2  =  Max[Abs[{max1, min1}]];
ge2 = 
        Plot[{max1, min1, (Subscript[e, 2])[x]}, {x, 0.0, 0.4}, 
          PlotRange -> {{0, 1}, Automatic}, 
          PlotStyle -> {Cyan, Cyan, Blue}];
Print[" f[x] = ", f[x]];
\(\(Print[\*""\!\(p\_2\)[x] = "", 
        p\_2[x]];\)\)
Print["The interval for interpolation is [0.0, 0.4]."];
\(\(Print[\*""Graph of the error (Subscript[e, 2])[x] = f[x]-(Subscript[p, 2])[x]""];\)\)
\(\(Print[\*""Extrema for \!\(e\_2\)[x] are "", 
        extrema\ ];\)\)
\(\(Print[\*""|\!\(e\_2\)[x]| ≤ "", 
        maxerror2\ ];\)\)

Use formula (ii).	\(\(\(|\)\(\(R\_2\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((Subscript[M, 3])/(9 (Sqrt[3]))) (h^3)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 2])]  =  [0.0,  0.4],  and find the error bound for this example.

h = (1/5);
(Subscript[M, 3]) =Derivative[3][f][(2/5)];
Plot[Derivative[3][f][x] ,{x, 0, 0.4},PlotRange -> {{0, 1}, {0, 1}},PlotStyle -> Green];
Print[" 
(f^((3)))[x]  = ",Derivative[3][f][x]]; 
\(Print[\*""|(f^((3)))(x)| ≤  (Subscript[M, 3]) = "", M\_3, " = ", 
        N[M\_3]];\)
Print["h = ", 
        h];
\(Print[\*""The remainder term (Subscript[R, 2])(x) has the form""];\)
\(Print[\*""|(Subscript[R, 2])(x)| ≤ "", \*""\!\(M\_3\/\(9\\ \@3\)\)(h^3) = "", \(M\_3\/\(9\ \@3\)\) h\^3, " = ", 
        N[\(M\_3\/\(9\ \@3\)\) h\^3]];\)

Thus,  \(\(\(|\)\(\(R\_2\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\(0.00019984958645269084`\)\)\)\)  is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 2])]  =  [0.0,  0.4],  which is a little bit larger than the maximum error  0.000116469.  After all, it is an error bound.

Solution 2 (c).

Investigate the error over the interval  [(Subscript[x, 0]), (Subscript[x, 3])]  =  [0.0,  0.6]  for the Newton interpolation polynomial  (Subscript[P, 3])[x],  of degree n = 3.

(Subscript[e, 3])[x_] = f[x] - (Subscript[p, 3])[x];
extrema1  =  
        FindMinimum[-(Subscript[e, 3])[x], {x, 0.3}];
max1 = -(Subscript[extrema1, [[1]]]);
extrema2  =  
        FindMinimum[(Subscript[e, 3])[x], {x, 0.1}];
min1 = 
        (Subscript[extrema2, [[1]]]);
extrema3  =  
        FindMinimum[(Subscript[e, 3])[x], {x, 0.5}];
min2 = 
        (Subscript[extrema3, [[1]]]);
extrema  =  {min1, max1, min2};
maxerror3  =  
        Max[Abs[{min1, max1, min2}]];
ge3 = 
        Plot[{min1, max1, min2, (Subscript[e, 3])[x]}, {x, 0.0, 0.6}, 
          PlotRange -> {{0, 1}, Automatic}, 
          PlotStyle -> {Cyan, Cyan, Cyan, Blue}];
Print[" f[x] = ", f[x]];
\(\(Print[\*""\!\(p\_3\)[x] = "", 
        p\_3[x]];\)\)
Print["The interval for interpolation is [0.0, 0.6]."];
\(\(Print[\*""Graph of the error (Subscript[e, 3])[x] = f[x]-(Subscript[p, 3])[x]""];\)\)
\(\(Print[\*""Extrema for \!\(e\_3\)[x] are "", 
        extrema\ ];\)\)
\(\(Print[\*""|\!\(e\_3\)[x]| ≤ "", 
        maxerror3\ ];\)\)

Use formula (iii).	\(\(\(|\)\(\(R\_3\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((Subscript[M, 4])/24) (h^4)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 3])]  =  [0.0,  0.6],  and find the error bound for this example.

h = (1/5);
(Subscript[M, 4]) =Derivative[4][f][0];
Plot[Derivative[4][f][x] ,{x, 0, 0.6},PlotRange -> {{0, 1}, {0, 1}},PlotStyle -> Green];
Print[" 
(f^((4)))[x]  = ",Derivative[4][f][x]]; 
Print["|
(f^((4)))[x]| ≤  (Subscript[M, 4]) = ",(Subscript[M, 4])," = ",N[(Subscript[M, 4])]];
Print["h = ", 
        h];
\(Print[\*""The remainder term (Subscript[R, 3])(x) has the form""];\)
\(Print[\*""|(Subscript[R, 3])(x)| ≤ "", \*""((Subscript[M, 4])/24)(h^4) = "", \(M\_4\/24\) h\^4, " = ", 
        N[\(M\_4\/24\) h\^4]];\)

Thus,  \(\(\(|\)\(\(R\_3\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \)\(0.00006666666666666667`\)\)\)\)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 3])]  =  [0.0,  0.6],  which is a little bit larger than the maximum error  0.0000642481.  After all, it is an error bound.

Solution 2 (d).

Investigate the error over the interval  [(Subscript[x, 0]), (Subscript[x, 4])]  =  [0.0,  0.8]  for the Newton interpolation polynomial  (Subscript[P, 4])[x],  of degree n = 4.

(Subscript[e, 4])[x_] = f[x] - (Subscript[p, 4])[x];
extrema1  =  
        FindMinimum[-(Subscript[e, 4])[x], {x, 0.3}];
max1 = -(Subscript[extrema1, [[1]]]);
extrema2  =  
        FindMinimum[-(Subscript[e, 4])[x], {x, 0.7}];
max2 = -(Subscript[extrema2, [[1]]]);
extrema3  =  
        FindMinimum[(Subscript[e, 4])[x], {x, 0.1}];
min1 = 
        (Subscript[extrema3, [[1]]]);
extrema4  =  
        FindMinimum[(Subscript[e, 4])[x], {x, 0.5}];
min2 = 
        (Subscript[extrema4, [[1]]]);
extrema  =  {min1, max1, min2, max2};
maxerror4  =  
        Max[Abs[{min1, max1, min2, max2}]];
ge4 = 
        Plot[{max1, max2, min1, min2, (Subscript[e, 4])[x]}, {x, 0.0, 0.8}, 
          PlotRange -> {{0, 1}, Automatic}, 
          PlotStyle -> {Cyan, Cyan, Cyan, Cyan, Blue}];
Print[" f[x] = ", f[x]];
\(\(Print[\*""\!\(p\_4\)[x] = "", 
        p\_4[x]];\)\)
Print["The interval for interpolation is [0.0, 0.8]."];
\(\(Print[\*""Graph of the error (Subscript[e, 4])[x] = f[x]-(Subscript[p, 4])[x]""];\)\)
\(\(Print[\*""Extrema for \!\(e\_4\)[x] are "", 
        extrema\ ];\)\)
\(\(Print[\*""|\!\(e\_4\)[x]| ≤ "", 
        maxerror4\ ];\)\)

Use formula (iv).	\(\(\(|\)\(\(R\_4\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((Sqrt[4750 + 290 (Sqrt[145])])/3000) (Subscript[M, 5]) (h^5)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 4])]  =  [0.0,  0.8],  and find the error bound for this example.

h = (1/5);
(Subscript[M, 5]) =Abs[Derivative[5][f][(4/5)]];
Plot[Abs[Derivative[5][f][x]] ,{x, 0, 0.8},PlotRange -> {{0, 1}, {0, 1}},PlotStyle -> Green];
Print["|
(f^((5)))[x]| ≤  (Subscript[M, 5]) = ",(Subscript[M, 5])," = ",N[(Subscript[M, 5])]];
Print["h = ", 
        h];
\(Print[\*""The remainder term (Subscript[R, 4])(x) has the form""];\)
\(Print[\*""|(Subscript[R, 4])(x)| ≤ "", \*""\!\(\@\(4750 + 290\\ \@145\)\/3000\) (Subscript[M, 5]) (h^5) = "", \@\(4750 + 290\ \@145\)\/3000\ M\_5\ h\^5, " = ", 
        N[\@\(4750 + 290\ \@145\)\/3000\ M\_5\ h\^5]];\)

Thus,  \(\(\(|\)\(\(R\_4\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \ \)\(6.94675×10\^\(-6\)\)\)\)\)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 4])]  =  [0.0,  0.8],  which is a little bit larger than the maximum error  4.23054×(10^(-6)).  After all, it is an error bound.

Solution 2 (e).

Investigate the error over the interval  [(Subscript[x, 0]), (Subscript[x, 5])]  =  [0.0,  1.0]  for the Newton interpolation polynomial  (Subscript[P, 5])[x],  of degree n = 5.

e5[x_] = f[x] - (Subscript[p, 5])[x];
extrema1  =  
        FindMinimum[-e5[x], {x, 0.1}];
max1 = -(Subscript[extrema1, [[1]]]);
extrema2  =  
        FindMinimum[-e5[x], {x, 0.5}];
max2 = -(Subscript[extrema2, [[1]]]);
extrema3  =  
        FindMinimum[-e5[x], {x, 0.9}];
max3 = -(Subscript[extrema3, [[1]]]);
extrema4  =  
        FindMinimum[e5[x], {x, 0.3}];
min1 = 
        (Subscript[extrema4, [[1]]]);
extrema5  =  
        FindMinimum[e5[x], {x, 0.7}];
min2 = 
        (Subscript[extrema5, [[1]]]);
extrema  =  {max1, min1, max2, min2, 
          max3};
maxerror5  =  
        Max[Abs[{max1, min1, max2, min2, 
              max3}]];
ge5 = 
        Plot[{max1, max2, max3, min1, min2, e5[x]}, {x, 0.0, 1.0}, 
          PlotRange -> {{0, 1}, Automatic}, 
          PlotStyle -> {Cyan, Cyan, Cyan, Cyan, Cyan, 
              Blue}];
Print[" f[x] = ", f[x]];
\(\(Print[\*""\!\(p\_5\)[x] = "", 
        p\_5[x]];\)\)
Print["The interval for interpolation is [0.0, 1.0]."];
\(\(Print[\*""Graph of the error (Subscript[e, 5])[x] = f[x]-(Subscript[p, 5])[x]""];\)\)
\(\(Print[\*""Extrema for \!\(e\_5\)[x] are "", 
        extrema\ ];\)\)
\(\(Print[\*""|\!\(e\_5\)[x]| ≤ "", 
        maxerror5\ ];\)\)

Use formula (v).	\(\(\(|\)\(\(R\_5\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \)\)\)\)((10 + 7 (Sqrt[7]))/1215) (Subscript[M, 6]) (h^6)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 5])]  =  [0.0,  1.0],  and find the error bound for this example.

h = (1/5);
(Subscript[M, 6]) =Abs[Derivative[6][f][0]];
Plot[Abs[Derivative[6][f][x]] ,{x, 0, 1},PlotRange -> {{0, 1}, {0, 1}},PlotStyle -> Green];
Print["|
(f^((6)))[x]| = ",Abs[Derivative[6][f][x]]]; 
Print["|
(f^((6)))[x]| ≤  (Subscript[M, 6]) = ",(Subscript[M, 6])];
Print["h = ", 
        h];
\(Print[\*""The remainder term (Subscript[R, 5])(x) has the form""];\)
\(Print[\*""|(Subscript[R, 5])(x)| ≤ "", \*""((10 + 7 (Sqrt[ 7]))/1215)(Subscript[M, 6]) (h^6) = "", \(\(10 + 7\ \@7\)\/1215\) \(M\_6\) 
          h\^6, " = ", 
        N[\(\(10 + 7\ \@7\)\/1215\) \(M\_6\) h\^6]];\)

Thus,  \(\(\(|\)\(\(R\_5\) \((x)\)\)\(|\)\(\ \)\(\(≤\)\(\ \ \ \)\(1.5023017179892483`*^-6\)\)\)\)   is valid for  x ∈ [(Subscript[x, 0]), (Subscript[x, 5])]  =  [0.0,  1.0],  which is a little bit larger than the maximum error  1.34999`*^-6.  After all, it is an error bound.

Example 3.  Summary of the maximum error and the error bounds over the intervals  [(Subscript[x, 0]), (Subscript[x, k])]  for the Newton polynomials in Example 2.
3 (a).  Find  (max_(o  ≤  x  ≤  0.2)) {|f[x] - 
              (Subscript[P, 1])[x]|}.  
3 (b).  Find  (max_(o  ≤  x  ≤  0.4)) {|f[x] - 
              (Subscript[P, 2])[x]|}.  
3 (c).  Find  (max_(o  ≤  x  ≤  0.6)) {|f[x] - 
              (Subscript[P, 3])[x]|}.  
3 (d).  Find  (max_(o  ≤  x  ≤  0.8)) {|f[x] - 
              (Subscript[P, 4])[x]|}.  
3 (e).  Find  (max_(o  ≤  x  ≤  1.0)) {|f[x] - 
              (Subscript[P, 5])[x]|}.

Solution 3.

Show[ge1]; Show[ge2]; Show[ge3]; Show[ge4]; 
    Show[ge5];
\(\(Print[\*""For (Subscript[p, 1])[x] over the interval [0,1] we have  |(Subscript[e, 1])[x]| ≤ "", 
        maxerror1\ ];\)\ \ \)
\(\(Print[\*""For (Subscript[p, 2])[x] over the interval [0,1] we have  |(Subscript[e, 2])[x]| ≤ "", 
        maxerror2\ ];\)\ \ \)
\(\(Print[\*""For (Subscript[p, 3])[x] over the interval [0,1] we have  |(Subscript[e, 3])[x]| ≤ "", 
        maxerror3\ ];\)\ \ \)
\(\(Print[\*""For (Subscript[p, 4])[x] over the interval [0,1] we have  |(Subscript[e, 4])[x]| ≤ "", 
        maxerror4\ ];\)\ \ \)
\(\(Print[\*""For (Subscript[p, 5])[x] over the interval [0,1] we have  |(Subscript[e, 5])[x]| ≤ "", maxerror5\ ];\)\ \ \)

\(\(Print[\*""The remainder term (Subscript[R, n])(x) has the form""];\)\)
\(\(Print[\*""|(Subscript[R, 2])(x)| ≤ "", \*""((Subscript[M, 2])/8)(h^2) = "", \(M\_2\/8\) 
          h\^2, \*""((1/5))^2"", " = ", 
        N[\(M\_2\/8\) \((1\/5)\)\^2]];\)\)
\(\(Print[\*""|(Subscript[R, 2])(x)| ≤ "", \*""\!\(M\_3\/\(9\\ \@3\)\) (h^3) = "", 
        M\_3\/\(9\ \@3\), \*""((1/5))^3"", " = ", 
        N[\(M\_3\/\(9\ \@3\)\) \((1\/5)\)\^3]];\)\)
\(\(Print[\*""|(Subscript[R, 3])(x)| ≤ "", \*""((Subscript[M, 4])/24) (h^4) = "", 
        M\_4\/24, \*""((1/5))^4"", " = ", 
        N[\(M\_4\/24\) \((1\/5)\)\^4]];\)\)
\(\(Print[\*""|(Subscript[R, 4])(x)| ≤ "", \*""\!\(\@\(4750 + 290\\ \@145\)\/3000\) (Subscript[M, 5]) (h^5) = "", \@\(4750 + 290\ \@145\)\/3000\ M\_5\ , \*""((1/5))^5"", " = ", 
        N[\@\(4750 + 290\ \@145\)\/3000\ M\_5\ \((1\/5)\)\^5]];\)\)
\(\(Print[\*""|(Subscript[R, 5])(x)| ≤ "", \*""\!\(\(10 + 7\\ \@7\)\/1215\)(Subscript[M, 6]) (h^6) = "", \(\(10 + 7\ \@7\)\/1215\) 
          M\_6, \*""((1/5))^6"", " = ", 
        N[\(\(10 + 7\ \@7\)\/1215\) \(M\_6\) \((1\/5)\)\^6]];\)\)

Print["The maximum error is less than the error bound."];
\(\(Print[\*""|\!\(e\_1\)[x]| ≤ "", 
        maxerror1\ , "   ≤ ", 
        N[\(M\_2\/8\) \((1\/5)\)\^2]];\)\ \ \)
\(\(Print[\*""|\!\(e\_2\)[x]| ≤ "", 
        maxerror2\ , "   ≤ ", 
        N[\(M\_3\/\(9\ \@3\)\) \((1\/5)\)\^3]];\)\ \ \)
\(\(Print[\*""|\!\(e\_3\)[x]| ≤ "", 
        maxerror3\ , "  ≤ ", 
        N[\(M\_4\/24\) \((1\/5)\)\^4]];\)\)
\(\(Print[\*""|\!\(e\_4\)[x]| ≤ "", 
        maxerror4, " ≤ ", 
        N[\@\(4750 + 290\ \@145\)\/3000\ M\_5\ \((1\/5)\)\^5]];\)\)
\(\(Print[\*""|\!\(e\_5\)[x]| ≤ "", 
        maxerror5, " ≤ ", 
        N[\(\(10 + 7\ \@7\)\/1215\) \(M\_6\) \((1\/5)\)\^6]];\)\)

Example 4.  Application to number theory. 
4 (a).  Find the formula for the sum of the first  n  integers.
4 (b).  Find the formula for the sum of the squares of the first  n integers.

Solution 4.

4 (a).  Find the formula for the sum of the first  n  integers.

xy = {{1, 1}, {2, 3}, {3, 6}, {4, 10}, {5, 15}} 
P[x_] = NewtonPoly[xy]; 
Print["P[x] = ", P[x]]; 
Print["P[x] = ", Together[P[x]] ];

Notice that the divided difference table has zeros for the 2nd and higher order differences.

TableForm[Dtable]

Thus the formula for the sum of the first n integers is:

Clear[n];
Print["P[n] = ", Together[P[n]] ];

Compare this formula with Mathematica's knowledge.

(∑_(i = 1)^n)i

4 (b).  Find the formula for the sum of the squares of the first  n  integers is:

xy = {{0, 0}, {1, 1}, {2, 5}, {3, 14}, {4, 30}, {5, 55}} 
P[x_] = NewtonPoly[xy]; 
Print["P[x] = ", P[x]]; 
Print["P[x] = ", Together[P[x]] ];

Notice that the divided difference table has zeros for the 3rd and higher order differences.

TableForm[Dtable]

Thus the formula for the sum of the squares of the first  n  integers is:

Clear[n];
Print["P[n] = ", Together[P[n]] ];

Factor[Together[ P[n] ] ]

Compare this formula with Mathematica's knowledge.

(∑_(i = 1)^n)(i^2)

Various Scenarios and Animations for Newton Polynomials.

\(\(NewtonPoly[
          XY_] := 
Module[{d, j, k, n, X, 
            Y}, 
X\  = \ \(Transpose[XY]\)\_\(\([[\)\(1\)\(]]\)\); 
Y\  = \ \(Transpose[XY]\)\_\(\([[\)\(2\)\(]]\)\); 
n\  = \ Length[XY] - 1; 
d\  = \ 
            Table["", {n + 1}, {n + 
                  1}\ ]; 
d\_\(\([[\)\(All, 1\)\(]]\)\)\  = \ 
            Y\_\(\([[\)\(All\)\(]]\)\); 
For[\ j = 1, \ 
            j ≤ 
              n, \ \(j++\), 
\(For[\ k = j, \ 
                k ≤ 
                  n, \ \(k++\), 
d\_\(\([[\)\(k + 1, j + 1\)\(]]\)\)\  = \ \(d\_\(\([[\)\(k + 1, j\)\(]]\)\)\  - \ d\_\(\([[\)\(k, j\)\(]]\)\)\)\/\(X\_\(\([[\)\(k + 1\)\(]]\)\)\  - \ X\_\(\([[\)\(k + 1 - j\)\(]]\)\)\)\ ];\)\ ]; 
For[\ k = 0, \ 
            k ≤ 
              n, \ \(k++\), 
p[k + 1, 
                x_]\  = \ ∏\+\(j = 1\)\%k\((x - 
                    X\_\(\([[\)\(j\)\(]]\)\))\)\ ]; 
Return[\ \ ∑\+\(k = 0\)\%n d\_\(\([[\)\(k + 1, k + 1\)\(]]\)\)\ p[k + 1, 
                  x]\ ];\ ];\)\)

Example 5.  Find the Newton polynomial approximation for  f[x] = (Sqrt[x]), on the interval [0, 2], [0, 8], and [0, 10].

Solution 5 (a).

f[x_] = (Sqrt[x]);  
x0 = 1;  
a = 0; 
b =  10.5; 
c = 0.0; 
d = 3.5; 
a0 = 0.0; 
b0 = 2.0; 
Needs["Graphics`Colors`"]; 
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
For[n = 2, 
      n ≤ 5, n++, 
aa = 
        0.0; 
bb = 2.0; 
XY = 
        N[Table[{aa + ((bb - aa)/n) k, 
              f[aa + ((bb - aa)/n) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY];    
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[0, 10, 2], Range[0, 4, 1]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; 
Print["\
f[x] = ", f[x], "\n", "P[x] = ", P[x], "\n", "P[x] = ", 
        Chop[Expand[P[x]]]];   ]

Animation 4.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Solution 5 (b).

f[x_] = (Sqrt[x]);  
x0 = 4;  
a = 0; 
b =  10.5; 
c = 0.0; 
d = 3.5; 
a0 = 0.0; 
b0 = 8.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 2, 
      n ≤ 5, \(n++\), 
aa = 
        0.0; 
bb = 8.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ \ 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[0, 10, 2], Range[0, 4, 1]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 5.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Solution 5 (c).

f[x_] = (Sqrt[x]);  
x0 = 5;  
a = 0; 
b =  10.5; 
c = 0.0; 
d = 3.5; 
a0 = 0.0; 
b0 = 10.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 2, 
      n ≤ 5, \(n++\), 
aa = 
        0.0; 
bb = 10.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ \ \ 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[0, 10, 2], Range[0, 4, 1]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 6.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 6.  Find the Newton polynomial approximation for  f[x] = (1/(Sqrt[1 - x])), on the interval [-1.5, 0.95].

Solution 6.

f[x_] = (1/(Sqrt[1 - x]));  
x0 = 0;  
a = -1.5; 
b =  1.5; 
c = 0.0; 
d = 8.5; 
a0 = -1.0; 
b0 = 1.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 2, 
      n ≤ 
        5, \(n++\), 
aa = \(-1.5\); 
bb = 
        0.95; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ 
graph1 = 
        Plot[f[x], {x, \(-1.5\),  .9999}, PlotStyle -> Magenta, 
          DisplayFunction -> Identity]; 
graph2 = 
        Plot[P[x], {x, \(-1.5\), 1.5}, PlotStyle -> Blue, 
          DisplayFunction -> Identity]; 
Show[graph1, 
        graph2, Cdot, Rends, PlotRange -> {{a, b}, {c, d}}, 
        Ticks -> {Range[\(-1.5\), 1.5, 0.5], Range[0, 5, 1]}, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 7.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 7.  Find the Newton polynomial approximation for  f[x] = (1/(1  +  10  (x^2))), on the interval [-1, 1].

Solution 7.

f[x_] = (1/(1 + 10  (x^2)));  
x0 = 0;  
a = -1.5; 
b =  1.5; 
c = 0.0; 
d = 1.25; 
a0 = -1.0; 
b0 = 1.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 3, 
      n ≤ 
        6, \(n++\), 
aa = \(-1.0\); 
bb = 
        1.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ 
graph = 
        Plot[{f[x], P[x]}, {x, \(-1.5\), 1.5}, 
          PlotStyle -> {Magenta, Blue}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, PlotRange -> {{a, b}, {c, d}}, 
        Ticks -> {Range[\(-1.5\), 1.5, 0.5], Range[0, 1, 0.5]}, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 8.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 8.  Find the Newton polynomial approximation for  f[x] = Log[x], on the interval [0.02, 2], and [0.02, 4].

Solution 8 (a).

f[x_] = Log[x];  
x0 = 1;  
a = -0.5; 
b =  4.1; 
c = -6.0; 
d = 3.0; 
a0 = 0; 
b0 = 2.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 2, 
      n ≤ 5, \(n++\), 
aa = 
        0.02; 
bb = 2.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ 
graph1 = 
        Plot[{f[x], P[x]}, {x, 0.001, b}, PlotStyle -> {Magenta, Blue}, 
          DisplayFunction -> Identity]; 
graph2 = 
        Plot[P[x], {x, \(-0.5\), b}, PlotStyle -> {Magenta, Blue}, 
          DisplayFunction -> Identity]; 
Show[graph1, 
        graph2, Cdot, Rends, PlotRange -> {{a, b}, {c, d}}, 
        Ticks -> {Range[\(-0.5\), 4.0, 0.5], Range[\(-5\), 0, 1]}, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 9.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Solution 8 (b).

f[x_] = Log[x];  
x0 = 2;  
a = -0.5; 
b =  4.1; 
c = -6.0; 
d = 3.0; 
a0 = 0; 
b0 = 4.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 2, 
      n ≤ 5, \(n++\), 
aa = 
        0.02; 
bb = 4.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ 
graph1 = 
        Plot[{f[x], P[x]}, {x, 0.001, b}, PlotStyle -> {Magenta, Blue}, 
          DisplayFunction -> Identity]; 
graph2 = 
        Plot[P[x], {x, \(-0.5\), b}, PlotStyle -> {Magenta, Blue}, 
          DisplayFunction -> Identity]; 
Show[graph1, 
        graph2, Cdot, Rends, PlotRange -> {{a, b}, {c, d}}, 
        Ticks -> {Range[\(-0.5\), 4.0, 0.5], Range[\(-6\), 0, 1]}, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 10.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 9.  Find the Newton polynomial approximation for  f[x] = Sin[x], on the interval [-10, 10].

Solution 9.

f[x_] = Sin[x];  
x0 = 0;  
a = -3 π; 
b = 3  π; 
c = -1.3; 
d =  1.3; 
a0 = -10; 
b0 = 10; 
Needs["Graphics`Colors`"];
Cdot = 
        Graphics[{{Red, PointSize[0.02], 
              Point[{x0, f[x0]}]}}]; 
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[i = 2, 
      i ≤ 6, \(i++\), 
n = 
        2  i + 4; 
XY = 
        N[Table[{a + \(\(b - a\)\/n\) k, f[a + \(\(b - a\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-3\) π, 3  π, π], 
              Range[\(-0.5\), 1, 0.5]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", 
        f[x], "
", "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 12.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 10.  Find the Newton polynomial approximation for  f[x] = Cos[x], on the interval [-10, 10].

Solution 10.

f[x_] = Cos[x];  
x0 = 0;  
a = -3 π; 
b = 3  π; 
c = -1.3; 
d =  1.3; 
a0 = -10; 
b0 = 10; 
Needs["Graphics`Colors`"];
Cdot = 
        Graphics[{{Red, PointSize[0.02], 
              Point[{x0, f[x0]}]}}]; 
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[i = 3, 
      i ≤ 6, \(i++\), 
n = 
        2  i; 
XY = 
        N[Table[{a + \(\(b - a\)\/n\) k, f[a + \(\(b - a\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-3\) π, 3  π, π], 
              Range[\(-0.5\), 1, 0.5]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 13.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 11.  Find the Newton polynomial approximation for  f[x] = Tan[x], on the interval [-1.45, 1.45].

Solution 11.

f[x_] = Tan[x];  
x0 = 0;  
a = -π; 
b = π; 
c = -10; 
d =  10; 
a0 = -1.45; 
b0 = 1.45; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 3, 
      n ≤ 
        6, \(n++\), 
aa = \(-1.45\); 
bb = 1.45; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ 
graph1 = 
        Plot[f[x], {x, a, \(-π\)\/2}, PlotStyle -> {Magenta}, 
          DisplayFunction -> Identity]; 
graph2 = 
        Plot[f[x], {x, \(-π\)\/2, π\/2}, PlotStyle -> {Magenta}, 
          DisplayFunction -> Identity]; 
graph3 = 
        Plot[f[x], {x, π\/2, b}, PlotStyle -> {Magenta}, 
          DisplayFunction -> Identity]; 
graph4 = 
        Plot[P[x], {x, \(-π\), π}, PlotStyle -> {Blue}, 
          DisplayFunction -> Identity]; 
Show[graph1, 
        graph2, graph3, graph4, Cdot, Rends, 
        PlotRange -> {{a, b}, {c, d}}, 
        Ticks -> {Range[\(-π\), π, π\/2], 
            Range[\(-6\), 10, 2]}, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 14.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 12.  Find the Newton polynomial approximation for  f[x] = (e^x), on the interval [-5, 5].

Solution 12.

f[x_] = (E^x);  
x0 = 0;  
a = -2.0; 
b =  3.0; 
c = -4.0; 
d = 20.5; 
a0 = -5.0; 
b0 = 5.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 2, 
      n ≤ 5, \(n++\), 
XY = 
        N[Table[{a + \(\(b - a\)\/n\) k, f[a + \(\(b - a\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; 
graph1 = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-2\), 3, 1], Range[0, 20, 5]}, 
          DisplayFunction -> Identity]; 
figure = 
        Show[graph1, Cdot, Rends, 
          DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 15.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 13.  Find the Newton polynomial approximation for  f[x] = ArcTan[x], on the interval [-2, 2].

Solution 13.

f[x_] = ArcTan[x];  
x0 = 0;  
a = -2; 
b =  2; 
c = -1.6; 
d =  1.6; 
a0 = -2; 
b0 = 2; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 3, 
      n ≤ 
        6, \(n++\), 
aa = \(-2.0\); 
bb = 
        2.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ \ 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-1.5\), 1.5, 0.5], 
              Range[\(-π\)\/4, π\/4, π\/4]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 18.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 14.  Find the Newton polynomial approximation for  f[x] = ArcSin[x], on the interval [-1, 1].

Solution 14.

f[x_] = ArcSin[x];  
x0 = 0;  
a = -1.5; 
b =  1.5; 
c = -2.2; 
d = 2.2; 
a0 = -1.0; 
b0 = 1.0; 
Needs["Graphics`Colors`"];
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 3, 
      n ≤ 
        6, \(n++\), 
aa = \(-1.0\); 
bb = 
        1.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; \ 
graph1 = 
        Plot[{f[x], P[x]}, {x, \(-0.9999\), 0.9999}, 
          PlotStyle -> Magenta, 
          DisplayFunction -> Identity]; 
graph2 = 
        Plot[P[x], {x, \(-1.5\), 1.5}, PlotStyle -> Blue, 
          DisplayFunction -> Identity]; 
Show[graph1, 
        graph2, Cdot, Rends, PlotRange -> {{a, b}, {c, d}}, 
        Ticks -> {Range[\(-1.5\), 1.5, 0.5], 
            Range[\(-π\)\/2, π\/2, π\/4]}, 
        DisplayFunction -> $DisplayFunction]; \ 
Show[
        graph, Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        P[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 19.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 15.  Find the Newton polynomial approximation for  f[x] = BesselJ[0, x], on the interval [-10, 10].

Solution 15.

f[x_] = BesselJ[0, x];  
x0 = 0;  
a = -10.2; 
b = 10.2; 
c = -0.8; 
d =  1.3; 
a0 = -11; 
b0 = 11; 
Needs["Graphics`Colors`"];
Cdot = 
        Graphics[{{Red, PointSize[0.02], 
              Point[{x0, f[x0]}]}}]; 
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[n = 5, 
      n ≤ 
        8, \(n++\), 
aa = \(-10\); 
bb = 
        10; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-10\), 10, 2], Range[\(-0.5\), 1, 0.5]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Show[graph, Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 20.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 16.  Find the Newton polynomial approximation for  f[x] = (1/(Sqrt[2  π])) (e^(-(x^2)/2)), on the interval [-3, 3].

Solution 16.

f[
          x_] = (1/(Sqrt[2  π])) (E^(-(x^2)/2));  
x0 = 0;  
a = -3.0; 
b =  3.0; 
c = -0.15; 
d = 0.5; 
a0 = -4.0; 
b0 = 4.0; 
Needs["Graphics`Colors`"];
Cdot = 
        Graphics[{{Red, PointSize[0.02], 
              Point[{x0, f[x0]}]}}]; 
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[i = 2, 
      i ≤ 5, \(i++\), 
n = 
        2  i; 
aa = \(-3.0\); 
bb = 
        3.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.02]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-3\), 3, 1], Range[0, 0.4, 0.1]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Show[graph, Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 22.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Example 17.  Find the Newton polynomial approximation for  f[x] = (1/2)  + (1/2) Erf[(x/(Sqrt[2]))], on the interval [-3, 3].

Solution 17.

f[x_] = (1/2)  + (1/2) Erf[(x/(Sqrt[2]))];  
x0 = 0;  
a = -3.0; 
b =  3.0; 
c = -0.25; 
d = 1.25; 
a0 = -4. ; 
b0 = 4.0; 
Needs["Graphics`Colors`"];
Cdot = 
        Graphics[{{Red, PointSize[0.02], 
              Point[{x0, f[x0]}]}}]; 
Rends = 
        Graphics[{{Green, Line[{{a0, c}, {a0, d}}]}, {Red, 
              Line[{{x0, c}, {x0, d}}]}, {Green, 
              Line[{{b0, c}, {b0, d}}]}}]; 
\(For[i = 2, 
      i ≤ 5, \(i++\), 
n = 
        2  i; 
aa = \(-3.0\); 
bb = 
        3.0; 
XY = 
        N[Table[{aa + \(\(bb - aa\)\/n\) k, 
              f[aa + \(\(bb - aa\)\/n\) k]}, {k, 0, 
              n}]]; 
Cdot = 
        ListPlot[XY, PlotStyle -> {Red, PointSize[0.011]}, 
          DisplayFunction -> Identity]; 
P[x_] = 
        NewtonPoly[XY]; 
graph = 
        Plot[{f[x], P[x]}, {x, a, b}, PlotStyle -> {Magenta, Blue}, 
          PlotRange -> {{a, b}, {c, d}}, 
          Ticks -> {Range[\(-3\), 3, 1], Range[0, 0.8, 0.25]}, 
          DisplayFunction -> Identity]; 
Show[graph, 
        Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ \ 
Show[graph, Cdot, Rends, 
        DisplayFunction -> $DisplayFunction]; \ 
Print["f[x] = ", f[x], "
", "P[x] = ", 
        Chop[Expand[P[x]]]];\ \ \ ]\ \ \)

Animation 23.  ( Newton Polynomial  Newton Polynomial ).  Internet hyperlink.

Animations (Newton Polynomial Approximation  Newton Polynomial Approximation).  Internet hyperlinks to animations.

Old Lab Project (Newton Interpolation Polynomial  Newton Interpolation Polynomial).  Internet hyperlinks to an old lab project.

Research Experience for Undergraduates

Newton Polynomials  Newton Polynomials  Internet hyperlinks to web sites and a bibliography of articles.