Thursday, January 3, 2013

Numerical Differentiation using Mathematica: Redux

While I was demonstrating this Mathematica script in class, I noticed there was an obvious error that had escaped into the wild. In particular, when I tried to get the center-difference approximation to the second derivative via AppDeriv[3, 2], I got the following:

Notice the "zero" truncation error for the second formula. That should have tipped me off (if I had done any testing, duh). The other errors are also incorrect. It can be fixed relatively easily by the following code.

AppDeriv[m_, n_] := Module[{points, IP, nthDeriv, e1},
  points = Table[{x0 + i h, Subscript[f, i]}, {i, 0, m - 1}];
  IP = InterpolatingPolynomial[points, x];
  nthDeriv = D[IP, {x, n}];
  e1 = D[F[X[x]]* Apply[Times,
            (x - Table[{x0 + i h}, {i, 0, m - 1}])], {x, n}];
  Formula = TableForm[
    Table[{Subsuperscript[f, i, n], nthDeriv /. x -> (x0 + i h), 
       O Superscript[h, Exponent[e1[[1]] /. 
       x -> (x0 + i h), h]]}, {i,0, m - 1}] // Simplify ]]

Here is my Mathematica screenshot with the new code. It at least seems to fix the problem above.

No comments: