Add files via upload

Author: Paulnkk

armijo.py

@@ -0,0 +1,15 @@
+import numpy as np
+from scipy import optimize as opt
+
+"""
+Armijo linesearch Algorithm
+"""
+
+def schrittweite(f, fd, xk, d, sigma, rho, gamma):
+
+    dnorm = (np.linalg.norm(d))**2
+    t = - gamma * (np.dot(fd(xk), d) / dnorm)
+    while f(xk + t * d) > f(xk) + t * sigma * np.dot(fd(xk), d):
+        t = rho * t
+
+    return t

bfgs.py

@@ -0,0 +1,57 @@
+import numpy as np
+from numpy.linalg import inv
+from scipy import optimize as opt
+import math
+import armijo
+
+class BFGS(object):
+    """
+    Constructor BFGS
+    """
+    def __init__ (self, f, fd, H, xk, eps):
+        self.fd = fd
+        self.H = H
+        self.xk = xk
+        self.eps = eps
+        self.f = f
+        return
+    """
+    BFGS-Method 
+    """
+
+    def work (self):
+        f = self.f
+        fd = self.fd
+        H = self.H
+        xk = self.xk
+        eps = self.eps
+        """
+        Initial Matrix for BFGS (Identitiy Matrix)
+        """
+        E =  np.array([   [1.,         0.],
+                            [0.,         1.] ])
+        xprev = xk
+        it = 0
+        maxit = 10000
+
+        while (np.linalg.norm(fd(xk)) > eps) and (it < maxit):
+            Hfd = inv(E)@fd(xprev)
+            xk = xprev - Hfd
+            sk = np.subtract(xk, xprev)
+            yk = np.subtract(fd(xk), fd(xprev))
+
+            b1 = (1 / np.dot(yk, sk))*(np.outer(yk, yk))
+            sub1b2 = np.outer(sk, sk)
+            Esk = E @ (sk)
+            sub2b2 = (1 / np.dot(sk, Esk))
+            sub3b2 = np.matmul(E, sub1b2)
+            sub4b2 = np.matmul(sub3b2, E)
+            b2 = sub2b2 * sub4b2
+            E1 = np.add(E, b1)
+            E = np.subtract(E1, b2)
+
+            xprev = xk
+            print("Log-Values(BFGS): ", math.log10(f(xk)))
+            it += 1
+
+        return xk, it

cg.py

@@ -0,0 +1,45 @@
+import numpy as np
+from numpy.linalg import inv
+from scipy import optimize as opt
+import math
+import armijo
+
+class CGMethod(object):
+    """
+    Constructor CG-Method
+    """
+    def __init__ (self, f, fd, H, xk, eps):
+        self.fd = fd
+        self.f = f
+        self.H = H
+        self.xk = xk
+        self.eps = eps
+        return
+    """
+    CG-Method
+    """
+    def work (self):
+        f = self.f
+        fd = self.fd
+        H = self.H
+        xk = self.xk
+        eps = self.eps
+        t = 1
+        it = 0
+        #maxit = 10000
+        dprev = - fd(xk)
+        xprev = self.xk
+        dk    = None
+
+        while(np.linalg.norm(fd(xk)) > eps): #and (it < maxit):
+            t = armijo.schrittweite(f, fd, xprev, dprev, sigma = 0.02, rho = 0.5, gamma = 2)
+            xk = xprev + t * dprev
+            normprev = (np.linalg.norm(fd(xprev)))**2
+            xprev = xk
+            normk = (np.linalg.norm(fd(xk)))**2
+            dk = - fd(xk) + (normk / normprev)*(dprev)
+            dprev = dk
+            print("Log-Values(CG): ", math.log10(f(xk)))
+            it += 1
+
+        return xk, it

gradv.py

@@ -0,0 +1,35 @@
+import numpy as np
+from scipy import optimize as opt
+import armijo
+import math
+
+
+
+class GradientMethod(object):
+    """
+    Constructor Gradient-Method
+    """
+    def __init__ (self, f, fd, H, xk, eps):
+        self.xk = xk
+        self.eps = eps
+        self.fd = fd
+        self.f = f
+        self.H = H
+        return
+    """
+    Gradient-Method
+    """
+
+    def work (self):
+        f = self.f
+        xk = self.xk
+        eps = self.eps
+        fd = self.fd
+        it = 0
+        maxit = 2000
+        while (np.linalg.norm(fd(xk)) > eps) and (it < maxit):
+            t = armijo.schrittweite(f, fd, xk, -fd(xk), sigma = 0.02, rho = 0.5, gamma = 0.0001)
+            xk = xk - t * fd(xk) 
+            print("Log-Values(Gradient): ", math.log10(f(xk)))
+            it += 1
+        return xk, it

newtonm.py

@@ -0,0 +1,36 @@
+import numpy as np
+from numpy.linalg import inv
+from scipy import optimize as opt
+import math
+
+
+class NewtonMethod(object):
+    """
+    Constructor Newton-Method
+    """
+    def __init__ (self, f, fd, H, xk, eps):
+        self.fd = fd
+        self.H = H
+        self.xk = xk
+        self.eps = eps
+        self.f = f
+        return
+    """
+    Newton-Method
+    """
+    def work (self):
+        f = self.f
+        fd = self.fd
+        H = self.H
+        xk = self.xk
+        eps = self.eps
+        it = 0
+        #maxit = 10000
+
+        while (np.linalg.norm(fd(xk)) > eps): #and (it < maxit):
+            Hfd = inv(H(xk))@(fd(xk))
+            xk = xk - Hfd
+            it += 1
+            print("Log-Values(Newton): ", math.log10(f(xk)))
+
+        return xk, it

