diff --git a/src/eko/io/bases.py b/src/eko/io/bases.py deleted file mode 100644 index 50ee74c64..000000000 --- a/src/eko/io/bases.py +++ /dev/null @@ -1,118 +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/manipulate.py b/src/eko/io/manipulate.py index d2bf0ff76..362c8d77d 100644 --- a/src/eko/io/manipulate.py +++ b/src/eko/io/manipulate.py @@ -1,8 +1,9 @@ """Manipulate output generate by EKO.""" +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 @@ -10,7 +11,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__) @@ -19,13 +20,13 @@ 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 +) -> npt.NDArray: """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 @@ -33,21 +34,23 @@ def rotation(new: Optional[Basis], old: Basis, 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): +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 @@ -62,81 +65,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( + targetrot = rotation( targetgrid, - eko.bases.targetgrid, + xgrid, check, lambda new, old: crot(new, old, interpdeg), ) - newinput, inputrot = rotation( + 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] - - 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 + operands = [elem.operator] + operands_errs = [elem.error] - 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) + 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 - 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" @@ -146,47 +134,38 @@ 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. +) -> Operator: + """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 inv_inputpids = np.zeros_like(inputpids) @@ -194,60 +173,47 @@ def flavor_reshape( 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) + return 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() - - -def to_evol(eko: EKO, source: bool = True, target: bool = False): +def to_evol(elem: Operator, source: bool = True, target: bool = False) -> Operator: """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 @@ -258,27 +224,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 + return flavor_reshape(elem, inputpids=inputpids, targetpids=targetpids) - eko.update() - -def to_uni_evol(eko: EKO, 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. - This assigns also the pids. The operation is in-place. - Parameters ---------- - eko : + elem : the operator to be rotated source : rotate on the input tensor @@ -289,14 +243,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/eko/io/metadata.py b/src/eko/io/metadata.py index 94cfad613..ce8b13998 100644 --- a/src/eko/io/metadata.py +++ b/src/eko/io/metadata.py @@ -9,7 +9,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 @@ -35,10 +35,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 076139d71..1520e1003 100644 --- a/src/eko/io/struct.py +++ b/src/eko/io/struct.py @@ -17,7 +17,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 @@ -106,20 +105,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 @@ -554,7 +548,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 a260a4cec..bf9f2083d 100644 --- a/src/ekobox/apply.py +++ b/src/ekobox/apply.py @@ -91,12 +91,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 @@ -107,9 +107,9 @@ def apply_pdf_flavor( error_final = np.einsum(CONTRACTION, elem.error, pdfs, optimize="optimal") else: error_final = None - out_grid[ep] = PdfResult(dict(zip(eko.bases.targetpids, pdf_final))) + out_grid[ep] = PdfResult(dict(zip(br.flavor_basis_pids, pdf_final))) 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)) qed = eko.theory_card.order[1] > 0 # rotate to evolution basis @@ -132,7 +132,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/src/ekomark/benchmark/runner.py b/src/ekomark/benchmark/runner.py index 8ac48f5d3..2b003e1b3 100644 --- a/src/ekomark/benchmark/runner.py +++ b/src/ekomark/benchmark/runner.py @@ -14,7 +14,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 @@ -117,22 +116,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 15fa19192..5e2bee420 100644 --- a/src/ekomark/plots.py +++ b/src/ekomark/plots.py @@ -2,7 +2,7 @@ import io import pprint -from typing import Dict, List, Tuple +from typing import Dict, List import matplotlib.pyplot as plt import numpy as np @@ -11,6 +11,7 @@ from matplotlib.colors import LogNorm from eko import basis_rotation as br +from eko.io import manipulate from eko.io.struct import EKO @@ -211,7 +212,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 @@ -226,15 +235,12 @@ def save_operators_to_pdf(path, theory, ops, me: EKO, skip_pdfs, change_lab=Fals 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 = list(me.bases.targetpids) - if np.allclose(ops_names, br.rotate_flavor_to_evolution): - ops_names = list(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}") @@ -247,6 +253,11 @@ 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 ep, op in me.items(): + 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) results = op.operator errors = op.error @@ -254,8 +265,8 @@ def save_operators_to_pdf(path, theory, ops, me: EKO, skip_pdfs, change_lab=Fals for label_out, res, res_err in zip(ops_names, results, errors): if label_out in skip_pdfs: continue - new_op: Dict[Tuple[int, ...], List[float]] = {} - new_op_err: Dict[Tuple[int, ...], List[float]] = {} + new_op: Dict[int, List[float]] = {} + new_op_err: Dict[int, List[float]] = {} # loop on xgrid point for j in range(len(me.xgrid)): # loop on pid in @@ -271,17 +282,10 @@ def save_operators_to_pdf(path, theory, ops, me: EKO, skip_pdfs, change_lab=Fals for label_in in ops_names: if label_in in skip_pdfs: continue - lab_in = label_in - lab_out = label_out - if change_lab: - 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] - lab_out = br.evol_basis[index_out] fig = None try: fig = plot_operator( - f"Operator ({lab_in};{lab_out}) µ_F^2 = {ep[0]} GeV^2, nf = {ep[1]}", + f"Operator ({label_in};{label_out}) µ_F^2 = {ep[0]} GeV^2, nf = {ep[1]}", new_op[label_in], new_op_err[label_in], ) 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_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 + ) 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)