PaperlayerView exact revision on arXiv

Source-linked AI summary

QCDNUM: Fast QCD Evolution and Convolution

M. Botje

arXiv:1005.1481v3hep-ph

TL;DR

QCDNUM addresses the need for numerical evolution of QCD parton densities and fragmentation functions and for convolution-based observable calculations. It uses discrete-grid spline evolution and weighted convolution sums, supporting NNLO un-polarised evolution and a broad set of schemes and inputs. The program evolves a full NNLO parton-density set over a large kinematic range in less than 10 ms on a 2 GHz processor.

  • Problem

    QCDNUM addresses numerical evolution of QCD parton densities and fragmentation functions and computation of observables from evolved densities and partonic cross-sections.

  • Method

    QCDNUM uses spline representations on discrete x and µ2 grids, forward-substitution evolution, and weighted convolution sums from user-provided kernels.

  • Results

    Less than 10 ms evolves the gluon density and 12 quark densities at NNLO over a large kinematic range on a 2 GHz Intel Core 2 Duo processor.

  • Takeaways & Limitations

    QCDNUM supports NNLO un-polarised evolution, NLO polarised or fragmentation-function evolution, alternative kernels, external density libraries, and multiple observable schemes.

  • Takeaways & Limitations

    Backward evolution of quadratic splines can sometimes lead to severe oscillations.

Abstract

from arXiv · show

The QCDNUM program numerically solves the evolution equations for parton densities and fragmentation functions in perturbative QCD. Un-polarised parton densities can be evolved up to next-to-next-to-leading order in powers of the strong coupling constant, while polarised densities or fragmentation functions can be evolved up to next-to-leading order. Other types of evolution can be accessed by feeding alternative sets of evolution kernels into the program. A versatile convolution engine provides tools to compute parton luminosities, cross-sections in hadron-hadron scattering, and deep inelastic structure functions in the zero-mass scheme or in generalised mass schemes. Input to these calculations are either the QCDNUM evolved densities, or those read in from an external parton density repository. Included in the software distribution are packages to calculate zero-mass structure functions in un-polarised deep inelastic scattering, and heavy flavour contributions to these structure functions in the fixed flavour number scheme.

PROGRAM SUMMARY

QCDNUM evolves parton densities and computes observables through spline-based numerical methods and fast convolutions. It provides broad QCD functionality, with speed and accuracy governed by the evolution-grid density.

  • PROGRAM SUMMARY: The software is distributed as a gzipped tar file under the GNU Public License and is implemented in Fortran-77 without external libraries beyond included add-on packages.
  • PROGRAM SUMMARY: QCDNUM evolves the strong coupling and parton densities up to NNLO and computes observables through Mellin convolution with partonic cross-sections.Its method represents densities with linear or quadratic splines and evolves spline coefficients by forward substitution.
  • PROGRAM SUMMARY: Accuracy and speed are determined by the density of the evolution grid.
  • PROGRAM SUMMARY: Less than 10 ms evolves the gluon density and 12 quark densities at NNLO over a large kinematic range on a 2 GHz Intel Core 2 Duo processor.

1 Introduction

QCDNUM17 upgrades QCDNUM16 with quadratic-spline evolution that improves accuracy and speed while supporting multiple evolution orders, schemes, and convolution-based observables. The program solves DGLAP evolution numerically on discrete grids and supports external parton-density inputs.

  • 1 Introduction: QCDNUM numerically solves DGLAP evolution equations on a discrete grid in x and µ2 using parton densities and αs as inputs.The renormalisation scale can be varied relative to the mass factorisation scale to study scale uncertainties.
  • 1 Introduction: QCDNUM17 is the NNLO upgrade of QCDNUM16, with quadratic-spline evolution yielding large gains in accuracy and speed.On a 2 GHz processor, it evolves the full parton-density set at NNLO over a large kinematic range in less than 10 ms.
  • 1 Introduction: Un-polarised densities evolve up to NNLO, while polarised densities and fragmentation functions evolve up to NLO; alternative evolution kernels are also supported.
  • 1 Introduction: Parton-density sets can be read from an external library instead of being evolved from the input scale.
  • 1 Introduction: QCDNUM provides convolution routines for hadron-hadron cross-sections, parton luminosities, and deep inelastic structure functions in zero-mass or generalised mass schemes.The zmstf and hqstf packages calculate zero-mass structure functions and fixed-flavour heavy-quark contributions.

