Polarized NNLO anomalous dimensions

benchmarks/DSSV_bench.py

@@ -0,0 +1,56 @@
+"""Benchmark DSSV pdf family.
+
+Note that the PDF set is private, but can be obtained from the
+authors upon request.
+"""
+import numpy as np
+from banana import register
+
+from eko import interpolation
+from ekomark.benchmark.runner import Runner
+
+register(__file__)
+
+
+class BenchmarkDSSV(Runner):
+    """Benchmark DSSV pdf family."""
+
+    external = "LHAPDF"
+
+    # Rotate to evolution basis
+    rotate_to_evolution_basis = True
+
+    def skip_pdfs(self, _theory):
+        return [-6, 6, 5, -5, 4, -4, 22, "ph", "T35", "V35", "T24", "V24", "T15", "V15"]
+
+    def benchmark(self, Q0=1.65, Q2grid=(100,)):
+        theory_card = {
+            "Qref": 1.0,
+            "mc": 1.4,
+            "mb": 4.75,
+            "mt": 175,
+            "kcThr": 1.0,
+            "kbThr": 1.0,
+            "ktThr": 1.0,
+            "alphas": 0.49127999999999999,
+            "FNS": "ZM-VFNS",
+            "ModEv": "EXA",
+            "Q0": Q0,
+            "PTO": 1,
+            "MaxNfPdf": 3,
+            "MaxNfAs": 3,
+        }
+
+        operator_card = {
+            "Q2grid": list(Q2grid),
+            "polarized": [True],
+            "interpolation_xgrid": interpolation.lambertgrid(50, 1e-5),
+        }
+        self.run([theory_card], [operator_card], ["DSSV_REP_LHAPDF6"])
+
+
+if __name__ == "__main__":
+    low2 = 5**2
+    high2 = 30**2
+    dssv = BenchmarkDSSV()
+    dssv.benchmark(Q0=np.sqrt(low2), Q2grid=[high2])

benchmarks/NNPDF_bench.py

@@ -4,6 +4,7 @@
 import numpy as np
 from banana import register
 
+from eko import interpolation
 from ekomark.benchmark.runner import Runner
 
 register(__file__)
@@ -90,14 +91,45 @@ def benchmark_nnlo(self, Q0=1.65, Q2grid=(100,)):
         self.run([theory_card], [operator_card], ["NNPDF40_nnlo_as_01180"])
 
 
+class BenchmarkNNPDFpol11(BenchmarkNNPDF):
+    """Benchmark NNPDFpol11"""
+
+    def benchmark(self, Q0=1.65, Q2grid=(100,)):
+        theory_card = {
+            "Qref": 91.2,
+            "mc": 1.41421,
+            "mb": 4.75,
+            "mt": 175,
+            "kcThr": 1.0,
+            "kbThr": 1.0,
+            "ktThr": 1.0,
+            "alphas": 0.119002,
+            "alphaqed": 0.007496,
+            "FNS": "ZM-VFNS",
+            "ModEv": "TRN",
+            "Q0": Q0,
+            "PTO": 1,
+        }
+
+        operator_card = {
+            **base_operator,
+            "Q2grid": list(Q2grid),
+            "polarized": [True],
+            "interpolation_xgrid": interpolation.lambertgrid(50, 1e-5),
+        }
+        self.run([theory_card], [operator_card], ["NNPDFpol11_100"])
+
+
 if __name__ == "__main__":
-    low2 = 4**2
+    low2 = 5**2
     high2 = 30**2
     # nn31 = BenchmarkNNPDF31()
     # # test forward
     # nn31.benchmark_nlo(Q0=np.sqrt(low2), Q2grid=[10])
     # # test backward
     # #nn31.benchmark_nlo(Q0=np.sqrt(high2), Q2grid=[low2])
-    nn40 = BenchmarkNNPDF40()
-    # nn40.benchmark_nnlo(Q2grid=[100])
-    nn40.benchmark_nnlo(Q0=np.sqrt(high2), Q2grid=[low2])
+    # nn40 = BenchmarkNNPDF40()
+    # # nn40.benchmark_nnlo(Q2grid=[100])
+    # nn40.benchmark_nnlo(Q0=np.sqrt(high2), Q2grid=[low2])
+    nnpol = BenchmarkNNPDFpol11()
+    nnpol.benchmark(Q0=np.sqrt(low2), Q2grid=[high2])

benchmarks/apfel_bench.py

@@ -158,6 +158,16 @@ def benchmark_plain(self, pto):
             cartesian_product(th), operators.build(operators.apfel_config), ["ToyLH"]
         )
 
+    def benchmark_plain_pol(self, pto):
+        """Plain configuration"""
+
+        th = self.ffns_theory.copy()
+        th.update({"PTO": [pto]})
+        th["ModEv"] = ["EXA"]  # TODO for the time one is sufficient
+        op = operators.apfel_config.copy()
+        op["polarized"] = [True]
+        self.run(cartesian_product(th), operators.build(op), ["ToyLH"])
+
     def benchmark_sv(self, pto, svmode):
         """Scale Variation"""
 
@@ -189,10 +199,11 @@ def benchmark_sv(self, pto, svmode):
 
 if __name__ == "__main__":
 
-    obj = BenchmarkVFNS()
-    # obj = BenchmarkFFNS()
+    # obj = BenchmarkVFNS()
+    obj = BenchmarkFFNS()
 
+    obj.benchmark_plain_pol(2)
     # obj.benchmark_plain(2)
-    obj.benchmark_sv(2, "exponentiated")
+    # obj.benchmark_sv(2, "exponentiated")
     # obj.benchmark_kthr(2)
     # obj.benchmark_msbar(2)

benchmarks/lha_paper_bench.py

@@ -3,9 +3,11 @@
 """
 import numpy as np
 from banana import register
+from banana.data import cartesian_product
 
 from eko.interpolation import lambertgrid
 from ekomark.benchmark.runner import Runner
+from ekomark.data import operators
 
 register(__file__)
 
@@ -219,12 +221,34 @@ def benchmark_sv(self, pto):
         self.run_lha([low, high])
 
 
+class BenchmarkFFNS_polarized(BenchmarkFFNS):
+    def run_lha(self, theory_updates):
+        """Enforce operator grid and PDF.
+
+        Parameters
+        ----------
+        theory_updates : list(dict)
+            theory updates
+        """
+        self.run(
+            theory_updates,
+            [
+                {
+                    "Q2grid": [1e4],
+                    "ev_op_iterations": 10,
+                    "interpolation_xgrid": lambertgrid(60).tolist(),
+                    "polarized": True,
+                }
+            ],
+            ["ToyLH_polarized"],
+        )
+
+
 if __name__ == "__main__":
     # Benchmark to LHA
-    obj = BenchmarkVFNS()
+    obj = BenchmarkFFNS_polarized()
     # obj = BenchmarkFFNS()
-
-    # obj.benchmark_plain(0)
+    # obj.benchmark_plain(1)
     obj.benchmark_sv(1)
 
     # # VFNS benchmarks with LHA settings

benchmarks/pegasus_bench.py

@@ -61,12 +61,12 @@ class BenchmarkVFNS(PegasusBenchmark):
         "nfref": 3,
         "nf0": 3,
     }
-    zm_theory = tolist(zm_theory)
+    vfns_theory = tolist(zm_theory)
 
     def benchmark_plain(self, pto):
         """Plain configuration"""
 
-        th = self.zm_theory.copy()
+        th = self.vfns_theory.copy()
         th.update(
             {
                 "PTO": [pto],
@@ -79,7 +79,7 @@ def benchmark_plain(self, pto):
     def benchmark_sv(self, pto, svmode):
         """Scale Variation"""
 
-        th = self.zm_theory.copy()
+        th = self.vfns_theory.copy()
         th.update(
             {
                 "PTO": [pto],
@@ -125,6 +125,15 @@ def benchmark_plain(self, pto):
             cartesian_product(th), operators.build(operators.pegasus_config), ["ToyLH"]
         )
 
+    def benchmark_plain_pol(self, pto):
+        """Plain polarized configuration"""
+
+        th = self.ffns_theory.copy()
+        th.update({"PTO": [pto]})
+        op = operators.pegasus_config.copy()
+        op["polarized"] = [True]
+        self.run(cartesian_product(th), operators.build(op), ["ToyLH_polarized"])
+
     def benchmark_sv(self, pto, svmode):
         """Scale Variation"""
 
@@ -143,9 +152,10 @@ def benchmark_sv(self, pto, svmode):
 
 if __name__ == "__main__":
 
-    obj = BenchmarkVFNS()
-    # obj = BenchmarkFFNS()
+    # obj = BenchmarkVFNS()
+    obj = BenchmarkFFNS()
+    obj.benchmark_plain_pol(1)
     # obj.benchmark_plain(1)
 
-    obj.benchmark_sv(2, "exponentiated")
+    # obj.benchmark_sv(2, "exponentiated")
     # vfns.benchmark_sv()

doc/source/development/Benchmarks.rst

@@ -45,6 +45,13 @@ List of bugs in :cite:`Dittmar:2005ed`
 - in table 15, part 1: :math:`xd_v(x=10^{-4}, \mu_F^2 = 10^4~\mathrm{GeV}^2) = 1.0699\cdot 10^{-4}` (wrong exponent) and
   :math:`xg(x=10^{-4}, \mu_F^2 = 10^4~\mathrm{GeV}^2) = 9.9694\cdot 10^{2}` (wrong exponent)
 
+- in table 16: the distribution :math:`L_{m}` is not vanishing despite the benchmark table is produced using the
+  toy PDF used, defined in :eqref:`4.57`, has :math:`\bar{u} = \bar{d}`. This means that :math:`T_3`
+  is not correct. The whole colum of :math:`L_{m}` should be set to 0, sine at |LO| it is not possible
+  to generate a difference between :math:`\bar{u}` and :math:`\bar{d}`
+
+- in table 17: same as before :math:`L_{m}` is bugged.
+
 LHAPDF
 ------
 

doc/source/refs.bib

@@ -776,3 +776,55 @@ @article{Albino:2000cp
     pages = "93--102",
     year = "2001"
 }
+
+@article{Gluck:1995yr,
+    author = "Gluck, M. and Reya, E. and Stratmann, M. and Vogelsang, W.",
+    title = "{Next-to-leading order radiative parton model analysis of polarized deep inelastic lepton - nucleon scattering}",
+    eprint = "hep-ph/9508347",
+    archivePrefix = "arXiv",
+    reportNumber = "DO-TH-95-13, RAL-TR-95-042",
+    doi = "10.1103/PhysRevD.53.4775",
+    journal = "Phys. Rev. D",
+    volume = "53",
+    pages = "4775--4786",
+    year = "1996"
+}
+
+@article{Floratos:1981hs,
+    author = "Floratos, E. G. and Kounnas, C. and Lacaze, R.",
+    title = "{Higher Order QCD Effects in Inclusive Annihilation and Deep Inelastic Scattering}",
+    reportNumber = "LPTENS-81-3",
+    doi = "10.1016/0550-3213(81)90434-X",
+    journal = "Nucl. Phys. B",
+    volume = "192",
+    pages = "417--462",
+    year = "1981"
+}
+
+@article{Moch:2015usa,
+    author = "Moch, S. and Vermaseren, J. A. M. and Vogt, A.",
+    title = "{On \ensuremath{\gamma}5 in higher-order QCD calculations and the NNLO evolution of the polarized valence distribution}",
+    eprint = "1506.04517",
+    archivePrefix = "arXiv",
+    primaryClass = "hep-ph",
+    reportNumber = "DESY-15-061, NIKHEF-15-018, LTH-1042",
+    doi = "10.1016/j.physletb.2015.07.027",
+    journal = "Phys. Lett. B",
+    volume = "748",
+    pages = "432--438",
+    year = "2015"
+}
+
+@article{Moch:2014sna,
+    author = "Moch, S. and Vermaseren, J. A. M. and Vogt, A.",
+    title = "{The Three-Loop Splitting Functions in QCD: The Helicity-Dependent Case}",
+    eprint = "1409.5131",
+    archivePrefix = "arXiv",
+    primaryClass = "hep-ph",
+    reportNumber = "DESY-14-157, NIKHEF-14-033, LTH-1023",
+    doi = "10.1016/j.nuclphysb.2014.10.016",
+    journal = "Nucl. Phys. B",
+    volume = "889",
+    pages = "351--400",
+    year = "2014"
+}

doc/source/theory/pQCD.rst

@@ -64,6 +64,45 @@ coupling :math:`a_s(\mu^2)` and are given by :cite:`Moch:2004pa,Vogt:2004mw`
 
 Note the additional minus in the definition of :math:`\gamma`.
 
+Polarized Splitting Functions
+-----------------------------
+
+Polarized Altarelli-Parisi splitting kernels are implemented up to |NNLO| and expanded in powers of the strong coupling as in the previous section.
+They are used to evolve longitudinally polarized parton distribution functions.
+Unlike in the unpolarized case, where the probability of the splitting describes the momentum of parent and daughter partons with averaged spins,
+the polarized splitting functions describe the parent and daughter momentums along with their spins
+and thus take into account positive or negative helicities.
+Throughout, the anomalous dimensions are defined as above and are represented with :math:`\gamma` and not :math:`\Delta\gamma` just like in the unpolarized case.
+
+The |LO| and |NLO| kernels are given in :cite:`Gluck:1995yr` and the |NNLO| in :cite:`Moch:2014sna` and :cite:`Moch:2015usa`.
+
+At |LO|, the non-singlet is the same in both the polarized and unpolarized case.
+Due to helicity conservation, the first moment of the anomalous dimension is :math:`\gamma^{(0)}_{qq} (N=1) = \gamma^{(0)}_{qg} (N=1) = 0`.
+
+At |NLO|, the singlet entry of the quark-quark anomalous dimension can be decomposed into the pure singlet
+(consisting of the flavour independent quark-quark and quark-antiquark anomalous dimensions) and the plus flavour asymmetry non-singlet:
+
+.. math ::
+    \gamma^{(1)}_{qq} =\gamma^{(1)}_{ps} + \gamma^{(1)}_{ns,+}
+
+The non-singlet sector in the polarized case swaps the plus and minus non-singlet relative to the unpolarized case.
+This is because the polarized non-singlet splitting functions are defined as the difference between the probability of the polarized parton splitting into daughter partons of the same flavour
+and daughters splitting into a different flavours and opposite helicity. The first moments of the anomalous dimensions are:
+
+.. math ::
+    \gamma^{(1)}_{ns,+} (N=1) &= 0 \\
+    \gamma^{(1)}_{qq} (N=1) &= 24 C_F T_R \\
+    \gamma^{(1)}_{qg} (N=1) &= 0  \\
+
+At |NNLO| the non-singlet is further decomposed into the helicity difference quark-antiquark anomalous dimension called the valence polarized non-singlet and defined as:
+
+.. math ::
+    \gamma^{(2)}_{ns,v} =\gamma^{(2)}_{ns,-} + \gamma^{(2)}_{ns,s}
+
+where :math:`\gamma^{(2)}_{ns,-}` is the minus flavour asymmetry non-singlet and :math:`\gamma^{(2)}_{ns,s}` the sea-like polarized non-singlet.
+The singlet entry :math:`\gamma^{(2)}_{qq}` is defined as above in the |NLO| case.
+
+
 Unified Splitting Functions
 ---------------------------
 
@@ -83,6 +122,7 @@ where :math:`a = \alpha/(4\pi)`.
 The expression of the pure |QED| and of the mixed |QED| :math:`\otimes` |QCD| splitting kernels are given in
 :cite:`deFlorian:2015ujt,deFlorian:2016gvk`
 
+
 Order specification
 -------------------