1 /* 2 * This file is part of the GROMACS molecular simulation package. 3 * 4 * Copyright (c) 2020, by the GROMACS development team, led by 5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, 6 * and including many others, as listed in the AUTHORS file in the 7 * top-level source directory and at http://www.gromacs.org. 8 * 9 * GROMACS is free software; you can redistribute it and/or 10 * modify it under the terms of the GNU Lesser General Public License 11 * as published by the Free Software Foundation; either version 2.1 12 * of the License, or (at your option) any later version. 13 * 14 * GROMACS is distributed in the hope that it will be useful, 15 * but WITHOUT ANY WARRANTY; without even the implied warranty of 16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 17 * Lesser General Public License for more details. 18 * 19 * You should have received a copy of the GNU Lesser General Public 20 * License along with GROMACS; if not, see 21 * http://www.gnu.org/licenses, or write to the Free Software Foundation, 22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. 23 * 24 * If you want to redistribute modifications to GROMACS, please 25 * consider that scientific software is very special. Version 26 * control is crucial - bugs must be traceable. We will be happy to 27 * consider code for inclusion in the official distribution, but 28 * derived work must not be called official GROMACS. Details are found 29 * in the README & COPYING files - if they are missing, get the 30 * official version at http://www.gromacs.org. 31 * 32 * To help us fund GROMACS development, we humbly ask that you cite 33 * the research papers on the package. Check out http://www.gromacs.org. 34 */ 35 /*! \libinternal \file 36 * 37 * \brief Declare classes to aid Nelder-Mead downhill simplex optimisation. 38 * 39 * \author Christian Blau <blau@kth.se> 40 * \ingroup module_math 41 */ 42 43 #ifndef GMX_MATH_NELDERMEAD_H 44 #define GMX_MATH_NELDERMEAD_H 45 46 #include <functional> 47 #include <list> 48 #include <vector> 49 50 #include "gromacs/utility/arrayref.h" 51 #include "gromacs/utility/real.h" 52 53 namespace gmx 54 { 55 56 /*! \internal 57 * \brief Tie together coordinate and function value at this coordinate. 58 */ 59 struct RealFunctionvalueAtCoordinate 60 { 61 //! Vertex coordinate 62 std::vector<real> coordinate_; 63 //! Function value at this coordinate 64 real value_; 65 }; 66 67 /*! \internal 68 * \brief The simplex for the Nelder-Mead algorithm. 69 * 70 * Contains the N+1 simplex N-dimensional coordinates and its function values. 71 * Allows for simplex manipulations as needed for the Nelder-Mead algorithm. 72 * 73 * \note Keeps the simplex sorted according to function values with the simplex 74 * at the lowest function value first. 75 */ 76 class NelderMeadSimplex 77 { 78 public: 79 /*! \brief Set up Nelder-Mead simplex from an initial guess. 80 * 81 * \note Triggers N+1 function evaluations at all simplex points. 82 * 83 * \param[in] f the function to be evaluated 84 * \param[in] initalGuess initial guess of the coordinates. 85 * 86 */ 87 NelderMeadSimplex(const std::function<real(ArrayRef<const real>)>& f, ArrayRef<const real> initalGuess); 88 89 //! Return the vertex with the lowest function value at any of the simplex vertices. 90 const RealFunctionvalueAtCoordinate& bestVertex() const; 91 92 //! Return the vertex of the simplex with the highest (worst) function value. 93 const RealFunctionvalueAtCoordinate& worstVertex() const; 94 95 //! Return the second largest function value at any of the simplex vertices. 96 real secondWorstValue() const; 97 98 //! Return the reflection point and the evaluated function value at this point. 99 RealFunctionvalueAtCoordinate 100 evaluateReflectionPoint(const std::function<real(ArrayRef<const real>)>& f) const; 101 102 //! Evaluate and return the expansion point and function value. 103 RealFunctionvalueAtCoordinate evaluateExpansionPoint(const std::function<real(ArrayRef<const real>)>& f) const; 104 105 //! Evaluate and return the contraction point and function value. 106 RealFunctionvalueAtCoordinate 107 evaluateContractionPoint(const std::function<real(ArrayRef<const real>)>& f) const; 108 109 /*! \brief Replace the simplex vertex with the largest function value. 110 * 111 * \param[in] newVertex to replace the worst vertex with 112 * \note keeps the simplex list sorted and reevaluates the reflection point 113 */ 114 void swapOutWorst(const RealFunctionvalueAtCoordinate& newVertex); 115 116 /*! \brief Shrink the simplex. 117 * 118 * All points move closer to the best point by a factor \f$\sigma\f$. 119 * 120 * Replace all point coordinates, except the best, with 121 * \f$x_i = x_{\mathrm{best}} + \sigma (x_i - x_{\mathrm{best}})\f$ 122 */ 123 void shrinkSimplexPointsExceptBest(const std::function<real(ArrayRef<const real>)>& f); 124 125 /*! \brief The oriented length of the vertex. 126 * 127 * The oriented length of the simplex is defined as the largest distance 128 * between the first simplex vertex coordinate (with the lowest, best function 129 * value) and any other simplex coordinate. 130 * 131 * The oriented length is used as a computationally fast and simple 132 * convergence criterion because it is proven that 133 * orientedLegnth < simplex_diameter < 2 * orientedLength 134 * 135 */ 136 real orientedLength() const; 137 138 private: 139 /*! \brief Update centroid and reflection point. 140 * 141 * The arithmetic mean of all vertex coordinates expect the one at the 142 * highest (worst) function value. 143 * 144 */ 145 void updateCentroidAndReflectionPoint(); 146 147 /*! \brief The points of the simplex with the function values. 148 * \note This list stays sorted according to function value during the 149 * life-time of this object. 150 */ 151 std::list<RealFunctionvalueAtCoordinate> simplex_; 152 153 //! The centroid of the simplex, skipping the worst point is updated once the simplex changes 154 std::vector<real> centroidWithoutWorstPoint_; 155 156 //! The reflection point and its function value is updated once the simplex changes 157 std::vector<real> reflectionPointCoordinates_; 158 }; 159 160 161 } // namespace gmx 162 163 #endif // GMX_MATH_NELDERMEAD_H 164