##
## xs = [x0,x1,x2,x3,x4]
##
- ## f'(x=x1) ~= [f(x2)-f(x1)]/(x2-x1)
+ ## f'(x1) ~= [f(x2)-f(x1)]/(x2-x1)
##
- ## This would return [0, -1/(x2-x1), 2/(x2-x1), 0, 0]. This aids
- ## the solution of linear systems.
+ ## This would return [0, -1/(x2-x1), 2/(x2-x1), 0, 0]. This aids the
+ ## solution of linear systems.
##
+ ##
+ ## INPUTS:
+ ##
+ ## * ``integer_order`` - The order of the derivative which we're
+ ## approximating.
+ ##
+ ## * ``xs`` - The vector of x-coordinates.
+ ##
+ ## * ``x`` - The point `x` at which you'd like to evaluate the
+ ## derivative of the specified `integer_order`. This should be an
+ ## element of `xs`.
+ ##
+ ##
+ ## OUTPUTS:
+ ##
+ ## * ``coefficients`` - The vector of coefficients, in order, of
+ ## f(x0), f(x1), ..., f(xn).
+ ##
+
+ if (integer_order < 0)
+ ## You have made a grave mistake.
+ coefficients = NA;
+ return;
+ end
+
if (integer_order == 0)
- df = x;
+ coefficients = x;
+ return;
+ end
+
+ if (length(xs) < 2)
+ ## You can't approximate a derivative of order greater than zero
+ ## with zero or one points!
+ coefficients = NA
return;
end
trailing_zeros = zeros(1, trailing_zero_count);
targets = xs(first_nonzero_idx : last_nonzero_idx);
- cs = divided_difference_coefficients(targets);
- coefficients = horzcat(leading_zeros, cs, trailing_zeros);
+ # The multiplier comes from the Taylor expansion.
+ multiplier = factorial(integer_order);
+ cs = divided_difference_coefficients(targets) * multiplier;
+
+ coefficients = horzcat(leading_zeros, cs, trailing_zeros);
end