author Michael Orlitzky Fri, 14 Sep 2012 23:45:14 +0000 (19:45 -0400) committer Michael Orlitzky Fri, 14 Sep 2012 23:45:14 +0000 (19:45 -0400)
 forward_euler.m [new file with mode: 0644] patch | blob

diff --git a/forward_euler.m b/forward_euler.m
new file mode 100644 (file)
index 0000000..c6ef217
--- /dev/null
@@ -0,0 +1,45 @@
+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