Merge pull request #274 from NNPDF/fix-apfelbench

Fix apfel benchmarks

benchmarks/CT14_bench.py

@@ -5,8 +5,6 @@
   whereas the lo family uses LO splitting functions with NLO alphas evolution
 """
 
-from math import nan
-
 from banana import register
 
 from ekomark.benchmark.runner import Runner
@@ -16,7 +14,7 @@
 
 base_theory = {
     "Qref": 91.1876,
-    "QrefQED": nan,
+    "Qedref": 0.0,
     "mc": 1.3,
     "mb": 4.75,
     "mt": 172,

benchmarks/CT18_bench.py

@@ -12,7 +12,7 @@
 
 base_theory = {
     "Qref": 91.1870,
-    "QrefQED": 91.1870,
+    "Qedref": 91.1870,
     "mc": 1.3,
     "mb": 4.75,
     "mt": 172.0,

benchmarks/HERA20_bench.py

@@ -2,7 +2,6 @@
 Benchmark HERAPDF2.0 pdf family
 
 """
-from math import nan
 
 from banana import register
 
@@ -14,7 +13,7 @@
 
 base_theory = {
     "Qref": 91.1876,
-    "QrefQED": nan,
+    "Qedref": 0.0,
     "mc": 1.43,
     "mb": 4.5,
     "mt": 173.0,

benchmarks/NNPDF_bench.py

@@ -72,7 +72,7 @@ def benchmark_nnlo(self, Q0=1.65, mugrid=(10,)):
             "QED": 2,
             "Q0": Q0,
         }
