--- /dev/null
+function df = forward_euler(integer_order, h, f, x)
+ ##
+ ## Use the forward Euler method to compute the derivative of `f` at
+ ## a point `x`.
+ ##
+ ## INPUTS:
+ ##
+ ## * ``integer_order`` - The order of the derivative.
+ ##
+ ## * ``h`` - The step size.
+ ##
+ ## * ``f`` - The function whose derivative we're computing.
+ ##
+ ## * ``x`` - The point at which to compute the derivative.
+ ##
+
+ if (integer_order == 0)
+ df = x;
+ return;
+ end
+
+ ## We need a few points around `x` to compute the derivative at `x`.
+ ## The number of points depends on the order.
+ if (even(integer_order))
+ offset_b = integer_order / 2;
+ offset_f = offset_b;
+ else
+ ## When the order is odd, we need one more "forward" point than we
+ ## do "backward" points.
+ offset_b = (integer_order - 1) / 2;
+ offset_f = offset_b + 1;
+ end
+
+ backward_xs = [x-(offset_b*h) : h : x];
+
+ ## We'll always have at least one forward point, so start this vector
+ ## from (x + h) and include `x` itself in the backward points.
+ forward_xs = [x+h : h : x+(offset_f*h)];
+
+ xs = horzcat(backward_xs, forward_xs);
+
+ ## Now that we have all of the points that we need in xs, we can use
+ ## the divided_difference function to compute the derivative.
+ df = divided_difference(f, xs);
+end