reinstate RelaxedNewton

nifty5/__init__.py

@@ -60,7 +60,8 @@ from .minimization.minimizer import Minimizer
 from .minimization.conjugate_gradient import ConjugateGradient
 from .minimization.nonlinear_cg import NonlinearCG
 from .minimization.descent_minimizers import (
-    DescentMinimizer, SteepestDescent, VL_BFGS, L_BFGS, NewtonCG)
+    DescentMinimizer, SteepestDescent, VL_BFGS, L_BFGS, RelaxedNewton,
+    NewtonCG)
 from .minimization.scipy_minimizer import (ScipyMinimizer, L_BFGS_B, ScipyCG)
 from .minimization.energy import Energy
 from .minimization.quadratic_energy import QuadraticEnergy

nifty5/minimization/descent_minimizers.py

@@ -135,6 +135,24 @@ class SteepestDescent(DescentMinimizer):
         return -energy.gradient
 
 
+class RelaxedNewton(DescentMinimizer):
+    """ Calculates the descent direction according to a Newton scheme.
+
+    The descent direction is determined by weighting the gradient at the
+    current parameter position with the inverse local metric.
+    """
+
+    def __init__(self, controller, line_searcher=None):
+        if line_searcher is None:
+            line_searcher = LineSearchStrongWolfe(
+                preferred_initial_step_size=1.)
+        super(RelaxedNewton, self).__init__(controller=controller,
+                                            line_searcher=line_searcher)
+
+    def get_descent_direction(self, energy):
+        return -energy.metric.inverse_times(energy.gradient)
+
+
 class NewtonCG(DescentMinimizer):
     """ Calculates the descent direction according to a Newton-CG scheme.
 

nifty5/minimization/energy.py

@@ -100,6 +100,15 @@ class Energy(NiftyMetaBase()):
             self._gradnorm = self.gradient.norm()
         return self._gradnorm
 
+    @property
+    def metric(self):
+        """
+        LinearOperator : implicitly defined metric.
+            A positive semi-definite operator or function describing the
+            metric of the potential at the given `position`.
+        """
+        raise NotImplementedError
+
     def apply_metric(self, x):
         """
         Parameters

nifty5/minimization/energy_adapter.py

@@ -37,5 +37,9 @@ class EnergyAdapter(Energy):
     def gradient(self):
         return self._grad
 
+    @property
+    def metric(self):
+        return self._metric
+
     def apply_metric(self, x):
         return self._metric(x)

nifty5/minimization/kl_energy.py

@@ -48,6 +48,10 @@ class KL_Energy(Energy):
     def apply_metric(self, x):
         return self._metric(x)
 
+    @property
+    def metric(self):
+        return self._metric
+
     @property
     def samples(self):
         return self._samples

nifty5/minimization/quadratic_energy.py

@@ -74,5 +74,9 @@ class QuadraticEnergy(Energy):
     def gradient(self):
         return self._grad
 
+    @property
+    def metric(self):
+        return self._A
+
     def apply_metric(self, x):
         return self._A(x)

test/test_minimization/test_minimizers.py

@@ -38,6 +38,7 @@ minimizers = ['ift.VL_BFGS(IC)',
               'ift.L_BFGS(IC)',
               'ift.NewtonCG(IC)']
 
+newton_minimizers = ['ift.RelaxedNewton(IC)']
 quadratic_only_minimizers = ['ift.ConjugateGradient(IC)',
                              'ift.ScipyCG(tol=1e-5, maxiter=300)']
 slow_minimizers = ['ift.SteepestDescent(IC)']
@@ -45,8 +46,8 @@ slow_minimizers = ['ift.SteepestDescent(IC)']
 
 class Test_Minimizers(unittest.TestCase):
 
-    @expand(product(minimizers + quadratic_only_minimizers + slow_minimizers,
-                    spaces))
+    @expand(product(minimizers + newton_minimizers +
+                    quadratic_only_minimizers + slow_minimizers, spaces))
     def test_quadratic_minimization(self, minimizer, space):
         np.random.seed(42)
         starting_point = ift.Field.from_random('normal', domain=space)*10
@@ -69,7 +70,7 @@ class Test_Minimizers(unittest.TestCase):
                         1./covariance_diagonal.local_data,
                         rtol=1e-3, atol=1e-3)
 
-    @expand(product(minimizers))
+    @expand(product(minimizers+newton_minimizers))
     def test_rosenbrock(self, minimizer):
         try:
             from scipy.optimize import rosen, rosen_der, rosen_hess_prod
@@ -93,6 +94,25 @@ class Test_Minimizers(unittest.TestCase):
                 out = ift.Field.from_global_data(space, rosen_der(inp))
                 return out
 
+            @property
+            def metric(self):
+                class RBCurv(ift.EndomorphicOperator):
+                    def __init__(self, loc):
+                        self._loc = loc.to_global_data_rw()
+                        self._capability = self.TIMES
+                        self._domain = space
+
+                    def apply(self, x, mode):
+                        self._check_input(x, mode)
+                        inp = x.to_global_data_rw()
+                        out = ift.Field.from_global_data(
+                            space, rosen_hess_prod(self._loc.copy(), inp))
+                        return out
+
+                t1 = ift.GradientNormController(tol_abs_gradnorm=1e-5,
+                                                iteration_limit=1000)
+                return ift.InversionEnabler(RBCurv(self._position), t1)
+
             def apply_metric(self, x):
                 inp = x.to_global_data_rw()
                 pos = self._position.to_global_data_rw()
@@ -149,7 +169,7 @@ class Test_Minimizers(unittest.TestCase):
         assert_allclose(energy.position.local_data, 0.,
                         atol=1e-3)
 
-    @expand(product(minimizers+slow_minimizers))
+    @expand(product(minimizers+newton_minimizers+slow_minimizers))
     def test_cosh(self, minimizer):
         space = ift.UnstructuredDomain((1,))
         starting_point = ift.Field.full(space, 3.)
@@ -168,6 +188,13 @@ class Test_Minimizers(unittest.TestCase):
                 x = self.position.to_global_data()[0]
                 return ift.Field.full(self.position.domain, np.sinh(x))
 
+            @property
+            def metric(self):
+                x = self.position.to_global_data()[0]
+                v = np.cosh(x)
+                return ift.DiagonalOperator(
+                    ift.Field.full(self.position.domain, v))
+
             def apply_metric(self, x):
                 p = self.position.to_global_data()[0]
                 v = np.cosh(p)