nonlinear_cg.py 3.69 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
#
Martin Reinecke's avatar
Martin Reinecke committed
14
# Copyright(C) 2013-2018 Max-Planck-Society
15 16 17 18 19 20
#
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
# and financially supported by the Studienstiftung des deutschen Volkes.

from __future__ import division
from .minimizer import Minimizer
Martin Reinecke's avatar
Martin Reinecke committed
21
from .line_search_strong_wolfe import LineSearchStrongWolfe
22 23 24


class NonlinearCG(Minimizer):
Martin Reinecke's avatar
Martin Reinecke committed
25 26 27
    """ Nonlinear Conjugate Gradient scheme according to Polak-Ribiere.

    Algorithm 5.4 from Nocedal & Wright.
Martin Reinecke's avatar
Martin Reinecke committed
28

29 30 31 32
    Parameters
    ----------
    controller : IterationController
        Object that decides when to terminate the minimization.
33
    beta_heuristics : str
Martin Reinecke's avatar
Martin Reinecke committed
34 35 36 37 38
        One of 'Polak-Ribiere', 'Fletcher-Reeves', 'Hestenes-Stiefel' or '5.49'

    Notes
    -----
    No restarting procedure has been implemented yet.
39 40 41 42 43 44 45

    References
    ----------
    Jorge Nocedal & Stephen Wright, "Numerical Optimization", Second Edition,
    2006, Springer-Verlag New York
    """

46
    def __init__(self, controller, beta_heuristics='Polak-Ribiere'):
Martin Reinecke's avatar
Martin Reinecke committed
47
        valid_beta_heuristics = ['Polak-Ribiere', 'Fletcher-Reeves',
48
                                 'Hestenes-Stiefel', "5.49"]
49
        if not (beta_heuristics in valid_beta_heuristics):
Martin Reinecke's avatar
Martin Reinecke committed
50
            raise ValueError("beta heuristics must be either 'Polak-Ribiere', "
Martin Reinecke's avatar
Martin Reinecke committed
51
                             "'Fletcher-Reeves', 'Hestenes-Stiefel, or '5.49'")
52
        self._beta_heuristic = beta_heuristics
53
        self._controller = controller
54
        self._line_searcher = LineSearchStrongWolfe(c2=0.1)
55 56 57 58 59 60 61 62 63 64 65 66 67

    def __call__(self, energy):
        controller = self._controller
        status = controller.start(energy)
        if status != controller.CONTINUE:
            return energy, status
        f_k_minus_1 = None

        p = -energy.gradient

        while True:
            grad_old = energy.gradient
            f_k = energy.value
Martin Reinecke's avatar
stage 1  
Martin Reinecke committed
68 69 70 71
            energy, success = self._line_searcher.perform_line_search(
                energy, p, f_k_minus_1)
            if not success:
                return energy, controller.ERROR
72 73 74 75 76
            f_k_minus_1 = f_k
            status = self._controller.check(energy)
            if status != controller.CONTINUE:
                return energy, status
            grad_new = energy.gradient
Martin Reinecke's avatar
Martin Reinecke committed
77

78
            if self._beta_heuristic == 'Hestenes-Stiefel':
79
                # Eq. (5.46) in Nocedal & Wright.
Martin Reinecke's avatar
Martin Reinecke committed
80 81 82
                beta = max(0.0, (grad_new.vdot(grad_new-grad_old) /
                                 (grad_new-grad_old).vdot(p)).real)
            elif self._beta_heuristic == 'Polak-Ribiere':
83
                # Eq. (5.44) in Nocedal & Wright. (with (5.45) additionally)
Martin Reinecke's avatar
Martin Reinecke committed
84 85
                beta = max(0.0, (grad_new.vdot(grad_new-grad_old) /
                                 (grad_old.vdot(grad_old))).real)
86 87
            elif self._beta_heuristic == 'Fletcher-Reeves':
                # Eq. (5.41a) in Nocedal & Wright.
88
                beta = (grad_new.vdot(grad_new)/(grad_old.vdot(grad_old))).real
89 90 91 92
            else:
                # Eq. (5.49) in Nocedal & Wright.
                beta = (grad_new.vdot(grad_new) /
                        ((grad_new-grad_old).vdot(p))).real
93
            p = beta*p - grad_new