cepgen is hosted by Hepforge, IPPP Durham
CepGen 1.2.5
Central exclusive processes event generator
Loading...
Searching...
No Matches
PlainIntegrator.cpp
Go to the documentation of this file.
1/*
2 * CepGen: a central exclusive processes event generator
3 * Copyright (C) 2013-2023 Laurent Forthomme
4 *
5 * This program is free software: you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation, either version 3 of the License, or
8 * any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program. If not, see <http://www.gnu.org/licenses/>.
17 */
18
19#include <gsl/gsl_monte_plain.h>
20
24
25namespace cepgen {
27 class PlainIntegrator final : public GSLIntegrator {
28 public:
29 explicit PlainIntegrator(const ParametersList& params)
30 : GSLIntegrator(params), ncvg_(steer<int>("numFunctionCalls")) {}
31
33 auto desc = GSLIntegrator::description();
34 desc.setDescription("Plain (trial/error) integrator");
35 desc.add<int>("numFunctionCalls", 50'000);
36 return desc;
37 }
38
39 Value integrate(Integrand& integrand) override {
40 setIntegrand(integrand);
41
42 //--- launch integration
43 std::unique_ptr<gsl_monte_plain_state, decltype(&gsl_monte_plain_free)> pln_state(
44 gsl_monte_plain_alloc(function_->dim), gsl_monte_plain_free);
45 double result, abserr;
46 if (int res = gsl_monte_plain_integrate(function_.get(),
47 &xlow_[0],
48 &xhigh_[0],
49 function_->dim,
50 ncvg_,
51 rnd_gen_->engine<gsl_rng>(),
52 pln_state.get(),
53 &result,
54 &abserr);
55 res != GSL_SUCCESS)
56 throw CG_FATAL("Integrator:integrate") << "Error while performing the integration!\n\t"
57 << "GSL error: " << gsl_strerror(res) << ".";
58
59 return Value{result, abserr};
60 }
61
62 private:
63 const int ncvg_;
64 };
65} // namespace cepgen
66REGISTER_INTEGRATOR("plain", PlainIntegrator);
#define CG_FATAL(mod)
Definition Exception.h:61
#define REGISTER_INTEGRATOR(name, obj)
Add a generic process definition to the list of handled processes.
std::vector< double > xhigh_
std::vector< double > xlow_
void setIntegrand(Integrand &)
std::unique_ptr< gsl_monte_function > function_
GSL structure storing the function to be integrated by this integrator instance (along with its param...
static ParametersDescription description()
An integrand wrapper placeholder.
Definition Integrand.h:27
const std::unique_ptr< utils::RandomGenerator > rnd_gen_
Definition Integrator.h:56
A description object for parameters collection.
Plain integration algorithm randomly sampling points in the phase space.
PlainIntegrator(const ParametersList &params)
Value integrate(Integrand &integrand) override
Perform the multidimensional Monte Carlo integration.
static ParametersDescription description()
T steer(const std::string &key) const
Retrieve a parameters as previously steered.
Definition Steerable.h:39
A scalar value with its uncertainty.
Definition Value.h:26
Common namespace for this Monte Carlo generator.