X-Git-Url: http://gitweb.michael.orlitzky.com/?p=octave.git;a=blobdiff_plain;f=central_difference.m;fp=central_difference.m;h=ad58fd44546fcd3031ef2d5e7514ecd5d0c6c08b;hp=0000000000000000000000000000000000000000;hb=9ecd26a43862541b1c341cb318472d56ad0b382b;hpb=26295a1c102d026b8c539d2b2ef35cf6be27272c

diff --git a/central_difference.m b/central_difference.m
new file mode 100644
index 0000000..ad58fd4
--- /dev/null
+++ b/central_difference.m
@@ -0,0 +1,51 @@
+function coefficients = central_difference(xs, x)
+  ##
+  ## The first order central difference at x1 is,
+  ##
+  ##   f'(x1) = (f(x2) - f(x0))/2
+  ##
+  ## where the index x1 is of course arbitrary but x2, x0 are adjacent
+  ## to x1. The coefficients we seek are the coefficients of f(xj) for
+  ## j = 1,...,N-2, where N is the length of ``xs``. We omit the first
+  ## and last coefficient because at x0 and xN, the previous/next
+  ## value is not available.
+  ##
+  ## This should probably take an 'order' parameter as well; see
+  ## forward_euler().
+  ##
+  ## INPUT:
+  ##
+  ##   * ``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`.
+  ##
+  ## OUTPUT:
+  ##
+  ##   * ``coefficients`` - The vector of coefficients, in order, of
+  ##   f(x0), f(x1), ..., f(xn).
+  ##
+
+  if (length(xs) < 3)
+    ## We need at least one point other than the first and last.
+    coefficients = NA;
+    return;
+  end
+
+  x_idx = find(xs == x);
+
+  if (x_idx == 1 || x_idx == length(xs))
+    ## You asked for the difference at the first or last element, which
+    ## we can't do.
+    coefficients = NA;
+    return;
+  end
+
+  ## Start with a vector of zeros.
+  coefficients = zeros(1, length(xs));
+
+  ## And fill in the two values that we know.
+  coefficients(x_idx - 1) = -1/2;
+  coefficients(x_idx + 1) = 1/2;
+end