2 QCD Evolution

QCDNUM evolves the strong coupling and parton densities across perturbative orders and flavour schemes by exploiting singlet/non-singlet structure. Its formalism includes threshold matching, perturbative splitting functions, and distinct evolution behaviour for non-singlet combinations.

  • Evolution of the Strong Coupling Constant: QCDNUM evolves the strong coupling αs and parton densities with flavour-dependent beta functions and splitting functions across perturbative orders.The number of active flavours changes at heavy-quark thresholds in the variable flavour number scheme.
  • Decomposition into Singlet and Non-singlets: The DGLAP system separates into a singlet equation coupled to the gluon and independently evolving non-singlet equations.This decomposition simplifies the coupled evolution of the original parton densities.
  • DGLAP Evolution Equations: QCDNUM uses perturbative expansions of splitting functions in αs, with the un-polarised kernels known through NNLO.The formalism separates the x dependence from the scale dependence in these expansions.
  • DGLAP Evolution Equations: At NNLO, all three non-singlet combinations evolve differently, whereas at LO they obey the same evolution equations.At NLO, q−ij and qv evolve identically but differently from q+ij.
  • Flavour Number Schemes: In the variable flavour number scheme, heavy-flavour distributions are generated at thresholds, with NNLO starting discontinuities exceeding 10% for charm at low x but remaining below 3% for bottom.The starting distributions are negative below x ≈ 10^-2.

3 Numerical Method

QCDNUM represents parton densities with polynomial B-splines on logarithmic grids, converting DGLAP evolution and convolutions into efficient algebraic operations. Quadratic interpolation improves the practical method, but grid choice and backward evolution impose accuracy and stability constraints.

  • Grid-based evolution: QCDNUM solves DGLAP equations on discrete x and µ2 grids, with convolution integrals evaluated as precomputed weighted sums.The operation count is quadratic in the x-grid size n and linear in the µ2-grid size m.
  • Spline representation: Parton densities are interpolated with linear or quadratic splines on an equidistant logarithmic x-grid, while µ2 interpolation is always quadratic.This transforms the evolution equations into a triangular set of linear equations for spline coefficients.
  • Spline representation: Increasing spline order adds only one free parameter rather than n because continuity constraints couple the piecewise polynomials.QCDNUM represents the resulting spline as a linear combination of B-splines.
  • Spline representation: Spline functions use B-spline bases whose order, grid distribution, and continuity constraints determine the interpolation basis.For qcdnum, the grid in y is equidistant, enabling translation invariance of the remaining basis functions.
  • Stability limits: Cubic and higher-order splines were abandoned because the translation-invariant cubic basis is numerically unstable and its corrective endpoint constraint does not fit the existing evolution algorithm.The proposed cubic basis imposes h(y0) = h′(y0) = h′′(y0) = 0 before the constraint is modified at the opposite endpoint.
  • Convolution integrals: Quadratic-spline convolutions require no more than O(5n) operations because only the first five convolution weights are non-zero.Limited B-spline support also restricts linear interpolation to three non-zero weight terms.
  • Stability limits: Backward evolution works well with linear splines but quadratic splines can produce severe oscillations.For linear interpolation, one iteration substantially improves large-x accuracy, while qcdnum’s default is one iteration, perhaps two.

4 The QCDNUM Program

QCDNUM provides a grid-based program for initializing, evolving, and accessing parton densities across fixed- and variable-flavour schemes. Its validation shows that multiple grids and quadratic splines substantially improve accuracy, while quadratic evolution also reduces computational cost.

  • Program setup: QCDNUM initializes x-µ2 grids, calculates convolution weights, sets perturbative and flavour-threshold parameters, and evolves input parton densities.Evolution parameters invalidate stored densities so subsequent results remain consistent with the selected αs, perturbative order, and scheme.
  • Parton-density evolution: The program evolves 2n_f + 1 momentum densities and internally decomposes input quark densities into singlet and non-singlet components.In the variable-flavour scheme, the input scale and heavy-flavour thresholds must satisfy the grid constraints described in the program documentation.
  • Input densities: Input quark flavours may be any linearly independent set, but sum rules must be imposed in the input parametrisation because QCDNUM lacks information below the lowest x-grid point.Once imposed at the input scale, the evolution conserves the sum rules.
  • Accessing results: The fvalxq and asfunc routines return interpolated flavour-basis densities and evolved αs values at requested x and µ2 scales.QCDNUM automatically converts its internal singlet/non-singlet representation to gluon, quark, or antiquark momentum densities.

