num = 32;
num2 = num (num - 1)/2;
step = 0.05;
ClearAll[\[Rho], Thiele]
\[Rho][x_List, y_List, i_Integer, n_Integer] := Module[{idx},
  If[n < 0,
   0
   ,
   If[n == 0,
    y[[i + 1]]
    ,
    idx = (num - 1 - n) (num - n)/2 + i + 1;
    If[r[[idx]] === Null,
     r[[idx]] = (x[[1 + i]] -
           x[[1 + i + n]])/(\[Rho][x, y, i, n - 1] - \[Rho][x, y,
            i + 1, n - 1]) + \[Rho][x, y, i + 1, n - 2];
     ];
    r[[idx]]
    ]
   ]
  ]
Thiele[x_List, y_List, xin_, n_Integer] := Module[{},
  If[n > num - 1,
   1
   ,
   \[Rho][x, y, 0, n] - \[Rho][x, y, 0, n - 2] + (xin - x[[n + 1]])/
     Thiele[x, y, xin, n + 1]
   ]
  ]
xval = Range[0, num - 1] step;
funcvals = Sin[xval];
r = ConstantArray[Null, num2];
6 Thiele[funcvals, xval, 0.5, 0]
funcvals = Cos[xval];
r = ConstantArray[Null, num2];
3 Thiele[funcvals, xval, 0.5, 0]
funcvals = Tan[xval];
r = ConstantArray[Null, num2];
4 Thiele[funcvals, xval, 1.0, 0]