From e07ebd54ed44418d2782268f4a88d1a49f1f0b34 Mon Sep 17 00:00:00 2001 From: niccolo Date: Tue, 19 Mar 2024 18:07:06 +0100 Subject: [PATCH] Clean un po' tutto --- inc/adani/ApproximateCoefficientFunction.h | 25 ++ inc/adani/AsymptoticCoefficientFunction.h | 4 + inc/adani/CoefficientFunction.h | 16 ++ inc/adani/Convolution.h | 40 ++- inc/adani/ExactCoefficientFunction.h | 37 +-- inc/adani/HighEnergyCoefficientFunction.h | 26 +- inc/adani/HighScaleCoefficientFunction.h | 4 + src/AsymptoticCoefficientFunction.cc | 19 +- src/CoefficientFunction.cc | 42 +++- src/Convolution.cc | 130 ++++++++-- src/ExactCoefficientFunction.cc | 273 ++++++--------------- src/HighEnergyCoefficientFunction.cc | 42 +++- src/HighScaleCoefficientFunction.cc | 21 +- 13 files changed, 416 insertions(+), 263 deletions(-) diff --git a/inc/adani/ApproximateCoefficientFunction.h b/inc/adani/ApproximateCoefficientFunction.h index 47921e9..bbc3163 100644 --- a/inc/adani/ApproximateCoefficientFunction.h +++ b/inc/adani/ApproximateCoefficientFunction.h @@ -22,6 +22,10 @@ #include "adani/AsymptoticCoefficientFunction.h" #include "adani/ExactCoefficientFunction.h" +//==========================================================================================// +// struct approximation_parameters: parameters of the damping functions +//------------------------------------------------------------------------------------------// + struct approximation_parameters { double A; double B; @@ -29,11 +33,19 @@ struct approximation_parameters { double D; }; +//==========================================================================================// +// struct variation_parameters: parameters of the variation of the damping functions +//------------------------------------------------------------------------------------------// + struct variation_parameters { double var; double fact; }; +//==========================================================================================// +// struct klmv_params: parameters for klmv approximation +//------------------------------------------------------------------------------------------// + struct klmv_params { double gamma; double C; @@ -42,6 +54,10 @@ struct klmv_params { double const_coeff; }; +//==========================================================================================// +// class AbstractApproximate +//------------------------------------------------------------------------------------------// + class AbstractApproximate : public CoefficientFunction { public: AbstractApproximate(const int& order, const char& kind, const char& channel, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000, const int& method_flag = 1, const int& MCcalls = 25000); @@ -56,6 +72,10 @@ class AbstractApproximate : public CoefficientFunction { ExactCoefficientFunction* muterms_; }; +//==========================================================================================// +// class ApproximateCoefficientFunction +//------------------------------------------------------------------------------------------// + class ApproximateCoefficientFunction : public AbstractApproximate { public: ApproximateCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL = true, const bool& exact_highscale = false, const bool& revised_approx_highscale = true, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000, const int& method_flag = 1, const int& MCcalls = 25000) ; @@ -75,6 +95,11 @@ class ApproximateCoefficientFunction : public AbstractApproximate { }; +//==========================================================================================// +// class ApproximateCoefficientFunctionKLMV: Approximate coefficient functions from [arXiv:1205.5727] +// klmv = Kawamura, Lo Presti, Moch, Vogt +//------------------------------------------------------------------------------------------// + class ApproximateCoefficientFunctionKLMV : public AbstractApproximate { public: ApproximateCoefficientFunctionKLMV(const int& order, const char& kind, const char& channel, const bool& revised_approx_highscale = true, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000, const int& method_flag = 1, const int& MCcalls = 25000) ; diff --git a/inc/adani/AsymptoticCoefficientFunction.h b/inc/adani/AsymptoticCoefficientFunction.h index eeb4c9c..081b716 100644 --- a/inc/adani/AsymptoticCoefficientFunction.h +++ b/inc/adani/AsymptoticCoefficientFunction.h @@ -21,6 +21,10 @@ #include "adani/HighEnergyCoefficientFunction.h" #include "adani/HighScaleCoefficientFunction.h" +//==========================================================================================// +// class AsymptoticCoefficientFunction +//------------------------------------------------------------------------------------------// + class AsymptoticCoefficientFunction : public AbstractHighEnergyCoefficientFunction { public: AsymptoticCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL = true, const bool& exact_highscale = false, const bool& revised_approx_highscale = true); diff --git a/inc/adani/CoefficientFunction.h b/inc/adani/CoefficientFunction.h index 57e9d35..77f9052 100644 --- a/inc/adani/CoefficientFunction.h +++ b/inc/adani/CoefficientFunction.h @@ -22,6 +22,22 @@ #include #include +//==========================================================================================// +// The notation used is the following: +// C^{(k)} = k-th order expansion in +// terms of a_s^{[nf]} +// +//------------------------------------------------------------------------------------------// +//==========================================================================================// +// Notation: +// m2Q2 = m^2/Q^2 +// m2mu2 = m^2/mu^2 +//------------------------------------------------------------------------------------------// + +//==========================================================================================// +// class CoefficientFunction: abstract base class for all the CoefficeintFunction +//------------------------------------------------------------------------------------------// + class CoefficientFunction { public: diff --git a/inc/adani/Convolution.h b/inc/adani/Convolution.h index 2c8becd..16b9602 100644 --- a/inc/adani/Convolution.h +++ b/inc/adani/Convolution.h @@ -9,8 +9,7 @@ * Author: Dribbla tutti, anche i cammelli del * deserto * - * In this file there are the convolutions between the - * functions defined in this library. + * In this file there are the convolutions between the coefficeint functions and splitting functions. * * For the convolution between a regular function and a * function containing regular, singular (plus distribution) @@ -31,13 +30,19 @@ #include #include +//==========================================================================================// +// class AbstractConvolution +//------------------------------------------------------------------------------------------// + class AbstractConvolution { public: AbstractConvolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000); virtual ~AbstractConvolution() = 0; + // result of the convolution double Convolute(double x, double m2Q2, int nf) const ; + // integrals of the reular, singular and local parts of the splittings virtual double RegularPart(double x, double m2Q2, int nf) const = 0; virtual double SingularPart(double x, double m2Q2, int nf) const = 0; virtual double LocalPart(double x, double m2Q2, int nf) const = 0; @@ -46,7 +51,6 @@ class AbstractConvolution { double GetAbserr() const {return abserr_;}; double GetRelerr() const {return relerr_;}; int GetDim() const {return dim_;} ; - CoefficientFunction* GetCoeffFunc() const {return coefffunc_;}; AbstractSplittingFunction* GetSplitFunc() const {return splitfunc_;}; @@ -66,6 +70,10 @@ class AbstractConvolution { }; +//==========================================================================================// +// class Convolution +//------------------------------------------------------------------------------------------// + class Convolution : public AbstractConvolution { public: Convolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000) : AbstractConvolution(coefffunc, splitfunc, abserr, relerr, dim) {} ; @@ -75,11 +83,15 @@ class Convolution : public AbstractConvolution { double SingularPart(double x, double m2Q2, int nf) const override; double LocalPart(double x, double m2Q2, int nf) const override; + // support function for the integral. it is static in order to be passed to gsl static double regular_integrand(double z, void *p) ; static double singular_integrand(double z, void *p) ; }; +//==========================================================================================// +// class ConvolutedCoefficientFunction: convolution between a CoefficientFunction and a SplittingFunction +//------------------------------------------------------------------------------------------// class ConvolutedCoefficientFunction : public CoefficientFunction { public: @@ -99,10 +111,14 @@ class ConvolutedCoefficientFunction : public CoefficientFunction { }; -class MonteCarloDoubleConvolution : public AbstractConvolution { +//==========================================================================================// +// class DoubleConvolution +//------------------------------------------------------------------------------------------// + +class DoubleConvolution : public AbstractConvolution { public: - MonteCarloDoubleConvolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000, const int& method_flag = 1, const int& MCcalls = 25000) ; - ~MonteCarloDoubleConvolution() override ; + DoubleConvolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr = 1e-3, const double& relerr = 1e-3, const int& dim = 1000, const int& method_flag = 1, const int& MCcalls = 25000) ; + ~DoubleConvolution() override ; double RegularPart(double x, double m2Q2, int nf) const override; double SingularPart(double x, double m2Q2, int nf) const override; @@ -116,6 +132,7 @@ class MonteCarloDoubleConvolution : public AbstractConvolution { void SetMethodFlag(const int& method_flag) ; void SetMCcalls(const int& MCcalls) ; + // support function for the integral. it is static in order to be passed to gsl static double regular1_integrand(double z[], size_t /*dim*/, void *p) ; static double regular2_integrand(double z[], size_t /*dim*/, void *p) ; static double regular3_integrand(double z, void *p) ; @@ -133,6 +150,10 @@ class MonteCarloDoubleConvolution : public AbstractConvolution { }; +//==========================================================================================// +// struct function_params to be passed to gsl +//------------------------------------------------------------------------------------------// + struct function_params { double x; double m2Q2; @@ -140,11 +161,4 @@ struct function_params { const AbstractConvolution* conv; }; -//==========================================================================================// -// Convolution between the first order massive gluon -// coefficient functions and the convolution between the -// splitting functions Pgg0 and Pgg0 using monte carlo -// mathods -//------------------------------------------------------------------------------------------// - #endif diff --git a/inc/adani/ExactCoefficientFunction.h b/inc/adani/ExactCoefficientFunction.h index 8fac33e..d4572ec 100644 --- a/inc/adani/ExactCoefficientFunction.h +++ b/inc/adani/ExactCoefficientFunction.h @@ -18,20 +18,13 @@ #define Exact_h #include "adani/CoefficientFunction.h" +#include "adani/SplittingFunction.h" #include "adani/Convolution.h" #include //==========================================================================================// -// The notation used is the following: -// C^{(k)} = k-th order expansion in -// terms of a_s^{[nf]} -// -//------------------------------------------------------------------------------------------// -//==========================================================================================// -// Notation: -// m2Q2 = m^2/Q^2 -// m2mu2 = m^2/mu^2 +// class ExactCoefficientFunction //------------------------------------------------------------------------------------------// class ExactCoefficientFunction : public CoefficientFunction { @@ -63,7 +56,6 @@ class ExactCoefficientFunction : public CoefficientFunction { void SetFunctions(); private: - int method_flag_ ; double abserr_ ; double relerr_ ; @@ -100,17 +92,6 @@ class ExactCoefficientFunction : public CoefficientFunction { double C2_g1(double x, double m2Q2, int /*nf*/) const; double CL_g1(double x, double m2Q2, int /*nf*/) const; - //==========================================================================================// - // Exact massive coefficient functions - // O(as^2) - //------------------------------------------------------------------------------------------// - - // double C2_g2(double x, double m2Q2, double m2mu2) const; - // double C2_ps2(double x, double m2Q2, double m2mu2) const; - - // double CL_g2(double x, double m2Q2, double m2mu2) const; - // double CL_ps2(double x, double m2Q2, double m2mu2) const; - //==========================================================================================// // Exact massive coefficient functions O(as^2): // mu-independent terms @@ -129,6 +110,11 @@ class ExactCoefficientFunction : public CoefficientFunction { double C_g21(double x, double m2Q2) const; double C_ps21(double x, double m2Q2) const; + //==========================================================================================// + // Exact massive coefficient functions O(as^2): + // mu-dependent terms + //------------------------------------------------------------------------------------------// + double C_ps2_MuDep(double x, double m2Q2, double m2mu2, int /*nf*/) const ; double C_g2_MuDep(double x, double m2Q2, double m2mu2, int /*nf*/) const ; @@ -146,12 +132,13 @@ class ExactCoefficientFunction : public CoefficientFunction { //------------------------------------------------------------------------------------------// double C_ps32(double x, double m2Q2, int nf) const; - - // These two expressions are integrated with montcarlo - // methods - double C_g32(double x, double m2Q2, int nf) const; + //==========================================================================================// + // Exact massive coefficient functions O(as^3): + // mu-dependent terms + //------------------------------------------------------------------------------------------// + double C_ps3_MuDep(double x, double m2Q2, double m2mu2, int nf) const ; double C_g3_MuDep(double x, double m2Q2, double m2mu2, int nf) const ; diff --git a/inc/adani/HighEnergyCoefficientFunction.h b/inc/adani/HighEnergyCoefficientFunction.h index 696cc8b..42d385d 100644 --- a/inc/adani/HighEnergyCoefficientFunction.h +++ b/inc/adani/HighEnergyCoefficientFunction.h @@ -29,6 +29,10 @@ // C_highenergy_highscale //------------------------------------------------------------------------------------------// +//==========================================================================================// +// class AbstractHighEnergyCoefficientFunction +//------------------------------------------------------------------------------------------// + class AbstractHighEnergyCoefficientFunction : public CoefficientFunction { public: AbstractHighEnergyCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL = true) ; @@ -53,6 +57,10 @@ class AbstractHighEnergyCoefficientFunction : public CoefficientFunction { double a_21(int nf) const; }; +//==========================================================================================// +// class HighEnergyCoefficientFunction +//------------------------------------------------------------------------------------------// + class HighEnergyCoefficientFunction : public AbstractHighEnergyCoefficientFunction { public: @@ -102,11 +110,19 @@ class HighEnergyCoefficientFunction : public AbstractHighEnergyCoefficientFuncti Value CL_g3_highenergy(double x, double m2Q2, double m2mu2, int nf) const ; Value CL_ps3_highenergy(double x, double m2Q2, double m2mu2, int nf) const ; + //==========================================================================================// + // Function needed to make LL_ and NLL_ point to a zero function + //------------------------------------------------------------------------------------------// + double ZeroFunction(double /*x*/, double /*m2Q2*/, double /*m2mu2*/) const {return 0.;}; Value ZeroFunctionBand(double /*x*/, double /*m2Q2*/, double /*m2mu2*/, int /*nf*/) const {return Value(0.);}; }; +//==========================================================================================// +// class HighEnergyHighScaleCoefficientFunction +//------------------------------------------------------------------------------------------// + class HighEnergyHighScaleCoefficientFunction : public AbstractHighEnergyCoefficientFunction { public: @@ -118,7 +134,6 @@ class HighEnergyHighScaleCoefficientFunction : public AbstractHighEnergyCoeffici void SetFunctions(); private: - double (HighEnergyHighScaleCoefficientFunction::*LL_)(double, double, double) const; Value (HighEnergyHighScaleCoefficientFunction::*NLL_)(double, double, double, int) const; @@ -161,11 +176,19 @@ class HighEnergyHighScaleCoefficientFunction : public AbstractHighEnergyCoeffici Value CL_ps3_highenergy_highscale(double x, double m2Q2, double m2mu2, int nf) const; + //==========================================================================================// + // Function needed to make LL_ and NLL_ point to a zero function + //------------------------------------------------------------------------------------------// + double ZeroFunction(double /*x*/, double /*m2Q2*/, double /*m2mu2*/) const {return 0.;}; Value ZeroFunctionBand(double /*x*/, double /*m2Q2*/, double /*m2mu2*/, int /*nf*/) const {return Value(0.);}; }; +//==========================================================================================// +// class PowerTermsCoefficientFunction +//------------------------------------------------------------------------------------------// + class PowerTermsCoefficientFunction : public AbstractHighEnergyCoefficientFunction { public: @@ -175,7 +198,6 @@ class PowerTermsCoefficientFunction : public AbstractHighEnergyCoefficientFuncti Value fxBand(double x, double m2Q2, double m2mu2, int nf) const override ; private: - HighEnergyCoefficientFunction *highenergy_ ; HighEnergyHighScaleCoefficientFunction *highenergyhighscale_ ; diff --git a/inc/adani/HighScaleCoefficientFunction.h b/inc/adani/HighScaleCoefficientFunction.h index 6f5bbf7..2ba3540 100644 --- a/inc/adani/HighScaleCoefficientFunction.h +++ b/inc/adani/HighScaleCoefficientFunction.h @@ -30,6 +30,10 @@ // \alpha_s^{[nf+1]} //------------------------------------------------------------------------------------------// +//==========================================================================================// +// class HighScaleCoefficientFunction +//------------------------------------------------------------------------------------------// + class HighScaleCoefficientFunction : public CoefficientFunction { public: HighScaleCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& exact, const bool& revised_approx) ; diff --git a/src/AsymptoticCoefficientFunction.cc b/src/AsymptoticCoefficientFunction.cc index 100a355..c5d9351 100644 --- a/src/AsymptoticCoefficientFunction.cc +++ b/src/AsymptoticCoefficientFunction.cc @@ -1,24 +1,39 @@ #include "adani/AsymptoticCoefficientFunction.h" -#include "adani/HighEnergyCoefficientFunction.h" -#include "adani/HighScaleCoefficientFunction.h" + #include +//==========================================================================================// +// ExactCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + AsymptoticCoefficientFunction::AsymptoticCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL, const bool& exact_highscale, const bool& revised_approx_highscale) : AbstractHighEnergyCoefficientFunction(order, kind, channel, NLL) { highscale_ = new HighScaleCoefficientFunction(GetOrder(), GetKind(), GetChannel(), exact_highscale, revised_approx_highscale); powerterms_ = new PowerTermsCoefficientFunction(GetOrder(), GetKind(), GetChannel(), GetNLL()); } +//==========================================================================================// +// ExactCoefficientFunction: destructor +//------------------------------------------------------------------------------------------// + AsymptoticCoefficientFunction::~AsymptoticCoefficientFunction() { delete highscale_; delete powerterms_; } +//==========================================================================================// +// ExactCoefficientFunction: central value of fx +//------------------------------------------------------------------------------------------// + double AsymptoticCoefficientFunction::fx(double x, double m2Q2, double m2mu2, int nf) const { return highscale_->fx(x, m2Q2, m2mu2, nf) + powerterms_->fx(x, m2Q2, m2mu2, nf) ; } +//==========================================================================================// +// ExactCoefficientFunction: band of fx +//------------------------------------------------------------------------------------------// + Value AsymptoticCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { return highscale_->fxBand(x, m2Q2, m2mu2, nf) + powerterms_->fxBand(x, m2Q2, m2mu2, nf) ; diff --git a/src/CoefficientFunction.cc b/src/CoefficientFunction.cc index 8f75334..7006b07 100644 --- a/src/CoefficientFunction.cc +++ b/src/CoefficientFunction.cc @@ -3,6 +3,10 @@ using std::cout ; using std::endl ; +//==========================================================================================// +// CoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + CoefficientFunction::CoefficientFunction(const int& order, const char& kind, const char& channel) { SetOrder(order); @@ -11,8 +15,16 @@ CoefficientFunction::CoefficientFunction(const int& order, const char& kind, con } +//==========================================================================================// +// CoefficientFunction: destructor +//------------------------------------------------------------------------------------------// + CoefficientFunction::~CoefficientFunction() {}; +//==========================================================================================// +// CoefficientFunction: set method for order: order = 1, 2, 3 +//------------------------------------------------------------------------------------------// + void CoefficientFunction::SetOrder(const int& order) { // check order if (order < 1 || order > 3) { @@ -22,6 +34,10 @@ void CoefficientFunction::SetOrder(const int& order) { order_ = order ; } +//==========================================================================================// +// CoefficientFunction: set method for kind: kind = '2', 'L' +//------------------------------------------------------------------------------------------// + void CoefficientFunction::SetKind(const char& kind) { // check kind if (kind != '2' && kind !='L') { @@ -31,6 +47,10 @@ void CoefficientFunction::SetKind(const char& kind) { kind_ = kind ; } +//==========================================================================================// +// CoefficientFunction: set method for channel: channel = 'g', 'q' +//------------------------------------------------------------------------------------------// + void CoefficientFunction::SetChannel(const char& channel) { //check channel if (channel != 'g' && channel != 'q') { @@ -44,22 +64,42 @@ void CoefficientFunction::SetChannel(const char& channel) { channel_ = channel ; } +//==========================================================================================// +// CoefficientFunction: central value of fx +//------------------------------------------------------------------------------------------// + double CoefficientFunction::fx(double x, double m2Q2, double m2mu2, int nf) const { return fxBand(x, m2Q2, m2mu2, nf).GetCentral(); } +//==========================================================================================// +// CoefficientFunction: central value of mu independent terms +//------------------------------------------------------------------------------------------// + double CoefficientFunction::MuIndependentTerms(double x, double m2Q2, int nf) const { return fx(x, m2Q2, 1., nf); } +//==========================================================================================// +// CoefficientFunction: central value mu dependent terms +//------------------------------------------------------------------------------------------// + double CoefficientFunction::MuDependentTerms(double x, double m2Q2, double m2mu2, int nf) const { return fx(x, m2Q2, m2mu2, nf) - MuIndependentTerms(x, m2Q2, nf); } +//==========================================================================================// +// CoefficientFunction: band for mu independent terms +//------------------------------------------------------------------------------------------// + Value CoefficientFunction::MuIndependentTermsBand(double x, double m2Q2, int nf) const { return fxBand(x, m2Q2, 1., nf);; } +//==========================================================================================// +// CoefficientFunction: band for mu dependent terms +//------------------------------------------------------------------------------------------// + Value CoefficientFunction::MuDependentTermsBand(double x, double m2Q2, double m2mu2, int nf) const { double central = fxBand(x, m2Q2, m2mu2, nf).GetCentral() - MuIndependentTermsBand(x, m2Q2, nf).GetCentral(); @@ -67,4 +107,4 @@ Value CoefficientFunction::MuDependentTermsBand(double x, double m2Q2, double m2 double lower = fxBand(x, m2Q2, m2mu2, nf).GetLower() - MuIndependentTermsBand(x, m2Q2, nf).GetLower(); return Value(central, higher, lower); -} \ No newline at end of file +} diff --git a/src/Convolution.cc b/src/Convolution.cc index 068a5a6..65b0ac3 100644 --- a/src/Convolution.cc +++ b/src/Convolution.cc @@ -1,10 +1,4 @@ #include "adani/Convolution.h" -#include "adani/Constants.h" -#include "adani/ExactCoefficientFunction.h" -#include "adani/MasslessCoefficientFunction.h" -#include "adani/MatchingCondition.h" -#include "adani/SpecialFunctions.h" -#include "adani/SplittingFunction.h" #include #include @@ -15,6 +9,10 @@ using std::endl; // TODO : in all the numerical convolutions the gsl default error is first // switched off and then swithced on again: see if it can be done once for all +//==========================================================================================// +// AbstractConvolution: constructor +//------------------------------------------------------------------------------------------// + AbstractConvolution::AbstractConvolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr, const double& relerr, const int& dim) { abserr_ = abserr; relerr_ = relerr; @@ -24,8 +22,16 @@ AbstractConvolution::AbstractConvolution(CoefficientFunction* coefffunc, Abstrac splitfunc_ = splitfunc; } +//==========================================================================================// +// AbstractConvolution: destructor +//------------------------------------------------------------------------------------------// + AbstractConvolution::~AbstractConvolution() {}; +//==========================================================================================// +// AbstractConvolution: set method for abserr +//------------------------------------------------------------------------------------------// + void AbstractConvolution::SetAbserr(const double& abserr) { // check abserr if (abserr <= 0) { @@ -35,6 +41,10 @@ void AbstractConvolution::SetAbserr(const double& abserr) { abserr_ = abserr; } +//==========================================================================================// +// AbstractConvolution: set method for relerr +//------------------------------------------------------------------------------------------// + void AbstractConvolution::SetRelerr(const double& relerr) { // check relerr if (relerr <= 0) { @@ -44,6 +54,10 @@ void AbstractConvolution::SetRelerr(const double& relerr) { relerr_ = relerr; } +//==========================================================================================// +// AbstractConvolution: set method for dim +//------------------------------------------------------------------------------------------// + void AbstractConvolution::SetDim(const int& dim) { // check dim if (dim <= 0) { @@ -53,10 +67,18 @@ void AbstractConvolution::SetDim(const int& dim) { dim_ = dim; } +//==========================================================================================// +// AbstractConvolution: convolute splitting function with coefficient function +//------------------------------------------------------------------------------------------// + double AbstractConvolution::Convolute(double x, double m2Q2, int nf) const { return RegularPart(x, m2Q2, nf) + SingularPart(x, m2Q2, nf) + LocalPart(x, m2Q2, nf) ; } +//==========================================================================================// +// AbstractConvolution: integrand of the regular part +//------------------------------------------------------------------------------------------// + double Convolution::regular_integrand(double z, void *p) { struct function_params *params = (struct function_params *)p; @@ -71,6 +93,10 @@ double Convolution::regular_integrand(double z, void *p) { return cf -> MuIndependentTerms(z, m2Q2, nf) * sf -> Regular(x / z, nf) / z ; } +//==========================================================================================// +// AbstractConvolution: integrand of the singular part +//------------------------------------------------------------------------------------------// + double Convolution::singular_integrand(double z, void *p) { struct function_params *params = (struct function_params *)p; @@ -95,6 +121,10 @@ double Convolution::singular_integrand(double z, void *p) { // will not be revealed) //------------------------------------------------------------------------------------------// +//==========================================================================================// +// AbstractConvolution: regular part +//------------------------------------------------------------------------------------------// + double Convolution::RegularPart(double x, double m2Q2, int nf) const { double x_max = 1. / (1. + 4 * m2Q2); @@ -127,6 +157,10 @@ double Convolution::RegularPart(double x, double m2Q2, int nf) const { return regular; } +//==========================================================================================// +// AbstractConvolution: singular part +//------------------------------------------------------------------------------------------// + double Convolution::SingularPart(double x, double m2Q2, int nf) const { double x_max = 1. / (1. + 4 * m2Q2); @@ -160,6 +194,10 @@ double Convolution::SingularPart(double x, double m2Q2, int nf) const { return singular; } +//==========================================================================================// +// AbstractConvolution: local part +//------------------------------------------------------------------------------------------// + double Convolution::LocalPart(double x, double m2Q2, int nf) const { double x_max = 1. / (1. + 4 * m2Q2); @@ -168,7 +206,11 @@ double Convolution::LocalPart(double x, double m2Q2, int nf) const { } -MonteCarloDoubleConvolution::MonteCarloDoubleConvolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr, const double& relerr, const int& dim, const int& method_flag, const int& MCcalls) : AbstractConvolution(coefffunc, splitfunc, abserr, relerr, dim) { +//==========================================================================================// +// DoubleConvolution: constructor +//------------------------------------------------------------------------------------------// + +DoubleConvolution::DoubleConvolution(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double& abserr, const double& relerr, const int& dim, const int& method_flag, const int& MCcalls) : AbstractConvolution(coefffunc, splitfunc, abserr, relerr, dim) { SetMethodFlag(method_flag); SetMCcalls(MCcalls); @@ -182,13 +224,21 @@ MonteCarloDoubleConvolution::MonteCarloDoubleConvolution(CoefficientFunction* co } } -MonteCarloDoubleConvolution::~MonteCarloDoubleConvolution() { +//==========================================================================================// +// DoubleConvolution: destructor +//------------------------------------------------------------------------------------------// + +DoubleConvolution::~DoubleConvolution() { delete convolution_; delete conv_coeff_; } -void MonteCarloDoubleConvolution::SetMethodFlag(const int& method_flag) { +//==========================================================================================// +// DoubleConvolution: set method for method_flag +//------------------------------------------------------------------------------------------// + +void DoubleConvolution::SetMethodFlag(const int& method_flag) { // check dim if (method_flag != 0 && method_flag != 1) { cout << "Error: method_flag must be 0 or 1. Got " << method_flag << endl; @@ -197,7 +247,11 @@ void MonteCarloDoubleConvolution::SetMethodFlag(const int& method_flag) { method_flag_ = method_flag; } -void MonteCarloDoubleConvolution::SetMCcalls(const int& MCcalls) { +//==========================================================================================// +// DoubleConvolution: set method for MCcalls +//------------------------------------------------------------------------------------------// + +void DoubleConvolution::SetMCcalls(const int& MCcalls) { // check dim if (MCcalls <= 0) { cout << "Error: MCcalls must be positive. Got " << MCcalls << endl; @@ -206,7 +260,11 @@ void MonteCarloDoubleConvolution::SetMCcalls(const int& MCcalls) { MCcalls_ = MCcalls; } -double MonteCarloDoubleConvolution::regular1_integrand(double z[], size_t /*dim*/, void *p) { +//==========================================================================================// +// DoubleConvolution: integrand of the first regular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::regular1_integrand(double z[], size_t /*dim*/, void *p) { struct function_params *params = (struct function_params *)p; @@ -227,7 +285,11 @@ double MonteCarloDoubleConvolution::regular1_integrand(double z[], size_t /*dim* } } -double MonteCarloDoubleConvolution::regular2_integrand(double z[], size_t /*dim*/, void *p) { +//==========================================================================================// +// DoubleConvolution: integrand of the second regular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::regular2_integrand(double z[], size_t /*dim*/, void *p) { struct function_params *params = (struct function_params *)p; double m2Q2 = (params->m2Q2); @@ -247,7 +309,11 @@ double MonteCarloDoubleConvolution::regular2_integrand(double z[], size_t /*dim* } } -double MonteCarloDoubleConvolution::regular3_integrand(double z, void *p) { +//==========================================================================================// +// DoubleConvolution: integrand of the third regular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::regular3_integrand(double z, void *p) { struct function_params *params = (struct function_params *)p; @@ -264,7 +330,11 @@ double MonteCarloDoubleConvolution::regular3_integrand(double z, void *p) { * splitfunc -> SingularIntegrated(z / x_max, nf); } -double MonteCarloDoubleConvolution::RegularPart(double x, double m2Q2, int nf) const { +//==========================================================================================// +// DoubleConvolution: regular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::RegularPart(double x, double m2Q2, int nf) const { if (method_flag_ == 0) return convolution_ -> RegularPart(x, m2Q2, nf); @@ -327,12 +397,10 @@ double MonteCarloDoubleConvolution::RegularPart(double x, double m2Q2, int nf) c } //==========================================================================================// -// Convolution between the first order massive gluon coefficient functions for -// F2 and the convolution between the splitting functions Pgg0 and Pgg0 using -// monte carlo methods +// DoubleConvolution: integrand of the first singular part //------------------------------------------------------------------------------------------// -double MonteCarloDoubleConvolution::singular1_integrand(double z[], size_t /*dim*/, void *p) { +double DoubleConvolution::singular1_integrand(double z[], size_t /*dim*/, void *p) { struct function_params *params = (struct function_params *)p; @@ -355,7 +423,11 @@ double MonteCarloDoubleConvolution::singular1_integrand(double z[], size_t /*dim return splitfunc -> Singular(z1, nf) * (tmp - splitfunc -> Regular(x / z2, nf)) * coefffunc -> MuIndependentTerms(z2, m2Q2, nf) / z2; } -double MonteCarloDoubleConvolution::singular2_integrand(double z[], size_t /*dim*/, void *p) { +//==========================================================================================// +// DoubleConvolution: integrand of the second singular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::singular2_integrand(double z[], size_t /*dim*/, void *p) { struct function_params *params = (struct function_params *)p; @@ -381,7 +453,11 @@ double MonteCarloDoubleConvolution::singular2_integrand(double z[], size_t /*dim * (tmp - (coefffunc -> MuIndependentTerms(x / z2, m2Q2, nf) / z2 - coefffunc -> MuIndependentTerms(x, m2Q2, nf))); } -double MonteCarloDoubleConvolution::singular3_integrand(double z, void *p) { +//==========================================================================================// +// DoubleConvolution: integrand of the third singular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::singular3_integrand(double z, void *p) { struct function_params *params = (struct function_params *)p; @@ -401,7 +477,11 @@ double MonteCarloDoubleConvolution::singular3_integrand(double z, void *p) { ); } -double MonteCarloDoubleConvolution::SingularPart(double x, double m2Q2, int nf) const { +//==========================================================================================// +// DoubleConvolution: singular part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::SingularPart(double x, double m2Q2, int nf) const { if (method_flag_ == 0) return convolution_ -> SingularPart(x, m2Q2, nf); @@ -466,7 +546,11 @@ double MonteCarloDoubleConvolution::SingularPart(double x, double m2Q2, int nf) return singular1 + singular2 + singular3 + singular4; } -double MonteCarloDoubleConvolution::LocalPart(double x, double m2Q2, int nf) const { +//==========================================================================================// +// DoubleConvolution: local part +//------------------------------------------------------------------------------------------// + +double DoubleConvolution::LocalPart(double x, double m2Q2, int nf) const { if (method_flag_ == 0) return convolution_ -> LocalPart(x, m2Q2, nf); @@ -474,5 +558,3 @@ double MonteCarloDoubleConvolution::LocalPart(double x, double m2Q2, int nf) con return convolution_ -> Convolute(x, m2Q2, nf) * (splitfunc_ -> Local(nf) - splitfunc_ -> SingularIntegrated(x / x_max, nf)); } - -//------------------------------------------------------------------------------------------// diff --git a/src/ExactCoefficientFunction.cc b/src/ExactCoefficientFunction.cc index 1f940a0..9025930 100644 --- a/src/ExactCoefficientFunction.cc +++ b/src/ExactCoefficientFunction.cc @@ -1,8 +1,6 @@ #include "adani/ExactCoefficientFunction.h" #include "adani/Constants.h" -#include "adani/MatchingCondition.h" #include "adani/SpecialFunctions.h" -#include "adani/SplittingFunction.h" #include #include @@ -10,6 +8,10 @@ using std::cout ; using std::endl ; +//==========================================================================================// +// ExactCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + ExactCoefficientFunction::ExactCoefficientFunction(const int& order, const char& kind, const char& channel, const double& abserr, const double& relerr, const int& dim, const int& method_flag, const int& MCcalls) : CoefficientFunction(order, kind, channel) { SetAbserr(abserr); @@ -78,7 +80,7 @@ ExactCoefficientFunction::ExactCoefficientFunction(const int& order, const char& convolutions_lmu1_.push_back( new Convolution(gluon_nlo_, Pgg0_, abserr, relerr, dim) ); convolutions_lmu1_.push_back( new Convolution(gluon_nlo_, delta_) ); - convolutions_lmu2_.push_back( new MonteCarloDoubleConvolution(gluon_lo_, Pgg0_, abserr, relerr, dim, method_flag, MCcalls) ); + convolutions_lmu2_.push_back( new DoubleConvolution(gluon_lo_, Pgg0_, abserr, relerr, dim, method_flag, MCcalls) ); convolutions_lmu2_.push_back( new Convolution(gluon_lo_, Pgq0Pqg0_, abserr, relerr, dim) ); convolutions_lmu2_.push_back( new Convolution(gluon_lo_, Pgg0_, abserr, relerr, dim) ); convolutions_lmu2_.push_back( new Convolution(gluon_lo_, delta_) ); @@ -89,6 +91,10 @@ ExactCoefficientFunction::ExactCoefficientFunction(const int& order, const char& } +//==========================================================================================// +// ExactCoefficientFunction: destructor +//------------------------------------------------------------------------------------------// + ExactCoefficientFunction::~ExactCoefficientFunction() { delete gluon_lo_; @@ -117,22 +123,43 @@ ExactCoefficientFunction::~ExactCoefficientFunction() { } +//==========================================================================================// +// ExactCoefficientFunction: central value of the exact coefficient function +//------------------------------------------------------------------------------------------// + double ExactCoefficientFunction::fx(double x, double m2Q2, double m2mu2, int nf) const { return MuIndependentTerms(x, m2Q2, nf) + MuDependentTerms(x, m2Q2, m2mu2, nf); } +//==========================================================================================// +// ExactCoefficientFunction: mu-independent terms +//------------------------------------------------------------------------------------------// + double ExactCoefficientFunction::MuIndependentTerms(double x, double m2Q2, int nf) const { return (this->*mu_indep_)(x, m2Q2, nf); } +//==========================================================================================// +// ExactCoefficientFunction: mu dependent terms +//------------------------------------------------------------------------------------------// + double ExactCoefficientFunction::MuDependentTerms(double x, double m2Q2, double m2mu2, int nf) const { return (this->*mu_dep_)(x, m2Q2, m2mu2, nf); } +//==========================================================================================// +// ExactCoefficientFunction: implemented only because it is pure virtual in base class. +// Returning three identical values +//------------------------------------------------------------------------------------------// + Value ExactCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { return Value(fx(x, m2Q2, m2mu2, nf)); } +//==========================================================================================// +// ExactCoefficientFunction: function that sets the pointer for mu_indep_ and mu_dep_ +//------------------------------------------------------------------------------------------// + void ExactCoefficientFunction::SetFunctions() { if (GetOrder() == 1) { if (GetKind() == '2') { @@ -167,6 +194,10 @@ void ExactCoefficientFunction::SetFunctions() { } } +//==========================================================================================// +// ExactCoefficientFunction: set method for method_flag +//------------------------------------------------------------------------------------------// + void ExactCoefficientFunction::SetMethodFlag(const int& method_flag) { // check method_flag if (method_flag != 0 && method_flag != 1) { @@ -177,6 +208,10 @@ void ExactCoefficientFunction::SetMethodFlag(const int& method_flag) { } +//==========================================================================================// +// ExactCoefficientFunction: set method for abserr +//------------------------------------------------------------------------------------------// + void ExactCoefficientFunction::SetAbserr(const double& abserr) { // check abserr if (abserr <= 0) { @@ -186,6 +221,10 @@ void ExactCoefficientFunction::SetAbserr(const double& abserr) { abserr_ = abserr; } +//==========================================================================================// +// ExactCoefficientFunction: set method for relerr +//------------------------------------------------------------------------------------------// + void ExactCoefficientFunction::SetRelerr(const double& relerr) { // check relerr if (relerr <= 0) { @@ -195,6 +234,10 @@ void ExactCoefficientFunction::SetRelerr(const double& relerr) { relerr_ = relerr; } +//==========================================================================================// +// ExactCoefficientFunction: set method for MCcalls +//------------------------------------------------------------------------------------------// + void ExactCoefficientFunction::SetMCcalls(const int& MCcalls) { // check MCcalls if (MCcalls <= 0) { @@ -204,6 +247,10 @@ void ExactCoefficientFunction::SetMCcalls(const int& MCcalls) { MCcalls_ = MCcalls; } +//==========================================================================================// +// ExactCoefficientFunction: set method for dim +//------------------------------------------------------------------------------------------// + void ExactCoefficientFunction::SetDim(const int& dim) { // check dim if (dim <= 0) { @@ -261,84 +308,12 @@ double ExactCoefficientFunction::CL_g1(double x, double m2Q2, int /*nf*/) const //------------------------------------------------------------------------------------------// //==========================================================================================// -// Exact massive gluon coefficient functions for F2 at O(as^2) -// -// Exact (but numerical) result from [arXiv:hep-ph/9411431]. -// Taken from the Fortran code 'src/hqcoef.f' -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::C2_g2(double x, double m2Q2, double m2mu2) const { - -// double xi = 1. / m2Q2; -// double eta = 0.25 * xi * (1 - x) / x - 1.; - -// if (eta > 1e6 || eta < 1e-6 || xi < 1e-3 || xi > 1e5) -// return 0.; - -// return C2_g20(x, m2Q2) + C2_g21(x, m2Q2) * log(1. / m2mu2); -// } - -//==========================================================================================// -// Exact massive quark coefficient functions for F2 at O(as) -// -// Exact (but numerical) result from [arXiv:hep-ph/9411431]. -// Taken from the Fortran code 'src/hqcoef.f' -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::C2_ps2(double x, double m2Q2, double m2mu2) const { - -// double xi = 1. / m2Q2; -// double eta = 0.25 * xi * (1 - x) / x - 1.; - -// if (eta > 1e6 || eta < 1e-6 || xi < 1e-3 || xi > 1e5) -// return 0.; - -// return C2_ps20(x, m2Q2) + C2_ps21(x, m2Q2) * log(1. / m2mu2); -// } - -//==========================================================================================// -// Exact massive gluon coefficient functions for FL at O(as) +// mu independent part of the exact massive quark coefficient functions for F2 at O(as^2) // // Exact (but numerical) result from [arXiv:hep-ph/9411431]. // Taken from the Fortran code 'src/hqcoef.f' //------------------------------------------------------------------------------------------// -// double ExactCoefficientFunction::CL_g2(double x, double m2Q2, double m2mu2) const { - -// double xi = 1. / m2Q2; -// double eta = 0.25 * xi * (1 - x) / x - 1.; - -// if (eta > 1e6 || eta < 1e-6 || xi < 1e-3 || xi > 1e5) -// return 0.; - -// return CL_g20(x, m2Q2) + CL_g21(x, m2Q2) * log(1. / m2mu2); -// } - -//==========================================================================================// -// Exact massive quarkk coefficient functions for FL at O(as) -// -// Exact (but numerical) result from [arXiv:hep-ph/9411431]. -// Taken from the Fortran code 'src/hqcoef.f' -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::CL_ps2(double x, double m2Q2, double m2mu2) const { - -// double xi = 1. / m2Q2; -// double eta = 0.25 * xi * (1 - x) / x - 1.; - -// if (eta > 1e6 || eta < 1e-6 || xi < 1e-3 || xi > 1e5) -// return 0.; - -// return CL_ps20(x, m2Q2) + CL_ps21(x, m2Q2) * log(1. / m2mu2); -// } - -//==========================================================================================// -// Exact massive quark coefficient functions for F2 at O(as^2): -// mu independent term. -// -// Eq. (4.4) of Ref. [arXiv:1205.5727] -//------------------------------------------------------------------------------------------// - double ExactCoefficientFunction::C2_ps20(double x, double m2Q2, int /*nf*/) const { double xi = 1. / m2Q2; @@ -351,7 +326,7 @@ double ExactCoefficientFunction::C2_ps20(double x, double m2Q2, int /*nf*/) cons } //==========================================================================================// -// Exact massive quark coefficient functions for F2 at O(as^2): +// Exact massive quark coefficient functions at O(as^2): // Term proportional to log(mu^2/m^2) // // Eq. (4.1) of Ref. [arXiv:1205.5727] @@ -369,10 +344,10 @@ double ExactCoefficientFunction::C_ps21(double x, double m2Q2) const { } //==========================================================================================// -// Exact massive quark coefficient functions for FL at O(as^2): -// mu independent term. +// mu independent part of the exact massive quark coefficient functions for FL at O(as) // -// Eq. (4.4) of Ref. [arXiv:1205.5727] +// Exact (but numerical) result from [arXiv:hep-ph/9411431]. +// Taken from the Fortran code 'src/hqcoef.f' //------------------------------------------------------------------------------------------// double ExactCoefficientFunction::CL_ps20(double x, double m2Q2, int /*nf*/) const { @@ -387,19 +362,10 @@ double ExactCoefficientFunction::CL_ps20(double x, double m2Q2, int /*nf*/) cons } //==========================================================================================// -// Exact massive quark coefficient functions for FL at O(as^2): -// Term proportional to log(mu^2/m^2) -// -// Eq. (4.1) of Ref. [arXiv:1205.5727] for FL -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::CL_ps21(double x, double m2Q2) const { return -CL_g1_x_Pgq0(x, m2Q2); } - -//==========================================================================================// -// Exact massive gluon coefficient functions for F2 at O(as^2): -// mu independent term. +// mu independent part of the exact massive gluon coefficient functions for F2 at O(as^2) // -// Eq. (4.4) of Ref. [arXiv:1205.5727] +// Exact (but numerical) result from [arXiv:hep-ph/9411431]. +// Taken from the Fortran code 'src/hqcoef.f' //------------------------------------------------------------------------------------------// double ExactCoefficientFunction::C2_g20(double x, double m2Q2, int /*nf*/) const { @@ -414,7 +380,7 @@ double ExactCoefficientFunction::C2_g20(double x, double m2Q2, int /*nf*/) const } //==========================================================================================// -// Exact massive gluon coefficient functions for F2 at O(as^2): +// Exact massive gluon coefficient functions at O(as^2): // Term proportional to log(mu^2/m^2) // // Eq. (4.4) of Ref. [arXiv:1205.5727] @@ -429,10 +395,10 @@ double ExactCoefficientFunction::C_g21(double x, double m2Q2) const { } //==========================================================================================// -// Exact massive gluon coefficient functions for FL at O(as^2): -// mu independent term. +// mu independent part of the exact massive gluon coefficient functions for FL at O(as^2) // -// Eq. (4.4) of Ref. [arXiv:1205.5727] +// Exact (but numerical) result from [arXiv:hep-ph/9411431]. +// Taken from the Fortran code 'src/hqcoef.f' //------------------------------------------------------------------------------------------// double ExactCoefficientFunction::CL_g20(double x, double m2Q2, int /*nf*/) const { @@ -447,32 +413,26 @@ double ExactCoefficientFunction::CL_g20(double x, double m2Q2, int /*nf*/) const } //==========================================================================================// -// Exact massive gluon coefficient functions for FL at O(as^2): -// Term proportional to log(mu^2/m^2) -// -// Eq. (4.4) of Ref. [arXiv:1205.5727] for FL +// ExactCoefficientFunction: mu dependent part of the exact massive gluon coefficient function at O(as^2): //------------------------------------------------------------------------------------------// -// double ExactCoefficientFunction::CL_g21(double x, double m2Q2) const { - -// int nf = 1; -// // Put nf to 1 since the nf contribution cancels for any value of nf - -// return -(CL_g1_x_Pgg0(x, m2Q2, nf) - CL_g1(x, m2Q2) * beta(0, nf)); -// } double ExactCoefficientFunction::C_g2_MuDep(double x, double m2Q2, double m2mu2, int /*nf*/) const { return C_g21(x, m2Q2) * log(1./m2mu2) ; } +//==========================================================================================// +// ExactCoefficientFunction: mu dependent part of the exact massive quark coefficient function at O(as^2): +//------------------------------------------------------------------------------------------// + double ExactCoefficientFunction::C_ps2_MuDep(double x, double m2Q2, double m2mu2, int /*nf*/) const { return C_ps21(x, m2Q2) * log(1./m2mu2) ; } //==========================================================================================// -// Exact massive quark coefficient functions for F2 at O(as^3): +// Exact massive quark coefficient functions at O(as^3): // Term proportional to log(mu^2/m^2) // // Eq. (4.2) of Ref. [arXiv:1205.5727] @@ -487,22 +447,7 @@ double ExactCoefficientFunction::C_ps31(double x, double m2Q2, int nf) const { } //==========================================================================================// -// Exact massive quark coefficient functions for FL at O(as^3): -// Term proportional to log(mu^2/m^2) -// -// Eq. (4.2) of Ref. [arXiv:1205.5727] for FL -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::CL_ps31(double x, double m2Q2, int nf) const { - -// return -( -// CL_g1_x_Pgq1(x, m2Q2, nf) + CL_g20_x_Pgq0(x, m2Q2) -// + CL_ps20_x_Pqq0(x, m2Q2, nf) - 2. * beta(0, nf) * CL_ps20(x, m2Q2) -// ); -// } - -//==========================================================================================// -// Exact massive quark coefficient functions for F2 at O(as^3): +// Exact massive quark coefficient functions at O(as^3): // Term proportional to log(mu^2/m^2)^2 // // Eq. (4.3) of Ref. [arXiv:1205.5727] @@ -514,21 +459,6 @@ double ExactCoefficientFunction::C_ps32(double x, double m2Q2, int nf) const { - 3./2 * beta(0, nf) * convolutions_lmu2_[2] -> Convolute(x, m2Q2, nf); } -//==========================================================================================// -// Exact massive quark coefficient functions for FL at O(as^3): -// Term proportional to log(mu^2/m^2)^2 -// -// Eq. (4.3) of Ref. [arXiv:1205.5727] for FL -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::CL_ps32(double x, double m2Q2, int nf) const { - -// return 0.5 -// * (CL_g1_x_Pgg0_x_Pgq0(x, m2Q2, nf) -// + CL_g1_x_Pqq0_x_Pgq0(x, m2Q2, nf)) -// - 3. / 2 * beta(0, nf) * CL_g1_x_Pgq0(x, m2Q2); -// } - //==========================================================================================// // Exact massive gluon coefficient functions for F2 at O(as^3): // Term proportional to log(mu^2/m^2) @@ -546,23 +476,7 @@ double ExactCoefficientFunction::C_g31(double x, double m2Q2, int nf) const { } //==========================================================================================// -// Exact massive gluon coefficient functions for FL at O(as^3): -// Term proportional to log(mu^2/m^2) -// -// Eq. (4.5) of Ref. [arXiv:1205.5727] for FL -//------------------------------------------------------------------------------------------// - -// double ExactCoefficientFunction::CL_g31(double x, double m2Q2, int nf) const { - -// return -( -// CL_g1_x_Pgg1(x, m2Q2, nf) - beta(1, nf) * CL_g1(x, m2Q2) -// + CL_ps20_x_Pqg0(x, m2Q2, nf) + CL_g20_x_Pgg0(x, m2Q2, nf) -// - 2 * beta(0, nf) * CL_g20(x, m2Q2) -// ); -// } - -//==========================================================================================// -// Exact massive gluon coefficient functions for F2 at O(as^3): +// Exact massive gluon coefficient functions at O(as^3): // Term proportional to log(mu^2/m^2)^2 // // Eq. (4.6) of Ref. [arXiv:1205.5727] @@ -570,54 +484,17 @@ double ExactCoefficientFunction::C_g31(double x, double m2Q2, int nf) const { double ExactCoefficientFunction::C_g32(double x, double m2Q2, int nf) const { - // double C2_g1xPgg0xPgg0; - double beta0 = beta(0, nf); - // if (method_flag_ == 0) - // C2_g1xPgg0xPgg0 = C2_g1_x_Pgg0_x_Pgg0(x, m2Q2, nf); - // else if (method_flag_ == 1) - // C2_g1xPgg0xPgg0 = C2_g1_x_Pgg0_x_Pgg0_MC(x, m2Q2, nf); - // else { - // cout << "C2_g32: Choose either method_flag = 0 or method_flag = 1" - // << endl; - // exit(-1); - // } - return 0.5 * (convolutions_lmu2_[0] -> Convolute(x, m2Q2, nf) + convolutions_lmu2_[1] -> Convolute(x, m2Q2, nf)) - 3./2 * beta0 * convolutions_lmu2_[2] -> Convolute(x, m2Q2, nf) + beta0 * beta0 * convolutions_lmu2_[3] -> Convolute(x, m2Q2, nf); } //==========================================================================================// -// Exact massive gluon coefficient functions for FL at O(as^3): -// Term proportional to log(mu^2/m^2)^2 -// -// Eq. (4.6) of Ref. [arXiv:1205.5727] for FL +// ExactCoefficientFunction: mu dependent part of the exact massive gluon coefficient function at O(as^3): //------------------------------------------------------------------------------------------// -// double ExactCoefficientFunction::CL_g32(double x, double m2Q2, int nf) const { - -// double CL_g1xPgg0xPgg0; - -// double beta0 = beta(0, nf); - -// if (method_flag_ == 0) -// CL_g1xPgg0xPgg0 = CL_g1_x_Pgg0_x_Pgg0(x, m2Q2, nf); -// else if (method_flag_ == 1) -// CL_g1xPgg0xPgg0 = CL_g1_x_Pgg0_x_Pgg0_MC(x, m2Q2, nf); -// else { -// std::cout << "C2_g32: Choose either method_flag = 0 or method_flag = 1" -// << std::endl; -// exit(-1); -// } - -// return 0.5 * CL_g1xPgg0xPgg0 + 0.5 * CL_g1_x_Pqg0_x_Pgq0(x, m2Q2, nf) -// - 3. / 2 * beta0 * CL_g1_x_Pgg0(x, m2Q2, nf) -// + beta0 * beta0 * CL_g1(x, m2Q2); -// } - - double ExactCoefficientFunction::C_g3_MuDep(double x, double m2Q2, double m2mu2, int nf) const { double lmu = log(1. / m2mu2) ; @@ -625,6 +502,10 @@ double ExactCoefficientFunction::C_g3_MuDep(double x, double m2Q2, double m2mu2, return C_g31(x, m2Q2, nf) * lmu + C_g32(x, m2Q2, nf) * lmu * lmu ; } +//==========================================================================================// +// ExactCoefficientFunction: mu dependent part of the exact massive quark coefficient function at O(as^3): +//------------------------------------------------------------------------------------------// + double ExactCoefficientFunction::C_ps3_MuDep(double x, double m2Q2, double m2mu2, int nf) const { double lmu = log(1. / m2mu2) ; @@ -632,6 +513,10 @@ double ExactCoefficientFunction::C_ps3_MuDep(double x, double m2Q2, double m2mu2 return C_ps31(x, m2Q2, nf) * lmu + C_ps32(x, m2Q2, nf) * lmu * lmu ; } +//==========================================================================================// +// ExactCoefficientFunction: print warning for mu independent part of O(as^3) +//------------------------------------------------------------------------------------------// + double ExactCoefficientFunction::WarningFunc(double /*x*/, double /*m2Q2*/, int /*nf*/) const { cout << "Error: mu-independent terms of the exact coefficient function at O(a_s^3) are not known!" << endl; exit(-1); diff --git a/src/HighEnergyCoefficientFunction.cc b/src/HighEnergyCoefficientFunction.cc index 53b52cf..ca47a87 100644 --- a/src/HighEnergyCoefficientFunction.cc +++ b/src/HighEnergyCoefficientFunction.cc @@ -7,20 +7,36 @@ using std::cout; using std::endl; +//==========================================================================================// +// AbstractHighEnergyCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + AbstractHighEnergyCoefficientFunction::AbstractHighEnergyCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL) : CoefficientFunction(order, kind, channel) { SetNLL(NLL); } +//==========================================================================================// +// HighEnergyCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + HighEnergyCoefficientFunction::HighEnergyCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL) : AbstractHighEnergyCoefficientFunction(order, kind, channel, NLL) { SetFunctions(); } +//==========================================================================================// +// HighEnergyCoefficientFunction: band of the high energy coefficient function +//------------------------------------------------------------------------------------------// + Value HighEnergyCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { return (this->*LL_)(x, m2Q2, m2mu2) + (this->*NLL_)(x, m2Q2, m2mu2, nf) ; } +//==========================================================================================// +// HighEnergyCoefficientFunction: function that sets the pointer for LL_ and NLL_ +//------------------------------------------------------------------------------------------// + void HighEnergyCoefficientFunction::SetFunctions() { if (GetOrder() == 1) { @@ -58,14 +74,26 @@ void HighEnergyCoefficientFunction::SetFunctions() { if (!GetNLL()) NLL_ = &HighEnergyCoefficientFunction::ZeroFunctionBand; } +//==========================================================================================// +// HighEnergyHighScaleCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + HighEnergyHighScaleCoefficientFunction::HighEnergyHighScaleCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL) : AbstractHighEnergyCoefficientFunction(order, kind, channel, NLL) { SetFunctions(); } +//==========================================================================================// +// HighEnergyHighScaleCoefficientFunction: band of the high scale limit of the high energy +//------------------------------------------------------------------------------------------// + Value HighEnergyHighScaleCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { return (this->*LL_)(x, m2Q2, m2mu2) + (this->*NLL_)(x, m2Q2, m2mu2, nf) ; } +//==========================================================================================// +// HighEnergyHighScaleCoefficientFunction: function that sets the pointer for LL_ and NLL_ +//------------------------------------------------------------------------------------------// + void HighEnergyHighScaleCoefficientFunction::SetFunctions() { if (GetOrder() == 1) { @@ -103,18 +131,30 @@ void HighEnergyHighScaleCoefficientFunction::SetFunctions() { if (!GetNLL()) NLL_ = &HighEnergyHighScaleCoefficientFunction::ZeroFunctionBand; } +//==========================================================================================// +// HighEnergyHighScaleCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + PowerTermsCoefficientFunction::PowerTermsCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& NLL) : AbstractHighEnergyCoefficientFunction(order, kind, channel, NLL) { highenergy_ = new HighEnergyCoefficientFunction(GetOrder(), GetKind(), GetChannel(), GetNLL()); highenergyhighscale_ = new HighEnergyHighScaleCoefficientFunction(GetOrder(), GetKind(), GetChannel(), GetNLL()); } +//==========================================================================================// +// HighEnergyHighScaleCoefficientFunction: destructor +//------------------------------------------------------------------------------------------// + PowerTermsCoefficientFunction::~PowerTermsCoefficientFunction() { delete highenergy_; delete highenergyhighscale_; } +//==========================================================================================// +// PowerTermsCoefficientFunction: band of the power terms +//------------------------------------------------------------------------------------------// + Value PowerTermsCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { - // TODO: in this way the error is very small: should take all the compbinations + // TODO: in this way the error is very small: should take all the combinations double central = (highenergy_->fx(x, m2Q2, m2mu2, nf)) - (highenergyhighscale_->fx(x, m2Q2, m2mu2, nf)); double higher = (highenergy_->fxBand(x, m2Q2, m2mu2, nf)).GetHigher() - (highenergyhighscale_->fxBand(x, m2Q2, m2mu2, nf)).GetHigher(); double lower = (highenergy_->fxBand(x, m2Q2, m2mu2, nf)).GetLower() - (highenergyhighscale_->fxBand(x, m2Q2, m2mu2, nf)).GetLower(); diff --git a/src/HighScaleCoefficientFunction.cc b/src/HighScaleCoefficientFunction.cc index 1b43255..28afd6e 100644 --- a/src/HighScaleCoefficientFunction.cc +++ b/src/HighScaleCoefficientFunction.cc @@ -1,6 +1,5 @@ #include "adani/HighScaleCoefficientFunction.h" #include "adani/Constants.h" -#include "adani/MasslessCoefficientFunction.h" #include "adani/SpecialFunctions.h" #include @@ -9,6 +8,10 @@ using std::cout; using std::endl; +//==========================================================================================// +// HighScaleCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + HighScaleCoefficientFunction::HighScaleCoefficientFunction(const int& order, const char& kind, const char& channel, const bool& exact, const bool& revised_approx) : CoefficientFunction(order, kind, channel) { massless_lo_ = nullptr; massless_nlo_ = nullptr; @@ -29,6 +32,10 @@ HighScaleCoefficientFunction::HighScaleCoefficientFunction(const int& order, con } +//==========================================================================================// +// HighScaleCoefficientFunction: destructor +//------------------------------------------------------------------------------------------// + HighScaleCoefficientFunction::~HighScaleCoefficientFunction() { delete massless_lo_; delete massless_nlo_; @@ -36,14 +43,26 @@ HighScaleCoefficientFunction::~HighScaleCoefficientFunction() { delete a_muindep_; } +//==========================================================================================// +// HighScaleCoefficientFunction: central value of the highscale coefficient function +//------------------------------------------------------------------------------------------// + double HighScaleCoefficientFunction::fx(double x, double m2Q2, double m2mu2, int nf) const { return fxBand(x, m2Q2, m2mu2, nf).GetCentral(); } +//==========================================================================================// +// HighScaleCoefficientFunction: band of the highscale coefficient function +//------------------------------------------------------------------------------------------// + Value HighScaleCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { return (this->*fx_)(x, m2Q2, m2mu2, nf); } +//==========================================================================================// +// HighScaleCoefficientFunction: function that sets the pointer for fxBand +//------------------------------------------------------------------------------------------// + void HighScaleCoefficientFunction::SetFunctions() { if (GetOrder() == 1) {