I am teaching myself Mathematica over summer, and coded up a simple module that computes numerical differentiation formulas automatically. Essentially, the module enables one to compute arbitrary differentiation rules like those explained here.
This module spits out the different formulae for a m-point approximation to the n-th derivative, and the leading error term.
Here is the (updated) module:
This seems to work alright for me (click to enlarge).
This module spits out the different formulae for a m-point approximation to the n-th derivative, and the leading error term.
Here is the (updated) module:
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 = (1/Factorial[m])*
D[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),
Superscript[f, m] e1 /. x -> (x0 + i h)}, {i, 0, m - 1}] //
Simplify ]]
points = Table[{x0 + i h, Subscript[f, i]}, {i, 0, m - 1}];
IP = InterpolatingPolynomial[points, x];
nthDeriv = D[IP, {x, n}];
e1 = (1/Factorial[m])*
D[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),
Superscript[f, m] e1 /. x -> (x0 + i h)}, {i, 0, m - 1}] //
Simplify ]]
This seems to work alright for me (click to enlarge).
AppDeriv[2, 1]
ReplyDeleteInterpolatingPolynomial::list: List expected at position {At Line = 3, the input was:,AppDeriv[2,1],1} in {At Line = 3, the input was:,AppDeriv[2,1],InterpolatingPolynomial[points,x]}. >>
You might want to learn programming a little more before posting such code.
You're right. I think I fixed it now. As bonus, it also gives the leading error term.
ReplyDelete