Example on the usage of the adaptive 1D integration algorithm of MathMore it calculates the numerically cumulative integral of a distribution (like in this case the BreitWigner) to execute the macro type it (you need to compile with AClic)
root[0] .x mathmoreIntegration.C+
This tutorials require having libMathMore built with ROOT.
To build mathmore you need to have a version of GSL >= 1.8 installed in your system The ROOT configure will automatically find GSL if the script gsl-config (from GSL) is in your PATH,. otherwise you need to configure root with the options –gsl-incdir and –gsl-libdir.
Processing /mnt/build/workspace/root-makedoc-v614/rootspi/rdoc/src/v6-14-00-patches/tutorials/math/mathmoreIntegration.C...
***************************************************************
Drawing cumulatives of BreitWigner in interval [ -2 , 2 ]
***************************************************************
***************************************************************
Test integration performances in interval [ -2 , 2 ]
Time using ROOT::Math::Integrator : 0.141767
Number of function calls = 69
42201.6649413923442
Time using TF1::Integral : 0.442968
Number of function calls = 91
42201.6649413923442
#include <iostream>
#include <limits>
int nc = 0;
double exactIntegral( double a, double b) {
}
double func( double x){
nc++;
}
double func2(const double *x, const double * = 0){
nc++;
}
void testIntegPerf(double x1, double x2, int n = 100000){
std::cout << "\n\n***************************************************************\n";
std::cout << "Test integration performances in interval [ " << x1 << " , " << x2 << " ]\n\n";
double dx = (x2-
x1)/
double(n);
double s1 = 0.0;
nc = 0;
for (
int i = 0; i <
n; ++i) {
double x = x1 + dx*i;
s1+= ig.Integral(x1,x);
}
std::cout <<
"Time using ROOT::Math::Integrator :\t" << timer.
RealTime() << std::endl;
std::cout << "Number of function calls = " << nc/n << std::endl;
int pr = std::cout.precision(18); std::cout << s1 << std::endl; std::cout.precision(pr);
TF1 *fBW =
new TF1(
"fBW",func2,x1, x2,0);
nc = 0;
double s2 = 0;
for (
int i = 0; i <
n; ++i) {
double x = x1 + dx*i;
}
std::cout <<
"Time using TF1::Integral :\t\t\t" << timer.
RealTime() << std::endl;
std::cout << "Number of function calls = " << nc/n << std::endl;
pr = std::cout.precision(18); std::cout << s1 << std::endl; std::cout.precision(pr);
}
void DrawCumulative(double x1, double x2, int n = 100){
std::cout << "\n\n***************************************************************\n";
std::cout << "Drawing cumulatives of BreitWigner in interval [ " << x1 << " , " << x2 << " ]\n\n";
double dx = (x2-
x1)/
double(n);
TH1D *cum0 =
new TH1D(
"cum0",
"", n, x1, x2);
for (
int i = 1; i <=
n; ++i) {
double x = x1 + dx*i;
}
TH1D *cum1 =
new TH1D(
"cum1",
"", n, x1, x2);
for (
int i = 1; i <=
n; ++i) {
double x = x1 + dx*i;
}
TF1 *fBW =
new TF1(
"fBW",
"TMath::BreitWigner(x, 0, 1)",x1, x2);
TH1D *cum2 =
new TH1D(
"cum2",
"", n, x1, x2);
for (
int i = 1; i <=
n; ++i) {
double x = x1 + dx*i;
}
TH1D *cum10 =
new TH1D(
"cum10",
"", n, x1, x2);
TH1D *cum20 =
new TH1D(
"cum23",
"", n, x1, x2);
for (
int i = 1; i <=
n; ++i) {
}
cum0->
SetTitle(
"BreitWigner - the cumulative");
l->
AddEntry(cum10,
"GSL integration - analytical ");
l->
AddEntry(cum20,
"TF1::Integral - analytical ");
std::cout << "\n***************************************************************\n";
}
void mathmoreIntegration(double a = -2, double b = 2)
{
DrawCumulative(a, b);
testIntegPerf(a, b);
}
- Authors
- M. Slawinska, L. Moneta
Definition in file mathmoreIntegration.C.