From: Bjørn Rustad <rustadbjornen@gmail.com>
Date: Wed, 28 Sep 2011 12:57:11 +0000 (+0200)
Subject: Add newton iteration
X-Git-Url: http://git.rustad.me/?a=commitdiff_plain;h=af6095cd790d8c6d2c86f0fbe1a0033b6c7728c3;p=nummat1

Add newton iteration
---

diff --git a/newton.py b/newton.py
new file mode 100644
index 0000000..8ea7538
--- /dev/null
+++ b/newton.py
@@ -0,0 +1,32 @@
+from function import function_g, gradient_g, hessian_g
+import generate
+import surf
+import numpy as np
+import matplotlib.pyplot as plt
+from numpy.linalg import norm, inv
+
+def newton_iter(x, grad_g, hess_g):
+    return x - inv(hess_g(x)).dot(grad_g(x))
+
+H = generate.spdmatrix(2, 3)
+b = generate.vector(2, -1, 1)
+c = generate.vector(2, -3, 3)
+
+alpha = 0.1
+xiter = np.array([1, 1])
+x0_norm = norm(gradient_g(xiter, H, b, c))
+
+tolerance = 0.0001
+
+relative_residual = norm(gradient_g(xiter, H, b, c)) / x0_norm
+relative_residuals = [relative_residual]
+
+while relative_residual > tolerance:
+    xiter = newton_iter(xiter, lambda x: gradient_g(x, H, b, c), lambda x:
+            hessian_g(x, H, b, c))
+    relative_residual = norm(gradient_g(xiter, H, b, c)) / x0_norm
+    relative_residuals.append(relative_residual)
+
+    print xiter
+
+print "Cirka null: ", gradient_g(xiter, H, b, c)