root_robbins_monro.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2017-2020, 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_ROOT_ROBBINS_MONRO_H
24 #define O2SCL_ROOT_ROBBINS_MONRO_H
25 
26 /** \file root_robbins_monro.h
27  \brief File defining \ref o2scl::root_robbins_monro
28 */
29 
30 #include <limits>
31 #include <o2scl/funct.h>
32 #include <o2scl/root.h>
33 
34 #ifndef DOXYGEN_NO_O2NS
35 namespace o2scl {
36 #endif
37 
38  /** \brief One-dimensional root-finding for noisy functions
39  */
40  template<class func_t=funct, class dfunc_t=func_t>
41  class root_robbins_monro : public root<func_t,dfunc_t> {
42 
43  public:
44 
46  coeff=1.0;
47  gamma=1.0;
48  }
49 
50  /// If true, apply Polyak-Ruppert averaging
52 
53  /// Coefficient for sequence (default 1.0)
54  double coeff;
55 
56  /// Exponent for sequence (default 1.0)
57  double gamma;
58 
59  /// Return the type, \c "root_robbins_monro".
60  virtual const char *type() { return "root_robbins_monro"; }
61 
62  /** \brief Solve \c func using \c x as an initial guess
63  */
64  virtual int solve(double &x, func_t &func) {
65  double x2;
66 
67  if (pr_averaging) {
68 
69  double avg=0.0;
70  for(size_t i=0;i<this->ntrial/2;i++) {
71  double fx=func(x);
72  if (fabs(fx)<this->tol_rel) {
73  return o2scl::success;
74  }
75  x2=x-fx*coeff/pow(((double)(i+1)),gamma);
76  avg+=(x2*((double)i)+avg)/((double)(i+1));
77  double favg=func(avg);
78  if (fabs(favg)<this->tol_rel) {
79  x=avg;
80  return o2scl::success;
81  }
82  x=x2;
83  }
84  O2SCL_CONV2_RET("Function root_robbins_monro::solve() exceeded ",
85  "maximum number of iterations.",o2scl::exc_emaxiter,
86  this->err_nonconv);
87  }
88 
89  for(size_t i=0;i<this->ntrial;i++) {
90  double fx=func(x);
91  if (fabs(fx)<this->tol_rel) {
92  return o2scl::success;
93  }
94  x2=x-f(x)*coeff/pow(((double)(i+1)),gamma);
95  x=x2;
96  }
97  fx=func(x);
98  if (fabs(fx)<this->tol_rel) {
99  return o2scl::success;
100  }
101  O2SCL_CONV2_RET("Function root_robbins_monro::solve() exceeded ",
102  "maximum number of iterations.",o2scl::exc_emaxiter,
103  this->err_nonconv);
104  return o2scl::exc_emaxiter;
105  }
106 
107  /// Solve \c func in region \f$ x_1<x<x_2 \f$ returning \f$ x_1 \f$.
108  virtual int solve_bkt(double &x1, double x2, func_t &f) {
109  x1=(x1+x2)/2.0;
110  return solve(x1,f);
111  }
112 
113  /** \brief Solve \c func using \c x as an initial
114  guess using derivatives \c df.
115  */
116  virtual int solve_de(double &x, func_t &func, dfunc_t &df) {
117  double x2;
118 
119  if (pr_averaging) {
120 
121  double avg=0.0;
122  for(size_t i=0;i<this->ntrial/2;i++) {
123  double fx=func(x);
124  double fpx=df(x);
125  if (fabs(fx)<this->tol_rel) {
126  return o2scl::success;
127  }
128  x2=x-fx/df*coeff/pow(((double)(i+1)),gamma);
129  avg+=(x2*((double)i)+avg)/((double)(i+1));
130  double favg=func(avg);
131  if (fabs(favg)<this->tol_rel) {
132  x=avg;
133  return o2scl::success;
134  }
135  x=x2;
136  }
137  O2SCL_CONV2_RET("Function root_robbins_monro::solve() exceeded ",
138  "maximum number of iterations.",o2scl::exc_emaxiter,
139  this->err_nonconv);
140  }
141 
142  for(size_t i=0;i<this->ntrial;i++) {
143  double fx=func(x);
144  double fpx=df(x);
145  if (fabs(fx)<this->tol_rel) {
146  return o2scl::success;
147  }
148  x2=x-fx/df*coeff/pow(((double)(i+1)),gamma);
149  x=x2;
150  }
151  fx=func(x);
152  if (fabs(fx)<this->tol_rel) {
153  return o2scl::success;
154  }
155  O2SCL_CONV2_RET("Function root_robbins_monro::solve() exceeded ",
156  "maximum number of iterations.",o2scl::exc_emaxiter,
157  this->err_nonconv);
158  return o2scl::exc_emaxiter;
159  }
160 
161  };
162 
163 #ifndef DOXYGEN_NO_O2NS
164 }
165 #endif
166 
167 #endif
o2scl::root< funct, funct >::tol_rel
double tol_rel
The maximum value of the functions for success (default )
Definition: root.h:66
O2SCL_CONV2_RET
#define O2SCL_CONV2_RET(d, d2, n, b)
Set an error and return the error value, two-string version.
Definition: err_hnd.h:303
o2scl::root_robbins_monro::type
virtual const char * type()
Return the type, "root_robbins_monro".
Definition: root_robbins_monro.h:60
o2scl::root_robbins_monro::solve_bkt
virtual int solve_bkt(double &x1, double x2, func_t &f)
Solve func in region returning .
Definition: root_robbins_monro.h:108
o2scl
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
Definition: anneal.h:42
o2scl::exc_emaxiter
@ exc_emaxiter
exceeded max number of iterations
Definition: err_hnd.h:73
o2scl::root
One-dimensional solver [abstract base].
Definition: root.h:48
o2scl_inte_qng_coeffs::x2
static const double x2[5]
Definition: inte_qng_gsl.h:66
o2scl::success
@ success
Success.
Definition: err_hnd.h:47
o2scl::root_robbins_monro::solve
virtual int solve(double &x, func_t &func)
Solve func using x as an initial guess.
Definition: root_robbins_monro.h:64
o2scl_inte_qng_coeffs::x1
static const double x1[5]
Definition: inte_qng_gsl.h:48
o2scl::root_robbins_monro
One-dimensional root-finding for noisy functions.
Definition: root_robbins_monro.h:41
o2scl::root_robbins_monro::coeff
double coeff
Coefficient for sequence (default 1.0)
Definition: root_robbins_monro.h:54
o2scl::root_robbins_monro::solve_de
virtual int solve_de(double &x, func_t &func, dfunc_t &df)
Solve func using x as an initial guess using derivatives df.
Definition: root_robbins_monro.h:116
o2scl::root_robbins_monro::pr_averaging
bool pr_averaging
If true, apply Polyak-Ruppert averaging.
Definition: root_robbins_monro.h:51
o2scl::root_robbins_monro::gamma
double gamma
Exponent for sequence (default 1.0)
Definition: root_robbins_monro.h:57
o2scl::root< funct, funct >::err_nonconv
bool err_nonconv
If true, call the error handler if the solver does not converge (default true)
Definition: root.h:82
o2scl::root< funct, funct >::ntrial
int ntrial
Maximum number of iterations (default 100)
Definition: root.h:77

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).