5 Subroutine Calls

QCDNUM’s subroutine interface initializes the program, configures numerical parameters and grids, and provides access to grid and PDF-related operations. Grid choices control calculation speed, accuracy, storage, and spline bases.

  • Initialization and parameters: QCINIT initializes qcdnum and defines the output stream before any other operation.
  • Initialization and parameters: SETINT|GETINT and SETVAL|GETVAL set or retrieve integer and floating-point parameters after initialization.
  • Initialization and parameters: The iter parameter controls backward-evolution iterations and is relevant only when quadratic splines are used.
  • Grid definition: Grid definitions determine calculation speed, accuracy, internal storage partitioning, and the B-spline bases.Changing either the x or µ2 grid invalidates the weight tables and PDF store.
  • Grid definition: GQMAKE generates a logarithmic µ2 grid from supplied endpoints, optional intermediate points, and relative regional densities.The µ2 grid need not be equidistant, and higher density can be assigned where αs changes rapidly.

6 Convolution Engine

The convolution engine computes convolutions for luminosities, scattering cross-sections, and structure functions using precomputed kernel weights and interpolated grid-point results. It supports zero-mass and generalized-mass calculations, with mass-dependent tables requiring additional scale dependence.

  • Scope and workflow: The engine calculates deep-inelastic structure functions, hadron-hadron cross-sections, and parton luminosities.The zmstf and hqstf packages use it for zero-mass and fixed-flavour-number heavy-quark structure-function contributions.
  • Convolution formalism: Convolution integrals combine parton densities and kernels, including single, nested, density-density, and density-density-kernel forms.Kernels must be defined consistently with convolution by a number density or transformed before input.
  • Scope and workflow: A standard calculation partitions stores, fills them with weights, evaluates a user function on grid points, and interpolates to requested x and µ2.The straightforward workflow is suitable for prototyping and debugging.
  • Scope and workflow: The fast calculation scheme provides at least an order-of-magnitude speed gain over the straightforward interpolation workflow.
  • Mass schemes: The rescaling variable χ = ax restricts the convolution range to x ≤ 1/a and reduces to χ = x in the zero-mass limit.Generalized-mass schemes use a ≥ 1 to incorporate heavy-quark kinematic constraints.
  • Mass schemes: Generalized-mass weight tables generally depend on x and µ2 because quark masses and the µ2–Q2 relation enter through χ.Such tables therefore use type-3 or type-4 storage.

A Singularities

The appendix decomposes splitting and coefficient kernels into regular, singular, product, and delta-function parts. This decomposition organizes convolution integrals and supports the treatment of singular terms.

  • Splitting functions: Time-like fragmentation evolution exchanges the relevant gluon-to-quark splitting functions, while the ‘+’ prescription defines singular terms.
  • Kernel decomposition: QCDNUM represents convolution kernels through regular, singular, product, and delta-function components, with absent terms omitted.
  • Integral formulation: The appendix formulates convolution integrals in logarithmic variables y = −ln(x), splitting them into contributions I1 through I4.The final integrals of I2 and I3 may remain in x = exp(−y) to avoid integration extending to infinity.

B Triangular Systems in the DGLAP Evolution

QCDNUM solves the discretized DGLAP systems as lower-triangular Toeplitz equations. Forward substitution handles non-singlet evolution and extends to coupled singlet-gluon systems with a modest additional cost.

  • Non-singlet systems: The non-singlet evolution matrix is lower triangular Toeplitz and is solved iteratively by forward substitution.Its entries are determined by storing the first column in a one-dimensional vector.
  • Non-singlet systems: A system of n equations requires n(n + 1)/2 operations because each iteration accumulates the relevant sums.
  • Singlet-gluon systems: The forward-substitution method extends to coupled singlet-gluon evolution through repeated 2 × 2 matrix solves.Its operation count is four times that of the non-singlet algorithm plus small overhead for the 2 × 2 solves.

C.1 General Formalism

