X-Git-Url: http://gitweb.michael.orlitzky.com/?p=octave.git;a=blobdiff_plain;f=newtons_method.m;fp=newtons_method.m;h=8bdd090945b49b1e99f5c7d9790f73ac218909da;hp=0000000000000000000000000000000000000000;hb=d7a73d33b0567c9abd72151b2972d7c0eb66e6a1;hpb=e99d78db1b964f4b3662b61012aa96d003722313

diff --git a/newtons_method.m b/newtons_method.m
new file mode 100644
index 0000000..8bdd090
--- /dev/null
+++ b/newtons_method.m
@@ -0,0 +1,43 @@
+function [root, iterations] = newtons_method(f, f_prime, epsilon, x0)
+  ## Find a root of the function `f` with initial guess x0.
+  ##
+  ## INPUTS:
+  ##
+  ##   * ``f`` - The function whose root we seek. Must return a column
+  ##     vector or scalar.
+  ##
+  ##   * ``f_prime`` - The derivative or Jacobian of ``f``.
+  ##
+  ##   * ``epsilon`` - We stop when the value of `f` becomes less
+  ##     than epsilon.
+  ##
+  ##   * ``x0`` - An initial guess. Either 1d or a column vector.
+  ##
+  ## OUTPUTS:
+  ##
+  ##   * ``root`` - The root that we found.
+  ##
+  ##   * ``iterations`` - The number of iterations that we performed
+  ##     during the search.
+  ##
+
+  if (columns(x0) > 1)
+    root = NA;
+    return;
+  end
+
+  iterations = 0;
+  xn = x0;
+  f_val = f(x0);
+
+  while (norm(f_val, Inf) > epsilon)
+    ## This uses the modified algorithm (2.11.5) in Atkinson.
+    ## Should work in 1d, too.    
+    delta = f_prime(xn) \ f_val;
+    xn = xn - delta;
+    f_val = f(xn);
+    iterations = iterations + 1;
+  end
+
+  root = xn;
+end