mcarlo_plain.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2017, 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 /* monte/plain.c
24  *
25  * Copyright (C) 1996, 1997, 1998, 1999, 2000 Michael Booth
26  *
27  * This program is free software; you can redistribute it and/or modify
28  * it under the terms of the GNU General Public License as published by
29  * the Free Software Foundation; either version 3 of the License, or (at
30  * your option) any later version.
31  *
32  * This program is distributed in the hope that it will be useful, but
33  * WITHOUT ANY WARRANTY; without even the implied warranty of
34  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
35  * General Public License for more details.
36  *
37  * You should have received a copy of the GNU General Public License
38  * along with this program; if not, write to the Free Software
39  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
40  * 02110-1301, USA.
41  */
42 #ifndef O2SCL_MCARLO_PLAIN_H
43 #define O2SCL_MCARLO_PLAIN_H
44 
45 /** \file mcarlo_plain.h
46  \brief File defining \ref o2scl::mcarlo_plain
47 */
48 #include <iostream>
49 
50 #include <boost/numeric/ublas/vector.hpp>
51 
52 #include <o2scl/multi_funct.h>
53 #include <o2scl/mcarlo.h>
54 #include <o2scl/rng_gsl.h>
55 
56 #ifndef DOXYGEN_NO_O2NS
57 namespace o2scl {
58 #endif
59 
60  /** \brief Multidimensional integration using plain Monte Carlo (GSL)
61  */
62  template<class func_t=multi_funct11,
64  class rng_t=rng_gsl>
65  class mcarlo_plain : public mcarlo<func_t,vec_t,rng_t> {
66 
67  public:
68 
69  virtual ~mcarlo_plain() {}
70 
71  /// Integrate function \c func from x=a to x=b.
72  virtual int minteg_err(func_t &func, size_t ndim, const vec_t &a,
73  const vec_t &b, double &res, double &err) {
74 
75  double r;
76 
77  vec_t cc(ndim);
78 
79  double vol=1.0, m=0.0, q=0.0;
80  for(size_t i=0;i<ndim;i++) vol*=b[i]-a[i];
81 
82  for(size_t n=0;n<this->n_points;n++) {
83  for(size_t i=0;i<ndim;i++) {
84  do {
85  r=this->rng_dist(this->rng);
86  } while (r==0.0);
87  cc[i]=a[i]+r*(b[i]-a[i]);
88  }
89  double fval;
90  fval=func(ndim,cc);
91  double d=fval-m;
92  m+=d/(n+1.0);
93  q+=d*d*(n/(n+1.0));
94  }
95 
96  res=vol*m;
97 
98  if (this->n_points<2) {
99  err=0.0;
100  } else {
101  err=vol*sqrt(q/(this->n_points*(this->n_points-1.0)));
102  }
103 
104  return 0;
105 
106  }
107 
108  /** \brief Integrate function \c func over the hypercube from
109  \f$ x_i=a_i \f$ to \f$ x_i=b_i \f$ for \f$ 0<i< \f$ ndim-1
110  */
111  virtual double minteg(func_t &func, size_t ndim, const vec_t &a,
112  const vec_t &b) {
113  double res;
114  minteg_err(func,ndim,a,b,res,this->interror);
115  return res;
116  }
117 
118  /// Return string denoting type ("mcarlo_plain")
119  virtual const char *type() { return "mcarlo_plain"; }
120 
121  };
122 
123 #ifndef DOXYGEN_NO_O2NS
124 }
125 #endif
126 
127 #endif
128 
rng_gsl_uniform_real rng_dist
The random number distribution.
Definition: mcarlo.h:67
The main O<span style=&#39;position: relative; top: 0.3em; font-size: 0.8em&#39;>2</span>scl O$_2$scl names...
Definition: anneal.h:42
virtual int minteg_err(func_t &func, size_t ndim, const vec_t &a, const vec_t &b, double &res, double &err)
Integrate function func from x=a to x=b.
Definition: mcarlo_plain.h:72
unsigned long n_points
Number of integration points (default 1000)
Definition: mcarlo.h:63
Monte-Carlo integration [abstract base].
Definition: mcarlo.h:51
virtual double minteg(func_t &func, size_t ndim, const vec_t &a, const vec_t &b)
Integrate function func over the hypercube from to for ndim-1.
Definition: mcarlo_plain.h:111
Multidimensional integration using plain Monte Carlo (GSL)
Definition: mcarlo_plain.h:65
double interror
The uncertainty for the last integration computation.
Definition: inte_multi.h:110
virtual const char * type()
Return string denoting type ("mcarlo_plain")
Definition: mcarlo_plain.h:119
std::function< double(size_t, const boost::numeric::ublas::vector< double > &)> multi_funct11
Multi-dimensional function typedef.
Definition: multi_funct.h:45
rng_t rng
The random number generator.
Definition: mcarlo.h:70

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