The formalism computes singlet/gluon and non-singlet structure-function contributions through Mellin convolutions of coefficient functions with evolved parton densities. Coefficients are expanded perturbatively through NNLO, with a distinct alternative expansion for F′L and a stated NNLO electroweak limitation.

  • Singlet/gluon sector: Singlet/gluon contributions to F2 and FL are convolutions of coefficient functions with singlet and gluon densities, while xF3 has no singlet/gluon contribution.The pure non-singlet nature of xF3 removes its singlet/gluon term.
  • Non-singlet sector: Non-singlet contributions to F2, FL, and xF3 are convolutions of non-singlet coefficient functions with the corresponding non-singlet densities.The non-singlet labels are +, −, and v.
  • Assembly: Structure functions combine singlet/gluon and non-singlet parts with the appropriate electroweak coupling combinations.The notation identifies F2, FL, and F3 with F2, FL, and xF3, respectively.
  • Perturbative expansion: Coefficient functions are expanded in powers of αs, with ℓ = (0, 1, 2) denoting LO, NLO, and NNLO and as = αs/2π.LO coefficient functions are zero or trivial delta functions.
  • Alternative longitudinal structure function: For F2, FL, and xF3, coefficients are calculated through C(2,2), whereas all coefficients are computed for the alternative F′L expansion.The alternative expansion gives a non-zero FL already at LO and uses the 3-loop coefficient function at NNLO.
  • Limitation: The NNLO 3-loop coefficient function used for F′L applies only to electromagnetic current exchange, so Z0 or W± contributions are unavailable at NNLO.This limitation is stated for the 3-loop calculation underlying the alternative longitudinal structure function.

C.2 Renormalisation and Factorisation Scale Dependence

QCDNUM handles scale dependence by expanding coefficient functions and replacing coupling powers with truncated Taylor series when renormalisation and factorisation scales differ. The formalism permits varying either scale relative to the other, but not both simultaneously.

  • Scale treatment: When the renormalisation scale differs from the factorisation scale, coefficient-function powers of as are replaced by a truncated Taylor series.The coefficient functions are organized using vector notation for non-singlet and singlet/gluon sectors.
  • Standard expansion: The expansion truncates at as in NLO and a2_s in NNLO for the standard coefficient-function expansion.This applies when the expansion in C.3 is used.
  • F′L expansion: For F′L, truncation occurs at as in LO, a2_s in NLO, and a3_s in NNLO, matching the splitting-function treatment.The alternative expansion therefore retains one additional power of as at each perturbative order.
  • Scale constraint: The standard coefficient-function formalism applies when the factorisation and renormalisation scales are equal.The scale-dependent expressions are constructed around this equality condition.
  • Scale constraint: Only one of the renormalisation scale or Q2 can be varied with respect to the other; both cannot be varied simultaneously.This restriction defines the available scale-uncertainty variations.

C.3 The zmstf Package

The zmstf add-on computes un-polarised deep-inelastic structure functions F2, FL, and xF3 by convolving parton densities with zero-mass coefficient functions. It provides routines for weight-table preparation, scale definition, density-set selection, and structure-function evaluation.

  • Package purpose: The zmstf package calculates F2, FL, and xF3 in un-polarised deep inelastic scattering using QCDNUM’s convolution engine.Its inputs are parton densities and zero-mass coefficient functions.
  • Weight tables: ZMFILLW fills weight tables for all flavour numbers 3 ≤ nf ≤ 6 and perturbative orders LO, NLO, and NNLO.It must be called after defining the x-µ2 grid and before the first structure-function calculation.
  • Weight-table storage: Weight tables can be dumped to disk and read back, but unformatted weight files are not exchangeable across machines.Read errors, version mismatches, and incompatible grid definitions are reported through ierr.
  • Scale definition: ZMDEFQ2 defines the relation between Q2 and the factorisation scale, with default parameters a = 1 and b = 0.The allowed ranges are 0.1 ≤ a ≤ 10 and −100 ≤ b ≤ 100, and Q2 variation requires equal renormalisation and factorisation scales in QCDNUM.
  • Density selection: ZSWITCH selects QCDNUM-evolved, custom-evolution, or external PDF sets, while polarised PDFs and fragmentation functions are invalid choices.The default is the un-polarised QCDNUM-evolved density set.
  • Structure-function evaluation: ZMSTFUN evaluates a structure function for a specified linear combination of parton densities at lists of x and Q2 points.Passing all interpolation points in one call is recommended because the internal loop is optimized for speed.

D Heavy Quark Structure Functions