-        theory_card.update({"ModEv": "iterate-exact", "FNS": "VFNS", "QrefQED": 91.2})
+        theory_card.update({"ModEv": "iterate-exact", "FNS": "VFNS", "Qedref": 91.2})
 
         self.skip_pdfs = lambda _theory: [
             -6,

benchmarks/apfel_bench.py

@@ -202,7 +202,7 @@ class BenchmarkFFNS_qed(ApfelBenchmark):
 
     ffns_theory = {
         "Qref": 91.1870,
-        "QrefQED": 91.1870,
+        "Qedref": 91.1870,
         "mc": 1.3,
         "mb": 4.75,
         "mt": 172.0,
@@ -281,7 +281,7 @@ class BenchmarkVFNS_qed(ApfelBenchmark):
 
     vfns_theory = {
         "Qref": 91.1870,
-        "QrefQED": 91.1870,
+        "Qedref": 91.1870,
         "mc": 1.3,
         "mb": 4.75,
         "mt": 172.0,

benchmarks/lha_paper_bench.py

@@ -2,7 +2,6 @@
 Benchmark to :cite:`Giele:2002hx` (LO + NLO) and :cite:`Dittmar:2005ed` (NNLO).
 """
 import os
-from math import nan
 
 import numpy as np
 import pytest
@@ -29,7 +28,7 @@
     "alphas": 0.35,  # Eq. (4.55) :cite:`Dittmar:2005ed`
     "alphaqed": 0.007496,
     "QED": 0,
-    "QrefQED": nan,
+    "Qedref": 0.0,
 }
 """Global theory settings"""
 

benchmarks/pegasus_bench.py

@@ -1,7 +1,6 @@
 """
 Benchmark to Pegasus :cite:`Vogt:2004ns`
 """
-from math import nan
 
 import numpy as np
 from banana import register
@@ -56,7 +55,7 @@ class BenchmarkVFNS(PegasusBenchmark):
         "kbThr": 1.0,
         "ktThr": 1.0,
         "Qref": np.sqrt(2.0),
-        "QrefQED": nan,
+        "Qedref": 0.0,
         "alphas": 0.35,
         "alphaqed": 0.007496,
         "QED": 0,
@@ -110,7 +109,7 @@ class BenchmarkFFNS(PegasusBenchmark):
         "kbThr": np.inf,
         "ktThr": np.inf,
         "Qref": np.sqrt(2.0),
-        "QrefQED": nan,
+        "Qedref": 0.0,
         "alphas": 0.35,
         "alphaqed": 0.007496,
         "Q0": np.sqrt(2.0),

benchmarks/sandbox.py

@@ -1,5 +1,4 @@
 # pylint: skip-file
-from math import nan
 
 import numpy as np
 from banana import register
@@ -27,7 +26,7 @@
 }
 nnpdf_base_theory = {
     "Qref": 91.2,
-    "QrefQED": nan,
+    "Qedref": 0.0,
     "mc": 1.51,
     "mb": 4.92,
     "mt": 172.5,

doc/source/code/IO.rst

@@ -69,7 +69,7 @@ and environment. The benchmark settings are available at :mod:`banana.data.theor
   * - ``alphaqed``
     - :py:obj:`float`
     - Reference value of the electromagnetic coupling :math:`\alpha_{em}`.
-  * - ``QrefQED``
+  * - ``Qedref``
     - :py:obj:`float`
     - Reference scale at which the ``alphaqed`` value is given (in GeV).
   * - ``HQ``

pyproject.toml

@@ -44,7 +44,7 @@ PyYAML = "^6.0"
 lz4 = "^4.0.2"
 numba = "^0.57.0"
 # ekomark
-banana-hep = { version = "^0.6.8", optional = true }
+banana-hep = { version = "^0.6.9", optional = true }
 sqlalchemy = { version = "^1.4.21", optional = true }
 pandas = { version = "^1.3.0", optional = true }
 matplotlib = { version = "^3.5.1", optional = true }

src/eko/io/runcards.py

@@ -173,11 +173,16 @@ def new_theory(self):
 
         new["order"] = [old["PTO"] + 1, old["QED"]]
         alphaem = self.fallback(old.get("alphaqed"), old.get("alphaem"), default=0.0)
+        if "Qedref" in old:
+            em_running = bool(np.isclose(old["Qedref"], old["Qref"]))
+        else:
+            em_running = False
         ms = self.heavies("m%s", old)
         ks = self.heavies("k%sThr", old)
         new["couplings"] = dict(
             alphas=old["alphas"],
             alphaem=alphaem,
+            em_running=em_running,
             scale=old["Qref"],
             num_flavs_ref=old["nfref"],
             max_num_flavs=old["MaxNfAs"],
@@ -220,7 +225,7 @@ def new_operator(self):
             mugrid = old["mugrid"]
         else:
             mu2grid = old["Q2grid"] if "Q2grid" in old else old["mu2grid"]
-            mugrid = np.sqrt(mu2grid)
+            mugrid = np.sqrt(mu2grid).tolist()
         new["mugrid"] = flavored_mugrid(
             mugrid,
             list(self.heavies("m%s", old_th)),

src/ekomark/benchmark/external/pegasus_utils.py

@@ -82,12 +82,12 @@ def compute_pegasus_data(theory, operators, skip_pdfs, rotate_to_evolution_basis
 
     # run pegaus
     out_tabs = {}
-    for mu2 in np.array(operators["mugrid"]) ** 2:
+    for mu in operators["mugrid"]:
         tab = {}
         for x in target_xgrid:
             # last two numbers are the min and max pid to calculate,
             # keep everthing for simplicity.
-            xf, _ = pegasus.xparton(x, mu2, -6, 6)
+            xf, _ = pegasus.xparton(x, mu**2, -6, 6)
             temp = dict(zip(labels, xf))
             for pid in labels:
                 if pid in skip_pdfs:
@@ -114,7 +114,7 @@ def compute_pegasus_data(theory, operators, skip_pdfs, rotate_to_evolution_basis
             evol_pdf = rotate_flavor_to_evolution @ pdfs
             tab = dict(zip(evol_basis, evol_pdf))
 
-        out_tabs[mu2] = tab
+        out_tabs[mu**2] = tab
 
     ref = {"target_xgrid": target_xgrid, "values": out_tabs}
 

src/ekomark/benchmark/runner.py

@@ -220,7 +220,6 @@ def log(self, theory, _, pdf, me, ext):
         # return a proper log table
         log_tabs = {}
         xgrid = ext["target_xgrid"]
-        q2s = list(ext["values"].keys())
 
         # LHA NNLO VFNS needs a special treatment
         # Valence contains only u and d
@@ -246,9 +245,8 @@ def log(self, theory, _, pdf, me, ext):
                 flavor_rotation=rotate_to_evolution,
                 qed=qed,
             )
-        for q2 in q2s:
+        for q2, ref_pdfs in ext["values"].items():
             log_tab = dfdict.DFdict()
-            ref_pdfs = ext["values"][q2]
             # find the first hit, regardless of the flavor, since others can not deal with them appearing in parallel
             eps = [ep for ep in pdf_grid.keys() if np.isclose(ep[0], q2)]
             if len(eps) != 1:

tests/eko/io/test_runcards.py

@@ -47,6 +47,7 @@ def test_runcards_opcard():
     np.testing.assert_allclose(no.mu0, tc["Q0"])
     np.testing.assert_allclose(no.mu20, tc["Q0"] ** 2.0)
     assert len(no.pids) == 14
+    check_dumpable(no)
     del oc["mugrid"]
     # or mu2grid
     oc["mu2grid"] = [9.0, 30.0, 32.0]
@@ -57,6 +58,7 @@ def test_runcards_opcard():
     assert no.evolgrid[0][-1] == 4
     assert no.evolgrid[1][-1] == 5
     assert no.evolgrid[2][-1] == 5
+    check_dumpable(no)
     del oc["mu2grid"]
     # or Q2grid
     oc["Q2grid"] = [15.0, 130.0, 140.0, 1e5]
@@ -68,6 +70,15 @@ def test_runcards_opcard():
     assert no.evolgrid[1][-1] == 5
     assert no.evolgrid[2][-1] == 5
     assert no.evolgrid[3][-1] == 6
+    check_dumpable(no)
+
+
+def check_dumpable(no):
+    """Check we can write and read to yaml."""
+    so = StringIO()
+    yaml.dump(no.raw, so)
+    so.seek(0)
+    noo = yaml.safe_load(so)
 
 
 def test_runcards_ekomark():