Efficient interpolation and evolution of parton distribution - - PowerPoint PPT Presentation

Efficient interpolation and evolution

• f parton distribution functions.

Riccardo Nagar

Deutsches Elektronen-Synchrotron (DESY)

XXVII International Workshop on Deep Inelastic Scattering 8–12 April 2019, Turin, Italy

work in collaboration with Markus Diehl and Frank Tackmann (DESY)

1 / 15

Motivation

Why do we need efficient and precise PDF interpolation and evolution? Precision predictions ◮ current interpolation precision becomes insufficient at N3LO

[see Dulat et al. 1710.03016]

◮ fast on-the-fly Mellin convolutions with “complicated” kernels required for analytical higher-order resummation ◮ PDF derivatives appear in subleading power calculations Double parton scattering ◮ cross-section formula for DPS includes double parton distributions (DPDs) Fa1a2(x1, x2, ②, µ1, µ2) → 5D grids? optimistic est. is 1 TB! ◮ DPDs needed for DPS phenomenology (e.g. W +jets) at different factorization scales (µ1, µ2) and differential in transverse separation ② ◮ already implemented solutions [Gaunt, Stirling ’11] [Elias, Golec-Biernat, Sta´

sto ’18]

◮ our goal: “fast” on-the-fly evolution and interpolation for DPDs

2 / 15

PDF interpolation methods

The usual interpolation [LHAPDF, APFEL, HOPPET, QCDNUM, ...] ◮ one or more grids in log x in interval [xmin, 1] (usually O(100) points) ◮ one or more grids in log αs(µ) or log(µ/GeV) ◮ splines or polynomial interpolation on equispaced grids

(not always on LHAPDF!)

10-6 10-5 10-4 0.001 0.010 0.100 1 10-10 10-8 10-6 10-4 10-2

◮ LHAPDF: log cubic splines, continuous 1st derivative ◮ Mathematica: log cubic splines, continuous 2nd derivative

3 / 15

PDF interpolation methods

The usual interpolation [LHAPDF, APFEL, HOPPET, QCDNUM, ...] ◮ one or more grids in log x in interval [xmin, 1] (usually O(100) points) ◮ one or more grids in log αs(µ) or log(µ/GeV) ◮ splines or polynomial interpolation on equispaced grids

(not always on LHAPDF!)

10-6 10-5 10-4 0.001 0.010 0.100 1 10-10 10-8 10-6 10-4 10-2

◮ LHAPDF: log cubic splines, continuous 1st derivative ◮ Mathematica: log cubic splines, continuous 2nd derivative

3 / 15

Chebyshev interpolation

Our interpolation ◮ one or more grids in u = log x ◮ interpolating N-th order polynomial in the Chebyshev points ◮ Chebyshev (N + 1)-grid on [−1, 1]: ˜ uk = cos kπ N → shifted grids in u ◮ barycentric formula is simple, fast and numerically stable: f (u) ≃

• j

f (uj) bj(u), with bj(u) = (−1)jcj u − uj

• i

(−1)ici u − ui ,

• barycentric basis function

(c0,N= 1

2 , ci =1)

Some advantages ◮ higher accuracy with considerably smaller number of points ◮ built-in Clenshaw-Curtis integration on full grid ◮ high-precision differentiation

4 / 15

Comparing splines vs barycentric

MMHT grid: 64 points

10-6 10-5 10-4 0.001 0.010 0.100 1 10-15 10-11 10-7 10-3

◮ LHAPDF: 64 pts ◮ Mathematica: 64 pts ◮ Chebyshev: 63 pts Relatively small amount of points can reach remarkable accuracy.

5 / 15

Comparing splines vs barycentric

HERAPDF grid: 199 points

10-6 10-5 10-4 0.001 0.010 0.100 1 10-15 10-11 10-7 10-3

◮ LHAPDF: 199 pts ◮ Mathematica: 199 pts ◮ Chebyshev: 71 pts Zero-crossings do not degrade the accuracy too much.

5 / 15

Comparing splines vs barycentric

Comparing: ◮ single grid [10−6, 1] with 64 points ◮ double grid [10−6, 0.2] and [0.2, 1] with 32 points each.

10-6 10-5 10-4 0.001 0.010 0.100 1 10-15 10-11 10-7 10-3

Generally, better accuracy with (few) subintervals.

5 / 15

Mellin moments: a measure of accuracy

