1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
17 
18 
19 
20 #ifndef LIBMESH_NLOPT_OPTIMIZATION_SOLVER_H
21 #define LIBMESH_NLOPT_OPTIMIZATION_SOLVER_H
22 
23 #include "libmesh/libmesh_config.h"
24 
25 // Petsc include files.
26 #if defined(LIBMESH_HAVE_NLOPT) && !defined(LIBMESH_USE_COMPLEX_NUMBERS)
27 
28 // Local includes
29 #include "libmesh/optimization_solver.h"
30 
31 // NLopt include
32 #include "nlopt.h"
33 
34 // C++ includes
35 
36 namespace libMesh
37 {
38 
39 // NLopt callback to set the objective function.
40 double __libmesh_nlopt_objective(unsigned n,
41                                  const double * x,
42                                  double * gradient,
43                                  void * data);
44 
45 // NLopt callback to set equality constraints.
46 void __libmesh_nlopt_equality_constraints(unsigned m,
47                                           double * result,
48                                           unsigned n,
49                                           const double * x,
50                                           double * gradient,
51                                           void * data);
52 
53 // NLopt callback to set inequality constraints.
54 void __libmesh_nlopt_inequality_constraints(unsigned m,
55                                             double * result,
56                                             unsigned n,
57                                             const double * x,
58                                             double * gradient,
59                                             void * data);
60 
61 /**
62  * This class provides an interface to the NLopt optimization solvers.
63  * http://ab-initio.mit.edu/wiki/index.php/NLopt
64  *
65  * \author David Knezevic
66  * \author John Peterson
67  * \date 2015
68  */
69 template <typename T>
70 class NloptOptimizationSolver : public OptimizationSolver<T>
71 {
72 public:
73 
74   /**
75    * The type of system that we use in conjunction with this solver.
76    */
77   typedef OptimizationSystem sys_type;
78 
79   /**
80    * Constructor.
81    */
82   explicit
83   NloptOptimizationSolver (sys_type & system);
84 
85   /**
86    * Destructor.
87    */
88   ~NloptOptimizationSolver ();
89 
90   /**
91    * Release all memory and clear data structures.
92    */
93   virtual void clear () override;
94 
95   /**
96    * Initialize data structures if not done so already.
97    */
98   virtual void init () override;
99 
100   /**
101    * \returns The raw NLopt object.
102    */
get_nlopt_object()103   nlopt_opt get_nlopt_object() { this->init(); return _opt; }
104 
105   /**
106    * Call the NLopt solver.
107    */
108   virtual void solve () override;
109 
110   /**
111    * Prints a useful message about why the latest optimization solve
112    * con(di)verged.
113    */
114   virtual void print_converged_reason() override;
115 
116   /**
117    * \returns The currently-available (or most recently obtained, if the NLopt object has
118    * been destroyed) convergence reason.  Refer to NLopt docs for the meaning of different
119    * the value.
120    */
121   virtual int get_converged_reason() override;
122 
123   /**
124    * \returns A writable reference to the current iteration count
125    * which can be incremented in the objective function.
126    */
get_iteration_count()127   unsigned & get_iteration_count() { return _iteration_count; }
128 
129 protected:
130 
131   /**
132    * Optimization solver context
133    */
134   nlopt_opt _opt;
135 
136   /**
137    * Store the result (i.e. convergence/divergence) for the most recent NLopt solve.
138    */
139   nlopt_result _result;
140 
141   /**
142    * Stores the current iteration index (incremented at each call of __libmesh_nlopt_objective).
143    */
144   unsigned _iteration_count;
145 
146   /**
147    * NLopt requires us to specify a tolerance for the constraints.
148    */
149   double _constraints_tolerance;
150 
151 private:
152 
153   // Make NLopt callback functions friends
154   friend double __libmesh_nlopt_objective (unsigned n,
155                                            const double * x,
156                                            double * gradient,
157                                            void * data);
158 
159   friend void __libmesh_nlopt_equality_constraints(unsigned m,
160                                                    double * result,
161                                                    unsigned n,
162                                                    const double * x,
163                                                    double * gradient,
164                                                    void * data);
165 
166   friend void __libmesh_nlopt_inequality_constraints(unsigned m,
167                                                      double * result,
168                                                      unsigned n,
169                                                      const double * x,
170                                                      double * gradient,
171                                                      void * data);
172 
173   // Map between strings and NLopt algorithms for command line parsing.
174   static std::map<std::string, nlopt_algorithm> _nlopt_algorithms;
175 
176   // Static function used to initialize the _nlopt_algorithms map, see
177   // below for its use.  Naming scheme:
178   // G/L == global/local optimization
179   // N/D == no/yes gradient required
180   // See the full list of algorithms at:
181   // http://ab-initio.mit.edu/wiki/index.php/NLopt_Algorithms
build_map()182   static std::map<std::string, nlopt_algorithm> build_map()
183   {
184     std::map<std::string, nlopt_algorithm> ret;
185     ret["LD_SLSQP"]                   = NLOPT_LD_SLSQP;
186     ret["LD_MMA"]                     = NLOPT_LD_MMA;
187     ret["LD_CCSAQ"]                   = NLOPT_LD_CCSAQ;
188     ret["LD_LBFGS"]                   = NLOPT_LD_LBFGS;
189     ret["LD_LBFGS_NOCEDAL"]           = NLOPT_LD_LBFGS_NOCEDAL;
190     ret["LD_TNEWTON"]                 = NLOPT_LD_TNEWTON;
191     ret["LD_TNEWTON_RESTART"]         = NLOPT_LD_TNEWTON_RESTART;
192     ret["LD_TNEWTON_PRECOND"]         = NLOPT_LD_TNEWTON_PRECOND;
193     ret["LD_TNEWTON_PRECOND_RESTART"] = NLOPT_LD_TNEWTON_PRECOND_RESTART;
194     ret["LD_AUGLAG"]                  = NLOPT_LD_AUGLAG;
195     ret["LD_VAR1"]                    = NLOPT_LD_VAR1;
196     ret["LD_VAR2"]                    = NLOPT_LD_VAR2;
197     ret["LN_COBYLA"]                  = NLOPT_LN_COBYLA;
198     ret["LN_BOBYQA"]                  = NLOPT_LN_BOBYQA;
199     ret["LN_PRAXIS"]                  = NLOPT_LN_PRAXIS;
200     ret["LN_NELDERMEAD"]              = NLOPT_LN_NELDERMEAD;
201     ret["LN_SBPLX"]                   = NLOPT_LN_SBPLX;
202     ret["GN_ISRES"]                   = NLOPT_GN_ISRES;
203     return ret;
204   }
205 };
206 
207 
208 
209 // Call the class-static function to define the class-static member.
210 // Since it's a template class, you actually do this in the header,
211 // not the source file.
212 template <typename T>
213 std::map<std::string, nlopt_algorithm>
214 NloptOptimizationSolver<T>::_nlopt_algorithms = NloptOptimizationSolver<T>::build_map();
215 
216 } // namespace libMesh
217 
218 
219 #endif // #if defined(LIBMESH_HAVE_NLOPT) && !defined(LIBMESH_USE_COMPLEX_NUMBERS)
220 #endif // LIBMESH_NLOPT_OPTIMIZATION_SOLVER_H
221