yat/utility/MultiMinimizerDerivative.h

Code
Comments
Other
Rev Date Author Line
4176 01 Jun 22 peter 1 #ifndef theplu_yat_utility_multi_minimizer_derivative
4176 01 Jun 22 peter 2 #define theplu_yat_utility_multi_minimizer_derivative
4176 01 Jun 22 peter 3
4176 01 Jun 22 peter 4 // $Id$
4176 01 Jun 22 peter 5
4346 24 Apr 23 peter 6 /*
4346 24 Apr 23 peter 7   Copyright (C) 2022, 2023 Peter Johansson
4346 24 Apr 23 peter 8
4346 24 Apr 23 peter 9   This file is part of the yat library, https://dev.thep.lu.se/yat
4346 24 Apr 23 peter 10
4346 24 Apr 23 peter 11   The yat library is free software; you can redistribute it and/or
4346 24 Apr 23 peter 12   modify it under the terms of the GNU General Public License as
4346 24 Apr 23 peter 13   published by the Free Software Foundation; either version 3 of the
4346 24 Apr 23 peter 14   License, or (at your option) any later version.
4346 24 Apr 23 peter 15
4346 24 Apr 23 peter 16   The yat library is distributed in the hope that it will be useful,
4346 24 Apr 23 peter 17   but WITHOUT ANY WARRANTY; without even the implied warranty of
4346 24 Apr 23 peter 18   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
4346 24 Apr 23 peter 19   General Public License for more details.
4346 24 Apr 23 peter 20
4346 24 Apr 23 peter 21   You should have received a copy of the GNU General Public License
4346 24 Apr 23 peter 22   along with yat. If not, see <https://www.gnu.org/licenses/>.
4346 24 Apr 23 peter 23 */
4346 24 Apr 23 peter 24
4252 18 Nov 22 peter 25 #include "multivariable/df.h"
4252 18 Nov 22 peter 26 #include "multivariable/f.h"
4176 01 Jun 22 peter 27
4176 01 Jun 22 peter 28 #include <yat/utility/Vector.h>
4176 01 Jun 22 peter 29 #include <yat/utility/VectorMutable.h>
4176 01 Jun 22 peter 30 #include <yat/utility/version.h>
4176 01 Jun 22 peter 31
4176 01 Jun 22 peter 32 #include <gsl/gsl_multimin.h>
4176 01 Jun 22 peter 33
4176 01 Jun 22 peter 34 #include <cassert>
4176 01 Jun 22 peter 35 #include <memory>
4176 01 Jun 22 peter 36
4176 01 Jun 22 peter 37 namespace theplu {
4176 01 Jun 22 peter 38 namespace yat {
4176 01 Jun 22 peter 39 namespace utility {
4176 01 Jun 22 peter 40
4176 01 Jun 22 peter 41   /**
4176 01 Jun 22 peter 42      \brief Wrapper class around gsl_multimin_fdfminimizer in GSL.
4176 01 Jun 22 peter 43
4176 01 Jun 22 peter 44      Class is the abstract base class for several minimisation
4176 01 Jun 22 peter 45      algorithms using the derivative.
4176 01 Jun 22 peter 46   */
4176 01 Jun 22 peter 47   class MultiMinimizerDerivative
4176 01 Jun 22 peter 48   {
4176 01 Jun 22 peter 49   public:
4176 01 Jun 22 peter 50     /// Copy is not allowed
4176 01 Jun 22 peter 51     MultiMinimizerDerivative(const MultiMinimizerDerivative&) = delete;
4176 01 Jun 22 peter 52
4176 01 Jun 22 peter 53     /// Assignment not allowed
4176 01 Jun 22 peter 54     MultiMinimizerDerivative&
4176 01 Jun 22 peter 55     operator=(const MultiMinimizerDerivative&) = delete;
4176 01 Jun 22 peter 56
4176 01 Jun 22 peter 57     /**
4176 01 Jun 22 peter 58        \brief Number of epochs (iterations) used in last minimisation.
4176 01 Jun 22 peter 59      */
4176 01 Jun 22 peter 60     unsigned int epochs(void) const;
4176 01 Jun 22 peter 61
4176 01 Jun 22 peter 62     /**
4176 01 Jun 22 peter 63        The size of the first trial step is given by \c step_size
4176 01 Jun 22 peter 64
4176 01 Jun 22 peter 65        \see gsl_multimin_fdfminimizer_set
4176 01 Jun 22 peter 66      */
4176 01 Jun 22 peter 67     double step_size(void) const;
4176 01 Jun 22 peter 68
4176 01 Jun 22 peter 69     /**
4176 01 Jun 22 peter 70        The size of the first trial step is given by \c step_size
4176 01 Jun 22 peter 71
4176 01 Jun 22 peter 72        \see gsl_multimin_fdfminimizer_set
4176 01 Jun 22 peter 73      */
4176 01 Jun 22 peter 74     void step_size(double ss);
4176 01 Jun 22 peter 75
4176 01 Jun 22 peter 76     /**
4176 01 Jun 22 peter 77        The accuracy of the line minimization is specified by \c
4176 01 Jun 22 peter 78        tolerance
4176 01 Jun 22 peter 79
4176 01 Jun 22 peter 80        \see <a href="https://www.gnu.org/software/gsl/doc/html/multimin.html#initializing-the-multidimensional-minimizer">gsl_multimin_fdfminimizer_set</a>
4176 01 Jun 22 peter 81      */
4176 01 Jun 22 peter 82     double tolerance(void) const;
4176 01 Jun 22 peter 83
4176 01 Jun 22 peter 84     /**
4176 01 Jun 22 peter 85        The accuracy of the line minimization is specified by \c
4176 01 Jun 22 peter 86        tolerance
4176 01 Jun 22 peter 87
4176 01 Jun 22 peter 88        \see <a href="https://www.gnu.org/software/gsl/doc/html/multimin.html#initializing-the-multidimensional-minimizer">gsl_multimin_fdfminimizer_set</a>
4176 01 Jun 22 peter 89      */
4176 01 Jun 22 peter 90     void tolerance(double tol);
4176 01 Jun 22 peter 91
4252 18 Nov 22 peter 92     /// Abstract class defining the interface for functions defining
4252 18 Nov 22 peter 93     /// the stop criteria in the minimization
4252 18 Nov 22 peter 94     class Stopper
4252 18 Nov 22 peter 95     {
4252 18 Nov 22 peter 96     public:
4252 18 Nov 22 peter 97       /// return true is search should stop
4252 18 Nov 22 peter 98       virtual bool operator()(const gsl_multimin_fdfminimizer*)=0;
4252 18 Nov 22 peter 99     };
4252 18 Nov 22 peter 100
4252 18 Nov 22 peter 101
4252 18 Nov 22 peter 102     /// wrapper around
4252 18 Nov 22 peter 103     /// <a href="https://www.gnu.org/software/gsl/doc/html/multimin.html#stopping-criteria">gsl_multimin_test_gradient</a>
4252 18 Nov 22 peter 104     class Gradient : public Stopper
4252 18 Nov 22 peter 105     {
4252 18 Nov 22 peter 106     public:
4252 18 Nov 22 peter 107       /// \param epsabs tolerance
4252 18 Nov 22 peter 108       explicit Gradient(double epsabs);
4252 18 Nov 22 peter 109
4252 18 Nov 22 peter 110       /**
4252 18 Nov 22 peter 111          \return \c true if <a
4252 18 Nov 22 peter 112          href="https://www.gnu.org/software/gsl/doc/html/multimin.html#stopping-criteria">gsl_multimin_test_gradient</a>(gradient,
4252 18 Nov 22 peter 113          epsabs) does not return \c GSL_CONTINUE, where \c gradient is
4252 18 Nov 22 peter 114          defined by \c gsl_multimin_fdfminimizer_gradient and \c
4252 18 Nov 22 peter 115          epsabs is defined in constructor.
4252 18 Nov 22 peter 116        */
4252 18 Nov 22 peter 117       bool operator()(const gsl_multimin_fdfminimizer*);
4252 18 Nov 22 peter 118     private:
4252 18 Nov 22 peter 119       double epsabs_;
4252 18 Nov 22 peter 120     };
4252 18 Nov 22 peter 121
4252 18 Nov 22 peter 122
4176 01 Jun 22 peter 123     /**
4176 01 Jun 22 peter 124        Function finds an \c x that minimizes the function defined by
4176 01 Jun 22 peter 125        \c func. It calls gsl_multimin_fdfminimizer_iterate until
4176 01 Jun 22 peter 126        either \c GSL_ENOPROG is returned or the norm of the gradient
4176 01 Jun 22 peter 127        is smaller than \c epsabs, as tested by
4176 01 Jun 22 peter 128        gsl_multimin_test_gradient.
4176 01 Jun 22 peter 129
4176 01 Jun 22 peter 130        Type Requirements:
4176 01 Jun 22 peter 131        - \c FUNC must have an operator defining the function
4176 01 Jun 22 peter 132        double operator()(const VectorBase& x)
4176 01 Jun 22 peter 133        - \c FUNC must have an operator defining the gradient
4176 01 Jun 22 peter 134        void operator()(const VectorBase& x, VectorMutable&)
4176 01 Jun 22 peter 135     */
4176 01 Jun 22 peter 136     template<class FUNC>
4338 14 Apr 23 peter 137     void operator()(VectorMutable&, FUNC& func,
4252 18 Nov 22 peter 138                     MultiMinimizerDerivative::Stopper&& stopper);
4176 01 Jun 22 peter 139
4176 01 Jun 22 peter 140     /**
4176 01 Jun 22 peter 141        Same as
4338 14 Apr 23 peter 142        operator()(VectorMutable&, FUNC& func, double epsabs)
4176 01 Jun 22 peter 143        but iterate at maximum \c max_epochs iterations.
4176 01 Jun 22 peter 144      */
4176 01 Jun 22 peter 145     template<class FUNC>
4338 14 Apr 23 peter 146     void operator()(VectorMutable&, FUNC& func,
4252 18 Nov 22 peter 147                     MultiMinimizerDerivative::Stopper&& stopper,
4252 18 Nov 22 peter 148                     unsigned int max_epochs);
4176 01 Jun 22 peter 149   protected:
4176 01 Jun 22 peter 150     /**
4176 01 Jun 22 peter 151        \brief Constructor
4176 01 Jun 22 peter 152
4176 01 Jun 22 peter 153        \param t defines type of GSL minimizer
4176 01 Jun 22 peter 154        \param size Number of dimensions in the input space.
4176 01 Jun 22 peter 155      */
4176 01 Jun 22 peter 156     MultiMinimizerDerivative(const gsl_multimin_fdfminimizer_type* t,
4176 01 Jun 22 peter 157                              size_t size);
4176 01 Jun 22 peter 158   private:
4176 01 Jun 22 peter 159     const gsl_multimin_fdfminimizer_type* type_;
4176 01 Jun 22 peter 160     unsigned int epochs_;
4176 01 Jun 22 peter 161     size_t size_;
4176 01 Jun 22 peter 162     struct GslFree
4176 01 Jun 22 peter 163     {
4176 01 Jun 22 peter 164       void operator()(gsl_multimin_fdfminimizer*) const;
4176 01 Jun 22 peter 165     };
4176 01 Jun 22 peter 166     std::unique_ptr<gsl_multimin_fdfminimizer, GslFree> solver_;
4176 01 Jun 22 peter 167     double step_size_;
4176 01 Jun 22 peter 168     double tol_;
4176 01 Jun 22 peter 169   };
4176 01 Jun 22 peter 170
4176 01 Jun 22 peter 171
4176 01 Jun 22 peter 172   // template implementation
4176 01 Jun 22 peter 173
4176 01 Jun 22 peter 174   template<class FUNC>
4252 18 Nov 22 peter 175   void
4338 14 Apr 23 peter 176   MultiMinimizerDerivative::operator()(VectorMutable& x,
4252 18 Nov 22 peter 177                                        FUNC& func,
4252 18 Nov 22 peter 178                                        MultiMinimizerDerivative::Stopper&& stopper)
4176 01 Jun 22 peter 179   {
4176 01 Jun 22 peter 180     unsigned int max_epochs = std::numeric_limits<unsigned int>::max();
4252 18 Nov 22 peter 181     (*this)(x, func, std::move(stopper), max_epochs);
4176 01 Jun 22 peter 182   }
4176 01 Jun 22 peter 183
4176 01 Jun 22 peter 184
4176 01 Jun 22 peter 185   template<class FUNC>
4252 18 Nov 22 peter 186   void
4338 14 Apr 23 peter 187   MultiMinimizerDerivative::operator()(VectorMutable& x,
4252 18 Nov 22 peter 188                                        FUNC& func,
4252 18 Nov 22 peter 189                                        MultiMinimizerDerivative::Stopper&& stopper,
4252 18 Nov 22 peter 190                                        unsigned int max_epochs)
4176 01 Jun 22 peter 191   {
4176 01 Jun 22 peter 192     assert(size_ == x.size());
4176 01 Jun 22 peter 193     gsl_multimin_function_fdf gsl_func;
4176 01 Jun 22 peter 194     gsl_func.n = size_;
4252 18 Nov 22 peter 195     gsl_func.f = multivariable::f<FUNC>;
4252 18 Nov 22 peter 196     gsl_func.df = multivariable::df<FUNC>;
4252 18 Nov 22 peter 197     gsl_func.fdf = multivariable::fdf<FUNC>;
4176 01 Jun 22 peter 198     gsl_func.params = &func;
4176 01 Jun 22 peter 199
4176 01 Jun 22 peter 200     gsl_multimin_fdfminimizer_set(solver_.get(), &gsl_func, x.gsl_vector_p(),
4176 01 Jun 22 peter 201                                   step_size_, tol_);
4176 01 Jun 22 peter 202
4176 01 Jun 22 peter 203     int status = 0;
4176 01 Jun 22 peter 204     for (epochs_=0; epochs_<max_epochs; ++epochs_) {
4176 01 Jun 22 peter 205       status = gsl_multimin_fdfminimizer_iterate(solver_.get());
4176 01 Jun 22 peter 206       if (status) {
4176 01 Jun 22 peter 207         if (status == GSL_ENOPROG)
4176 01 Jun 22 peter 208           break;
4338 14 Apr 23 peter 209         throw GSL_error("MultiMinimizerDerivative", status);
4176 01 Jun 22 peter 210       }
4176 01 Jun 22 peter 211
4252 18 Nov 22 peter 212       if (stopper(solver_.get()))
4176 01 Jun 22 peter 213         break;
4252 18 Nov 22 peter 214       //if (gsl_multimin_test_gradient(solver_->gradient,epsabs) != GSL_CONTINUE)
4252 18 Nov 22 peter 215       //break;
4176 01 Jun 22 peter 216     }
4176 01 Jun 22 peter 217
4176 01 Jun 22 peter 218     // copy result to passed x
4338 14 Apr 23 peter 219     VectorConstView view(solver_->x);
4176 01 Jun 22 peter 220     x = view;
4176 01 Jun 22 peter 221   }
4176 01 Jun 22 peter 222
4176 01 Jun 22 peter 223 }}}
4176 01 Jun 22 peter 224 #endif

Generated on Sat Oct 14 02:24:56 2023 (UTC) by svndigest 0.11.2pre (dev build 1672)