Newton's divided differences iteration [duplicate]

Question

I am trying to write a program that forms the interpolation polynomial for a given function on a given interval for any number of data points n. I wish to write a formula that will compute all of the necessary divided differences. Here is part of my code,

g[x_] := 1/(1 + x^2) 
f = Table[g[x], {x, -5.0, 5.0, 1}]
x = Table[a, {a, -5, 5, 1}]
 
Table[(f[[i]] - f[[i - 1]])/(x[[i + j]] - x[[i - 1]]), {j, 0, 9}, {i, 2, 11 - j}]

For j=0 I get a list of the correct outputs I need, {0.020362, 0.0411765, 0.1, 0.3, 0.5, -0.5, -0.3, -0.1, -0.0411765,-0.020362} but for j=1, I must find a way to have the terms in the numerator replaced with these output values over the loop, and so forth for each j. So for j=1,i=2 I wish to write a formula that will calculate (0.0411765-0.020362)/(x[[3]] - x[[1]]) for example, if this is possible.

I am a Mathematica novice, so bear with me and I hope I have been as clear as possible. I would appreciate any help, tips, tricks or guidance if it seems my approach is not a good one.

Answer 1

Here is one possibility (although a bit a tricky one):

We first define a function that returns the Newton polynomial as a function :

getNewton[xs_, ys_] := Module[{dd, i, newtoncoef},
  dd[{i_Integer}] := ys[[i]] ;(*dd implements Newton recursion*)
  dd[i : {__} ] := (dd[ Rest@i] - dd[Most@i])/(xs[[Last@i]] - 
      xs[[First@i]]);
  
  newtoncoef = 
   dd[Range[#]] & /@ Range[Length[xs]];(*Newton coefficients*)
  Function[x, newtoncoef.FoldList[#1 (x - #2) &, 1, Most@xs]]
  ]

Now we try it out. We first create some data and then apply our function to the x values:

g[x_] := 1/(1 + x^2);
ys = Table[g[x], {x, -5.0, 5.0, 1}];
xs = Table[x, {x, -5, 5, 1}];
Print["xs=", xs]; Print["ys=", ys];


Print["Newton poly at data points= ", getNewton[xs, ys] /@ xs]