Merge pull request #6 from N3PDF/x-space

[CLEAN UP]: x-space resummed expressions
N3PDF · Jun 9, 2021 · c149031 · c149031
2 parents 517bec4 2c0c0b5
commit c149031
Show file tree

Hide file tree

Showing 10 changed files with 426 additions and 2 deletions.
diff --git a/include/CombinedRes.h b/include/CombinedRes.h
@@ -30,6  30,7 @@
 #include "./ComplexDefs.h"
 #include "./SmallptExp.h"
 #include "./ThresExp.h"
 #include "./IntMellin.h"
 #include "higgs-fo/partonic.h"
 
 class CombinedRes {
@@ -50,6  51,7 @@ class CombinedRes {
   // Init. resummation classes
   SmallptExp *SMALLPT;
   ThresExp *THRESHOLD;
   MellinTrans *MELLIN;
 
   // Init. exact FO class
   CrossHiggs *MELLINPARTONIC;

diff --git a/include/IntMellin.h b/include/IntMellin.h
@@ -0,0  1,48 @@
 /*
  * =====================================================================================
  *
  *       Filename:
  *
  *    Description:
  *
  *        Version:  1.0
  *        Created:
  *       Revision:  none
  *       Compiler:  gcc
  *
  *         Author:
  *   Organization:  N3PDF
  *
  * =====================================================================================
  */
 
 #pragma once
 
 #include <math.h>
 #include <stdio.h>
 #include <stdlib.h>
 
 #include <cmath>
 #include <iostream>
 #include <string>
 #include <vector>
 
 #include "./ComplexDefs.h"
 #include "./Integration.h"
 #include "./ThresXspace.h"
 
 class MellinTrans {
  public:
   MellinTrans(int order, int channel, std::string pdfname, void *params);
   virtual ~MellinTrans();
 
   double ExtractThresMom(double x, double N, double pt);
   std::complex<long double> xSpaceThres(std::complex<long double> N,
 										long double pt);
 
  private:
   int NF, ORD, CHANNEL;
   double MH2, MUR2, MUF2;
 
   ThresXspace xThres;
 };
diff --git a/include/Integration.h b/include/Integration.h
@@ -0,0  1,38 @@
 /*
  * =====================================================================================
  *
  *       Filename:
  *
  *    Description:
  *
  *        Version:  1.0
  *        Created:
  *       Revision:  none
  *       Compiler:  gcc
  *
  *         Author:
  *   Organization:  N3PDF
  *
  * =====================================================================================
  */
 
 #pragma once
 
 #include <math.h>
 #include <stdio.h>
 #include <stdlib.h>
 
 #include <cmath>
 #include <iostream>
 #include <string>
 #include <vector>
 
 #include "./ComplexDefs.h"
 #include "./SmallptExp.h"
 #include "./ThresExp.h"
 #include "higgs-fo/integration.h"
 
 namespace ExtrIntegration {
 double IntegrateOverx(int method, double(Func)(double, void *), void *pp,
                       double *error);
 }
diff --git a/include/ThresXspace.h b/include/ThresXspace.h
@@ -0,0  1,58 @@
 /*
  * =====================================================================================
  *
  *       Filename:
  *
  *    Description:
  *
  *        Version:  1.0
  *        Created:
  *       Revision:  none
  *       Compiler:  gcc
  *
  *         Author:
  *   Organization:  N3PDF
  *
  * =====================================================================================
  */
 
 #pragma once
 
 #include <gsl/gsl_math.h>
 #include <gsl/gsl_sf_psi.h>
 #include <math.h>
 #include <stdio.h>
 #include <stdlib.h>
 
 #include <cmath>
 #include <iostream>
 #include <string>
 #include <vector>
 
 #include "higgs-fo/params.h"
 
 class ThresXspace {
  public:
   ThresXspace(int order, int channel, void *params);
   virtual ~ThresXspace();
 
   double ThresXspaceExpr(double x, double N, double pt);
 
   // LO pt-distribution
   double LOgggH(double x, double pt);
   double LOgqqH(double x, double pt);
   double LOqqgH(double x, double pt);
 
  private:
   int NC, NF, CA, ORD, CHANNEL;
   double LF, LR, LQ;
   double MH2, MUF2, MUR2, SROOT;
   double CF, aass, SIGMA0;
 
   double Beta0;
 
   // Polygamma functions
   double psi0one = gsl_sf_psi_n(0, 1);
   double psi0half = gsl_sf_psi_n(0, 0.5);
   double psi1half = gsl_sf_psi_n(1, 0.5);
 };
diff --git a/include/meson.build b/include/meson.build
@@ -6,7  6,10 @@ install_headers(
   	'AnomalousDim.h',
     'HarmonicSum.h',
     'SmallptExp.h',
-    'CombinedRes.h'
     'CombinedRes.h',
     'Integration.h',
 	'ThresXspace.h',
     'IntMellin.h'
   ],
   subdir : 'higgs-n3lo'
 )
diff --git a/src/CombinedRes.cpp b/src/CombinedRes.cpp
@@ -35,10  35,12 @@ CombinedRes::CombinedRes(int order, int channel, std::string pdfname,
 
   SMALLPT = new SmallptExp(order, channel, params);
   THRESHOLD = new ThresExp(order, channel, params);
   MELLIN = new MellinTrans(order, channel, pdfname, params);
   MELLINPARTONIC = new CrossHiggs(EXACT_ORD, channel, pdfname, params);
 }
 
 CombinedRes::~CombinedRes() {
   delete MELLIN;
   delete MELLINPARTONIC;
   delete SMALLPT;
   delete THRESHOLD;
@@ -84,10  86,12 @@ std::complex<long double> CombinedRes::CombinedResExpr(
   // Compute approximation from resummations
   std::complex<long double> SptMellin = SMALLPT->SmallptExpExpr(N, pt);
   std::complex<long double> ThresMellin = THRESHOLD->ThresExpExpr(N, pt);
   std::complex<long double> xThresMellin = MELLIN->xSpaceThres(N, pt);
 
   std::complex<long double> ExactMellinCmpx(ExactMellin[0], 0.);
   mres = (1. - Matching(N, pt, scheme)) * SptMellin  
          Matching(N, pt, scheme) * ThresMellin;
          /* Matching(N, pt, scheme) * xThresMellin; */
 
   return ExactMellinCmpx   mres;
 }
diff --git a/src/IntMellin.cpp b/src/IntMellin.cpp
@@ -0,0  1,73 @@
 /*
  * =====================================================================================
  *
  *       Filename:
  *
  *    Description:
  *
  *        Version:  1.0
  *        Created:
  *       Revision:  none
  *       Compiler:  gcc
  *
  *         Author:
  *   Organization:  N3PDF
  *
  * =====================================================================================
  */
 
 #include "../include/IntMellin.h"
 
 #include <higgs-fo/params.h>
 
 MellinTrans::MellinTrans(int order, int channel, std::string pdfname,
                          void *params)
     : xThres(order, channel, params) {
   PhysParams param = *reinterpret_cast<PhysParams *>(params);
 
   MH2 = static_cast<long double>(std::pow(param.mh, 2));
   MUR2 = static_cast<long double>(std::pow(param.mur, 2));
   MUF2 = static_cast<long double>(std::pow(param.muf, 2));
 }
 
 
 MellinTrans::~MellinTrans() {}
 
 
 struct IntStruct {
   MellinTrans *TransMellin;
   double N;
   double pt;
 };
 
 
 double MellinTrans::ExtractThresMom(double x, double N, double pt) {
   return xThres.ThresXspaceExpr(x, N, pt);
 }
 
 
 double xThresIntegrand(double x, void *p) {
   IntStruct par = *reinterpret_cast<IntStruct *>(p);
   return par.TransMellin->ExtractThresMom(x, par.N, par.pt);
 }
 
 
 std::complex<long double> MellinTrans::xSpaceThres(
     std::complex<long double> N, long double pt) {
   int method = 0;
   double reslt, error;
   double nn = static_cast<double>(N.real());
 
   // pass parameters
   IntStruct finalparams;
   finalparams.N = nn;
   finalparams.pt = static_cast<double>(pt);
   finalparams.TransMellin = this;
 
   // Integration
   reslt = ExtrIntegration::IntegrateOverx(method, xThresIntegrand, &finalparams,
                                           &error);
 
   std::complex<long double> finreslt(static_cast<long double>(reslt), 0.);
   return finreslt;
 }
diff --git a/src/Integration.cpp b/src/Integration.cpp
@@ -0,0  1,59 @@
 /*
  * =====================================================================================
  *
  *       Filename:
  *
  *    Description:
  *
  *        Version:  1.0
  *        Created:
  *       Revision:  none
  *       Compiler:  gcc
  *
  *         Author:
  *   Organization:  N3PDF
  *
  * =====================================================================================
  */
 
 #include "../include/Integration.h"
 
 struct IntData {
   std::function<double(double, void *)> fun;
   void *param;
 };
 
 
 int xIntegrand(int *ndim, double *x, int *ncomp, double *y, void *p) {
   IntData inparam = *reinterpret_cast<IntData *>(p);
   y[0] = inparam.fun(x[0], inparam.param);
   return 0;
 }
 
 
 double ExtrIntegration::IntegrateOverx(int method, double(Func)(double, void *),
                                        void *pp, double *error) {
   int fail;
   double prec = 1e-8;
   double *res = NULL;
   double *err = NULL;
   double *prb = NULL;
 
   res = new double[1];
   err = new double[1];
   prb = new double[1];
 
   int ndim = 2;
   int ncomp = 1;
   int verbose = 0;
 
   IntData UserData;
   UserData.fun = Func;
   UserData.param = pp;
 
   integration::CubaIntegrator(method, xIntegrand, ndim, ncomp, prec, &res,
                               &fail, &err, &prb, &UserData, verbose);
 
   *error = err[0];
   return res[0];
 }