All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
mmin.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2014, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_MULTI_MIN_H
24 #define O2SCL_MULTI_MIN_H
25 
26 /** \file mmin.h
27  \brief Function object for gradients and \ref o2scl::mmin_base
28 */
29 
30 #include <o2scl/multi_funct.h>
31 #include <o2scl/mm_funct.h>
32 #include <o2scl/string_conv.h>
33 
34 #ifndef DOXYGEN_NO_O2NS
35 namespace o2scl {
36 #endif
37 
38  /// Array of multi-dimensional functions typedef
39  typedef std::function<int(size_t,boost::numeric::ublas::vector<double> &,
42 
43  /** \brief Class for automatically computing gradients [abstract base]
44 
45  Default template arguments
46  - \c func_t - (no default)
47  - \c vec_t - \ref boost::numeric::ublas::vector < double >
48 
49  \future Consider making an exact_grad class for computing exact
50  gradients.
51  */
52  template<class func_t, class vec_t=boost::numeric::ublas::vector<double> >
53  class gradient {
54 
55  public:
56 
57  // (Need to have empty default constructor since we
58  // have private copy constructor)
59  gradient() {}
60 
61  virtual ~gradient() {}
62 
63  /// Set the function to compute the gradient of
64  virtual int set_function(func_t &f) {
65  func=&f;
66  return 0;
67  }
68 
69  /** \brief Compute the gradient \c g at the point \c x
70  */
71  virtual int operator()(size_t nv, vec_t &x, vec_t &g)=0;
72 
73 #ifndef DOXYGEN_INTERNAL
74 
75  protected:
76 
77  /// A pointer to the user-specified function
78  func_t *func;
79 
80  private:
81 
82  gradient(const gradient &);
83  gradient& operator=(const gradient&);
84 
85 #endif
86 
87  };
88 
89  /** \brief Simple automatic computation of gradient by finite
90  differencing
91 
92  \comment
93  \endcomment
94  */
95  template<class func_t, class vec_t> class gradient_gsl :
96  public gradient<func_t,vec_t> {
97 
98  public:
99 
100  gradient_gsl() {
101  epsrel=1.0e-6;
102  epsmin=1.0e-15;
103  }
104 
105  virtual ~gradient_gsl() {}
106 
107  /** \brief The relative stepsize for finite-differencing
108  (default \f$ 10^{-6} \f$ )
109  */
110  double epsrel;
111 
112  /// The minimum stepsize (default \f$ 10^{-15} \f$)
113  double epsmin;
114 
115  /** \brief Compute the gradient \c g at the point \c x
116  */
117  virtual int operator()(size_t nv, vec_t &x, vec_t &g) {
118  double fv1, fv2, h;
119 
120  fv1=(*this->func)(nv,x);
121 
122  for(size_t i=0;i<nv;i++) {
123 
124  h=epsrel*fabs(x[i]);
125  if (fabs(h)<=epsmin) h=epsrel;
126 
127  x[i]+=h;
128  fv2=(*this->func)(nv,x);
129  x[i]-=h;
130  g[i]=(fv2-fv1)/h;
131 
132  }
133 
134  return 0;
135  }
136 
137  };
138 
139  /** \brief Multidimensional minimization [abstract base]
140 
141  <b>The template parameters:</b>
142  The template parameter \c func_t specifies the function to
143  min and should be a class containing a definition
144  \code
145  func_t::operator()(size_t nv, const vec_t &x, double &f);
146  \endcode
147  where \c f is the value of the function at \c x ,
148  where \c x is a array-like class defining \c operator[] of size \c nv.
149  The parameter \c dfunc_t (if used) should provide the gradient with
150  \code
151  func_t::operator()(size_t nv, vec_t &x, vec_t &g);
152  \endcode
153  where \c g is the gradient of the function at \c x.
154 
155  Verbose I/O is sent through \c std::cout and \c std::cin by
156  default, but this can be modified using \ref
157  set_verbose_stream(). Note that this function
158  stores pointers to the user-specified output streams,
159  and these pointers are not copied in child copy
160  constructors.
161  */
162  template<class func_t=multi_funct11, class dfunc_t=func_t,
164 
165 #ifndef DOXYGEN_INTERNAL
166 
167  protected:
168 
169  /// Stream for verbose output
170  std::ostream *outs;
171 
172  /// Stream for verbose input
173  std::istream *ins;
174 
175 #endif
176 
177  public:
178 
179  mmin_base() {
180  verbose=0;
181  ntrial=100;
182  tol_rel=1.0e-4;
183  tol_abs=1.0e-4;
184  last_ntrial=0;
185  err_nonconv=true;
186  outs=&std::cout;
187  ins=&std::cin;
188  }
189 
190  virtual ~mmin_base() {}
191 
192  /// Output control
193  int verbose;
194 
195  /// Maximum number of iterations
196  int ntrial;
197 
198  /// Function value tolerance
199  double tol_rel;
200 
201  /// The independent variable tolerance
202  double tol_abs;
203 
204  /// The number of iterations for in the most recent minimization
206 
207  /// If true, call the error handler if the routine does not "converge"
209 
210  /** \brief Set streams for verbose I/O
211 
212  Note that this function stores pointers to the user-specified
213  output streams, and these pointers are not copied in child copy
214  constructors.
215  */
216  int set_verbose_stream(std::ostream &out, std::istream &in) {
217  outs=&out;
218  ins=&in;
219  return 0;
220  }
221 
222  /** \brief Calculate the minimum \c min of \c func w.r.t. the
223  array \c x of size \c nvar.
224  */
225  virtual int mmin(size_t nvar, vec_t &x, double &fmin,
226  func_t &func)=0;
227 
228  /** \brief Calculate the minimum \c min of \c func
229  w.r.t. the array \c x of size \c nvar with gradient
230  \c dfunc
231  */
232  virtual int mmin_de(size_t nvar, vec_t &x, double &fmin,
233  func_t &func, dfunc_t &dfunc)
234  {
235  return mmin(nvar,x,fmin,func);
236  }
237 
238  /** \brief Print out iteration information.
239 
240  Depending on the value of the variable verbose, this prints out
241  the iteration information. If verbose=0, then no information is
242  printed, while if verbose>1, then after each iteration, the
243  present values of x and y are output to std::cout along with the
244  iteration number. If verbose>=2 then each iteration waits for a
245  character.
246  */
247  template<class vec2_t>
248  int print_iter(size_t nv, vec2_t &x, double y, int iter,
249  double value, double limit, std::string comment)
250  {
251 
252  if (verbose<=0) return 0;
253 
254  int i;
255  char ch;
256 
257  (*outs) << comment << " Iteration: " << iter << std::endl;
258  {
259  (*outs) << "x: " << std::endl;
260  for(i=0;i<((int)nv);i++) (*outs) << x[i] << " ";
261  (*outs) << std::endl;
262  }
263  (*outs) << "y: " << y << " Val: " << value << " Lim: "
264  << limit << std::endl;
265  if (verbose>1) {
266  (*outs) << "Press a key and type enter to continue. ";
267  (*ins) >> ch;
268  }
269 
270  return 0;
271  }
272 
273  /// Return string denoting type ("mmin_base")
274  const char *type() { return "mmin_base"; }
275 
276  /** \brief Copy constructor
277  */
280  this->verbose=mb.verbose;
281  this->ntrial=mb.ntrial;
282  this->tol_rel=mb.tol_rel;
283  this->tol_abs=mb.tol_abs;
284  this->last_ntrial=mb.last_ntrial;
285  this->err_nonconv=mb.err_nonconv;
286  }
287 
288  /** \brief Copy constructor from operator=
289  */
292 
293  if (this != &mb) {
294  this->verbose=mb.verbose;
295  this->ntrial=mb.ntrial;
296  this->tol_rel=mb.tol_rel;
297  this->tol_abs=mb.tol_abs;
298  this->last_ntrial=mb.last_ntrial;
299  this->err_nonconv=mb.err_nonconv;
300  }
301 
302  return *this;
303  }
304 
305  };
306 
307 #ifndef DOXYGEN_NO_O2NS
308 }
309 #endif
310 
311 #endif
312 
Class for automatically computing gradients [abstract base].
Definition: mmin.h:53
double epsmin
The minimum stepsize (default )
Definition: mmin.h:113
std::istream * ins
Stream for verbose input.
Definition: mmin.h:173
int set_verbose_stream(std::ostream &out, std::istream &in)
Set streams for verbose I/O.
Definition: mmin.h:216
std::ostream * outs
Stream for verbose output.
Definition: mmin.h:170
const char * type()
Return string denoting type ("mmin_base")
Definition: mmin.h:274
virtual int mmin(size_t nvar, vec_t &x, double &fmin, func_t &func)=0
Calculate the minimum min of func w.r.t. the array x of size nvar.
int verbose
Output control.
Definition: mmin.h:193
virtual int mmin_de(size_t nvar, vec_t &x, double &fmin, func_t &func, dfunc_t &dfunc)
Calculate the minimum min of func w.r.t. the array x of size nvar with gradient dfunc.
Definition: mmin.h:232
virtual int set_function(func_t &f)
Set the function to compute the gradient of.
Definition: mmin.h:64
bool err_nonconv
If true, call the error handler if the routine does not "converge".
Definition: mmin.h:208
Multidimensional minimization [abstract base].
Definition: mmin.h:163
int print_iter(size_t nv, vec2_t &x, double y, int iter, double value, double limit, std::string comment)
Print out iteration information.
Definition: mmin.h:248
double epsrel
The relative stepsize for finite-differencing (default )
Definition: mmin.h:110
Simple automatic computation of gradient by finite differencing.
Definition: mmin.h:95
int last_ntrial
The number of iterations for in the most recent minimization.
Definition: mmin.h:205
virtual int operator()(size_t nv, vec_t &x, vec_t &g)
Compute the gradient g at the point x.
Definition: mmin.h:117
func_t * func
A pointer to the user-specified function.
Definition: mmin.h:78
std::function< int(size_t, boost::numeric::ublas::vector< double > &, boost::numeric::ublas::vector< double > &)> grad_funct11
Array of multi-dimensional functions typedef.
Definition: mmin.h:41
double tol_abs
The independent variable tolerance.
Definition: mmin.h:202
double tol_rel
Function value tolerance.
Definition: mmin.h:199
virtual int operator()(size_t nv, vec_t &x, vec_t &g)=0
Compute the gradient g at the point x.
std::function< double(size_t, const boost::numeric::ublas::vector< double > &)> multi_funct11
Multi-dimensional function typedef.
Definition: multi_funct.h:45
int ntrial
Maximum number of iterations.
Definition: mmin.h:196

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).
Hosted at Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads..