FK tables and the theory pipeline

From PineAPPL grids through EKO to the linear/quadratic contractions fits evaluate.

Here is the number that dictates the architecture of every modern PDF fit: our grid-model run made 137 million likelihood evaluations (PPDF-8). Each one needs NNLO QCD predictions for 648 data points. Computed honestly — solve DGLAP, then do the convolution integrals point by point — a single pass costs minutes; 137 million passes cost centuries. The fix is the FK table: precompute everything that doesn't depend on the fitted parameters, once, so that "NNLO QCD" becomes a single cached tensor contraction. This chapter derives that contraction, tours the pipeline that builds it, and retells our most expensive lesson about it.

What exactly is being cached

A DIS structure function at the data's scale $Q$ is a factorized convolution (Chapter 2.2), $F_I = \sum_j C_{Ij} \otimes f_j(Q)$, and the PDF at $Q$ is itself a linear map of the PDF at the fitting scale $Q_0$ — DGLAP evolution (Chapter 2.4) is a linear operator on $f$. The fitted parameters $\theta$ enter only through $f(\cdot\,, Q_0; \theta)$. Everything else — coefficient functions at NNLO, the evolution operator, the $x$-integrals — is fixed once the theory settings are fixed. So cache it all.

— one small tensor contraction per likelihood call. For a parametrized model like MSHT20 there is not even interpolation error at fit time: the model is simply evaluated at the nodes $x_k$ (the interpolation already happened, inside the FK integrals). The grid is the standard NNPDF one: $\sim$50 nodes spanning $x = 10^{-9}$ to 1, log-dense at small $x$.

EKO, and theory IDs as frozen assumption sets

The evolution operator $E(Q_I \leftarrow Q_0)$ inside the bracket is produced by EKO (Evolution Kernel Operator), which solves DGLAP once per target scale and stores the result as a matrix over (flavour ⊗ grid) — exactly the right shape to be fused with the coefficient-function grids (PineAPPL) into FK tables. The whole bundle is versioned by a theory ID. Ours, 40000000, pins: perturbative order NNLO; $\alpha_s(M_Z) = 0.118$ with NNLO running; heavy-quark masses $m_c = 1.51$, $m_b = 4.92$, $m_t = 172.5$ GeV; the FONLL general-mass scheme for crossing thresholds; parametrization scale $Q_0 = 1.65$ GeV; scales $\mu_R = \mu_F = Q$. One integer = one complete, frozen set of theory assumptions. That is what makes our model comparisons honest: PPDF-7, 9 and 12 share the ID, so their evidence differences are pure parametrization statements — change the ID and you change the question. (It also makes results reproducible years later: the ID resolves to immutable tables, not to "whatever the code did that week".)

Two of those pinned choices deserve a second look because they surfaced elsewhere in this study. The FONLL scheme and $m_c = 1.51$ GeV are why §3.1's MSHT20 charm block is all zeros: with perturbative charm, the charm PDF is generated by evolution at the threshold — it lives inside $E(Q \leftarrow Q_0)$, not in the fitted parameters. And $Q_0 = 1.65$ GeV is the scale at which everything in §3.1 and §3.2 is defined: the Chebyshev shapes, the sum-rule integrals, the prior boxes — change $Q_0$ and every fitted number changes meaning. Finally, note what the caching buys methodologically: the same frozen tables serve the Bayesian sampler, the gradient descent of the Hessian and Monte-Carlo methods (JAX differentiates through the contraction for free), and every closure test — so our three-method comparisons (PPDF-15) differ in inference only, never in theory.

Hadronic data: the quadratic version

DIS is linear because one proton is probed. Drell–Yan, $W/Z$, jets, top — two protons, two PDFs — give the double convolution, which discretizes the same way into a rank-5 object:

Cost jumps from $14 \times 50 = 700$ multiplies per point to $(14 \times 50)^2 \approx 5 \times 10^5$ — and the table itself from kilobytes to megabytes per point. This is why our SLAC+BCDMS fits run happily on a laptop while "add LHC data" (Maria's next step) is a GPU-box conversation, and why quadratic FK is the dominant cost driver in global fits.

The flavour dimension — and the masking hazard

The index $j$ runs over the 14-row evolution basis: $\gamma$ (photon), $\Sigma$, $g$, $V$, $V_3$, $V_8$, $V_{15}$, $V_{24}$, $V_{35}$, $T_3$, $T_8$, $T_{15}$, $T_{24}$, $T_{35}$. Every model must fill all 14 rows at the grid nodes — and every row it fills with zeros is a physics assumption, silently contracted against real FK columns. The Les Houches toy (§3.1) states its assumptions in this language: $T_3 = 0$, $T_8 = 0.4\,\Sigma$, all higher $V$'s equal to $V$.

Now the lesson we paid twelve hours and 137 million evaluations to learn (PPDF-8). Our first grid model fitted only $[\Sigma, g, V]$ and let the framework's flavour_mapping zero the other eleven rows — which masks the corresponding FK columns. But $T_3 = (u + \bar u) - (d + \bar d)$ is precisely the combination that drives the proton–deuteron $F_2$ difference, and half our 648 points are deuteron data. The model was structurally incapable of fitting the dataset, no matter its 24 parameters: $\chi^2/N = 8.02$, evidence worse than the 13-parameter toy by $\Delta \log Z \approx -808$. Nothing crashed; no warning fired; the contraction happily multiplied real data sensitivity into hard zeros. Restoring the full flavour set (PPDF-9, $[\Sigma, g, V, V_3, T_3, T_8]$) dropped $\chi^2/N$ to 5.84 immediately. The moral, stated technically: the FK contraction is total — audit which rows your model populates against which columns your data actually reads.

Provenance: from HepData to the numbers the sampler sees

The data side has its own pipeline, and knowing it is part of trusting a fit. Published measurements are deposited on HepData; the NNPDF commondata layer re-encodes each dataset as versioned YAML: central values, kinematics ($x, Q^2, y$ per point), and — most importantly for §3.4 — the uncertainty breakdown, each source tagged with a treatment (ADDitive or MULTiplicative) and a correlation type (UNCORR, CORR, or a named systematic shared across points and even across datasets). Datasets carry variants encoding curatorial choices; our runcards use legacy_dw for all four sets (SLAC p/d, BCDMS p/d), the variant that carries the legacy deuteron treatment with de-weighted nuclear-correction systematics — deuteron targets are not free neutrons, and that model dependence must live somewhere in the error budget. Then kinematic cuts ($Q^2 \ge 3.49$ GeV$^2$, $W^2 \ge 12.5$ GeV$^2$ — stay where Chapter 2's theory applies) reduce the published points to our fitted 648: SLAC p 33, SLAC d 34, BCDMS p 333, BCDMS d 248. Only after cuts are the FK tables assembled, one row per surviving point. The cuts are theory protection, not data hygiene: below $W^2 \approx 12.5$ GeV$^2$ the resonance region and higher-twist $\mathcal{O}(\Lambda^2/Q^2)$ corrections invalidate leading-twist factorization — precisely the failure modes Chapter 2.2 flagged. One habit worth copying: we counted the 648 through the framework's own API rather than assuming it from the papers (PPDF-1) — cuts are code, and code is what you verify.