◮ truncated Mellin moment of a PDF defined as Mj(f , x0) = 1

x0

dz zj−1 f (z) tends to the full Mellin moment as x0 → 0 ◮ global measure of interpolation accuracy

40 60 80 100 120 10-15 10-12 10-9 10-6 10-3 40 60 80 100 120 10-15 10-12 10-9 10-6 10-3 6 / 15

Mellin convolution

◮ Mellin convolution with a PDF (K ⊗ f )(x) = 1

x

dz z K(z) f x z

• ,

where K(z) is a kernel and f (x) the parton distribution ◮ use barycentric formula f (x) =

n fn bn(x),

(K ⊗ f )(xm) = Kmn fn, with Kmn = 1

xm

dz z K(z) bn xm z

• ◮ obtain (K ⊗ f ) at any x via interpolation

Benefits

• 1. precision: pre-compute kernel matrix Kij once at desired precision
• 2. efficiency: apply same matrix to any discretized PDF-like distribution

7 / 15

Mellin convolution accuracy

Comparing Chebyshev (63 pts.) vs LHAPDF (MMHT grid, 64 pts.)

10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 8 / 15

Discretized DGLAP equations

◮ DGLAP evolution equations d d log µ fa(x, µ) = (P(αs) ⊗ fa) (x, µ) ◮ define ˜ fa(x, µ) = x fa(x, µ) d d log µ ˜ f (x, µ) = 1

x

dz P(z, αs(µ)) ˜ f (x/z, µ) ◮ discretize splitting kernel d d log µ ˜ fm = Pmn ˜ fn with Pmn = 1

xm

dz P(z) bn(xm/z) ◮ solve the linear system of differential equations

9 / 15

Runge–Kutta methods

We solve the system of homogeneous differential equations using an explicit Runge–Kutta routine. Various Runge–Kutta implementations Accuracy in orders p of the step-size h: O(hp) ◮ classic RK4: “Fiat”, 4 function calls, 4th order ◮ Cash–Karp: “Alfa Romeo”, 6 function calls, 5th order ◮ Dormand–Prince: “Lamborghini”, 8 function calls, 6th order

• ■

■ ■ ■ ■ ■ ■ ■ ■ ■ ■ ■ ■ ■ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆ ◆

• ■

◆ 10 100 1000 104 10-16 10-14 10-12 10-10 10-8 10-6

10 / 15

Comparison with benchmark

[hep-ph/0204316, hep-ph/0511119]

We agree with all the digits shown in the benchmark tables. ◮ benchmark evolution using HOPPET (G. Salam) and PEGASUS (A. Vogt) ◮ benchmark values given with 5 significant digits (some points with 4) ◮ evolution from µ0 = √ 2 GeV to µ = 100 GeV in variable flavour number scheme from Nf = 3 to Nf = 5, with matching at µ = mc and µ = mb ◮ HOPPET result obtained with a total of 1,170 pts in x ∈ [10−8, 1] and 220 pts in µ2 ∈ [2, 106] GeV2 Absolute difference of our NNLO evolution vs. the benchmark results

x x uv x dv x L− 2x L+ x s+ x c+ x s− x g 10−7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 10−6 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 . . . . . . . . . . . . . . . . . . . . . . . . . . . 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.9 0.0 0.0 4e-13 1e-12 3e-13 1e-13 1e-14 0.0

• ur settings: N = 70 pts in x and h = 0.02

with L± = ¯ u ± ¯ d, q± = q ± ¯ q

11 / 15

PDF evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta relative discrepancy between h = 0.02 and h = 0.004 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ Evolution to µ = 1.01 mb (just after matching from nF = 3 to 5) ◮ Inaccuracy grows in valence distributions at low x

10-5 0.001 0.100 10-13 10-10 10-7 10-5 0.001 0.100 10-13 10-10 10-7 12 / 15

PDF evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta relative discrepancy between h = 0.02 and h = 0.004 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ Evolution to µ = 1.01 mb (just after matching from nF = 3 to 5) ◮ Inaccuracy grows in valence distributions at low x

10-5 0.001 0.100 10-14 10-11 10-8 0.2 0.4 0.6 0.8 1.0 10-14 10-11 10-8 12 / 15

PDF evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta relative discrepancy between h = 0.02 and h = 0.004 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ Evolution to µ = 10 TeV ◮ Numerical errors still under control; Runge-Kutta error negligible

