cepgen is hosted by Hepforge, IPPP Durham
CepGen 1.2.5
Central exclusive processes event generator
Loading...
Searching...
No Matches
GSLIntegrator.cpp
Go to the documentation of this file.
1/*
2 * CepGen: a central exclusive processes event generator
3 * Copyright (C) 2013-2024 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
24
25namespace cepgen {
27 CG_DEBUG("Integrator:build") << "Random numbers generator: " << gsl_rng_name(rnd_gen_->engine<gsl_rng>()) << ".";
28 }
29
31 //--- specify the integrand through the GSL wrapper
32 funct_ = [&](double* x, size_t ndim, void*) -> double { return integrand.eval(std::vector<double>(x, x + ndim)); };
33 function_ = utils::GSLMonteFunctionWrapper<decltype(funct_)>::build(funct_, integrand.size());
34 if (!function_)
35 throw CG_FATAL("GSLIntegrator:setIntegrand") << "Integrand was not properly set.";
36 if (function_->dim <= 0)
37 throw CG_FATAL("GSLIntegrator:setIntegrand") << "Invalid phase space dimension: " << function_->dim << ".";
38
39 CG_DEBUG("GSLIntegrator:setIntegrand") << "Number of integration dimensions: " << function_->dim << ".";
40
41 checkLimits(integrand); // check the integration bounds
42 }
43
44 void GSLIntegrator::setLimits(const std::vector<Limits>& lims) {
46 xlow_.clear();
47 xhigh_.clear();
48 for (const auto& lim : limits_) {
49 xlow_.emplace_back(lim.min());
50 xhigh_.emplace_back(lim.max());
51 }
52 }
53
55 auto desc = Integrator::description();
56 desc.add<ParametersDescription>("randomGenerator", RandomGeneratorFactory::get().describeParameters("gsl"));
57 return desc;
58 }
59} // namespace cepgen
#define CG_FATAL(mod)
Definition Exception.h:61
#define CG_DEBUG(mod)
Definition Message.h:220
std::vector< double > xhigh_
std::vector< double > xlow_
std::function< double(double *, size_t, void *)> funct_
A functor wrapping GSL's function footprint.
void setIntegrand(Integrand &)
GSLIntegrator(const ParametersList &)
std::unique_ptr< gsl_monte_function > function_
GSL structure storing the function to be integrated by this integrator instance (along with its param...
void setLimits(const std::vector< Limits > &) override
Specify the variables limits on integration.
static ParametersDescription description()
An integrand wrapper placeholder.
Definition Integrand.h:27
virtual double eval(const std::vector< double > &)=0
Compute the integrand for a given coordinates set.
virtual size_t size() const =0
Phase space dimension.
Monte-Carlo integration algorithm.
Definition Integrator.h:28
void checkLimits(const Integrand &)
Ensure the integration bounds are properly set.
virtual void setLimits(const std::vector< Limits > &limits)
Specify the variables limits on integration.
Definition Integrator.h:38
std::vector< Limits > limits_
List of per-variable integration limits.
Definition Integrator.h:58
const std::unique_ptr< utils::RandomGenerator > rnd_gen_
Definition Integrator.h:56
static ParametersDescription description()
A description object for parameters collection.
GSL wrapper to define a functor as an integrable functional.
Common namespace for this Monte Carlo generator.