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
19
#include "
CepGen/Core/Exception.h
"
20
#include "
CepGen/Core/ParametersList.h
"
21
#include "
CepGen/Integration/GSLIntegrator.h
"
22
#include "
CepGen/Integration/Integrand.h
"
23
#include "
CepGen/Modules/RandomGeneratorFactory.h
"
24
25
namespace
cepgen
{
26
GSLIntegrator::GSLIntegrator
(
const
ParametersList
& params) :
Integrator
(params) {
27
CG_DEBUG
(
"Integrator:build"
) <<
"Random numbers generator: "
<< gsl_rng_name(
rnd_gen_
->engine<gsl_rng>()) <<
"."
;
28
}
29
30
void
GSLIntegrator::setIntegrand
(
Integrand
& integrand) {
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) {
45
Integrator::setLimits
(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
54
ParametersDescription
GSLIntegrator::description
() {
55
auto
desc =
Integrator::description
();
56
desc.add<
ParametersDescription
>(
"randomGenerator"
, RandomGeneratorFactory::get().describeParameters(
"gsl"
));
57
return
desc;
58
}
59
}
// namespace cepgen
Exception.h
CG_FATAL
#define CG_FATAL(mod)
Definition
Exception.h:61
GSLIntegrator.h
Integrand.h
CG_DEBUG
#define CG_DEBUG(mod)
Definition
Message.h:220
ParametersList.h
RandomGeneratorFactory.h
cepgen::GSLIntegrator::xhigh_
std::vector< double > xhigh_
Definition
GSLIntegrator.h:43
cepgen::GSLIntegrator::xlow_
std::vector< double > xlow_
Definition
GSLIntegrator.h:43
cepgen::GSLIntegrator::funct_
std::function< double(double *, size_t, void *)> funct_
A functor wrapping GSL's function footprint.
Definition
GSLIntegrator.h:39
cepgen::GSLIntegrator::setIntegrand
void setIntegrand(Integrand &)
Definition
GSLIntegrator.cpp:30
cepgen::GSLIntegrator::GSLIntegrator
GSLIntegrator(const ParametersList &)
Definition
GSLIntegrator.cpp:26
cepgen::GSLIntegrator::function_
std::unique_ptr< gsl_monte_function > function_
GSL structure storing the function to be integrated by this integrator instance (along with its param...
Definition
GSLIntegrator.h:42
cepgen::GSLIntegrator::setLimits
void setLimits(const std::vector< Limits > &) override
Specify the variables limits on integration.
Definition
GSLIntegrator.cpp:44
cepgen::GSLIntegrator::description
static ParametersDescription description()
Definition
GSLIntegrator.cpp:54
cepgen::Integrand
An integrand wrapper placeholder.
Definition
Integrand.h:27
cepgen::Integrand::eval
virtual double eval(const std::vector< double > &)=0
Compute the integrand for a given coordinates set.
cepgen::Integrand::size
virtual size_t size() const =0
Phase space dimension.
cepgen::Integrator
Monte-Carlo integration algorithm.
Definition
Integrator.h:28
cepgen::Integrator::checkLimits
void checkLimits(const Integrand &)
Ensure the integration bounds are properly set.
Definition
Integrator.cpp:32
cepgen::Integrator::setLimits
virtual void setLimits(const std::vector< Limits > &limits)
Specify the variables limits on integration.
Definition
Integrator.h:38
cepgen::Integrator::limits_
std::vector< Limits > limits_
List of per-variable integration limits.
Definition
Integrator.h:58
cepgen::Integrator::rnd_gen_
const std::unique_ptr< utils::RandomGenerator > rnd_gen_
Definition
Integrator.h:56
cepgen::Integrator::description
static ParametersDescription description()
Definition
Integrator.cpp:71
cepgen::ParametersDescription
A description object for parameters collection.
Definition
ParametersDescription.h:26
cepgen::ParametersList
Definition
ParametersList.h:52
cepgen::utils::GSLMonteFunctionWrapper
GSL wrapper to define a functor as an integrable functional.
Definition
GSLFunctionsWrappers.h:67
cepgen
Common namespace for this Monte Carlo generator.
Definition
CommandLineHandler.cpp:36
CepGen
Integration
GSLIntegrator.cpp
Generated on Mon Jul 29 2024 for CepGen by
1.9.7