10-5 0.001 0.100 10-13 10-10 10-7 10-5 0.001 0.100 10-13 10-10 10-7 12 / 15

PDF evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta relative discrepancy between h = 0.02 and h = 0.004 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ Evolution to µ = 10 TeV ◮ Numerical errors still under control; Runge-Kutta error negligible

10-5 0.001 0.100 10-14 10-11 10-8 0.2 0.4 0.6 0.8 1.0 10-14 10-11 10-8 12 / 15

Extension to DPDs

p p

H1 H2

◮ DPD cross section is factorized as dσDPS ∝ dσ(1)

a1b1(µ1) ⊗ dσ(2) a2b2(µ2)

⊗

• d2② Fa1a2(②, µ1, µ2) ⊗ Fb1b2(②, µ1, µ2)

◮ DGLAP equations extended to DPDs are separate evolution equations w.r.t. µ1 or µ2 [Diehl, Ostermeier, Sch¨

afer ’11]

d d log µi Fa1a2(x1, x2, y, µ1, µ2) =

• P(i)(αs) ⊗

i Fa1a2

• (x1, x2, y, µ1, µ2)

EXPERT NOTE: homogeneous equations only valid for y-dependent DPDs!

◮ defining Fa1a2(x1, x2) = 0 for x1 + x2 > 1, obtain µ1-evolution independent from x2 and viceversa ◮ same discretization as done for PDFs → matrices Pmn are the same

13 / 15

DPD evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta error est. comparing DOPRI6 h = 0.05 vs. h = 0.01 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ Representative sample of Fgg, Fu ¯

u, Fug, Fgu, Fsg at µ1 = µ2 = 1.01 mb

◮ Notice degradation for x1 + x2 close to 1

10-5 0.001 0.100 10-10 10-8 10-6 10-4 10-5 0.001 0.100 10-10 10-8 10-6 10-4 14 / 15

DPD evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta error est. comparing DOPRI6 h = 0.05 vs. h = 0.01 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ Representative sample of Fgg, Fu ¯

u, Fug, Fgu, Fsg at µ1 = µ2 = 1.01 mb

◮ Notice degradation for x1 + x2 close to 1

10-5 0.001 0.100 10-10 10-8 10-6 10-4 10-6 10-5 10-4 0.001 0.010 0.100 10-10 10-8 10-6 10-4 14 / 15

DPD evolution accuracy

Settings ◮ NNLO VFN evolution ◮ Runge–Kutta error est. comparing DOPRI6 h = 0.05 vs. h = 0.01 ◮ discretization error estimated comparing N = 70 and N = 106 grids ◮ evolution to µ1 = µ2 = 100 GeV ◮ Runge-Kutta error orders of magnitude below discretization error

10-5 0.001 0.100 10-12 10-9 10-6 10-5 0.001 0.100 10-12 10-9 10-6 14 / 15

Conclusions

This method is implemented in a C++14 library:

ChiliPDF

(Chebyshev-based Interpolation Library for PDFs) ◮ our results for PDFs reach higher numerical accuracy for a considerably smaller amount of points ◮ up to NNLO DGLAP evolution with O(α2

s) flavour matching for PDFs and

DPDs ◮ on-the-fly evolution: on Intel Core i5 @ 3.20 GHz, evolving a full PDF set from 1 GeV to 100 GeV

◮ FFN LO, h = 0.01, N = 70 → t ≈ 40 ms (+ 380 ms init.) ◮ VFN NNLO, h = 0.01, N = 70 → t ≈ 80 ms (+ 10 s init.) ǫ ≤ O(10−8) for x < 0.8 ◮ FFN NNLO full DPD set, h = 0.05, N = 70 → t ≈ 5 s (+ 11 s init.) ǫ ≤ O(10−4) for x1 + x2 < 0.8

Work in progress: ◮ include polarization ◮ full y-dependence in DPDs ◮ optimize code ...make the library publicly available

15 / 15

Backup: Mellin convolution with GSL integration

Comparing the accuracy of Chebyshev against LHAPDF for some kernels:

10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 16 / 15

Backup: Mellin convolution on HERAPDF grid

Taking the HERAPDF grid (199 points) vs Chebyshev with 71 points:

10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 10-5 10-4 0.001 0.010 0.100 1 10-14 10-11 10-8 10-5 0.01 10 17 / 15