The hqstf add-on supplies heavy-quark contributions to F2 and FL in NLO charged-lepton–proton scattering through photon–gluon fusion and quark subprocesses. Its calculation uses mass-dependent convolutions and imposes heavy-quark production kinematics and scale-setting restrictions.

  • Heavy-flavour setup: The heavy-quark calculation treats the heavy flavour as produced in the hard scattering rather than as a constituent of the incoming proton.Quarks below the heavy-quark pole-mass threshold are treated as massless light flavours.
  • Perturbative contributions: A heavy-flavour contribution to F2 or FL combines gluon, singlet, and charge-weighted light-quark distributions with heavy-quark coefficient functions.The first term is the LO photon–gluon fusion process γ∗g → h̄h; NLO terms include γ∗g → h̄hg and γ∗q → h̄hq.
  • Convolution and scales: The heavy-quark convolution uses the number-density form with zf(z, µ2), coefficient argument x/z, and mass-dependent variables.The factorisation and renormalisation scales are usually set to µ2 = Q2 or µ2 = Q2 + 4m2_h.
  • Production kinematics: Heavy-quark production is restricted by W2 = M2 + Q2(1 − x)/x ≥ M2 + 4m2_h, which sets the lower integration limit.The kinematic threshold ensures sufficient centre-of-mass energy to produce the heavy-quark pair.
  • Scale dependence: Only one of µ2_R or Q2 may be varied relative to µ2_F at a time, because their dependence cannot be factorised in the heavy-quark coefficient functions.Each relation between Q2 and µ2 therefore requires its own weight tables.
  • QCDNUM integration: The hqstf add-on makes the heavy-quark calculation available in qcdnum17, as in qcdnum16.The convolution is adapted to QCDNUM through a modified kernel because its integral differs from the general convolution form.

D.1 The hqstf Package

The hqstf package computes heavy-flavour contributions to F2 and FL at NLO from FFNS-evolved parton densities, using weight tables and routines that enforce flavour-scheme checks.

  • Package purpose: The hqstf package calculates heavy-flavour contributions to F2 or FL up to NLO from parton densities evolved in the FFNS with nf light flavours.Its routines accept evolved densities and produce the requested structure-function values.
  • Initialization: HQFILLW fills the weight tables and must be called before other hqstf operations.The routine accepts a structure-function selection and associated mass and coupling parameters.
  • Coefficient functions: The LO and first two NLO terms couple the virtual photon to the heavy quark, whereas the last NLO term uses a charge-weighted light-quark distribution sum.Some coefficient functions are supplied as interpolation tables because they are too complex for analytical form.
  • Structure-function calculation: HQSTFUN calculates FL or F2 for selected charm, bottom, or top contributions from input x and Q2 arrays.It returns an output array of structure functions and can check whether inputs lie within grid boundaries.
  • Consistency checks: For icbt = (1,2,3) = (c,b,t), the routine checks that densities were evolved in the FFNS with nf = (3,4,5) flavours.A preceding minus sign disables the flavour-number check while retaining the FFNS check.
  • Example application: In charged lepton–proton scattering, zmstf and hqstfun can calculate the d,u,s, charm, and total F2 contributions while neglecting bottom and top.The example uses densities evolved with nf = 3 flavours and combines the resulting arrays.

Index

The index maps QCDNUM concepts, routines, schemes, and implementation topics to their locations, including convolution, evolution, structure functions, and flavour schemes.

  • Numerical methods: Implementation topics include B-splines, spline interpolation and oscillation, interpolation meshes, Mellin convolution, and Gauss quadrature.The index also lists beta functions, the scale parameter Λ, and the number of active flavours nf.
  • Core capabilities: The index covers QCDNUM convolution integrals, convolution weights, and custom evolution.These entries point to the program’s convolution engine and alternative evolution capabilities.
  • Evolution: It identifies DGLAP evolution equations, alpha_s and PDF discontinuities, heavy-quark PDF evolution, and forward or backward evolution.The entries also reference evolution of alpha_s and the e± basis PDFs.
  • Flavour schemes and scales: The index distinguishes fixed-flavour and variable-flavour number schemes and lists flavour-threshold topics.It also references factorisation- and renormalisation-scale dependence.
  • Program documentation: The index records general program topics such as execution speed, example jobs, subroutine lists, PDFs, sum-rule integrals, and the QCDNUM web site.It also includes null values, qcdnum.inc parameters, and program hang-ups.
  • Structure functions: Structure-function entries include zero-mass structure functions, heavy-flavour contributions, and Wilson coefficients.The index separately references the FL structure function and rescaling variable χ.
Loading 1005.1481v3…