From 7564d21791ff7a91d9333d993110c6e306c78f21 Mon Sep 17 00:00:00 2001 From: niccolo Date: Wed, 20 Mar 2024 11:18:39 +0100 Subject: [PATCH] Move some stuff --- inc/adani/Convolution.h | 26 ++++++---------- inc/adani/HighScaleSplitLogs.h | 5 +++ src/Convolution.cc | 56 +++++++++++++++++++++++++++++++--- src/HighScaleSplitLogs.cc | 29 +++++++++++++++++- tests/test_mc_integral.py | 4 +-- 5 files changed, 96 insertions(+), 24 deletions(-) diff --git a/inc/adani/Convolution.h b/inc/adani/Convolution.h index 05bb216..28ee1ed 100644 --- a/inc/adani/Convolution.h +++ b/inc/adani/Convolution.h @@ -105,30 +105,24 @@ class Convolution : public AbstractConvolution { class ConvolutedCoefficientFunction : public CoefficientFunction { public: - ConvolutedCoefficientFunction(Convolution *conv) - : CoefficientFunction(conv->GetCoeffFunc()) { - conv_ = conv; - }; - ~ConvolutedCoefficientFunction() override{}; + ConvolutedCoefficientFunction(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double &abserr = 1e-3, const double &relerr = 1e-3, const int &dim = 1000); + ~ConvolutedCoefficientFunction() override; - double MuIndependentTerms(double x, double m2Q2, int nf) const override { - return conv_->Convolute(x, m2Q2, nf); - }; + // get method + Convolution* GetConv() const {return conv_;}; + + double MuIndependentTerms(double x, double m2Q2, int nf) const override; double MuDependentTerms( double /*x*/, double /*m2Q2*/, double /*m2mu2*/, int /*nf*/ ) const override { return 0.; }; - double fx(double x, double m2Q2, double m2mu2, int nf) const override { - return MuIndependentTerms(x, m2Q2, nf) - + MuDependentTerms(x, m2Q2, m2mu2, nf); - } - Value fxBand(double x, double m2Q2, double m2mu2, int nf) const override { - return Value(fx(x, m2Q2, m2mu2, nf)); - } + double fx(double x, double m2Q2, double m2mu2, int nf) const override; + Value fxBand(double x, double m2Q2, double m2mu2, int nf) const override; private: - Convolution *conv_; + Convolution* conv_; + }; //==========================================================================================// diff --git a/inc/adani/HighScaleSplitLogs.h b/inc/adani/HighScaleSplitLogs.h index de0f92f..82cf1d3 100644 --- a/inc/adani/HighScaleSplitLogs.h +++ b/inc/adani/HighScaleSplitLogs.h @@ -29,6 +29,10 @@ // (The expressions do not contain the log) //------------------------------------------------------------------------------------------// +//==========================================================================================// +// class HighScaleSplitLogs +//------------------------------------------------------------------------------------------// + class HighScaleSplitLogs : public CoefficientFunction { public: HighScaleSplitLogs( @@ -46,6 +50,7 @@ class HighScaleSplitLogs : public CoefficientFunction { const override; Value fxBand(double x, double m2Q2, int nf) const; + // division of the total result in the log terms double LL(double x, int nf) const { return (this->*LL_)(x, nf); } double NLL(double x, int nf) const { return (this->*NLL_)(x, nf); } double N2LL(double x, int nf) const { return (this->*N2LL_)(x, nf); } diff --git a/src/Convolution.cc b/src/Convolution.cc index dc10f8b..93a73b2 100644 --- a/src/Convolution.cc +++ b/src/Convolution.cc @@ -98,7 +98,7 @@ double Convolution::regular_integrand(double z, void *p) { } //==========================================================================================// -// AbstractConvolution: integrand of the singular part +// Convolution: integrand of the singular part //------------------------------------------------------------------------------------------// double Convolution::singular_integrand(double z, void *p) { @@ -127,7 +127,7 @@ double Convolution::singular_integrand(double z, void *p) { //------------------------------------------------------------------------------------------// //==========================================================================================// -// AbstractConvolution: regular part +// Convolution: regular part //------------------------------------------------------------------------------------------// double Convolution::RegularPart(double x, double m2Q2, int nf) const { @@ -163,7 +163,7 @@ double Convolution::RegularPart(double x, double m2Q2, int nf) const { } //==========================================================================================// -// AbstractConvolution: singular part +// Convolution: singular part //------------------------------------------------------------------------------------------// double Convolution::SingularPart(double x, double m2Q2, int nf) const { @@ -200,7 +200,7 @@ double Convolution::SingularPart(double x, double m2Q2, int nf) const { } //==========================================================================================// -// AbstractConvolution: local part +// Convolution: local part //------------------------------------------------------------------------------------------// double Convolution::LocalPart(double x, double m2Q2, int nf) const { @@ -212,6 +212,52 @@ double Convolution::LocalPart(double x, double m2Q2, int nf) const { - splitfunc_->SingularIntegrated(x / x_max, nf)); } +//==========================================================================================// +// ConvolutedCoefficientFunction: constructor +//------------------------------------------------------------------------------------------// + +ConvolutedCoefficientFunction::ConvolutedCoefficientFunction(CoefficientFunction* coefffunc, AbstractSplittingFunction* splitfunc, const double &abserr, const double &relerr, const int &dim) : CoefficientFunction(coefffunc) { + + conv_ = new Convolution(coefffunc, splitfunc, abserr, relerr, dim); +} + +//==========================================================================================// +// ConvolutedCoefficientFunction: destructor +//------------------------------------------------------------------------------------------// + +ConvolutedCoefficientFunction::~ConvolutedCoefficientFunction() { + + delete conv_; +} + +//==========================================================================================// +// ConvolutedCoefficientFunction: mu independent terms +//------------------------------------------------------------------------------------------// + +double ConvolutedCoefficientFunction::MuIndependentTerms(double x, double m2Q2, int nf) const { + + return conv_->Convolute(x, m2Q2, nf); +} + +//==========================================================================================// +// ConvolutedCoefficientFunction: central value of fx +//------------------------------------------------------------------------------------------// + +double ConvolutedCoefficientFunction::fx(double x, double m2Q2, double m2mu2, int nf) const { + + return MuIndependentTerms(x, m2Q2, nf) + + MuDependentTerms(x, m2Q2, m2mu2, nf); +} + +//==========================================================================================// +// ConvolutedCoefficientFunction: band of fx +//------------------------------------------------------------------------------------------// + +Value ConvolutedCoefficientFunction::fxBand(double x, double m2Q2, double m2mu2, int nf) const { + + return Value(fx(x, m2Q2, m2mu2, nf)); +} + //==========================================================================================// // DoubleConvolution: constructor //------------------------------------------------------------------------------------------// @@ -232,7 +278,7 @@ DoubleConvolution::DoubleConvolution( conv_coeff_ = nullptr; } else { conv_coeff_ = new ConvolutedCoefficientFunction( - new Convolution(coefffunc, splitfunc, abserr, relerr, dim) + coefffunc, splitfunc, abserr, relerr, dim ); convolution_ = new Convolution(conv_coeff_, splitfunc, abserr, relerr, dim); diff --git a/src/HighScaleSplitLogs.cc b/src/HighScaleSplitLogs.cc index 587715f..07fa64f 100644 --- a/src/HighScaleSplitLogs.cc +++ b/src/HighScaleSplitLogs.cc @@ -1,6 +1,5 @@ #include "adani/HighScaleSplitLogs.h" #include "adani/Constants.h" -#include "adani/Convolution.h" #include "adani/SpecialFunctions.h" #include @@ -9,6 +8,10 @@ using std::cout; using std::endl; +//==========================================================================================// +// HighScaleSplitLogs: constructor +//------------------------------------------------------------------------------------------// + HighScaleSplitLogs::HighScaleSplitLogs( const int &order, const char &kind, const char &channel, const bool &exact, const bool &revised_approx @@ -39,12 +42,20 @@ HighScaleSplitLogs::HighScaleSplitLogs( SetFunctions(); } +//==========================================================================================// +// HighScaleSplitLogs: destructor +//------------------------------------------------------------------------------------------// + HighScaleSplitLogs::~HighScaleSplitLogs() { delete massless_; delete massless_lo_; delete a_muindep_; } +//==========================================================================================// +// HighScaleSplitLogs: warning function that prevents from calling fx with Q != mu. +//------------------------------------------------------------------------------------------// + double HighScaleSplitLogs:: fx(double /*x*/, double /*m2Q2*/, double /*m2mu2*/, int /*nf*/) const { cout << "Error: HighScaleSplitLogs is implemented only in the case mu=Q!" @@ -53,10 +64,18 @@ double HighScaleSplitLogs:: exit(-1); } +//==========================================================================================// +// HighScaleSplitLogs: central value of the band +//------------------------------------------------------------------------------------------// + double HighScaleSplitLogs::fx(double x, double m2Q2, int nf) const { return fxBand(x, m2Q2, nf).GetCentral(); } +//==========================================================================================// +// HighScaleSplitLogs: warning function that prevents from calling fxBand with Q != mu. +//------------------------------------------------------------------------------------------// + Value HighScaleSplitLogs:: fxBand(double /*x*/, double /*m2Q2*/, double /*m2mu2*/, int /*nf*/) const { cout << "Error: HighScaleSplitLogs is implemented only in the case mu=Q!" @@ -66,6 +85,10 @@ Value HighScaleSplitLogs:: exit(-1); } +//==========================================================================================// +// HighScaleSplitLogs: band of the total result +//------------------------------------------------------------------------------------------// + Value HighScaleSplitLogs::fxBand(double x, double m2Q2, int nf) const { double Log = log(m2Q2); double Log2 = Log * Log; @@ -74,6 +97,10 @@ Value HighScaleSplitLogs::fxBand(double x, double m2Q2, int nf) const { + N3LL(x, nf); } +//==========================================================================================// +// HighScaleSplitLogs: function that sets the pointers to the correct function +//------------------------------------------------------------------------------------------// + void HighScaleSplitLogs::SetFunctions() { if (GetKind() == '2') { diff --git a/tests/test_mc_integral.py b/tests/test_mc_integral.py index 01e455c..a42f912 100644 --- a/tests/test_mc_integral.py +++ b/tests/test_mc_integral.py @@ -7,11 +7,11 @@ def test_mc_integrals(): for kind in ['2', 'L']: massive_mc = ad.ExactCoefficientFunction(3, kind, 'g', 1e-3, 1e-3, 1000, 1, 1000000) massive_nomc = ad.ExactCoefficientFunction(3, kind, 'g', 1e-3, 1e-3, 1000, 0, 25000) - for xi in np.geomspace(1e-2, 1e4, 4, endpoint=True): + for xi in np.geomspace(1e-2, 1e4, 5, endpoint=True): m2Q2 = 1/xi xmax = 1. / (1. + 4 * m2Q2) for x in np.geomspace(1e-5, xmax, 5, endpoint=False): - for nf in range(7): + for nf in range(1, 6 + 1): res1 = massive_mc.MuDependentTerms(x, m2Q2, m2Q2, nf) res2 = massive_nomc.MuDependentTerms(x, m2Q2, m2Q2, nf) np.testing.assert_allclose(res1, res2, rtol=5e-3)