From b9cfc9e5914287b550fb75e01443197dba3b870b Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Mon, 24 Jul 2023 18:11:01 +0200 Subject: [PATCH 01/20] Change manipulate act on Operator --- src/eko/io/manipulate.py | 212 +++++++++--------------- src/ekomark/benchmark/runner.py | 13 +- src/ekomark/plots.py | 30 ++-- tests/eko/io/test_manipulate.py | 274 ++++++++++---------------------- 4 files changed, 180 insertions(+), 349 deletions(-) diff --git a/src/eko/io/manipulate.py b/src/eko/io/manipulate.py index e9cc89538..1372e9e14 100644 --- a/src/eko/io/manipulate.py +++ b/src/eko/io/manipulate.py @@ -1,7 +1,8 @@ -"""Manipulate output generate by EKO.""" +"""Manipulate operators.""" +import copy import logging import warnings -from typing import Callable, Optional, Union +from typing import Callable, Optional import numpy as np import numpy.typing as npt @@ -9,7 +10,7 @@ from .. import basis_rotation as br from .. import interpolation from ..interpolation import XGrid -from .struct import EKO, Operator +from .struct import Operator logger = logging.getLogger(__name__) @@ -18,10 +19,8 @@ SIMGRID_ROTATION = "ij,ajbk,kl->aibl" """Simultaneous grid rotation contraction indices.""" -Basis = Union[XGrid, npt.NDArray] - -def rotation(new: Optional[Basis], old: Basis, check: Callable, compute: Callable): +def rotation(new: Optional[XGrid], old: XGrid, check: Callable, compute: Callable): """Define grid rotation. This function returns the new grid to be assigned and the rotation computed, @@ -61,81 +60,66 @@ def xgrid_compute_rotation(new: XGrid, old: XGrid, interpdeg: int, swap: bool = def xgrid_reshape( - eko: EKO, + elem: Operator, + xgrid: XGrid, + interpdeg: int, targetgrid: Optional[XGrid] = None, inputgrid: Optional[XGrid] = None, -): - """Reinterpolate operators on output and/or input grids. +) -> Operator: + """Reinterpolate the operator on output and/or input grid(s). Target corresponds to the output PDF. - - The operation is in-place. - """ - eko.assert_permissions(write=True) - # calling with no arguments is an error if targetgrid is None and inputgrid is None: raise ValueError("Nor inputgrid nor targetgrid was given") - interpdeg = eko.operator_card.configs.interpolation_polynomial_degree check = xgrid_check crot = xgrid_compute_rotation # construct matrices newtarget, targetrot = rotation( targetgrid, - eko.bases.targetgrid, + xgrid, check, lambda new, old: crot(new, old, interpdeg), ) newinput, inputrot = rotation( inputgrid, - eko.bases.inputgrid, + xgrid, check, lambda new, old: crot(new, old, interpdeg, swap=True), ) - # after the checks: if there is still nothing to do, skip if targetrot is None and inputrot is None: logger.debug("Nothing done.") - return - # if no rotation is done, the grids are not modified - if targetrot is not None: - eko.bases.targetgrid = newtarget - if inputrot is not None: - eko.bases.inputgrid = newinput + return copy.deepcopy(elem) # build new grid - for ep, elem in eko.items(): - assert elem is not None - - operands = [elem.operator] - operands_errs = [elem.error] + operands = [elem.operator] + operands_errs = [elem.error] - if targetrot is not None and inputrot is None: - contraction = TARGETGRID_ROTATION - elif inputrot is not None and targetrot is None: - contraction = INPUTGRID_ROTATION - else: - contraction = SIMGRID_ROTATION + if targetrot is not None and inputrot is None: + contraction = TARGETGRID_ROTATION + elif inputrot is not None and targetrot is None: + contraction = INPUTGRID_ROTATION + else: + contraction = SIMGRID_ROTATION - if targetrot is not None: - operands.insert(0, targetrot) - operands_errs.insert(0, targetrot) - if inputrot is not None: - operands.append(inputrot) - operands_errs.append(inputrot) - - new_operator = np.einsum(contraction, *operands, optimize="optimal") - if elem.error is not None: - new_error = np.einsum(contraction, *operands_errs, optimize="optimal") - else: - new_error = None + if targetrot is not None: + operands.insert(0, targetrot) + operands_errs.insert(0, targetrot) + if inputrot is not None: + operands.append(inputrot) + operands_errs.append(inputrot) - eko[ep] = Operator(operator=new_operator, error=new_error) + new_operator = np.einsum(contraction, *operands, optimize="optimal") + if elem.error is not None: + new_error = np.einsum(contraction, *operands_errs, optimize="optimal") + else: + new_error = None - eko.update() + return Operator(operator=new_operator, error=new_error) TARGETPIDS_ROTATION = "ca,ajbk->cjbk" @@ -145,107 +129,85 @@ def xgrid_reshape( def flavor_reshape( - eko: EKO, + elem: Operator, targetpids: Optional[npt.NDArray] = None, inputpids: Optional[npt.NDArray] = None, - update: bool = True, ): - """Change the operators to have in the output targetpids and/or in the input inputpids. - - The operation is in-place. + """Change the operator to have in the output targetpids and/or in the input inputpids. Parameters ---------- - eko : + elem : the operator to be rotated targetpids : target rotation specified in the flavor basis inputpids : input rotation specified in the flavor basis - update : - update :class:`~eko.io.struct.EKO` metadata after writing """ - eko.assert_permissions(write=True) - # calling with no arguments is an error if targetpids is None and inputpids is None: raise ValueError("Nor inputpids nor targetpids was given") # now check to the current status if targetpids is not None and np.allclose( - targetpids, np.eye(len(eko.bases.targetpids)) + targetpids, np.eye(elem.operator.shape[0]) ): targetpids = None warnings.warn("The new targetpids is close to current basis") - if inputpids is not None and np.allclose( - inputpids, np.eye(len(eko.bases.inputpids)) - ): + if inputpids is not None and np.allclose(inputpids, np.eye(elem.operator.shape[2])): inputpids = None warnings.warn("The new inputpids is close to current basis") # after the checks: if there is still nothing to do, skip if targetpids is None and inputpids is None: logger.debug("Nothing done.") - return + return copy.deepcopy(elem) # flip input around if inputpids is not None: inv_inputpids = np.linalg.inv(inputpids) # build new grid - for q2, elem in eko.items(): - ops = elem.operator - errs = elem.error - if targetpids is not None and inputpids is None: - ops = np.einsum(TARGETPIDS_ROTATION, targetpids, ops, optimize="optimal") - errs = ( - np.einsum(TARGETPIDS_ROTATION, targetpids, errs, optimize="optimal") - if errs is not None - else None - ) - elif inputpids is not None and targetpids is None: - ops = np.einsum(INPUTPIDS_ROTATION, ops, inv_inputpids, optimize="optimal") - errs = ( - np.einsum(INPUTPIDS_ROTATION, errs, inv_inputpids, optimize="optimal") - if errs is not None - else None - ) - else: - ops = np.einsum( - SIMPIDS_ROTATION, targetpids, ops, inv_inputpids, optimize="optimal" - ) - errs = ( - np.einsum( - SIMPIDS_ROTATION, - targetpids, - errs, - inv_inputpids, - optimize="optimal", - ) - if errs is not None - else None + ops = elem.operator + errs = elem.error + if targetpids is not None and inputpids is None: + ops = np.einsum(TARGETPIDS_ROTATION, targetpids, ops, optimize="optimal") + errs = ( + np.einsum(TARGETPIDS_ROTATION, targetpids, errs, optimize="optimal") + if errs is not None + else None + ) + elif inputpids is not None and targetpids is None: + ops = np.einsum(INPUTPIDS_ROTATION, ops, inv_inputpids, optimize="optimal") + errs = ( + np.einsum(INPUTPIDS_ROTATION, errs, inv_inputpids, optimize="optimal") + if errs is not None + else None + ) + else: + ops = np.einsum( + SIMPIDS_ROTATION, targetpids, ops, inv_inputpids, optimize="optimal" + ) + errs = ( + np.einsum( + SIMPIDS_ROTATION, + targetpids, + errs, + inv_inputpids, + optimize="optimal", ) + if errs is not None + else None + ) - eko[q2] = Operator(operator=ops, error=errs) - - # drop PIDs - keeping them int nevertheless - # there is no meaningful way to set them in general, after rotation - if inputpids is not None: - eko.bases.inputpids = np.array([0] * len(eko.bases.inputpids)) - if targetpids is not None: - eko.bases.targetpids = np.array([0] * len(eko.bases.targetpids)) - - if update: - eko.update() + return Operator(operator=ops, error=errs) -def to_evol(eko: EKO, source: bool = True, target: bool = False): +def to_evol(elem: Operator, source: bool = True, target: bool = False): """Rotate the operator into evolution basis. - This assigns also the pids. The operation is in-place. - Parameters ---------- - eko : + elem : the operator to be rotated source : rotate on the input tensor @@ -256,27 +218,15 @@ def to_evol(eko: EKO, source: bool = True, target: bool = False): # rotate inputpids = br.rotate_flavor_to_evolution if source else None targetpids = br.rotate_flavor_to_evolution if target else None - # prevent metadata update, since flavor_reshape has not enough information - # to determine inpupids and targetpids, and they will be updated after the - # call - flavor_reshape(eko, inputpids=inputpids, targetpids=targetpids, update=False) - # assign pids - if source: - eko.bases.inputpids = inputpids - if target: - eko.bases.targetpids = targetpids - - eko.update() + return flavor_reshape(elem, inputpids=inputpids, targetpids=targetpids) -def to_uni_evol(eko: EKO, source: bool = True, target: bool = False): +def to_uni_evol(elem: Operator, source: bool = True, target: bool = False): """Rotate the operator into evolution basis. - This assigns also the pids. The operation is in-place. - Parameters ---------- - eko : + elem : the operator to be rotated source : rotate on the input tensor @@ -287,14 +237,4 @@ def to_uni_evol(eko: EKO, source: bool = True, target: bool = False): # rotate inputpids = br.rotate_flavor_to_unified_evolution if source else None targetpids = br.rotate_flavor_to_unified_evolution if target else None - # prevent metadata update, since flavor_reshape has not enough information - # to determine inpupids and targetpids, and they will be updated after the - # call - flavor_reshape(eko, inputpids=inputpids, targetpids=targetpids, update=False) - # assign pids - if source: - eko.bases.inputpids = inputpids - if target: - eko.bases.targetpids = targetpids - - eko.update() + return flavor_reshape(elem, inputpids=inputpids, targetpids=targetpids) diff --git a/src/ekomark/benchmark/runner.py b/src/ekomark/benchmark/runner.py index 8c07e7153..b1e9fc579 100644 --- a/src/ekomark/benchmark/runner.py +++ b/src/ekomark/benchmark/runner.py @@ -13,7 +13,6 @@ import eko from eko import EKO from eko import basis_rotation as br -from eko.io import manipulate from ekobox import apply from .. import pdfname @@ -116,22 +115,14 @@ def run_me(self, theory, ocard, _pdf): os.makedirs(output_path) # rotating to evolution basis if requested with EKO.edit(path) as out_copy: - change_lab = False - if self.rotate_to_evolution_basis: - qed = theory["QED"] > 0 - if not qed: - manipulate.to_evol(out_copy, source=True, target=True) - else: - manipulate.to_uni_evol(out_copy, source=True, target=True) - change_lab = True - save_operators_to_pdf( output_path, theory, ocard, out_copy, self.skip_pdfs(theory), - change_lab, + self.rotate_to_evolution_basis, + theory["QED"] > 0, ) else: # else we always rerun diff --git a/src/ekomark/plots.py b/src/ekomark/plots.py index 92a155de3..5c300f996 100644 --- a/src/ekomark/plots.py +++ b/src/ekomark/plots.py @@ -9,6 +9,7 @@ from matplotlib.colors import LogNorm from eko import basis_rotation as br +from eko.io import manipulate from eko.io.struct import EKO @@ -209,7 +210,15 @@ def plot_operator(var_name, op, op_err, log_operator=True, abs_operator=True): return fig -def save_operators_to_pdf(path, theory, ops, me: EKO, skip_pdfs, change_lab=False): +def save_operators_to_pdf( + path, + theory, + ops, + me: EKO, + skip_pdfs, + rotate_to_evolution_basis: bool, + qed: bool = False, +): """Output all operator heatmaps to PDF. Parameters @@ -228,11 +237,7 @@ def save_operators_to_pdf(path, theory, ops, me: EKO, skip_pdfs, change_lab=Fals set whether to rename the labels """ - ops_names = list(me.bases.targetpids) - if np.allclose(ops_names, br.rotate_flavor_to_evolution): - ops_names = br.evol_basis_pids - else: - raise ValueError("Can not reconstruct PDF names") + ops_names = br.evol_basis_pids if rotate_to_evolution_basis else br.evol_basis_pids ops_id = f"o{ops['hash'][:6]}_t{theory['hash'][:6]}" path = f"{path}/{ops_id}.pdf" print(f"Plotting operators plots to {path}") @@ -244,9 +249,16 @@ def save_operators_to_pdf(path, theory, ops, me: EKO, skip_pdfs, change_lab=Fals # plot the operators # it's necessary to reshuffle the eko output - for mu2 in me.mu2grid: - results = me[mu2].operator - errors = me[mu2].error + for mu2, op in me.items(): + change_lab = False + if rotate_to_evolution_basis: + if not qed: + op = manipulate.to_evol(op, source=True, target=True) + else: + op = manipulate.to_uni_evol(op, source=True, target=True) + change_lab = True + results = op.operator + errors = op.error # loop on pids for label_out, res, res_err in zip(ops_names, results, errors): diff --git a/tests/eko/io/test_manipulate.py b/tests/eko/io/test_manipulate.py index 3431d93ce..56ceafc0f 100644 --- a/tests/eko/io/test_manipulate.py +++ b/tests/eko/io/test_manipulate.py @@ -1,5 +1,3 @@ -import pathlib - import numpy as np import pytest @@ -21,235 +19,125 @@ def chk_keys(a, b): class TestManipulate: - def test_xgrid_reshape(self, eko_factory: EKOFactory, tmp_path: pathlib.Path): + def test_xgrid_reshape(self): # create object - muout = 10.0 - mu2out = muout**2 - epout = (mu2out, 5) + interpdeg = 1 xg = interpolation.XGrid(np.geomspace(1e-5, 1.0, 21)) - eko_factory.operator.mugrid = [(muout, 5)] - eko_factory.operator.xgrid = xg - o1 = eko_factory.get() + xgp = interpolation.XGrid(np.geomspace(1e-5, 1.0, 11)) lpids = 2 - o1[epout] = eko.io.Operator( + o1 = eko.io.Operator( operator=eko_identity([1, lpids, len(xg), lpids, len(xg)])[0] ) - xgp = interpolation.XGrid(np.geomspace(1e-5, 1.0, 11)) # only target - otpath = tmp_path / "ot.tar" - o1.deepcopy(otpath) - with EKO.edit(otpath) as ot: - manipulate.xgrid_reshape(ot, xgp) - chk_keys(o1.raw, ot.raw) - assert ot[epout].operator.shape == (lpids, len(xgp), lpids, len(xg)) - ottpath = tmp_path / "ott.tar" - o1.deepcopy(ottpath) - with EKO.edit(ottpath) as ott: - with pytest.warns(Warning): - manipulate.xgrid_reshape(ott, xg) - chk_keys(o1.raw, ott.raw) - np.testing.assert_allclose(ott[epout].operator, o1[epout].operator) + ot = manipulate.xgrid_reshape(o1, xg, interpdeg, xgp) + assert ot.operator.shape == (lpids, len(xgp), lpids, len(xg)) + ott = manipulate.xgrid_reshape(ot, xgp, interpdeg, xg) + # when blowing up again a line 0 ... 0 0 1 0 0 ... 0 becomes + # 0 ... 0 0.5 0 0.5 0 ... 0 instead + np.testing.assert_allclose( + np.sum(ott.operator, axis=3), np.sum(o1.operator, axis=3) + ) # only input - oipath = tmp_path / "oi.tar" - o1.deepcopy(oipath) - with EKO.edit(oipath) as oi: - manipulate.xgrid_reshape(oi, inputgrid=xgp) - assert oi[epout].operator.shape == (lpids, len(xg), lpids, len(xgp)) - chk_keys(o1.raw, oi.raw) - oiipath = tmp_path / "oii.tar" - o1.deepcopy(oiipath) - with EKO.edit(oiipath) as oii: - with pytest.warns(Warning): - manipulate.xgrid_reshape(oii, inputgrid=xg) - chk_keys(o1.raw, oii.raw) - np.testing.assert_allclose(oii[epout].operator, o1[epout].operator) + oi = manipulate.xgrid_reshape(o1, xg, interpdeg, inputgrid=xgp) + assert oi.operator.shape == (lpids, len(xg), lpids, len(xgp)) + oii = manipulate.xgrid_reshape(oi, xgp, interpdeg, inputgrid=xg) + np.testing.assert_allclose( + np.sum(oii.operator, axis=3), np.sum(o1.operator, axis=3) + ) + with pytest.warns(Warning): + oiii = manipulate.xgrid_reshape(oii, xg, interpdeg, inputgrid=xg) + np.testing.assert_allclose(oiii.operator, oii.operator) # both - oitpath = tmp_path / "oit.tar" - o1.deepcopy(oitpath) - with EKO.edit(oitpath) as oit: - manipulate.xgrid_reshape(oit, xgp, xgp) - chk_keys(o1.raw, oit.raw) - op = eko_identity([1, 2, len(xgp), 2, len(xgp)]) - np.testing.assert_allclose(oit[epout].operator, op[0], atol=1e-10) - + oit = manipulate.xgrid_reshape(o1, xg, interpdeg, xgp, xgp) + op = eko_identity([1, 2, len(xgp), 2, len(xgp)]) + np.testing.assert_allclose(oit.operator, op[0], atol=1e-10) # op error handling - ep2 = (25, 5) - o1[ep2] = eko.io.Operator( + o1e = eko.io.Operator( operator=eko_identity([1, lpids, len(xg), lpids, len(xg)])[0], error=0.1 * eko_identity([1, lpids, len(xg), lpids, len(xg)])[0], ) - ot2path = tmp_path / "ot2.tar" - o1.deepcopy(ot2path) - with EKO.edit(ot2path) as ot2: - manipulate.xgrid_reshape(ot2, xgp) - chk_keys(o1.raw, ot2.raw) - assert ot2[ep2].operator.shape == (lpids, len(xgp), lpids, len(xg)) - assert ot2[epout].error is None - assert ot2[ep2].error is not None + assert ot.error is None + assert oi.error is None + ot2 = manipulate.xgrid_reshape(o1e, xg, interpdeg, xgp) + assert ot2.error is not None # Python error - with pytest.raises(ValueError): - manipulate.xgrid_reshape(o1) + with pytest.raises(ValueError, match="Nor inputgrid nor targetgrid"): + manipulate.xgrid_reshape(o1, xg, interpdeg) - def test_reshape_io(self, eko_factory: EKOFactory, tmp_path): - eko_factory.path = tmp_path / "eko.tar" - eko_factory.operator.configs.interpolation_polynomial_degree = 1 - # create object - o1 = eko_factory.get() - lpids = len(o1.bases.pids) - path_copy = tmp_path / "eko_copy.tar" - o1.deepcopy(path_copy) - newxgrid = interpolation.XGrid([0.1, 1.0]) - inputpids = np.eye(lpids) - inputpids[:2, :2] = np.array([[1, -1], [1, 1]]) - with EKO.edit(path_copy) as o2: - manipulate.xgrid_reshape(o2, newxgrid, newxgrid) - manipulate.flavor_reshape(o2, inputpids=inputpids) - # reload - with EKO.read(path_copy) as o3: - chk_keys(o1.raw, o3.raw) - np.testing.assert_allclose(o3.bases.inputgrid.raw, newxgrid.raw) - np.testing.assert_allclose(o3.bases.targetgrid.raw, newxgrid.raw) - # since we use a general rotation, the inputpids are erased, - # leaving just as many zeros as PIDs, as placeholders for missing - # values - np.testing.assert_allclose(o3.bases.inputpids, [0] * len(o3.bases.pids)) - # these has to be unchanged - np.testing.assert_allclose(o3.bases.targetpids, o3.bases.pids) - - def test_flavor_reshape(self, eko_factory: EKOFactory, tmp_path: pathlib.Path): + def test_flavor_reshape(self): # create object xg = interpolation.XGrid(np.geomspace(1e-5, 1.0, 21)) - muout = 10.0 - mu2out = muout**2 - epout = (mu2out, 5) - eko_factory.operator.xgrid = xg - eko_factory.operator.mugrid = [(muout, 5)] - o1 = eko_factory.get() - lpids = len(o1.bases.pids) + lpids = len(br.flavor_basis_pids) lx = len(xg) - o1[epout] = eko.io.Operator( + o1 = eko.io.Operator( operator=eko_identity([1, lpids, lx, lpids, lx])[0], error=None, ) - # only target - target_r = np.eye(lpids) - target_r[:2, :2] = np.array([[1, -1], [1, 1]]) - tpath = tmp_path / "ot.tar" - ttpath = tmp_path / "ott.tar" - o1.deepcopy(tpath) - with EKO.edit(tpath) as ot: - manipulate.flavor_reshape(ot, target_r) - chk_keys(o1.raw, ot.raw) - assert ot[epout].operator.shape == (lpids, len(xg), lpids, len(xg)) - ot.deepcopy(ttpath) - with EKO.edit(ttpath) as ott: - manipulate.flavor_reshape(ott, np.linalg.inv(target_r)) - np.testing.assert_allclose(ott[epout].operator, o1[epout].operator) - with pytest.warns(Warning): - manipulate.flavor_reshape(ott, np.eye(lpids)) - chk_keys(o1.raw, ott.raw) - np.testing.assert_allclose(ott[epout].operator, o1[epout].operator) # only input input_r = np.eye(lpids) input_r[:2, :2] = np.array([[1, -1], [1, 1]]) - ipath = tmp_path / "oi.tar" - iipath = tmp_path / "oii.tar" - o1.deepcopy(ipath) - with EKO.edit(ipath) as oi: - manipulate.flavor_reshape(oi, inputpids=input_r) - chk_keys(o1.raw, oi.raw) - assert oi[epout].operator.shape == (lpids, len(xg), lpids, len(xg)) - oi.deepcopy(iipath) - with EKO.edit(iipath) as oii: - manipulate.flavor_reshape(oii, inputpids=np.linalg.inv(input_r)) - np.testing.assert_allclose(oii[epout].operator, o1[epout].operator) - with pytest.warns(Warning): - manipulate.flavor_reshape(oii, inputpids=np.eye(lpids)) - chk_keys(o1.raw, oii.raw) - np.testing.assert_allclose(oii[epout].operator, o1[epout].operator) + oi = manipulate.flavor_reshape(o1, inputpids=input_r) + assert oi.operator.shape == (lpids, len(xg), lpids, len(xg)) + oii = manipulate.flavor_reshape(oi, inputpids=np.linalg.inv(input_r)) + np.testing.assert_allclose(oii.operator, o1.operator) + with pytest.warns(Warning): + oiii = manipulate.flavor_reshape(oii, inputpids=np.eye(lpids)) + np.testing.assert_allclose(oiii.operator, oii.operator) + + # only target + target_r = np.eye(lpids) + target_r[:2, :2] = np.array([[1, -1], [1, 1]]) + ot = manipulate.flavor_reshape(o1, target_r) + assert ot.operator.shape == (lpids, len(xg), lpids, len(xg)) + ott = manipulate.flavor_reshape(ot, np.linalg.inv(target_r)) + np.testing.assert_allclose(ott.operator, o1.operator) + with pytest.warns(Warning): + ottt = manipulate.flavor_reshape(ott, np.eye(lpids)) + np.testing.assert_allclose(ottt.operator, ott.operator) # both - itpath = tmp_path / "oit.tar" - o1.deepcopy(itpath) - with EKO.edit(itpath) as oit: - manipulate.flavor_reshape(oit, target_r, input_r) - chk_keys(o1.raw, oit.raw) - op = eko_identity([1, lpids, len(xg), lpids, len(xg)]).copy() - np.testing.assert_allclose(oit[epout].operator, op[0], atol=1e-10) + oit = manipulate.flavor_reshape(o1, target_r, input_r) + op = eko_identity([1, lpids, len(xg), lpids, len(xg)]).copy() + np.testing.assert_allclose(oit.operator, op[0], atol=1e-10) # error - fpath = tmp_path / "fail.tar" - o1.deepcopy(fpath) - with pytest.raises(ValueError): - with EKO.edit(fpath) as of: - manipulate.flavor_reshape(of) + with pytest.raises(ValueError, match="Nor inputpids nor targetpids"): + manipulate.flavor_reshape(o1) - def test_to_evol(self, eko_factory: EKOFactory, tmp_path): + def test_to_evol(self): self._test_to_all_evol( - eko_factory, - tmp_path, manipulate.to_evol, br.rotate_flavor_to_evolution, - br.flavor_basis_pids, ) - def test_to_uni_evol(self, eko_factory: EKOFactory, tmp_path): + def test_to_uni_evol(self): self._test_to_all_evol( - eko_factory, - tmp_path, manipulate.to_uni_evol, br.rotate_flavor_to_unified_evolution, - br.flavor_basis_pids, ) - def _test_to_all_evol( - self, eko_factory: EKOFactory, tmp_path, to_evol_fnc, rot_matrix, pids - ): - xgrid = interpolation.XGrid([0.5, 1.0]) - mu_out = 2.0 - mu2_out = mu_out**2 - nfout = 4 - epout = (mu2_out, nfout) - eko_factory.operator.mu0 = float(np.sqrt(1.0)) - eko_factory.operator.mugrid = [(mu_out, nfout)] - eko_factory.operator.xgrid = xgrid - eko_factory.operator.configs.interpolation_polynomial_degree = 1 - eko_factory.operator.configs.interpolation_is_log = False - eko_factory.operator.configs.ev_op_max_order = (2, 0) - eko_factory.operator.configs.ev_op_iterations = 1 - eko_factory.operator.configs.inversion_method = runcards.InversionMethod.EXACT - o00 = eko_factory.get() - o01_path = tmp_path / "o01.tar" - o00.deepcopy(o01_path) - with EKO.edit(o01_path) as o01: - to_evol_fnc(o01) - o10_path = tmp_path / "o10.tar" - o00.deepcopy(o10_path) - with EKO.edit(o10_path) as o10: - to_evol_fnc(o10, False, True) - o11_path = tmp_path / "o11.tar" - o00.deepcopy(o11_path) - with EKO.edit(o11_path) as o11: - to_evol_fnc(o11, True, True) - chk_keys(o00.raw, o11.raw) - - with EKO.edit(o01_path) as o01: - with EKO.edit(o10_path) as o10: - with EKO.read(o11_path) as o11: - # check the input rotated one - np.testing.assert_allclose(o01.bases.inputpids, rot_matrix) - np.testing.assert_allclose(o01.bases.targetpids, pids) - # rotate also target - to_evol_fnc(o01, False, True) - np.testing.assert_allclose(o01[epout].operator, o11[epout].operator) - chk_keys(o00.raw, o01.raw) - # check the target rotated one - np.testing.assert_allclose(o10.bases.inputpids, pids) - np.testing.assert_allclose(o10.bases.targetpids, rot_matrix) - # rotate also input - to_evol_fnc(o10) - np.testing.assert_allclose(o10[epout].operator, o11[epout].operator) - chk_keys(o00.raw, o10.raw) + def _test_to_all_evol(self, to_evol_fnc, rot_matrix): + # create object + xg = interpolation.XGrid(np.geomspace(1e-5, 1.0, 21)) + lpids = len(br.flavor_basis_pids) + lx = len(xg) + o = eko.io.Operator( + operator=eko_identity([1, lpids, lx, lpids, lx])[0], + error=None, + ) + + # do it once + o01 = to_evol_fnc(o, True, False) + o10 = to_evol_fnc(o, False, True) + o11 = to_evol_fnc(o, True, True) + + # do also the other one + np.testing.assert_allclose( + to_evol_fnc(o01, False, True).operator, o11.operator, atol=1e-15 + ) + np.testing.assert_allclose( + to_evol_fnc(o10, True, False).operator, o11.operator, atol=1e-15 + ) From 4699064f9e4c5e368e148f1afb9614e0e9908689 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Mon, 31 Jul 2023 11:26:54 +0200 Subject: [PATCH 02/20] Address review comments --- src/eko/io/manipulate.py | 12 ++++++------ src/ekomark/plots.py | 11 +++++------ 2 files changed, 11 insertions(+), 12 deletions(-) diff --git a/src/eko/io/manipulate.py b/src/eko/io/manipulate.py index 1372e9e14..72f3a641b 100644 --- a/src/eko/io/manipulate.py +++ b/src/eko/io/manipulate.py @@ -23,7 +23,7 @@ def rotation(new: Optional[XGrid], old: XGrid, check: Callable, compute: Callable): """Define grid rotation. - This function returns the new grid to be assigned and the rotation computed, + This function returns the necessary rotation, if the checks for a non-trivial new grid are passed. However, the check and the computation are delegated respectively to the @@ -31,13 +31,13 @@ def rotation(new: Optional[XGrid], old: XGrid, check: Callable, compute: Callabl """ if new is None: - return old, None + return None if check(new, old): warnings.warn("The new grid is close to the current one") - return old, None + return None - return new, compute(new, old) + return compute(new, old) def xgrid_check(new: Optional[XGrid], old: XGrid): @@ -78,13 +78,13 @@ def xgrid_reshape( crot = xgrid_compute_rotation # construct matrices - newtarget, targetrot = rotation( + targetrot = rotation( targetgrid, xgrid, check, lambda new, old: crot(new, old, interpdeg), ) - newinput, inputrot = rotation( + inputrot = rotation( inputgrid, xgrid, check, diff --git a/src/ekomark/plots.py b/src/ekomark/plots.py index 5c300f996..2e74586c6 100644 --- a/src/ekomark/plots.py +++ b/src/ekomark/plots.py @@ -233,9 +233,10 @@ def save_operators_to_pdf( DGLAP result skip_pdfs : list PDF to skip - change_lab : bool - set whether to rename the labels - + rotate_to_evolution_basis : bool + plot operators in evolution basis + qed : bool + plot operators in unified evolution basis, if the former is active """ ops_names = br.evol_basis_pids if rotate_to_evolution_basis else br.evol_basis_pids ops_id = f"o{ops['hash'][:6]}_t{theory['hash'][:6]}" @@ -250,13 +251,11 @@ def save_operators_to_pdf( # plot the operators # it's necessary to reshuffle the eko output for mu2, op in me.items(): - change_lab = False if rotate_to_evolution_basis: if not qed: op = manipulate.to_evol(op, source=True, target=True) else: op = manipulate.to_uni_evol(op, source=True, target=True) - change_lab = True results = op.operator errors = op.error @@ -283,7 +282,7 @@ def save_operators_to_pdf( continue lab_in = label_in lab_out = label_out - if change_lab: + if rotate_to_evolution_basis: index_in = br.evol_basis_pids.index(label_in) index_out = br.evol_basis_pids.index(label_out) lab_in = br.evol_basis[index_in] From 0ee166af669a0fbd50e48c4415ee4bd07cdd814d Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Mon, 31 Jul 2023 11:51:39 +0200 Subject: [PATCH 03/20] Add missing typ hints --- src/eko/io/manipulate.py | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/eko/io/manipulate.py b/src/eko/io/manipulate.py index 72f3a641b..c54883fdd 100644 --- a/src/eko/io/manipulate.py +++ b/src/eko/io/manipulate.py @@ -20,7 +20,9 @@ """Simultaneous grid rotation contraction indices.""" -def rotation(new: Optional[XGrid], old: XGrid, check: Callable, compute: Callable): +def rotation( + new: Optional[XGrid], old: XGrid, check: Callable, compute: Callable +) -> npt.NDArray: """Define grid rotation. This function returns the necessary rotation, @@ -40,12 +42,14 @@ def rotation(new: Optional[XGrid], old: XGrid, check: Callable, compute: Callabl return compute(new, old) -def xgrid_check(new: Optional[XGrid], old: XGrid): +def xgrid_check(new: Optional[XGrid], old: XGrid) -> bool: """Check validity of new xgrid.""" return new is not None and len(new) == len(old) and np.allclose(new.raw, old.raw) -def xgrid_compute_rotation(new: XGrid, old: XGrid, interpdeg: int, swap: bool = False): +def xgrid_compute_rotation( + new: XGrid, old: XGrid, interpdeg: int, swap: bool = False +) -> npt.NDArray: """Compute rotation from old to new xgrid. By default, the roation is computed for a target xgrid. Whether the function @@ -132,7 +136,7 @@ def flavor_reshape( elem: Operator, targetpids: Optional[npt.NDArray] = None, inputpids: Optional[npt.NDArray] = None, -): +) -> Operator: """Change the operator to have in the output targetpids and/or in the input inputpids. Parameters @@ -202,7 +206,7 @@ def flavor_reshape( return Operator(operator=ops, error=errs) -def to_evol(elem: Operator, source: bool = True, target: bool = False): +def to_evol(elem: Operator, source: bool = True, target: bool = False) -> Operator: """Rotate the operator into evolution basis. Parameters @@ -221,7 +225,7 @@ def to_evol(elem: Operator, source: bool = True, target: bool = False): return flavor_reshape(elem, inputpids=inputpids, targetpids=targetpids) -def to_uni_evol(elem: Operator, source: bool = True, target: bool = False): +def to_uni_evol(elem: Operator, source: bool = True, target: bool = False) -> Operator: """Rotate the operator into evolution basis. Parameters From 3e2d12a4ff0a1013ce33cd07b193a2c52c9e4939 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Tue, 1 Aug 2023 17:06:34 +0200 Subject: [PATCH 04/20] Drop Bases --- src/eko/io/bases.py | 117 ---------------------------------- src/eko/io/metadata.py | 8 +-- src/eko/io/struct.py | 12 +--- src/ekobox/apply.py | 12 ++-- tests/eko/io/test_bases.py | 83 ------------------------ tests/eko/io/test_metadata.py | 6 +- tests/eko/io/test_runcards.py | 1 - tests/eko/runner/__init__.py | 10 +-- tests/eko/runner/conftest.py | 3 +- tests/ekobox/test_apply.py | 11 ++-- 10 files changed, 27 insertions(+), 236 deletions(-) delete mode 100644 src/eko/io/bases.py delete mode 100644 tests/eko/io/test_bases.py diff --git a/src/eko/io/bases.py b/src/eko/io/bases.py deleted file mode 100644 index 97a28dc36..000000000 --- a/src/eko/io/bases.py +++ /dev/null @@ -1,117 +0,0 @@ -"""Operators bases.""" -from dataclasses import dataclass, fields -from typing import Optional - -import numpy as np -import numpy.typing as npt - -from .. import basis_rotation as br -from .. import interpolation -from .dictlike import DictLike - - -@dataclass -class Bases(DictLike): - """Rotations related configurations. - - Here "Rotation" is intended in a broad sense: it includes both rotations in - flavor space (labeled with suffix `pids`) and in :math:`x`-space (labeled - with suffix `grid`). - Rotations in :math:`x`-space correspond to reinterpolate the result on a - different basis of polynomials. - - """ - - xgrid: interpolation.XGrid - """Internal momentum fraction grid.""" - _targetgrid: Optional[interpolation.XGrid] = None - _inputgrid: Optional[interpolation.XGrid] = None - _targetpids: Optional[npt.NDArray] = None - _inputpids: Optional[npt.NDArray] = None - - def __post_init__(self): - """Adjust types when loaded from serialized object.""" - for attr in ("xgrid", "_inputgrid", "_targetgrid"): - value = getattr(self, attr) - if value is None: - continue - if isinstance(value, (np.ndarray, list)): - setattr(self, attr, interpolation.XGrid(value)) - elif not isinstance(value, interpolation.XGrid): - setattr(self, attr, interpolation.XGrid.load(value)) - - @property - def pids(self): - """Internal flavor basis, used for computation.""" - return np.array(br.flavor_basis_pids) - - @property - def inputpids(self) -> npt.NDArray: - """Provide pids expected on the input PDF.""" - if self._inputpids is None: - return self.pids - return self._inputpids - - @inputpids.setter - def inputpids(self, value): - self._inputpids = value - - @property - def targetpids(self) -> npt.NDArray: - """Provide pids corresponding to the output PDF.""" - if self._targetpids is None: - return self.pids - return self._targetpids - - @targetpids.setter - def targetpids(self, value): - self._targetpids = value - - @property - def inputgrid(self) -> interpolation.XGrid: - """Provide :math:`x`-grid expected on the input PDF.""" - if self._inputgrid is None: - return self.xgrid - return self._inputgrid - - @inputgrid.setter - def inputgrid(self, value: interpolation.XGrid): - self._inputgrid = value - - @property - def targetgrid(self) -> interpolation.XGrid: - """Provide :math:`x`-grid corresponding to the output PDF.""" - if self._targetgrid is None: - return self.xgrid - return self._targetgrid - - @targetgrid.setter - def targetgrid(self, value: interpolation.XGrid): - self._targetgrid = value - - @classmethod - def from_dict(cls, dictionary: dict): - """Deserialize rotation. - - Load from full state, but with public names. - - """ - d = dictionary.copy() - for f in fields(cls): - if f.name.startswith("_"): - d[f.name] = d.pop(f.name[1:]) - return cls._from_dict(d) - - @property - def raw(self): - """Serialize rotation. - - Pass through interfaces, access internal values but with a public name. - - """ - d = self._raw() - for key in d.copy(): - if key.startswith("_"): - d[key[1:]] = d.pop(key) - - return d diff --git a/src/eko/io/metadata.py b/src/eko/io/metadata.py index 415b2dfa5..00f189ab4 100644 --- a/src/eko/io/metadata.py +++ b/src/eko/io/metadata.py @@ -8,7 +8,7 @@ import yaml from .. import version as vmod -from .bases import Bases +from ..interpolation import XGrid from .dictlike import DictLike from .paths import InternalPaths from .types import EvolutionPoint as EPoint @@ -34,10 +34,8 @@ class Metadata(DictLike): origin: EPoint """Inital scale.""" - bases: Bases - """Manipulation information, describing the current status of the EKO (e.g. - `inputgrid` and `targetgrid`). - """ + xgrid: XGrid + """Interpolation grid""" # tagging information _path: Optional[pathlib.Path] = None """Path to temporary dir.""" diff --git a/src/eko/io/struct.py b/src/eko/io/struct.py index fecba3c57..dbbb231f9 100644 --- a/src/eko/io/struct.py +++ b/src/eko/io/struct.py @@ -16,7 +16,6 @@ from .. import interpolation from . import exceptions, raw from .access import AccessConfigs -from .bases import Bases from .inventory import Inventory from .items import Evolution, Matching, Operator, Recipe, Target from .metadata import Metadata @@ -105,20 +104,15 @@ def paths(self) -> InternalPaths: """Accessor for internal paths.""" return InternalPaths(self.metadata.path) - @property - def bases(self) -> Bases: - """Bases information.""" - return self.metadata.bases - @property def xgrid(self) -> interpolation.XGrid: """Momentum fraction internal grid.""" - return self.bases.xgrid + return self.metadata.xgrid @xgrid.setter def xgrid(self, value: interpolation.XGrid): """Set `xgrid` value.""" - self.bases.xgrid = value + self.metadata.xgrid = value self.update() @property @@ -551,7 +545,7 @@ def build(self) -> EKO: metadata = Metadata( _path=self.path, origin=(self.operator.mu20, self.theory.heavy.num_flavs_init), - bases=Bases(xgrid=self.operator.xgrid), + xgrid=self.operator.xgrid, ) InternalPaths(self.path).bootstrap( theory=self.theory.raw, diff --git a/src/ekobox/apply.py b/src/ekobox/apply.py index c392ffe78..5722a6ad6 100644 --- a/src/ekobox/apply.py +++ b/src/ekobox/apply.py @@ -74,12 +74,12 @@ def apply_pdf_flavor( output PDFs and their associated errors for the computed mu2grid """ # create pdfs - pdfs = np.zeros((len(eko.bases.inputpids), len(eko.bases.inputgrid))) - for j, pid in enumerate(eko.bases.inputpids): + pdfs = np.zeros((len(br.flavor_basis_pids), len(eko.xgrid))) + for j, pid in enumerate(br.flavor_basis_pids): if not lhapdf_like.hasFlavor(pid): continue pdfs[j] = np.array( - [lhapdf_like.xfxQ2(pid, x, eko.mu20) / x for x in eko.bases.inputgrid.raw] + [lhapdf_like.xfxQ2(pid, x, eko.mu20) / x for x in eko.xgrid.raw] ) # build output @@ -91,11 +91,11 @@ def apply_pdf_flavor( else: error_final = None out_grid[ep] = { - "pdfs": dict(zip(eko.bases.targetpids, pdf_final)), + "pdfs": dict(zip(br.flavor_basis_pids, pdf_final)), "errors": None, } if error_final is not None: - out_grid[ep]["errors"] = dict(zip(eko.bases.targetpids, error_final)) + out_grid[ep]["errors"] = dict(zip(br.flavor_basis_pids, error_final)) # rotate to evolution basis if flavor_rotation is not None: @@ -117,7 +117,7 @@ def apply_pdf_flavor( # rotate/interpolate to target grid if targetgrid is not None: b = interpolation.InterpolatorDispatcher( - xgrid=eko.bases.targetgrid, + xgrid=eko.xgrid, polynomial_degree=eko.operator_card.configs.interpolation_polynomial_degree, mode_N=False, ) diff --git a/tests/eko/io/test_bases.py b/tests/eko/io/test_bases.py deleted file mode 100644 index a0165243f..000000000 --- a/tests/eko/io/test_bases.py +++ /dev/null @@ -1,83 +0,0 @@ -from dataclasses import fields - -import numpy as np - -from eko import basis_rotation as br -from eko import interpolation -from eko.io.bases import Bases - - -class TestBases: - XGRID_TEST = [1e-3, 1e-2, 1e-1, 1.0] - - def test_serialization(self): - rot = Bases(interpolation.XGrid(self.XGRID_TEST)) - - d = rot.raw - rot1 = rot.from_dict(d) - - for f in fields(Bases): - assert getattr(rot, f.name) == getattr(rot1, f.name) - - assert d["targetgrid"] is None - assert "_targetgrid" not in d - - def test_pids(self): - rot = Bases(interpolation.XGrid(self.XGRID_TEST)) - - # no check on correctness of value set - rot.inputpids = [0, 1] - # but the internal grid is unmodified - assert len(rot.pids) == 14 - # and fallback implemented for unset external bases - assert np.all(rot.targetpids == rot.pids) - - def test_grids(self): - rot = Bases(interpolation.XGrid(self.XGRID_TEST)) - - # no check on correctness of value set - rot.inputgrid = interpolation.XGrid([0.1, 1]) - # but the internal grid is unmodified - assert len(rot.xgrid) == len(self.XGRID_TEST) - # and fallback implemented for unset external grids - assert np.all(rot.targetgrid == rot.xgrid) - - def test_fallback(self): - xg = interpolation.XGrid([0.1, 1.0]) - r = Bases(xgrid=xg) - np.testing.assert_allclose(r.targetpids, r.pids) - np.testing.assert_allclose(r.inputpids, r.pids) - assert r.xgrid == xg - assert r.targetgrid == xg - assert r.inputgrid == xg - - def test_overwrite(self): - tpids = np.array([3, 4] + list(br.flavor_basis_pids[2:])) - ipids = np.array([5, 6] + list(br.flavor_basis_pids[2:])) - xg = interpolation.XGrid([0.1, 1.0]) - txg = interpolation.XGrid([0.2, 1.0]) - ixg = interpolation.XGrid([0.3, 1.0]) - r = Bases( - xgrid=xg, - _targetgrid=txg, - _inputgrid=ixg, - _targetpids=tpids, - _inputpids=ipids, - ) - np.testing.assert_allclose(r.targetpids, tpids) - np.testing.assert_allclose(r.inputpids, ipids) - assert r.xgrid == xg - assert r.targetgrid == txg - assert r.inputgrid == ixg - - def test_init(self): - xg = interpolation.XGrid([0.1, 1.0]) - txg = np.array([0.2, 1.0]) - ixg = {"grid": [0.3, 1.0], "log": True} - r = Bases(xgrid=xg, _targetgrid=txg, _inputgrid=ixg) - assert isinstance(r.xgrid, interpolation.XGrid) - assert isinstance(r.targetgrid, interpolation.XGrid) - assert isinstance(r.inputgrid, interpolation.XGrid) - assert r.xgrid == xg - assert r.targetgrid == interpolation.XGrid(txg) - assert r.inputgrid == interpolation.XGrid.load(ixg) diff --git a/tests/eko/io/test_metadata.py b/tests/eko/io/test_metadata.py index a5651e0c8..fcfa8dd09 100644 --- a/tests/eko/io/test_metadata.py +++ b/tests/eko/io/test_metadata.py @@ -3,11 +3,11 @@ import pytest import yaml -from eko.io import bases, metadata, paths +from eko.io import metadata, paths def test_metadata(tmp_path, caplog): - m = metadata.Metadata((1.0, 3), bases.Bases([0.1, 1.0])) + m = metadata.Metadata(origin=(1.0, 3), xgrid=[0.1, 1.0]) # errors with caplog.at_level(logging.INFO): m.update() @@ -26,7 +26,7 @@ def test_metadata(tmp_path, caplog): m.version = "0.0.0-a1~really1.0.0" m.update() # if I read back the thing should be what I set - mn = metadata.Metadata((1.0, 3), bases.Bases([0.1, 1.0])) + mn = metadata.Metadata(origin=(1.0, 3), xgrid=[0.1, 1.0]) mm = metadata.Metadata.load(tmp_path) assert m.path == tmp_path assert mm.version != mn.version diff --git a/tests/eko/io/test_runcards.py b/tests/eko/io/test_runcards.py index 6af37edcb..4b4996336 100644 --- a/tests/eko/io/test_runcards.py +++ b/tests/eko/io/test_runcards.py @@ -7,7 +7,6 @@ from banana.data.theories import default_card as theory_card from eko.io import runcards as rc -from eko.io.bases import Bases from ekomark.data.operators import default_card as operator_card diff --git a/tests/eko/runner/__init__.py b/tests/eko/runner/__init__.py index c8ec5b854..0e3ea68c0 100644 --- a/tests/eko/runner/__init__.py +++ b/tests/eko/runner/__init__.py @@ -1,19 +1,21 @@ import numpy as np +from eko import basis_rotation as br + def check_shapes(o, txs, ixs, theory_card, operators_card): - tpids = len(o.bases.targetpids) - ipids = len(o.bases.inputpids) + tpids = len(br.flavor_basis_pids) + ipids = len(br.flavor_basis_pids) op_shape = (tpids, len(txs), ipids, len(ixs)) # check output = input np.testing.assert_allclose(o.xgrid.raw, operators_card.xgrid.raw) # targetgrid and inputgrid in the opcard are now ignored, we are testing this np.testing.assert_allclose( - o.bases.targetgrid.raw, + o.xgrid.raw, txs.raw, ) - np.testing.assert_allclose(o.bases.inputgrid.raw, ixs.raw) + np.testing.assert_allclose(o.xgrid.raw, ixs.raw) np.testing.assert_allclose(o.mu20, operators_card.mu20) # check available operators ~o.operators diff --git a/tests/eko/runner/conftest.py b/tests/eko/runner/conftest.py index a1516b2a7..a2966c5e9 100644 --- a/tests/eko/runner/conftest.py +++ b/tests/eko/runner/conftest.py @@ -4,6 +4,7 @@ import pytest from eko import EKO +from eko import basis_rotation as br from eko.io.items import Operator from eko.io.runcards import OperatorCard, TheoryCard from eko.runner import commons, recipes @@ -21,7 +22,7 @@ def neweko(theory_card: TheoryCard, operator_card: OperatorCard, tmp_path: Path) @pytest.fixture def identity(neweko: EKO): xs = len(neweko.xgrid.raw) - flavs = len(neweko.bases.pids) + flavs = len(br.flavor_basis_pids) return Operator(operator=np.eye(xs * flavs).reshape((xs, flavs, xs, flavs))) diff --git a/tests/ekobox/test_apply.py b/tests/ekobox/test_apply.py index 8fa23548f..70eb74adc 100644 --- a/tests/ekobox/test_apply.py +++ b/tests/ekobox/test_apply.py @@ -1,5 +1,6 @@ import numpy as np +from eko import basis_rotation as br from ekobox import apply from tests.conftest import EKOFactory @@ -12,13 +13,13 @@ def test_apply(self, eko_factory: EKOFactory, fake_pdf): pdf_grid = apply.apply_pdf(eko, fake_pdf) assert len(pdf_grid) == len(eko.evolgrid) pdfs = pdf_grid[ep_out]["pdfs"] - assert list(pdfs.keys()) == list(eko.bases.targetpids) + assert list(pdfs.keys()) == list(br.flavor_basis_pids) # rotate to target_grid target_grid = [0.75] pdf_grid = apply.apply_pdf(eko, fake_pdf, target_grid) assert len(pdf_grid) == 1 pdfs = pdf_grid[ep_out]["pdfs"] - assert list(pdfs.keys()) == list(eko.bases.targetpids) + assert list(pdfs.keys()) == list(br.flavor_basis_pids) def test_apply_flavor(self, eko_factory: EKOFactory, fake_pdf, monkeypatch): eko = eko_factory.get() @@ -27,12 +28,8 @@ def test_apply_flavor(self, eko_factory: EKOFactory, fake_pdf, monkeypatch): monkeypatch.setattr( "eko.basis_rotation.rotate_flavor_to_evolution", np.ones((14, 14)) ) - monkeypatch.setattr( - "eko.basis_rotation.flavor_basis_pids", - eko.bases.targetpids, - ) fake_evol_basis = tuple( - ["a", "b"] + [str(x) for x in range(len(eko.bases.pids) - 2)] + ["a", "b"] + [str(x) for x in range(len(br.flavor_basis_pids) - 2)] ) monkeypatch.setattr("eko.basis_rotation.evol_basis", fake_evol_basis) pdf_grid = apply.apply_pdf(eko, fake_pdf, rotate_to_evolution_basis=True) From 78e533cee78e30e441973fdc219d4be8f27ac4eb Mon Sep 17 00:00:00 2001 From: juacrumar Date: Tue, 16 Jul 2024 09:41:34 +0200 Subject: [PATCH 05/20] use nf to generate alpha_s --- src/ekobox/info_file.py | 24 +++++++++--------------- 1 file changed, 9 insertions(+), 15 deletions(-) diff --git a/src/ekobox/info_file.py b/src/ekobox/info_file.py index 26e3cb170..fb744ce4b 100644 --- a/src/ekobox/info_file.py +++ b/src/ekobox/info_file.py @@ -68,19 +68,13 @@ def build( hqm_scheme=theory_card.heavy.masses_scheme, thresholds_ratios=np.power(list(iter(theory_card.heavy.matching_ratios)), 2.0), ) - alphas_values = np.array( - [ - 4.0 - * np.pi - * sc.a_s( - muf2, - ) - for muf2 in operators_card.mu2grid - ], - dtype=float, - ) - template_info["AlphaS_Vals"] = alphas_values.tolist() - template_info["AlphaS_Qs"] = np.array( - [mu for mu, _ in operators_card.mugrid] - ).tolist() + + alphas_values = [] + alphas_qs = [] + for mu, nf in operators_card.mugrid: + alphas_values.append(float(4.0 * np.pi * sc.a_s(mu * mu, nf_to=nf))) + alphas_qs.append(mu) + + template_info["AlphaS_Vals"] = alphas_values + template_info["AlphaS_Qs"] = alphas_qs return template_info From cb6953df9adeef319c6f23f6aedcc21993ad601b Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Tue, 16 Jul 2024 14:15:58 +0300 Subject: [PATCH 06/20] Check sorting of alphas in LHAPDF info --- src/ekobox/evol_pdf.py | 10 +----- src/ekobox/info_file.py | 62 +++++++++++++++++++++++++++------- src/ekobox/utils.py | 9 +++++ tests/ekobox/test_info_file.py | 40 ++++++++++++++++++++++ 4 files changed, 100 insertions(+), 21 deletions(-) diff --git a/src/ekobox/evol_pdf.py b/src/ekobox/evol_pdf.py index ddd80f75a..0045cdd28 100644 --- a/src/ekobox/evol_pdf.py +++ b/src/ekobox/evol_pdf.py @@ -1,13 +1,13 @@ """Tools to evolve actual PDFs.""" import pathlib -from collections import defaultdict from eko import basis_rotation as br from eko.io import EKO from eko.runner import managed from . import apply, genpdf, info_file +from .utils import regroup_evolgrid DEFAULT_NAME = "eko.tar" @@ -95,14 +95,6 @@ def evolve_pdfs( genpdf.install_pdf(name) -def regroup_evolgrid(evolgrid: list): - """Split evolution points by nf and sort by scale.""" - by_nf = defaultdict(list) - for q, nf in sorted(evolgrid, key=lambda ep: ep[1]): - by_nf[nf].append(q) - return {nf: sorted(qs) for nf, qs in by_nf.items()} - - def collect_blocks(evolved_PDF: dict, q2block_per_nf: dict, xgrid: list): """Collect all LHAPDF blocks for a given replica. diff --git a/src/ekobox/info_file.py b/src/ekobox/info_file.py index fb744ce4b..c9ac88b0c 100644 --- a/src/ekobox/info_file.py +++ b/src/ekobox/info_file.py @@ -9,6 +9,7 @@ from eko.io.runcards import OperatorCard, TheoryCard from .genpdf import load +from .utils import regroup_evolgrid def build( @@ -16,19 +17,19 @@ def build( operators_card: OperatorCard, num_members: int, info_update: dict, -): - """Generate a lhapdf info file from theory and operators card. +) -> dict: + """Generate a lhapdf info file. Parameters ---------- - theory_card : dict + theory_card : theory card - operators_card : dict - operators_card - num_members : int + operators_card : + operators card + num_members : number of pdf set members - info_update : dict - info to update + info_update : + additional info to update Returns ------- @@ -56,8 +57,45 @@ def build( template_info["MCharm"] = theory_card.heavy.masses.c.value template_info["MBottom"] = theory_card.heavy.masses.b.value template_info["MTop"] = theory_card.heavy.masses.t.value + # dump alphas + template_info.update(build_alphas(theory_card, operators_card)) + return template_info + + +def build_alphas( + theory_card: TheoryCard, + operators_card: OperatorCard, +) -> dict: + """Generate a couplings section of lhapdf info file. + + Parameters + ---------- + theory_card : dict + theory card + operators_card : dict + operators card + + Returns + ------- + dict + info file section in lhapdf format + """ + # start with meta stuff + template_info = {} template_info["AlphaS_MZ"] = theory_card.couplings.alphas template_info["AlphaS_OrderQCD"] = theory_card.order[0] - 1 + + # check we have disjoint scale ranges + evolgrid = regroup_evolgrid(operators_card.mugrid) + nfs = list(evolgrid.keys()) + for j in range(len(nfs) - 1): + # equal points are allowed by LHAPDF + if evolgrid[nfs[j]][-1] > evolgrid[nfs[j + 1]][0]: + raise ValueError( + f"Last scale point for nf={nfs[j]} is bigger than first in nf={nfs[j+1]}" + ) + + # add actual values evmod = couplings.couplings_mod_ev(operators_card.configs.evolution_method) quark_masses = [(x.value) ** 2 for x in theory_card.heavy.masses] sc = couplings.Couplings( @@ -68,12 +106,12 @@ def build( hqm_scheme=theory_card.heavy.masses_scheme, thresholds_ratios=np.power(list(iter(theory_card.heavy.matching_ratios)), 2.0), ) - alphas_values = [] alphas_qs = [] - for mu, nf in operators_card.mugrid: - alphas_values.append(float(4.0 * np.pi * sc.a_s(mu * mu, nf_to=nf))) - alphas_qs.append(mu) + for nf, mus in evolgrid.items(): + for mu in mus: + alphas_values.append(float(4.0 * np.pi * sc.a_s(mu * mu, nf_to=nf))) + alphas_qs.append(mu) template_info["AlphaS_Vals"] = alphas_values template_info["AlphaS_Qs"] = alphas_qs diff --git a/src/ekobox/utils.py b/src/ekobox/utils.py index 0047d0872..091cde15e 100644 --- a/src/ekobox/utils.py +++ b/src/ekobox/utils.py @@ -1,6 +1,7 @@ """Generic utilities to work with EKOs.""" import os +from collections import defaultdict from typing import Optional import numpy as np @@ -78,3 +79,11 @@ def ekos_product( if path is not None: final_eko.close() + + +def regroup_evolgrid(evolgrid: list): + """Split evolution points by nf and sort by scale.""" + by_nf = defaultdict(list) + for q, nf in sorted(evolgrid, key=lambda ep: ep[1]): + by_nf[nf].append(q) + return {nf: sorted(qs) for nf, qs in by_nf.items()} diff --git a/tests/ekobox/test_info_file.py b/tests/ekobox/test_info_file.py index 64af8eb93..f7267fbd4 100644 --- a/tests/ekobox/test_info_file.py +++ b/tests/ekobox/test_info_file.py @@ -1,6 +1,7 @@ import math import numpy as np +import pytest from ekobox import cards, info_file @@ -23,3 +24,42 @@ def test_build(): np.testing.assert_allclose(info["QMin"], math.sqrt(op.mu2grid[0]), rtol=1e-5) assert info["XMin"] == op.xgrid.raw[0] assert info["XMax"] == op.xgrid.raw[-1] == 1.0 + + +def test_build_alphas_good(): + """Good configurations.""" + theory = cards.example.theory() + theory.order = (2, 0) + theory.couplings.alphas = 0.2 + op = cards.example.operator() + op.mu0 = 1.0 + # base case + op.mugrid = [(10.0, 5), (100.0, 5)] + info = info_file.build_alphas(theory, op) + assert len(info["AlphaS_Vals"]) == 2 + np.testing.assert_allclose(info["AlphaS_Qs"], [10.0, 100.0]) + # several nf + op.mugrid = [(1.0, 3), (5.0, 3), (5.0, 4), (10.0, 5), (10.0, 5), (100.0, 5)] + info = info_file.build_alphas(theory, op) + assert len(info["AlphaS_Vals"]) == 6 + np.testing.assert_allclose(info["AlphaS_Qs"], [1.0, 5.0, 5.0, 10.0, 10.0, 100.0]) + # several nf with gap + op.mugrid = [(1.0, 3), (10.0, 3), (10.0, 5), (100.0, 5)] + info = info_file.build_alphas(theory, op) + assert len(info["AlphaS_Vals"]) == 4 + np.testing.assert_allclose(info["AlphaS_Qs"], [1.0, 10.0, 10.0, 100.0]) + + +def test_build_alphas_bad(): + """Bad configurations.""" + theory = cards.example.theory() + theory.order = (2, 0) + theory.couplings.alphas = 0.2 + op = cards.example.operator() + op.mu0 = 1.0 + op.mugrid = [(5.0, 3), (15.0, 4), (10.0, 5), (100.0, 5)] + with pytest.raises(ValueError, match="is bigger"): + info_file.build_alphas(theory, op) + op.mugrid = [(5.0, 3), (10.0, 3), (10.0, 4), (15.0, 4), (10.0, 5)] + with pytest.raises(ValueError, match="is bigger"): + info_file.build_alphas(theory, op) From d5f63fe7ac995fd9d0d38dba48f2a7969de01c28 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Tue, 16 Jul 2024 14:53:40 +0300 Subject: [PATCH 07/20] Lift scale check for LHAPDF --- src/ekobox/evol_pdf.py | 19 ++++++++++++++++--- src/ekobox/info_file.py | 13 ++----------- tests/ekobox/test_evol_pdf.py | 16 ++++++++++++++++ tests/ekobox/test_info_file.py | 20 ++------------------ 4 files changed, 36 insertions(+), 32 deletions(-) diff --git a/src/ekobox/evol_pdf.py b/src/ekobox/evol_pdf.py index 0045cdd28..abdb2ec6a 100644 --- a/src/ekobox/evol_pdf.py +++ b/src/ekobox/evol_pdf.py @@ -2,6 +2,8 @@ import pathlib +import numpy as np + from eko import basis_rotation as br from eko.io import EKO from eko.runner import managed @@ -46,6 +48,18 @@ def evolve_pdfs( info_update : dict dict of info to add or update to default info file """ + # separate by nf the evolgrid (and order per nf/q) + q2block_per_nf = regroup_evolgrid(operators_card.mugrid) + + # check we have disjoint scale ranges + nfs = list(q2block_per_nf.keys()) + for j in range(len(nfs) - 1): + # equal points are allowed by LHAPDF + if q2block_per_nf[nfs[j]][-1] > q2block_per_nf[nfs[j + 1]][0]: + raise ValueError( + f"Last scale point for nf={nfs[j]} is bigger than first in nf={nfs[j+1]}" + ) + # update op and th cards if path is not None: eko_path = pathlib.Path(path) @@ -79,9 +93,8 @@ def evolve_pdfs( info_update=info_update, ) - # separate by nf the evolgrid (and order per nf/q) - q2block_per_nf = regroup_evolgrid(eko_output.evolgrid) - + # in the eko scales are squared + q2block_per_nf = {nf: np.power(q2s, 2) for nf, q2s in q2block_per_nf.items()} # write all replicas all_member_blocks = [] for evolved_PDF in evolved_PDF_list: diff --git a/src/ekobox/info_file.py b/src/ekobox/info_file.py index c9ac88b0c..fc5c17f1d 100644 --- a/src/ekobox/info_file.py +++ b/src/ekobox/info_file.py @@ -84,18 +84,8 @@ def build_alphas( template_info = {} template_info["AlphaS_MZ"] = theory_card.couplings.alphas template_info["AlphaS_OrderQCD"] = theory_card.order[0] - 1 - - # check we have disjoint scale ranges + # prepare evolgrid = regroup_evolgrid(operators_card.mugrid) - nfs = list(evolgrid.keys()) - for j in range(len(nfs) - 1): - # equal points are allowed by LHAPDF - if evolgrid[nfs[j]][-1] > evolgrid[nfs[j + 1]][0]: - raise ValueError( - f"Last scale point for nf={nfs[j]} is bigger than first in nf={nfs[j+1]}" - ) - - # add actual values evmod = couplings.couplings_mod_ev(operators_card.configs.evolution_method) quark_masses = [(x.value) ** 2 for x in theory_card.heavy.masses] sc = couplings.Couplings( @@ -106,6 +96,7 @@ def build_alphas( hqm_scheme=theory_card.heavy.masses_scheme, thresholds_ratios=np.power(list(iter(theory_card.heavy.matching_ratios)), 2.0), ) + # add actual values alphas_values = [] alphas_qs = [] for nf, mus in evolgrid.items(): diff --git a/tests/ekobox/test_evol_pdf.py b/tests/ekobox/test_evol_pdf.py index 54bce6623..141787b49 100644 --- a/tests/ekobox/test_evol_pdf.py +++ b/tests/ekobox/test_evol_pdf.py @@ -1,4 +1,5 @@ import numpy as np +import pytest from banana import toy import eko @@ -50,6 +51,21 @@ def test_evolve_pdfs_dump_path(fake_lhapdf, cd): assert p.exists() +def test_evolve_pdfs_bad_scales(fake_lhapdf, cd): + """Bad scale configurations.""" + theory, op = init_cards() + n = "test_evolve_pdfs_bad_scales" + op = cards.example.operator() + op.mugrid = [(5.0, 3), (15.0, 4), (10.0, 5), (100.0, 5)] + with pytest.raises(ValueError, match="is bigger"): + with cd(fake_lhapdf): + ev_p.evolve_pdfs([toy.mkPDF("", 0)], theory, op, name=n, path=fake_lhapdf) + op.mugrid = [(5.0, 3), (10.0, 3), (10.0, 4), (15.0, 4), (10.0, 5)] + with pytest.raises(ValueError, match="is bigger"): + with cd(fake_lhapdf): + ev_p.evolve_pdfs([toy.mkPDF("", 0)], theory, op, name=n, path=fake_lhapdf) + + def test_evolve_pdfs_dump_file(fake_lhapdf, cd): theory, op = init_cards() n = "test_evolve_pdfs_dump_file" diff --git a/tests/ekobox/test_info_file.py b/tests/ekobox/test_info_file.py index f7267fbd4..acf3f7272 100644 --- a/tests/ekobox/test_info_file.py +++ b/tests/ekobox/test_info_file.py @@ -1,7 +1,6 @@ import math import numpy as np -import pytest from ekobox import cards, info_file @@ -34,12 +33,12 @@ def test_build_alphas_good(): op = cards.example.operator() op.mu0 = 1.0 # base case - op.mugrid = [(10.0, 5), (100.0, 5)] + op.mugrid = [(100.0, 5), (10.0, 5)] info = info_file.build_alphas(theory, op) assert len(info["AlphaS_Vals"]) == 2 np.testing.assert_allclose(info["AlphaS_Qs"], [10.0, 100.0]) # several nf - op.mugrid = [(1.0, 3), (5.0, 3), (5.0, 4), (10.0, 5), (10.0, 5), (100.0, 5)] + op.mugrid = [(5.0, 4), (10.0, 5), (1.0, 3), (5.0, 3), (10.0, 5), (100.0, 5)] info = info_file.build_alphas(theory, op) assert len(info["AlphaS_Vals"]) == 6 np.testing.assert_allclose(info["AlphaS_Qs"], [1.0, 5.0, 5.0, 10.0, 10.0, 100.0]) @@ -48,18 +47,3 @@ def test_build_alphas_good(): info = info_file.build_alphas(theory, op) assert len(info["AlphaS_Vals"]) == 4 np.testing.assert_allclose(info["AlphaS_Qs"], [1.0, 10.0, 10.0, 100.0]) - - -def test_build_alphas_bad(): - """Bad configurations.""" - theory = cards.example.theory() - theory.order = (2, 0) - theory.couplings.alphas = 0.2 - op = cards.example.operator() - op.mu0 = 1.0 - op.mugrid = [(5.0, 3), (15.0, 4), (10.0, 5), (100.0, 5)] - with pytest.raises(ValueError, match="is bigger"): - info_file.build_alphas(theory, op) - op.mugrid = [(5.0, 3), (10.0, 3), (10.0, 4), (15.0, 4), (10.0, 5)] - with pytest.raises(ValueError, match="is bigger"): - info_file.build_alphas(theory, op) From ac8dbe859d97747b4c3c144f991822987a85d418 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Thu, 18 Jul 2024 17:15:17 +0300 Subject: [PATCH 08/20] rust: Allow abbrev --- crates/doc-header.html | 65 ++++++++++++++++++++++++++++++++++ crates/eko/Cargo.toml | 2 +- crates/eko/doc-header.html | 1 + crates/eko/katex-header.html | 1 - crates/ekore/Cargo.toml | 2 +- crates/ekore/doc-header.html | 1 + crates/ekore/katex-header.html | 1 - crates/katex-header.html | 19 ---------- crates/parse-abbrev.py | 20 +++++++++++ 9 files changed, 89 insertions(+), 23 deletions(-) create mode 100644 crates/doc-header.html create mode 120000 crates/eko/doc-header.html delete mode 120000 crates/eko/katex-header.html create mode 120000 crates/ekore/doc-header.html delete mode 120000 crates/ekore/katex-header.html delete mode 100644 crates/katex-header.html create mode 100644 crates/parse-abbrev.py diff --git a/crates/doc-header.html b/crates/doc-header.html new file mode 100644 index 000000000..3b2f9651b --- /dev/null +++ b/crates/doc-header.html @@ -0,0 +1,65 @@ + + + + diff --git a/crates/eko/Cargo.toml b/crates/eko/Cargo.toml index 0fc9b2ee7..0fffdea1a 100644 --- a/crates/eko/Cargo.toml +++ b/crates/eko/Cargo.toml @@ -13,7 +13,7 @@ rust-version.workspace = true version.workspace = true [package.metadata.docs.rs] -rustdoc-args = ["--html-in-header", "katex-header.html"] +rustdoc-args = ["--html-in-header", "doc-header.html"] [lib] name = "ekors" diff --git a/crates/eko/doc-header.html b/crates/eko/doc-header.html new file mode 120000 index 000000000..22282661d --- /dev/null +++ b/crates/eko/doc-header.html @@ -0,0 +1 @@ +../doc-header.html \ No newline at end of file diff --git a/crates/eko/katex-header.html b/crates/eko/katex-header.html deleted file mode 120000 index 810e68dcb..000000000 --- a/crates/eko/katex-header.html +++ /dev/null @@ -1 +0,0 @@ -../katex-header.html \ No newline at end of file diff --git a/crates/ekore/Cargo.toml b/crates/ekore/Cargo.toml index 86e9784eb..419bc126a 100644 --- a/crates/ekore/Cargo.toml +++ b/crates/ekore/Cargo.toml @@ -13,7 +13,7 @@ rust-version.workspace = true version.workspace = true [package.metadata.docs.rs] -rustdoc-args = ["--html-in-header", "katex-header.html"] +rustdoc-args = ["--html-in-header", "doc-header.html"] # See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html diff --git a/crates/ekore/doc-header.html b/crates/ekore/doc-header.html new file mode 120000 index 000000000..22282661d --- /dev/null +++ b/crates/ekore/doc-header.html @@ -0,0 +1 @@ +../doc-header.html \ No newline at end of file diff --git a/crates/ekore/katex-header.html b/crates/ekore/katex-header.html deleted file mode 120000 index 810e68dcb..000000000 --- a/crates/ekore/katex-header.html +++ /dev/null @@ -1 +0,0 @@ -../katex-header.html \ No newline at end of file diff --git a/crates/katex-header.html b/crates/katex-header.html deleted file mode 100644 index ff018b9c3..000000000 --- a/crates/katex-header.html +++ /dev/null @@ -1,19 +0,0 @@ - - - - diff --git a/crates/parse-abbrev.py b/crates/parse-abbrev.py new file mode 100644 index 000000000..b928ed03c --- /dev/null +++ b/crates/parse-abbrev.py @@ -0,0 +1,20 @@ +"""Parse abbreviations from sphinx to Rust.""" + +import pathlib +import re + +SRC = ( + pathlib.Path(__file__).parents[1] + / "doc" + / "source" + / "shared" + / "abbreviations.rst" +) + +cnt = SRC.read_text("utf-8") +for el in cnt.split(".."): + test = re.match(r"\s*(.+)\s+replace::\n\s+:abbr:`(.+?)\((.+)\)`", el.strip()) + if test is None: + continue + # Print to terminal - the user can dump to the relevant file + print(f'"{test[2].strip()}": "{test[3].strip()}",') From 6b004dba87d494415a71a78535c65fbe32b7c142 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Fri, 19 Jul 2024 14:37:54 +0300 Subject: [PATCH 09/20] Make add SV to ModeMixin name --- src/eko/evolution_operator/__init__.py | 2 +- src/eko/evolution_operator/__init__.py.patch | 8 ++++---- src/eko/evolution_operator/grid.py | 2 +- src/eko/scale_variations/__init__.py | 2 +- 4 files changed, 7 insertions(+), 7 deletions(-) diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py index 1c759c5c5..16aa4bb5b 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -600,7 +600,7 @@ def quad_ker_qed( return ker -class Operator(sv.ModeMixin): +class Operator(sv.ScaleVariationModeMixin): """Internal representation of a single EKO. The actual matrices are computed upon calling :meth:`compute`. diff --git a/src/eko/evolution_operator/__init__.py.patch b/src/eko/evolution_operator/__init__.py.patch index a3bc24ec1..58c3f74fd 100644 --- a/src/eko/evolution_operator/__init__.py.patch +++ b/src/eko/evolution_operator/__init__.py.patch @@ -534,10 +534,10 @@ index 1c759c5c..5eb394d0 100644 - return ker - - - class Operator(sv.ModeMixin): + class Operator(sv.ScaleVariationModeMixin): """Internal representation of a single EKO. -@@ -780,50 +282,6 @@ class Operator(sv.ModeMixin): +@@ -780,50 +282,6 @@ class Operator(sv.ScaleVariationModeMixin): labels.extend(br.singlet_unified_labels) return labels @@ -588,7 +588,7 @@ index 1c759c5c..5eb394d0 100644 def initialize_op_members(self): """Init all operators with the identity or zeros.""" eye = OpMember( -@@ -846,10 +304,7 @@ class Operator(sv.ModeMixin): +@@ -846,10 +304,7 @@ class Operator(sv.ScaleVariationModeMixin): else: self.op_members[n] = zero.copy() @@ -600,7 +600,7 @@ index 1c759c5c..5eb394d0 100644 """Run the integration for each grid point. Parameters -@@ -864,18 +319,56 @@ class Operator(sv.ModeMixin): +@@ -864,18 +319,56 @@ class Operator(sv.ScaleVariationModeMixin): """ column = [] k, logx = log_grid diff --git a/src/eko/evolution_operator/grid.py b/src/eko/evolution_operator/grid.py index e64885418..2539c3fc7 100644 --- a/src/eko/evolution_operator/grid.py +++ b/src/eko/evolution_operator/grid.py @@ -29,7 +29,7 @@ """In particular, only the ``operator`` and ``error`` fields are expected.""" -class OperatorGrid(sv.ModeMixin): +class OperatorGrid(sv.ScaleVariationModeMixin): """Collection of evolution operators for several scales. The operator grid is the driver class of the evolution. diff --git a/src/eko/scale_variations/__init__.py b/src/eko/scale_variations/__init__.py index 4af80592b..9c43f1f9c 100644 --- a/src/eko/scale_variations/__init__.py +++ b/src/eko/scale_variations/__init__.py @@ -36,7 +36,7 @@ def sv_mode(s): return Modes.unvaried -class ModeMixin: +class ScaleVariationModeMixin: """Mixin to cast scale variation mode.""" @property From 7a6c1cc2c42bd63a42ad2e4bae4d362cd4daafb7 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Fri, 19 Jul 2024 15:15:12 +0300 Subject: [PATCH 10/20] Introduce EvoMethods --- src/eko/evolution_operator/__init__.py | 6 +++++ src/eko/kernels/__init__.py | 32 ++++++++++++++++++++++++++ src/eko/scale_variations/__init__.py | 13 ++++++----- tests/eko/kernels/test_init.py | 21 +++++++++++++++++ 4 files changed, 66 insertions(+), 6 deletions(-) create mode 100644 tests/eko/kernels/test_init.py diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py index 16aa4bb5b..8cd5829bd 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -20,6 +20,7 @@ from .. import basis_rotation as br from .. import interpolation, mellin from .. import scale_variations as sv +from ..kernels import ev_method from ..kernels import non_singlet as ns from ..kernels import non_singlet_qed as qed_ns from ..kernels import singlet as s @@ -780,6 +781,11 @@ def labels(self): labels.extend(br.singlet_unified_labels) return labels + @property + def ev_method(self): + """Return the evolution method.""" + return ev_method(self.config["method"]) + def quad_ker(self, label, logx, areas): """Return partially initialized integrand function. diff --git a/src/eko/kernels/__init__.py b/src/eko/kernels/__init__.py index 7640727b5..a10998b23 100644 --- a/src/eko/kernels/__init__.py +++ b/src/eko/kernels/__init__.py @@ -1 +1,33 @@ """The solutions to the |DGLAP| equations.""" + +import enum + + +class EvoMethods(enum.IntEnum): + """Enumerate evolution methods.""" + + ITERATE_EXACT = enum.auto() + ITERATE_EXPANDED = enum.auto() + PERTURBATIVE_EXACT = enum.auto() + PERTURBATIVE_EXPANDED = enum.auto() + TRUNCATED = enum.auto() + ORDERED_TRUNCATED = enum.auto() + DECOMPOSE_EXACT = enum.auto() + DECOMPOSE_EXPANDED = enum.auto() + + +def ev_method(s: EvoMethods) -> EvoMethods: + """Return the evolution methods. + + Parameters + ---------- + s : + string representation + + Returns + ------- + i : + int representation + + """ + return EvoMethods[s.value.upper().replace("-", "_")] diff --git a/src/eko/scale_variations/__init__.py b/src/eko/scale_variations/__init__.py index 9c43f1f9c..21d56a503 100644 --- a/src/eko/scale_variations/__init__.py +++ b/src/eko/scale_variations/__init__.py @@ -6,29 +6,30 @@ import enum +from ..io.types import ScaleVariationsMethod from . import expanded, exponentiated class Modes(enum.IntEnum): - """Enumerate scale Variation modes.""" + """Enumerate scale variation modes.""" unvaried = enum.auto() exponentiated = enum.auto() expanded = enum.auto() -def sv_mode(s): +def sv_mode(s: ScaleVariationsMethod) -> Modes: """Return the scale variation mode. Parameters ---------- - s : str + s : string representation Returns ------- - enum.IntEnum - enum representation + i : + int representation """ if s is not None: @@ -40,6 +41,6 @@ class ScaleVariationModeMixin: """Mixin to cast scale variation mode.""" @property - def sv_mode(self): + def sv_mode(self) -> Modes: """Return the scale variation mode.""" return sv_mode(self.config["ModSV"]) diff --git a/tests/eko/kernels/test_init.py b/tests/eko/kernels/test_init.py new file mode 100644 index 000000000..3717de275 --- /dev/null +++ b/tests/eko/kernels/test_init.py @@ -0,0 +1,21 @@ +from eko.io.types import EvolutionMethod +from eko.kernels import EvoMethods, ev_method + + +def test_ev_method(): + methods = { + "iterate-expanded": EvoMethods.ITERATE_EXPANDED, + "decompose-expanded": EvoMethods.DECOMPOSE_EXPANDED, + "perturbative-expanded": EvoMethods.PERTURBATIVE_EXPANDED, + "truncated": EvoMethods.TRUNCATED, + "ordered-truncated": EvoMethods.ORDERED_TRUNCATED, + "iterate-exact": EvoMethods.ITERATE_EXACT, + "decompose-exact": EvoMethods.DECOMPOSE_EXACT, + "perturbative-exact": EvoMethods.PERTURBATIVE_EXACT, + } + assert len(methods.keys()) == len(EvolutionMethod) + assert len(methods.keys()) == len(EvoMethods) + for s, i in methods.items(): + j = ev_method(EvolutionMethod(s)) + assert j == i + assert isinstance(j, int) From 5636b4a7e2e14f83ca0f5e34a24f71b53042e026 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Fri, 19 Jul 2024 15:45:25 +0300 Subject: [PATCH 11/20] Use EvoMethods in kernels --- src/eko/kernels/non_singlet.py | 27 ++++++------ src/eko/kernels/non_singlet_qed.py | 4 +- src/eko/kernels/singlet.py | 19 +++++---- src/eko/kernels/singlet_qed.py | 7 ++-- src/eko/kernels/valence_qed.py | 5 ++- tests/eko/kernels/test_kernels_QEDns.py | 3 +- tests/eko/kernels/test_kernels_QEDsinglet.py | 12 +++++- tests/eko/kernels/test_kernels_QEDvalence.py | 16 +++++-- tests/eko/kernels/test_ns.py | 32 ++++++-------- tests/eko/kernels/test_s.py | 44 ++++++++------------ 10 files changed, 87 insertions(+), 82 deletions(-) diff --git a/src/eko/kernels/non_singlet.py b/src/eko/kernels/non_singlet.py index bd5bd10f0..15b913438 100644 --- a/src/eko/kernels/non_singlet.py +++ b/src/eko/kernels/non_singlet.py @@ -4,6 +4,7 @@ import numpy as np from .. import beta +from . import EvoMethods from . import as4_evolution_integrals as as4_ei from . import evolution_integrals as ei @@ -355,7 +356,7 @@ def dispatcher( ---------- order : tuple(int,int) perturbation order - method : str + method : int method gamma_ns : numpy.ndarray non-singlet anomalous dimensions @@ -378,38 +379,38 @@ def dispatcher( # use always exact in LO if order[0] == 1: return lo_exact(gamma_ns, a1, a0, betalist) - if method == "ordered-truncated": + if method is EvoMethods.ORDERED_TRUNCATED: return eko_ordered_truncated( gamma_ns, a1, a0, betalist, order, ev_op_iterations ) - if method == "truncated": + if method is EvoMethods.TRUNCATED: return eko_truncated(gamma_ns, a1, a0, betalist, order, ev_op_iterations) # NLO if order[0] == 2: if method in [ - "iterate-expanded", - "decompose-expanded", - "perturbative-expanded", + EvoMethods.ITERATE_EXPANDED, + EvoMethods.DECOMPOSE_EXPANDED, + EvoMethods.PERTURBATIVE_EXPANDED, ]: return nlo_expanded(gamma_ns, a1, a0, betalist) - # if method in ["iterate-exact", "decompose-exact", "perturbative-exact"]: + # exact is left return nlo_exact(gamma_ns, a1, a0, betalist) # NNLO if order[0] == 3: if method in [ - "iterate-expanded", - "decompose-expanded", - "perturbative-expanded", + EvoMethods.ITERATE_EXPANDED, + EvoMethods.DECOMPOSE_EXPANDED, + EvoMethods.PERTURBATIVE_EXPANDED, ]: return nnlo_expanded(gamma_ns, a1, a0, betalist) return nnlo_exact(gamma_ns, a1, a0, betalist) # N3LO if order[0] == 4: if method in [ - "iterate-expanded", - "decompose-expanded", - "perturbative-expanded", + EvoMethods.ITERATE_EXPANDED, + EvoMethods.DECOMPOSE_EXPANDED, + EvoMethods.PERTURBATIVE_EXPANDED, ]: return n3lo_expanded(gamma_ns, a1, a0, betalist) return n3lo_exact(gamma_ns, a1, a0, betalist) diff --git a/src/eko/kernels/non_singlet_qed.py b/src/eko/kernels/non_singlet_qed.py index 7a2de816a..d81d5076b 100644 --- a/src/eko/kernels/non_singlet_qed.py +++ b/src/eko/kernels/non_singlet_qed.py @@ -146,7 +146,7 @@ def as4_exact(gamma_ns, a1, a0, beta): @nb.njit(cache=True) def dispatcher( order, - method, + _method, gamma_ns, as_list, aem_half, @@ -164,7 +164,7 @@ def dispatcher( ---------- order : tuple(int,int) perturbation order - method : str + method : int method gamma_ns : numpy.ndarray non-singlet anomalous dimensions diff --git a/src/eko/kernels/singlet.py b/src/eko/kernels/singlet.py index 5983c0c67..e637af046 100644 --- a/src/eko/kernels/singlet.py +++ b/src/eko/kernels/singlet.py @@ -6,6 +6,7 @@ from ekore import anomalous_dimensions as ad from .. import beta +from . import EvoMethods from . import as4_evolution_integrals as as4_ei from . import evolution_integrals as ei @@ -579,7 +580,7 @@ def dispatcher( # pylint: disable=too-many-return-statements ---------- order : tuple(int,int) perturbative order - method : str + method : int method gamma_singlet : numpy.ndarray singlet anomalous dimensions matrices @@ -609,9 +610,9 @@ def dispatcher( # pylint: disable=too-many-return-statements return lo_exact(gamma_singlet, a1, a0, betalist) # Common method for NLO and NNLO - if method in ["iterate-exact", "iterate-expanded"]: + if method in [EvoMethods.ITERATE_EXACT, EvoMethods.ITERATE_EXPANDED]: return eko_iterate(gamma_singlet, a1, a0, betalist, order, ev_op_iterations) - if method == "perturbative-exact": + if method is EvoMethods.PERTURBATIVE_EXACT: return eko_perturbative( gamma_singlet, a1, @@ -622,7 +623,7 @@ def dispatcher( # pylint: disable=too-many-return-statements ev_op_max_order, True, ) - if method == "perturbative-expanded": + if method is EvoMethods.PERTURBATIVE_EXPANDED: return eko_perturbative( gamma_singlet, a1, @@ -633,19 +634,19 @@ def dispatcher( # pylint: disable=too-many-return-statements ev_op_max_order, False, ) - if method in ["truncated", "ordered-truncated"]: + if method in [EvoMethods.TRUNCATED, EvoMethods.ORDERED_TRUNCATED]: return eko_truncated(gamma_singlet, a1, a0, betalist, order, ev_op_iterations) # These methods are scattered for nlo and nnlo - if method == "decompose-exact": + if method is EvoMethods.DECOMPOSE_EXACT: if order[0] == 2: return nlo_decompose_exact(gamma_singlet, a1, a0, betalist) - elif order[0] == 3: + if order[0] == 3: return nnlo_decompose_exact(gamma_singlet, a1, a0, betalist) return n3lo_decompose_exact(gamma_singlet, a1, a0, nf) - if method == "decompose-expanded": + if method is EvoMethods.DECOMPOSE_EXPANDED: if order[0] == 2: return nlo_decompose_expanded(gamma_singlet, a1, a0, betalist) - elif order[0] == 3: + if order[0] == 3: return nnlo_decompose_expanded(gamma_singlet, a1, a0, betalist) return n3lo_decompose_expanded(gamma_singlet, a1, a0, nf) raise NotImplementedError("Selected method is not implemented") diff --git a/src/eko/kernels/singlet_qed.py b/src/eko/kernels/singlet_qed.py index c13dee95c..2f2ba8eee 100644 --- a/src/eko/kernels/singlet_qed.py +++ b/src/eko/kernels/singlet_qed.py @@ -6,6 +6,7 @@ from ekore import anomalous_dimensions as ad from .. import beta +from . import EvoMethods @nb.njit(cache=True) @@ -66,7 +67,7 @@ def dispatcher( a_half, nf, ev_op_iterations, - ev_op_max_order, + _ev_op_max_order, ): """Determine used kernel and call it. @@ -74,7 +75,7 @@ def dispatcher( ---------- order : tuple(int,int) perturbative order - method : str + method : int method gamma_singlet : numpy.ndarray singlet anomalous dimensions matrices @@ -96,7 +97,7 @@ def dispatcher( e_s : numpy.ndarray singlet EKO """ - if method in ["iterate-exact", "iterate-expanded"]: + if method in [EvoMethods.ITERATE_EXACT, EvoMethods.ITERATE_EXPANDED]: return eko_iterate( gamma_singlet, as_list, a_half, nf, order, ev_op_iterations, 4 ) diff --git a/src/eko/kernels/valence_qed.py b/src/eko/kernels/valence_qed.py index 0d5b747c3..074dd95ed 100644 --- a/src/eko/kernels/valence_qed.py +++ b/src/eko/kernels/valence_qed.py @@ -2,6 +2,7 @@ import numba as nb +from . import EvoMethods from .singlet_qed import eko_iterate @@ -14,7 +15,7 @@ def dispatcher( a_half, nf, ev_op_iterations, - ev_op_max_order, + _ev_op_max_order, ): """ Determine used kernel and call it. @@ -45,7 +46,7 @@ def dispatcher( e_v : numpy.ndarray singlet EKO """ - if method in ["iterate-exact", "iterate-expanded"]: + if method in [EvoMethods.ITERATE_EXACT, EvoMethods.ITERATE_EXPANDED]: return eko_iterate( gamma_valence, as_list, a_half, nf, order, ev_op_iterations, 2 ) diff --git a/tests/eko/kernels/test_kernels_QEDns.py b/tests/eko/kernels/test_kernels_QEDns.py index bb18c3e17..9fce3097b 100644 --- a/tests/eko/kernels/test_kernels_QEDns.py +++ b/tests/eko/kernels/test_kernels_QEDns.py @@ -3,6 +3,7 @@ from eko import basis_rotation as br from eko.couplings import Couplings +from eko.kernels import EvoMethods from eko.kernels import non_singlet_qed as ns from eko.quantities.couplings import CouplingEvolutionMethod, CouplingsInfo from eko.quantities.heavy_quarks import QuarkMassScheme @@ -14,7 +15,7 @@ # "perturbative-expanded", # "truncated", # "ordered-truncated", - "iterate-exact", + EvoMethods.ITERATE_EXACT, # "decompose-exact", # "perturbative-exact", ] diff --git a/tests/eko/kernels/test_kernels_QEDsinglet.py b/tests/eko/kernels/test_kernels_QEDsinglet.py index b953c053c..2fee4e17f 100644 --- a/tests/eko/kernels/test_kernels_QEDsinglet.py +++ b/tests/eko/kernels/test_kernels_QEDsinglet.py @@ -1,6 +1,7 @@ import numpy as np import pytest +from eko.kernels import EvoMethods from eko.kernels import singlet_qed as s from ekore.anomalous_dimensions.unpolarized import space_like as ad @@ -10,7 +11,7 @@ # "perturbative-expanded", # "truncated", # "ordered-truncated", - "iterate-exact", + EvoMethods.ITERATE_EXACT, # "decompose-exact", # "perturbative-exact", ] @@ -125,5 +126,12 @@ def test_zero_true_gamma(monkeypatch): def test_error(): with pytest.raises(NotImplementedError): s.dispatcher( - (3, 2), "AAA", np.random.rand(4, 3, 2, 2), [0.2, 0.1], [0.01], 3, 10, 10 + (3, 2), + "iterate-exact", + np.random.rand(4, 3, 2, 2), + [0.2, 0.1], + [0.01], + 3, + 10, + 10, ) diff --git a/tests/eko/kernels/test_kernels_QEDvalence.py b/tests/eko/kernels/test_kernels_QEDvalence.py index 29f0f3d7a..4b12a492b 100644 --- a/tests/eko/kernels/test_kernels_QEDvalence.py +++ b/tests/eko/kernels/test_kernels_QEDvalence.py @@ -1,6 +1,7 @@ import numpy as np import pytest +from eko.kernels import EvoMethods from eko.kernels import valence_qed as val from ekore import anomalous_dimensions @@ -10,13 +11,13 @@ # "perturbative-expanded", # "truncated", # "ordered-truncated", - "iterate-exact", + EvoMethods.ITERATE_EXACT, # "decompose-exact", # "perturbative-exact", ] -def test_zero(monkeypatch): +def test_zero(): """No evolution results in exp(0)""" nf = 3 ev_op_iterations = 2 @@ -57,7 +58,7 @@ def test_zero(monkeypatch): ) -def test_zero_true_gamma(monkeypatch): +def test_zero_true_gamma(): """No evolution results in exp(0)""" nf = 3 ev_op_iterations = 2 @@ -101,5 +102,12 @@ def test_zero_true_gamma(monkeypatch): def test_error(): with pytest.raises(NotImplementedError): val.dispatcher( - (3, 2), "AAA", np.random.rand(4, 3, 2, 2), [0.2, 0.1], [0.01], 3, 10, 10 + (3, 2), + "iterate-exact", + np.random.rand(4, 3, 2, 2), + [0.2, 0.1], + [0.01], + 3, + 10, + 10, ) diff --git a/tests/eko/kernels/test_ns.py b/tests/eko/kernels/test_ns.py index f78282142..70e7faf8d 100644 --- a/tests/eko/kernels/test_ns.py +++ b/tests/eko/kernels/test_ns.py @@ -5,19 +5,9 @@ import ekore.anomalous_dimensions.unpolarized.space_like as ad from eko import beta +from eko.kernels import EvoMethods from eko.kernels import non_singlet as ns -methods = [ - "iterate-expanded", - "decompose-expanded", - "perturbative-expanded", - "truncated", - "ordered-truncated", - "iterate-exact", - "decompose-exact", - "perturbative-exact", -] - def test_zero(): """No evolution results in exp(0)""" @@ -25,7 +15,7 @@ def test_zero(): ev_op_iterations = 2 gamma_ns = np.array([1 + 0.0j, 1 + 0j, 1 + 0j, 1 + 0j]) for order in [1, 2, 3, 4]: - for method in methods: + for method in EvoMethods: np.testing.assert_allclose( ns.dispatcher( (order, 0), method, gamma_ns, 1.0, 1.0, nf, ev_op_iterations @@ -54,7 +44,7 @@ def test_ode_lo(): a0 = 0.3 for a1 in [0.1, 0.2]: r = a1 * gamma_ns / (beta.beta_qcd((2, 0), nf) * a1**2) - for method in methods: + for method in EvoMethods: rhs = r * ns.dispatcher( (1, 0), method, gamma_ns, a1, a0, nf, ev_op_iterations ) @@ -91,7 +81,7 @@ def test_ode_nlo(): r = (a1 * gamma_ns[0] + a1**2 * gamma_ns[1]) / ( beta.beta_qcd((2, 0), nf) * a1**2 + beta.beta_qcd((3, 0), nf) * a1**3 ) - for method in ["iterate-exact"]: + for method in [EvoMethods.ITERATE_EXACT]: rhs = r * ns.dispatcher( (2, 0), method, gamma_ns, a1, a0, nf, ev_op_iterations ) @@ -130,7 +120,7 @@ def test_ode_nnlo(): + beta.beta_qcd((3, 0), nf) * a1**2 + beta.beta_qcd((4, 0), nf) * a1**3 ) - for method in ["iterate-exact"]: + for method in [EvoMethods.ITERATE_EXACT]: rhs = r * ns.dispatcher( (3, 0), method, gamma_ns, a1, a0, nf, ev_op_iterations ) @@ -172,7 +162,7 @@ def test_ode_n3lo(): + beta.beta_qcd((4, 0), nf) * a1**3 + beta.beta_qcd((5, 0), nf) * a1**4 ) - for method in ["iterate-exact"]: + for method in [EvoMethods.ITERATE_EXACT]: rhs = r * ns.dispatcher( (4, 0), method, gamma_ns, a1, a0, nf, ev_op_iterations ) @@ -202,7 +192,9 @@ def test_ode_n3lo(): def test_error(monkeypatch): monkeypatch.setattr("eko.beta.beta_qcd", lambda *_args: 1.0) with pytest.raises(NotImplementedError, match="order is not implemented"): - ns.dispatcher((5, 0), "iterate-exact", np.random.rand(3) + 0j, 0.2, 0.1, 3, 10) + ns.dispatcher( + (5, 0), EvoMethods.ITERATE_EXACT, np.random.rand(3) + 0j, 0.2, 0.1, 3, 10 + ) with pytest.raises(NotImplementedError): ad.gamma_ns((2, 0), 10202, 1, (0, 0, 0, 0, 0, 0, 0), 3) @@ -216,7 +208,7 @@ def test_gamma_usage(): gamma_ns = np.full(4, np.nan) for order in range(1, 5): gamma_ns[order - 1] = np.random.rand() - for method in methods: + for method in EvoMethods: r = ns.dispatcher( (order, 0), method, gamma_ns, a1, a0, nf, ev_op_iterations ) @@ -225,8 +217,8 @@ def test_gamma_usage(): for order in range(1, 5): gamma_ns = np.random.rand(order) gamma_ns[order - 1] = np.nan - for method in methods: - if method == "ordered-truncated": + for method in EvoMethods: + if method is EvoMethods.ORDERED_TRUNCATED: # we are actually dividing by a np.nan, # since the sum of U vec is nan warnings.simplefilter("ignore", RuntimeWarning) diff --git a/tests/eko/kernels/test_s.py b/tests/eko/kernels/test_s.py index 5380ba319..ff4fdfc0d 100644 --- a/tests/eko/kernels/test_s.py +++ b/tests/eko/kernels/test_s.py @@ -3,20 +3,10 @@ import numpy as np import pytest +from eko.kernels import EvoMethods from eko.kernels import singlet as s from ekore import anomalous_dimensions as ad -methods = [ - "iterate-expanded", - "decompose-expanded", - "perturbative-expanded", - "truncated", - "ordered-truncated", - "iterate-exact", - "decompose-exact", - "perturbative-exact", -] - def test_zero_lo(monkeypatch): """No evolution results in exp(0)""" @@ -35,7 +25,7 @@ def test_zero_lo(monkeypatch): np.array([[0, 0], [0, 1]]), ), ) - for method in methods: + for method in EvoMethods: np.testing.assert_allclose( s.dispatcher( (1, 0), method, gamma_s, 1, 1, nf, ev_op_iterations, ev_op_max_order @@ -75,8 +65,8 @@ def test_zero_nlo_decompose(monkeypatch): ), ) for method in [ - "decompose-expanded", - "decompose-exact", + EvoMethods.DECOMPOSE_EXPANDED, + EvoMethods.DECOMPOSE_EXACT, ]: np.testing.assert_allclose( s.dispatcher( @@ -117,8 +107,8 @@ def test_zero_nnlo_decompose(monkeypatch): ), ) for method in [ - "decompose-expanded", - "decompose-exact", + EvoMethods.DECOMPOSE_EXPANDED, + EvoMethods.DECOMPOSE_EXACT, ]: np.testing.assert_allclose( s.dispatcher( @@ -159,8 +149,8 @@ def test_zero_n3lo_decompose(monkeypatch): ), ) for method in [ - "decompose-expanded", - "decompose-exact", + EvoMethods.DECOMPOSE_EXPANDED, + EvoMethods.DECOMPOSE_EXACT, ]: np.testing.assert_allclose( s.dispatcher( @@ -200,7 +190,7 @@ def test_similarity(): ]: ref = s.dispatcher( order, - "decompose-exact", + EvoMethods.DECOMPOSE_EXACT, gamma_s, a1, a0, @@ -208,7 +198,7 @@ def test_similarity(): ev_op_iterations, ev_op_max_order, ) - for method in methods: + for method in EvoMethods: np.testing.assert_allclose( s.dispatcher( order, @@ -227,7 +217,9 @@ def test_similarity(): def test_error(): with pytest.raises(NotImplementedError): - s.dispatcher((4, 0), "AAA", np.random.rand(3, 2, 2), 0.2, 0.1, 3, 10, 10) + s.dispatcher( + (4, 0), "iterate-exact", np.random.rand(3, 2, 2), 0.2, 0.1, 3, 10, 10 + ) def mk_almost_diag_matrix(n, max_ang=np.pi / 8.0): @@ -247,7 +239,7 @@ def test_gamma_usage(): gamma_s = np.full((4, 2, 2), np.nan) for order in range(1, 5): gamma_s[order - 1] = mk_almost_diag_matrix(1) - for method in methods: + for method in EvoMethods: r = s.dispatcher( (order, 0), method, @@ -263,8 +255,8 @@ def test_gamma_usage(): for order in range(1, 5): gamma_s = mk_almost_diag_matrix(4) gamma_s[order - 1] = np.full((2, 2), np.nan) - for method in methods: - if "iterate" in method: + for method in EvoMethods: + if method in [EvoMethods.ITERATE_EXACT, EvoMethods.ITERATE_EXPANDED]: # we are actually dividing by the determinant of # matrix full of np.nan warnings.simplefilter("ignore", RuntimeWarning) @@ -287,8 +279,8 @@ def test_singlet_back(): nf = 4 a1 = 3.0 a0 = 4.0 - s10 = s.dispatcher(order, "iterate-exact", gamma_s, a1, a0, nf, 15, 1) + s10 = s.dispatcher(order, EvoMethods.ITERATE_EXACT, gamma_s, a1, a0, nf, 15, 1) np.testing.assert_allclose( np.linalg.inv(s10), - s.dispatcher(order, "iterate-exact", gamma_s, a0, a1, nf, 15, 1), + s.dispatcher(order, EvoMethods.ITERATE_EXACT, gamma_s, a0, a1, nf, 15, 1), ) From db2ef3ed63894e806b78e4c0ff41281179a3ac64 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Fri, 19 Jul 2024 16:03:51 +0300 Subject: [PATCH 12/20] Use EvoMethods in ev_op --- src/eko/evolution_operator/__init__.py | 5 ++- tests/eko/evolution_operator/test_init.py | 43 ++++++++++++++------- tests/eko/scale_variations/test_diff.py | 4 +- tests/eko/scale_variations/test_expanded.py | 4 +- 4 files changed, 36 insertions(+), 20 deletions(-) diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py index 8cd5829bd..ec56b6db6 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -20,6 +20,7 @@ from .. import basis_rotation as br from .. import interpolation, mellin from .. import scale_variations as sv +from ..io.types import EvolutionMethod from ..kernels import ev_method from ..kernels import non_singlet as ns from ..kernels import non_singlet_qed as qed_ns @@ -784,7 +785,7 @@ def labels(self): @property def ev_method(self): """Return the evolution method.""" - return ev_method(self.config["method"]) + return ev_method(EvolutionMethod(self.config["method"])) def quad_ker(self, label, logx, areas): """Return partially initialized integrand function. @@ -809,7 +810,7 @@ def quad_ker(self, label, logx, areas): order=self.order, mode0=label[0], mode1=label[1], - method=self.config["method"], + method=self.ev_method, is_log=self.int_disp.log, logx=logx, areas=areas, diff --git a/tests/eko/evolution_operator/test_init.py b/tests/eko/evolution_operator/test_init.py index 597f645f8..23395f66b 100644 --- a/tests/eko/evolution_operator/test_init.py +++ b/tests/eko/evolution_operator/test_init.py @@ -11,6 +11,7 @@ from eko.evolution_operator import Operator, quad_ker from eko.interpolation import InterpolatorDispatcher from eko.io.runcards import OperatorCard, ScaleVariationsMethod, TheoryCard +from eko.kernels import EvoMethods from eko.kernels import non_singlet as ns from eko.kernels import non_singlet_qed as qed_ns from eko.kernels import singlet as s @@ -25,7 +26,7 @@ def test_quad_ker_errors(): order=(1, 0), mode0=mode0, mode1=0, - method="", + method="iterate-exact", is_log=True, logx=np.log(0.1), areas=[[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]], @@ -60,15 +61,29 @@ def test_quad_ker(monkeypatch): monkeypatch.setattr(qed_ns, "dispatcher", lambda *args: 1.0) monkeypatch.setattr(s, "dispatcher", lambda *args: np.identity(2)) params = [ - ((1, 0), br.non_singlet_pids_map["ns+"], 0, "", 0.0, 0.0), - ((1, 0), br.non_singlet_pids_map["ns+"], 0, "", 0.123, 1.0), - ((3, 1), br.non_singlet_pids_map["ns+u"], 0, "", 0.0, 0.0), - ((1, 0), 100, 100, "", 0.123, 1.0), - ((1, 0), 100, 21, "", 0.0, 0.0), - ((1, 1), 100, 100, "iterate-exact", 0.123, 1.0), - ((1, 1), 100, 21, "iterate-exact", 0.123, 0.0), - ((1, 1), 10200, 10200, "iterate-exact", 0.123, 1.0), - ((1, 1), 10200, 10204, "iterate-exact", 0.123, 0.0), + ((1, 0), br.non_singlet_pids_map["ns+"], 0, EvoMethods.ITERATE_EXACT, 0.0, 0.0), + ( + (1, 0), + br.non_singlet_pids_map["ns+"], + 0, + EvoMethods.ITERATE_EXACT, + 0.123, + 1.0, + ), + ( + (3, 1), + br.non_singlet_pids_map["ns+u"], + 0, + EvoMethods.ITERATE_EXACT, + 0.0, + 0.0, + ), + ((1, 0), 100, 100, EvoMethods.ITERATE_EXACT, 0.123, 1.0), + ((1, 0), 100, 21, EvoMethods.ITERATE_EXACT, 0.0, 0.0), + ((1, 1), 100, 100, EvoMethods.ITERATE_EXACT, 0.123, 1.0), + ((1, 1), 100, 21, EvoMethods.ITERATE_EXACT, 0.123, 0.0), + ((1, 1), 10200, 10200, EvoMethods.ITERATE_EXACT, 0.123, 1.0), + ((1, 1), 10200, 10204, EvoMethods.ITERATE_EXACT, 0.123, 0.0), ] for order, mode0, mode1, method, logx, res in params: for is_log in [True, False]: @@ -107,7 +122,7 @@ def test_quad_ker(monkeypatch): order=(1, 0), mode0=label[0], mode1=label[1], - method="", + method=EvoMethods.ITERATE_EXACT, is_log=True, logx=0.123, areas=np.zeros(3), @@ -143,7 +158,7 @@ def test_quad_ker(monkeypatch): order=(1, 1), mode0=label[0], mode1=label[1], - method="iterate-exact", + method=EvoMethods.ITERATE_EXACT, is_log=True, logx=0.123, areas=np.zeros(3), @@ -171,7 +186,7 @@ def test_quad_ker(monkeypatch): order=(1, 0), mode0=br.non_singlet_pids_map["ns+"], mode1=0, - method="", + method=EvoMethods.ITERATE_EXACT, is_log=True, logx=0.0, areas=np.zeros(3), @@ -436,7 +451,7 @@ def quad_ker_pegasus( order = (2, 0) mode0 = br.non_singlet_pids_map["ns+"] mode1 = 0 - method = "" + method = EvoMethods.ITERATE_EXACT logxs = np.log(int_disp.xgrid.raw) a1 = 1 a0 = 2 diff --git a/tests/eko/scale_variations/test_diff.py b/tests/eko/scale_variations/test_diff.py index c0dd67236..8ade8bcac 100644 --- a/tests/eko/scale_variations/test_diff.py +++ b/tests/eko/scale_variations/test_diff.py @@ -11,7 +11,7 @@ from eko import basis_rotation as br from eko.beta import beta_qcd_as2, beta_qcd_as3 from eko.couplings import CouplingEvolutionMethod, Couplings, CouplingsInfo -from eko.kernels import non_singlet, singlet +from eko.kernels import EvoMethods, non_singlet, singlet from eko.quantities.heavy_quarks import QuarkMassScheme from eko.scale_variations import expanded, exponentiated from ekore.anomalous_dimensions.unpolarized.space_like import gamma_ns, gamma_singlet @@ -19,7 +19,7 @@ NF = 4 Q02 = 1.65**2 Q12 = 100**2 -EV_METHOD = "truncated" +EV_METHOD = EvoMethods.TRUNCATED def compute_a_s(q2, order): diff --git a/tests/eko/scale_variations/test_expanded.py b/tests/eko/scale_variations/test_expanded.py index 6af22dd51..bfacc1c96 100644 --- a/tests/eko/scale_variations/test_expanded.py +++ b/tests/eko/scale_variations/test_expanded.py @@ -2,7 +2,7 @@ from eko import basis_rotation as br from eko.couplings import CouplingEvolutionMethod, Couplings, CouplingsInfo -from eko.kernels import non_singlet, singlet +from eko.kernels import EvoMethods, non_singlet, singlet from eko.quantities.heavy_quarks import QuarkMassScheme from eko.scale_variations import Modes, expanded from ekore.anomalous_dimensions.unpolarized.space_like import gamma_ns, gamma_singlet @@ -10,7 +10,7 @@ NF = 4 Q02 = 1.65**2 Q12 = 100**2 -EV_METHOD = "truncated" +EV_METHOD = EvoMethods.TRUNCATED def compute_a_s(q2, order): From 1b1e4ce4b651033fc1c23c89cb79f42a03019c66 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Fri, 19 Jul 2024 16:41:12 +0300 Subject: [PATCH 13/20] Introduce BackwardsMethods --- .../operator_matrix_element.py | 50 ++++++++++++++++--- src/eko/kernels/__init__.py | 6 ++- tests/eko/evolution_operator/test_ome.py | 19 +++---- 3 files changed, 57 insertions(+), 18 deletions(-) diff --git a/src/eko/evolution_operator/operator_matrix_element.py b/src/eko/evolution_operator/operator_matrix_element.py index 202fd18d9..b4b6731b4 100644 --- a/src/eko/evolution_operator/operator_matrix_element.py +++ b/src/eko/evolution_operator/operator_matrix_element.py @@ -1,6 +1,7 @@ """The |OME| for the non-trivial matching conditions in the |VFNS| evolution.""" import copy +import enum import functools import logging @@ -20,6 +21,33 @@ logger = logging.getLogger(__name__) +class BackwardsMethods(enum.IntEnum): + """Enumerate backward methods.""" + + FORWARD = enum.auto() + EXACT = enum.auto() + EXPANDED = enum.auto() + + +def bw_method(s: InversionMethod) -> BackwardsMethods: + """Return the backward method. + + Parameters + ---------- + s : + string representation + + Returns + ------- + i : + int representation + + """ + if s is not None: + return BackwardsMethods[s.value.upper()] + return BackwardsMethods.FORWARD + + @nb.njit(cache=True) def build_ome(A, matching_order, a_s, backward_method): r"""Construct the matching expansion in :math:`a_s` with the appropriate method. @@ -32,7 +60,7 @@ def build_ome(A, matching_order, a_s, backward_method): perturbation matching order a_s : float strong coupling, needed only for the exact inverse - backward_method : InversionMethod or None + backward_method : BackwardsMethods empty or method for inverting the matching condition (exact or expanded) Returns @@ -51,7 +79,7 @@ def build_ome(A, matching_order, a_s, backward_method): ome = np.eye(len(A[0]), dtype=np.complex_) A = A[:, :, :] A = np.ascontiguousarray(A) - if backward_method is InversionMethod.EXPANDED: + if backward_method is BackwardsMethods.EXPANDED: # expended inverse if matching_order[0] >= 1: ome -= a_s * A[0] @@ -68,7 +96,7 @@ def build_ome(A, matching_order, a_s, backward_method): if matching_order[0] >= 3: ome += a_s**3 * A[2] # need inverse exact ? so add the missing pieces - if backward_method is InversionMethod.EXACT: + if backward_method is BackwardsMethods.EXACT: ome = np.linalg.inv(ome) return ome @@ -216,7 +244,9 @@ class OperatorMatrixElement(Operator): def __init__(self, config, managers, nf, q2, is_backward, L, is_msbar): super().__init__(config, managers, Segment(q2, q2, nf)) - self.backward_method = config["backward_inversion"] if is_backward else None + self.backward_method = bw_method( + config["backward_inversion"] if is_backward else None + ) if is_backward: self.is_intrinsic = True else: @@ -241,7 +271,10 @@ def labels(self): logger.warning("%s: skipping non-singlet sector", self.log_label) else: labels.append((200, 200)) - if self.is_intrinsic or self.backward_method is not None: + if ( + self.is_intrinsic + or self.backward_method is not BackwardsMethods.FORWARD + ): # intrinsic labels, which are not zero at NLO labels.append((br.matching_hminus_pid, br.matching_hminus_pid)) # These contributions are always 0 for the moment @@ -257,7 +290,10 @@ def labels(self): (br.matching_hplus_pid, 100), ] ) - if self.is_intrinsic or self.backward_method is not None: + if ( + self.is_intrinsic + or self.backward_method is not BackwardsMethods.FORWARD + ): labels.extend( [ (21, br.matching_hplus_pid), @@ -329,7 +365,7 @@ def compute(self): self.log_label, self.order[0], self.order[1], - self.backward_method, + BackwardsMethods(self.backward_method).name, ) self.integrate() diff --git a/src/eko/kernels/__init__.py b/src/eko/kernels/__init__.py index a10998b23..594fbf1fa 100644 --- a/src/eko/kernels/__init__.py +++ b/src/eko/kernels/__init__.py @@ -2,6 +2,8 @@ import enum +from ..io.types import EvolutionMethod + class EvoMethods(enum.IntEnum): """Enumerate evolution methods.""" @@ -16,8 +18,8 @@ class EvoMethods(enum.IntEnum): DECOMPOSE_EXPANDED = enum.auto() -def ev_method(s: EvoMethods) -> EvoMethods: - """Return the evolution methods. +def ev_method(s: EvolutionMethod) -> EvoMethods: + """Return the evolution method. Parameters ---------- diff --git a/tests/eko/evolution_operator/test_ome.py b/tests/eko/evolution_operator/test_ome.py index b844290c9..584635967 100644 --- a/tests/eko/evolution_operator/test_ome.py +++ b/tests/eko/evolution_operator/test_ome.py @@ -8,6 +8,7 @@ from eko import interpolation, mellin from eko import scale_variations as sv from eko.evolution_operator.operator_matrix_element import ( + BackwardsMethods, OperatorMatrixElement, build_ome, quad_ker, @@ -33,7 +34,7 @@ def test_build_ome_as(): aS = A_singlet((o, 0), N, nf, L, is_msbar) for a in [aNS, aS]: - for method in [None, InversionMethod.EXPANDED, InversionMethod.EXACT]: + for method in BackwardsMethods: dim = len(a[0]) if o != 1: assert len(a) == o @@ -53,7 +54,7 @@ def test_build_ome_nlo(): aNSi = A_non_singlet((1, 0), N, nf, L) aSi = A_singlet((1, 0), N, nf, L, is_msbar) for a in [aNSi, aSi]: - for method in [None, InversionMethod.EXPANDED, InversionMethod.EXACT]: + for method in BackwardsMethods: dim = len(a[0]) # hh assert a[0, -1, -1] != 0.0 @@ -86,7 +87,7 @@ def test_quad_ker_errors(): is_log=True, logx=0.123, areas=[[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]], - backward_method=None, + backward_method=BackwardsMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -119,7 +120,7 @@ def test_quad_ker(monkeypatch): is_log=is_log, logx=0.123, areas=np.zeros(3), - backward_method=None, + backward_method=BackwardsMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -138,7 +139,7 @@ def test_quad_ker(monkeypatch): is_log=is_log, logx=0.123, areas=np.zeros(3), - backward_method=None, + backward_method=BackwardsMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -157,7 +158,7 @@ def test_quad_ker(monkeypatch): is_log=is_log, logx=0.0, areas=np.zeros(3), - backward_method=None, + backward_method=BackwardsMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -191,7 +192,7 @@ def test_quad_ker(monkeypatch): is_log=True, logx=0.123, areas=np.zeros(3), - backward_method=InversionMethod.EXPANDED, + backward_method=BackwardsMethods.EXPANDED, a_s=0.0, nf=3, L=0.0, @@ -228,7 +229,7 @@ def test_quad_ker(monkeypatch): is_log=True, logx=0.123, areas=np.zeros(3), - backward_method=InversionMethod.EXACT, + backward_method=BackwardsMethods.EXACT, a_s=0.0, nf=3, L=0.0, @@ -252,7 +253,7 @@ def test_quad_ker(monkeypatch): is_log=True, logx=0.0, areas=np.array([0.01, 0.1, 1.0]), - backward_method=None, + backward_method=BackwardsMethods.FORWARD, a_s=0.0, nf=3, L=0.0, From c88627f71da3df6780692cbcc1bd04c37d0776ae Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Mon, 22 Jul 2024 13:23:45 +0300 Subject: [PATCH 14/20] Rename matching method --- .../operator_matrix_element.py | 34 ++++++++----------- src/ekomark/data/operators.py | 2 +- tests/eko/evolution_operator/test_ome.py | 20 +++++------ 3 files changed, 25 insertions(+), 31 deletions(-) diff --git a/src/eko/evolution_operator/operator_matrix_element.py b/src/eko/evolution_operator/operator_matrix_element.py index b4b6731b4..cf6cbcb97 100644 --- a/src/eko/evolution_operator/operator_matrix_element.py +++ b/src/eko/evolution_operator/operator_matrix_element.py @@ -21,16 +21,16 @@ logger = logging.getLogger(__name__) -class BackwardsMethods(enum.IntEnum): +class MatchingMethods(enum.IntEnum): """Enumerate backward methods.""" FORWARD = enum.auto() - EXACT = enum.auto() - EXPANDED = enum.auto() + BACKWARD_EXACT = enum.auto() + BACKWARD_EXPANDED = enum.auto() -def bw_method(s: InversionMethod) -> BackwardsMethods: - """Return the backward method. +def matching_method(s: InversionMethod) -> MatchingMethods: + """Return the matching method. Parameters ---------- @@ -44,8 +44,8 @@ def bw_method(s: InversionMethod) -> BackwardsMethods: """ if s is not None: - return BackwardsMethods[s.value.upper()] - return BackwardsMethods.FORWARD + return MatchingMethods[s.value.upper()] + return MatchingMethods.FORWARD @nb.njit(cache=True) @@ -60,7 +60,7 @@ def build_ome(A, matching_order, a_s, backward_method): perturbation matching order a_s : float strong coupling, needed only for the exact inverse - backward_method : BackwardsMethods + backward_method : MatchingMethods empty or method for inverting the matching condition (exact or expanded) Returns @@ -79,7 +79,7 @@ def build_ome(A, matching_order, a_s, backward_method): ome = np.eye(len(A[0]), dtype=np.complex_) A = A[:, :, :] A = np.ascontiguousarray(A) - if backward_method is BackwardsMethods.EXPANDED: + if backward_method is MatchingMethods.BACKWARD_EXPANDED: # expended inverse if matching_order[0] >= 1: ome -= a_s * A[0] @@ -96,7 +96,7 @@ def build_ome(A, matching_order, a_s, backward_method): if matching_order[0] >= 3: ome += a_s**3 * A[2] # need inverse exact ? so add the missing pieces - if backward_method is BackwardsMethods.EXACT: + if backward_method is MatchingMethods.BACKWARD_EXACT: ome = np.linalg.inv(ome) return ome @@ -244,7 +244,7 @@ class OperatorMatrixElement(Operator): def __init__(self, config, managers, nf, q2, is_backward, L, is_msbar): super().__init__(config, managers, Segment(q2, q2, nf)) - self.backward_method = bw_method( + self.backward_method = matching_method( config["backward_inversion"] if is_backward else None ) if is_backward: @@ -271,10 +271,7 @@ def labels(self): logger.warning("%s: skipping non-singlet sector", self.log_label) else: labels.append((200, 200)) - if ( - self.is_intrinsic - or self.backward_method is not BackwardsMethods.FORWARD - ): + if self.is_intrinsic or self.backward_method is not MatchingMethods.FORWARD: # intrinsic labels, which are not zero at NLO labels.append((br.matching_hminus_pid, br.matching_hminus_pid)) # These contributions are always 0 for the moment @@ -290,10 +287,7 @@ def labels(self): (br.matching_hplus_pid, 100), ] ) - if ( - self.is_intrinsic - or self.backward_method is not BackwardsMethods.FORWARD - ): + if self.is_intrinsic or self.backward_method is not MatchingMethods.FORWARD: labels.extend( [ (21, br.matching_hplus_pid), @@ -365,7 +359,7 @@ def compute(self): self.log_label, self.order[0], self.order[1], - BackwardsMethods(self.backward_method).name, + MatchingMethods(self.backward_method).name, ) self.integrate() diff --git a/src/ekomark/data/operators.py b/src/ekomark/data/operators.py index ad928574d..03ad4d349 100644 --- a/src/ekomark/data/operators.py +++ b/src/ekomark/data/operators.py @@ -15,7 +15,7 @@ ev_op_max_order=10, ev_op_iterations=10, backward_inversion="expanded", - n_integration_cores=0, + n_integration_cores=1, debug_skip_non_singlet=False, debug_skip_singlet=False, mugrid=[10], diff --git a/tests/eko/evolution_operator/test_ome.py b/tests/eko/evolution_operator/test_ome.py index 584635967..553e5823c 100644 --- a/tests/eko/evolution_operator/test_ome.py +++ b/tests/eko/evolution_operator/test_ome.py @@ -8,7 +8,7 @@ from eko import interpolation, mellin from eko import scale_variations as sv from eko.evolution_operator.operator_matrix_element import ( - BackwardsMethods, + MatchingMethods, OperatorMatrixElement, build_ome, quad_ker, @@ -34,7 +34,7 @@ def test_build_ome_as(): aS = A_singlet((o, 0), N, nf, L, is_msbar) for a in [aNS, aS]: - for method in BackwardsMethods: + for method in MatchingMethods: dim = len(a[0]) if o != 1: assert len(a) == o @@ -54,7 +54,7 @@ def test_build_ome_nlo(): aNSi = A_non_singlet((1, 0), N, nf, L) aSi = A_singlet((1, 0), N, nf, L, is_msbar) for a in [aNSi, aSi]: - for method in BackwardsMethods: + for method in MatchingMethods: dim = len(a[0]) # hh assert a[0, -1, -1] != 0.0 @@ -87,7 +87,7 @@ def test_quad_ker_errors(): is_log=True, logx=0.123, areas=[[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]], - backward_method=BackwardsMethods.FORWARD, + backward_method=MatchingMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -120,7 +120,7 @@ def test_quad_ker(monkeypatch): is_log=is_log, logx=0.123, areas=np.zeros(3), - backward_method=BackwardsMethods.FORWARD, + backward_method=MatchingMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -139,7 +139,7 @@ def test_quad_ker(monkeypatch): is_log=is_log, logx=0.123, areas=np.zeros(3), - backward_method=BackwardsMethods.FORWARD, + backward_method=MatchingMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -158,7 +158,7 @@ def test_quad_ker(monkeypatch): is_log=is_log, logx=0.0, areas=np.zeros(3), - backward_method=BackwardsMethods.FORWARD, + backward_method=MatchingMethods.FORWARD, a_s=0.0, nf=3, L=0.0, @@ -192,7 +192,7 @@ def test_quad_ker(monkeypatch): is_log=True, logx=0.123, areas=np.zeros(3), - backward_method=BackwardsMethods.EXPANDED, + backward_method=MatchingMethods.BACKWARD_EXPANDED, a_s=0.0, nf=3, L=0.0, @@ -229,7 +229,7 @@ def test_quad_ker(monkeypatch): is_log=True, logx=0.123, areas=np.zeros(3), - backward_method=BackwardsMethods.EXACT, + backward_method=MatchingMethods.BACKWARD_EXACT, a_s=0.0, nf=3, L=0.0, @@ -253,7 +253,7 @@ def test_quad_ker(monkeypatch): is_log=True, logx=0.0, areas=np.array([0.01, 0.1, 1.0]), - backward_method=BackwardsMethods.FORWARD, + backward_method=MatchingMethods.FORWARD, a_s=0.0, nf=3, L=0.0, From d4a3f71889822f0311f342ec2d3ca97c049acf8f Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Mon, 22 Jul 2024 13:29:41 +0300 Subject: [PATCH 15/20] Fix MatchingMethods doc string --- src/eko/evolution_operator/operator_matrix_element.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/eko/evolution_operator/operator_matrix_element.py b/src/eko/evolution_operator/operator_matrix_element.py index cf6cbcb97..87c76183a 100644 --- a/src/eko/evolution_operator/operator_matrix_element.py +++ b/src/eko/evolution_operator/operator_matrix_element.py @@ -22,7 +22,7 @@ class MatchingMethods(enum.IntEnum): - """Enumerate backward methods.""" + """Enumerate matching methods.""" FORWARD = enum.auto() BACKWARD_EXACT = enum.auto() From 7f079bc6001ac3f38d0d55ebd0a3cd271da7e976 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Mon, 22 Jul 2024 15:12:05 +0300 Subject: [PATCH 16/20] Fix matching renaming --- src/eko/evolution_operator/operator_matrix_element.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/eko/evolution_operator/operator_matrix_element.py b/src/eko/evolution_operator/operator_matrix_element.py index 87c76183a..25e5df2dd 100644 --- a/src/eko/evolution_operator/operator_matrix_element.py +++ b/src/eko/evolution_operator/operator_matrix_element.py @@ -44,7 +44,7 @@ def matching_method(s: InversionMethod) -> MatchingMethods: """ if s is not None: - return MatchingMethods[s.value.upper()] + return MatchingMethods["BACKWARD_" + s.value.upper()] return MatchingMethods.FORWARD From e40acaec28c4212d5a3ffb809ece959d415c2016 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Wed, 24 Jul 2024 15:10:40 +0300 Subject: [PATCH 17/20] Allow only iterate-exact for QED --- src/eko/kernels/singlet_qed.py | 2 +- src/eko/kernels/valence_qed.py | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/eko/kernels/singlet_qed.py b/src/eko/kernels/singlet_qed.py index 2f2ba8eee..d63a2f1ac 100644 --- a/src/eko/kernels/singlet_qed.py +++ b/src/eko/kernels/singlet_qed.py @@ -97,7 +97,7 @@ def dispatcher( e_s : numpy.ndarray singlet EKO """ - if method in [EvoMethods.ITERATE_EXACT, EvoMethods.ITERATE_EXPANDED]: + if method is EvoMethods.ITERATE_EXACT: return eko_iterate( gamma_singlet, as_list, a_half, nf, order, ev_op_iterations, 4 ) diff --git a/src/eko/kernels/valence_qed.py b/src/eko/kernels/valence_qed.py index 074dd95ed..186c70219 100644 --- a/src/eko/kernels/valence_qed.py +++ b/src/eko/kernels/valence_qed.py @@ -24,7 +24,7 @@ def dispatcher( ---------- order : tuple(int,int) perturbative order - method : str + method : int method gamma_singlet : numpy.ndarray singlet anomalous dimensions matrices @@ -46,7 +46,7 @@ def dispatcher( e_v : numpy.ndarray singlet EKO """ - if method in [EvoMethods.ITERATE_EXACT, EvoMethods.ITERATE_EXPANDED]: + if method is EvoMethods.ITERATE_EXACT: return eko_iterate( gamma_valence, as_list, a_half, nf, order, ev_op_iterations, 2 ) From c92c3a24a0cd2cc440048f1028af9b98ab63b49a Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 29 Jul 2024 18:11:07 +0000 Subject: [PATCH 18/20] [pre-commit.ci] pre-commit autoupdate MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit updates: - [github.com/asottile/pyupgrade: v3.16.0 → v3.17.0](https://github.com/asottile/pyupgrade/compare/v3.16.0...v3.17.0) - [github.com/pre-commit/pre-commit: v3.7.1 → v3.8.0](https://github.com/pre-commit/pre-commit/compare/v3.7.1...v3.8.0) --- .pre-commit-config.yaml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 67e901972..d579fb856 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -32,7 +32,7 @@ repos: - id: isort args: ["--profile", "black"] - repo: https://github.com/asottile/pyupgrade - rev: v3.16.0 + rev: v3.17.0 hooks: - id: pyupgrade - repo: https://github.com/pycqa/pydocstyle @@ -53,6 +53,6 @@ repos: files: ^crates/.*\.rs$ args: [] - repo: https://github.com/pre-commit/pre-commit - rev: v3.7.1 + rev: v3.8.0 hooks: - id: validate_manifest From 48476f9f2a3ab8121ff694db7c413fbb3bdddf1c Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Thu, 8 Aug 2024 14:49:28 +0300 Subject: [PATCH 19/20] Fix Rust ev_op patch --- src/eko/evolution_operator/__init__.py.patch | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/src/eko/evolution_operator/__init__.py.patch b/src/eko/evolution_operator/__init__.py.patch index 58c3f74fd..d8258ee14 100644 --- a/src/eko/evolution_operator/__init__.py.patch +++ b/src/eko/evolution_operator/__init__.py.patch @@ -1,5 +1,5 @@ diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py -index 1c759c5c..5eb394d0 100644 +index ec56b6db..374d0d0b 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py @@ -3,15 +3,15 @@ r"""Contains the central operator classes. @@ -20,7 +20,7 @@ index 1c759c5c..5eb394d0 100644 import ekore.anomalous_dimensions.polarized.space_like as ad_ps import ekore.anomalous_dimensions.unpolarized.space_like as ad_us -@@ -27,92 +27,10 @@ from ..kernels import singlet_qed as qed_s +@@ -29,92 +29,10 @@ from ..kernels import singlet_qed as qed_s from ..kernels import valence_qed as qed_v from ..matchings import Segment, lepton_number from ..member import OpMember @@ -114,7 +114,7 @@ index 1c759c5c..5eb394d0 100644 spec = [ ("is_singlet", nb.boolean), ("is_QEDsinglet", nb.boolean), -@@ -184,422 +102,6 @@ class QuadKerBase: +@@ -186,422 +104,6 @@ class QuadKerBase: return self.path.prefactor * pj * self.path.jac @@ -537,9 +537,9 @@ index 1c759c5c..5eb394d0 100644 class Operator(sv.ScaleVariationModeMixin): """Internal representation of a single EKO. -@@ -780,50 +282,6 @@ class Operator(sv.ScaleVariationModeMixin): - labels.extend(br.singlet_unified_labels) - return labels +@@ -787,50 +289,6 @@ class Operator(sv.ScaleVariationModeMixin): + """Return the evolution method.""" + return ev_method(EvolutionMethod(self.config["method"])) - def quad_ker(self, label, logx, areas): - """Return partially initialized integrand function. @@ -564,7 +564,7 @@ index 1c759c5c..5eb394d0 100644 - order=self.order, - mode0=label[0], - mode1=label[1], -- method=self.config["method"], +- method=self.ev_method, - is_log=self.int_disp.log, - logx=logx, - areas=areas, @@ -588,7 +588,7 @@ index 1c759c5c..5eb394d0 100644 def initialize_op_members(self): """Init all operators with the identity or zeros.""" eye = OpMember( -@@ -846,10 +304,7 @@ class Operator(sv.ScaleVariationModeMixin): +@@ -853,10 +311,7 @@ class Operator(sv.ScaleVariationModeMixin): else: self.op_members[n] = zero.copy() @@ -600,7 +600,7 @@ index 1c759c5c..5eb394d0 100644 """Run the integration for each grid point. Parameters -@@ -864,18 +319,56 @@ class Operator(sv.ScaleVariationModeMixin): +@@ -871,18 +326,56 @@ class Operator(sv.ScaleVariationModeMixin): """ column = [] k, logx = log_grid From b651ebe9699dc2ffcacd561e05bb52026e4d1944 Mon Sep 17 00:00:00 2001 From: Felix Hekhorn Date: Thu, 8 Aug 2024 15:02:47 +0300 Subject: [PATCH 20/20] Fix Rust ev_op patch --- src/eko/evolution_operator/__init__.py.patch | 26 ++++++++++---------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/src/eko/evolution_operator/__init__.py.patch b/src/eko/evolution_operator/__init__.py.patch index 58c3f74fd..995d77b6b 100644 --- a/src/eko/evolution_operator/__init__.py.patch +++ b/src/eko/evolution_operator/__init__.py.patch @@ -1,8 +1,8 @@ diff --git a/src/eko/evolution_operator/__init__.py b/src/eko/evolution_operator/__init__.py -index 1c759c5c..5eb394d0 100644 +index fe07ade9..0f58c9e5 100644 --- a/src/eko/evolution_operator/__init__.py +++ b/src/eko/evolution_operator/__init__.py -@@ -3,15 +3,15 @@ r"""Contains the central operator classes. +@@ -3,16 +3,16 @@ r"""Contains the central operator classes. See :doc:`Operator overview `. """ @@ -11,6 +11,7 @@ index 1c759c5c..5eb394d0 100644 import os import time from multiprocessing import Pool + from typing import Dict, Tuple +import ekors import numba as nb @@ -20,7 +21,7 @@ index 1c759c5c..5eb394d0 100644 import ekore.anomalous_dimensions.polarized.space_like as ad_ps import ekore.anomalous_dimensions.unpolarized.space_like as ad_us -@@ -27,92 +27,10 @@ from ..kernels import singlet_qed as qed_s +@@ -30,92 +30,10 @@ from ..kernels import singlet_qed as qed_s from ..kernels import valence_qed as qed_v from ..matchings import Segment, lepton_number from ..member import OpMember @@ -114,7 +115,7 @@ index 1c759c5c..5eb394d0 100644 spec = [ ("is_singlet", nb.boolean), ("is_QEDsinglet", nb.boolean), -@@ -184,422 +102,6 @@ class QuadKerBase: +@@ -187,421 +105,6 @@ class QuadKerBase: return self.path.prefactor * pj * self.path.jac @@ -533,13 +534,12 @@ index 1c759c5c..5eb394d0 100644 - ) - return ker - -- - class Operator(sv.ScaleVariationModeMixin): - """Internal representation of a single EKO. -@@ -780,50 +282,6 @@ class Operator(sv.ScaleVariationModeMixin): - labels.extend(br.singlet_unified_labels) - return labels + OpMembers = Dict[OperatorLabel, OpMember] + """Map of all operators.""" +@@ -792,50 +295,6 @@ class Operator(sv.ScaleVariationModeMixin): + """Return the evolution method.""" + return ev_method(EvolutionMethod(self.config["method"])) - def quad_ker(self, label, logx, areas): - """Return partially initialized integrand function. @@ -564,7 +564,7 @@ index 1c759c5c..5eb394d0 100644 - order=self.order, - mode0=label[0], - mode1=label[1], -- method=self.config["method"], +- method=self.ev_method, - is_log=self.int_disp.log, - logx=logx, - areas=areas, @@ -588,7 +588,7 @@ index 1c759c5c..5eb394d0 100644 def initialize_op_members(self): """Init all operators with the identity or zeros.""" eye = OpMember( -@@ -846,10 +304,7 @@ class Operator(sv.ScaleVariationModeMixin): +@@ -858,10 +317,7 @@ class Operator(sv.ScaleVariationModeMixin): else: self.op_members[n] = zero.copy() @@ -600,7 +600,7 @@ index 1c759c5c..5eb394d0 100644 """Run the integration for each grid point. Parameters -@@ -864,18 +319,56 @@ class Operator(sv.ScaleVariationModeMixin): +@@ -876,18 +332,56 @@ class Operator(sv.ScaleVariationModeMixin): """ column = [] k, logx = log_grid