ros_test.py

@@ -0,0 +1,117 @@
+from gradv import GradientMethod
+import numpy as np
+import time
+from newtonm import NewtonMethod
+from numpy.linalg import inv
+from cg import CGMethod
+from bfgs import BFGS
+import armijo
+
+
+
+class Rosenb:
+
+    """
+    Select Solver (Gradient-Method, Newton-Method, CG-Method, BFGS-Method)
+    """
+
+
+    def algo(self, fo, alsolver):
+        x0 = fo.x0
+        f = fo.f
+        fd = fo.fd
+        H = fo.H
+        x = fo.x0
+        eps = 0.0001 
+        it = None
+        while(np.linalg.norm(fd(x)) > eps):
+            solver = None
+            if(alsolver == 'gradv'):
+                #Nutze Gradientenverfahren Loesung
+                solver = GradientMethod(f, fd, H, x, eps)
+            elif(alsolver == 'newtonm'):
+                solver = NewtonMethod(f, fd, H, x, eps)
+            elif(alsolver == 'cg'):
+                solver = CGMethod(f, fd, H, x, eps)
+            elif(alsolver == 'bfgs'):
+                solver = BFGS(f, fd, H, x, eps)
+            else:
+                print("ERROR, NO SOLVER SELECTED. EXITING")
+                exit()
+
+            x, it = solver.work()
+
+        return x, it
+
+"""
+Create Function that should be minimized
+"""
+class Function(object):
+    x0 = None
+    f = None
+    fd = None
+    H = None
+
+    def __init__(self, x0, f, fd, H):
+        self.x0 = x0
+        self.f = f
+        self.fd = fd
+        self.H = H
+    """
+    Rosenbrock function
+    """
+    def testFunctionRosenbrock():
+        f = lambda xy: (10*(xy[0] - xy[1]**2))**2 + (1-xy[0])**2
+        fd = lambda xy: np.array([202.*xy[0] - 200*xy[1]**2 - 2, -400*xy[1]*(xy[0] - xy[1]**2)])
+        H = lambda xy: np.array([   [202.,         -400.*xy[1]                              ],
+                                    [-400.*xy[1],   800.*xy[1]**2 - 400.*(xy[0] - xy[1]**2) ]
+                                ])
+        x0 = np.array([-2.,2.])
+        return Function(x0,f,fd,H)
+
+    """
+    Simple quadratic function ((x,y) -> (x + 3)^2 + y^2)
+    """
+    def testFunction2():
+        f = lambda xy: (xy[0] + 3)**2 + xy[1]**2
+        fd = lambda xy: np.array([2*xy[0] + 6, 2*xy[1]])
+        H = lambda xy: np.array([[2., 0.],[0., 2.]])
+        x0 = np.array([3.,0.])
+        return Function(x0,f,fd,H)
+
+
+print("STARTING TESTS")
+a = Rosenb()
+"""
+print("-----Testing quadratic(2) function-----")
+f = Function.testFunction2()
+x_gv, it_gv = a.algo(f,'gradv')
+x_nm, it_nm = a.algo(f, 'newtonm')
+x_cg, it_cg = a.algo(f, 'cg')
+#x_bfgs, it_bfgs = a.algo(f, 'bfgs')
+print("Result (Grad): ",x_gv, " steps: ", it_gv)
+print("Result (Newton): ",x_nm, " steps: ", it_nm)
+print("Result (CG): ",x_cg, "steps: ", it_cg)
+#print("Result (BFGS): ",x_bfgs, "steps: ", it_bfgs)
+assert np.allclose(x_gv, np.array([-3.,0.]), rtol=1e-04, atol=1e-05)
+assert np.allclose(x_nm, np.array([-3.,0.]), rtol=1e-04, atol=1e-05)
+assert np.allclose(x_cg, np.array([-3.,0.]), rtol=1e-04, atol=1e-05)
+#assert np.allclose(x_bfgs, np.array([-3.,0.]), rtol=1e-04, atol=1e-05)
+print(" + Pass")
+"""
+
+print("-----Testing Rosenbrock function-----")
+f = Function.testFunctionRosenbrock()
+x_gv, it_gv = a.algo(f,'gradv')
+x_nm, it_nm = a.algo(f, 'newtonm')
+x_cg, it_cg = a.algo(f, 'cg')
+x_bfgs, it_bfgs = a.algo(f, 'bfgs')
+print("Result (Grad): ",x_gv, " steps: ", it_gv)
+print("Result (Newton): ",x_nm, " steps: ", it_nm)
+print("Result (CG): ",x_cg, "steps: ", it_cg)
+print("Result (BFGS): ",x_bfgs, "steps: ", it_bfgs)
+assert np.allclose(x_gv, np.array([1.,1.]), rtol=1e-04, atol=1e-05)
+assert np.allclose(x_nm, np.array([1.,1.]), rtol=1e-04, atol=1e-05)
+assert np.allclose(x_cg, np.array([1.,1.]), rtol=1e-04, atol=1e-05)
+assert np.allclose(x_bfgs, np.array([1.,1.]), rtol=1e-04, atol=1e-05)
+print(" + Pass")