diff --git a/.github/actions/cache-test-data/action.yml b/.github/actions/cache-test-data/action.yml index f0629ee6..855d6704 100644 --- a/.github/actions/cache-test-data/action.yml +++ b/.github/actions/cache-test-data/action.yml @@ -10,10 +10,11 @@ runs: uses: actions/cache@v4 with: path: test-data - key: test-data-v24 + key: test-data-v30 - name: Download test data if cache miss if: steps.cache.outputs.cache-hit != 'true' run: | cd maintainer + # TODO: if this script fails to download some files, the cache is saved nevertheless ./download-test-data.sh shell: bash diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md index 666c1ec3..f6f36d88 100644 --- a/CONTRIBUTING.md +++ b/CONTRIBUTING.md @@ -56,6 +56,11 @@ increasing the MSRV make sure to set it everywhere to the same value: data with `curl`. To make Github refresh the cached test data when running the CI, increase the integer `XX` in the line `key: test-data-vXX` by one. +### Regression tests for Github Issues + +If you're writing a regression test for a Github Issue, name the test +`issue_XXX`, where `XXX` is the Github Issue number. + ## Git - When you commit, make sure the commit message is written properly. This diff --git a/Cargo.lock b/Cargo.lock index ea33ec9f..b9b51bf0 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -487,9 +487,6 @@ name = "float-cmp" version = "0.9.0" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "98de4bbd547a563b716d8dfa9aad1cb19bfab00f4fa09a6a4ed21dbcf44ce9c4" -dependencies = [ - "num-traits", -] [[package]] name = "form_urlencoded" @@ -712,22 +709,13 @@ version = "0.4.20" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "b5e6163cb8c49088c2c36f57875e58ccd8c87c7427f7fbd50ea6710b2f3f2e8f" -[[package]] -name = "lz4_flex" -version = "0.9.5" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "1a8cbbb2831780bc3b9c15a41f5b49222ef756b6730a95f3decfdd15903eb5a3" -dependencies = [ - "twox-hash 1.6.3", -] - [[package]] name = "lz4_flex" version = "0.11.6" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "373f5eceeeab7925e0c1098212f2fbc4d416adec9d35051a6ab251e824c1854a" dependencies = [ - "twox-hash 2.1.2", + "twox-hash", ] [[package]] @@ -930,27 +918,6 @@ dependencies = [ "sha2", ] -[[package]] -name = "pineappl" -version = "0.8.3" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "5041fcf611eb0c41f1f6562b498fabdd1319c8d4572fd137ac244ca4e73d999c" -dependencies = [ - "anyhow", - "arrayvec", - "bincode", - "bitflags 2.4.2", - "enum_dispatch", - "float-cmp", - "git-version", - "itertools", - "lz4_flex 0.9.5", - "ndarray", - "rustc-hash 1.1.0", - "serde", - "thiserror", -] - [[package]] name = "pineappl" version = "1.3.3" @@ -963,11 +930,11 @@ dependencies = [ "float-cmp", "git-version", "itertools", - "lz4_flex 0.11.6", + "lz4_flex", "managed-lhapdf", "ndarray", "num-complex", - "pineappl 0.8.3", + "pineappl_v0", "rand", "rand_pcg", "rayon", @@ -993,7 +960,7 @@ version = "1.3.3" dependencies = [ "itertools", "ndarray", - "pineappl 1.3.3", + "pineappl", ] [[package]] @@ -1012,11 +979,11 @@ dependencies = [ "float-cmp", "git-version", "itertools", - "lz4_flex 0.11.6", + "lz4_flex", "managed-lhapdf", "ndarray", "ndarray-npy", - "pineappl 1.3.3", + "pineappl", "pineappl_applgrid", "pineappl_fastnlo", "predicates", @@ -1044,10 +1011,21 @@ dependencies = [ "itertools", "ndarray", "numpy", - "pineappl 1.3.3", + "pineappl", "pyo3", ] +[[package]] +name = "pineappl_v0" +version = "1.3.3" +dependencies = [ + "bincode", + "enum_dispatch", + "ndarray", + "serde", + "thiserror", +] + [[package]] name = "pkg-config" version = "0.3.29" @@ -1448,12 +1426,6 @@ version = "1.13.2" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "3c5e1a9a646d36c3599cd173a41282daf47c44583ad367b8e6837255952e5c67" -[[package]] -name = "static_assertions" -version = "1.1.0" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "a2eb9349b6444b326872e140eb1cf5e7c522154d69e7a0ffb0fb81c06b37543f" - [[package]] name = "strsim" version = "0.11.1" @@ -1586,16 +1558,6 @@ dependencies = [ "winnow", ] -[[package]] -name = "twox-hash" -version = "1.6.3" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "97fee6b57c6a41524a810daee9286c02d7752c4253064d0b05472833a438f675" -dependencies = [ - "cfg-if", - "static_assertions", -] - [[package]] name = "twox-hash" version = "2.1.2" diff --git a/Cargo.toml b/Cargo.toml index e4fde002..dd002899 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -6,6 +6,7 @@ members = [ "pineappl_cli", "pineappl_fastnlo", "pineappl_py", + "pineappl_v0", "xtask", ] default-members = [ diff --git a/maintainer/download-test-data.sh b/maintainer/download-test-data.sh index 926abfed..63cab140 100755 --- a/maintainer/download-test-data.sh +++ b/maintainer/download-test-data.sh @@ -21,6 +21,10 @@ files=( 'https://data.nnpdf.science/dy_high_mass/NNPDF_DY_14TEV_BSM_AFB.pineappl.lz4' 'https://data.nnpdf.science/pineappl/test-data/ATLASWPT11-Wplus_tot.appl' 'https://data.nnpdf.science/pineappl/test-data/CMS_TTB_8TEV_2D_TTM_TRAP_TOT-opt.pineappl.lz4' + 'https://data.nnpdf.science/pineappl/test-data/DYAA_0.3.0_LagrangeSubgridV1.pineappl.lz4' + 'https://data.nnpdf.science/pineappl/test-data/DYAA_0.4.1_LagrangeSparseSubgridV1.pineappl.lz4' + 'https://data.nnpdf.science/pineappl/test-data/DYAA_0.4.1_LagrangeSubgridV1.pineappl.lz4' + 'https://data.nnpdf.science/pineappl/test-data/DYAA_0.4.1_LagrangeSubgridV2.pineappl.lz4' 'https://data.nnpdf.science/pineappl/test-data/CMS_TTB_8TEV_2D_TTM_TRAP_TOT.tar' 'https://data.nnpdf.science/pineappl/test-data/E906nlo_bin_00.pineappl.lz4' 'https://data.nnpdf.science/pineappl/test-data/E906nlo_bin_00.tar' diff --git a/maintainer/make-release.sh b/maintainer/make-release.sh index e3afe71f..372daa7f 100755 --- a/maintainer/make-release.sh +++ b/maintainer/make-release.sh @@ -70,6 +70,7 @@ for crate in . ${crates[@]}; do -e "s:^\(pineappl_applgrid = .*\)version = \".*\":\1version = \"=${version}\":" \ -e "s:^\(pineappl_cli = .*\)version = \".*\":\1version = \"=${version}\":" \ -e "s:^\(pineappl_fastnlo = .*\)version = \".*\":\1version = \"=${version}\":" \ + -e "s:^\(pineappl_v0 = .*\)version = \".*\":\1version = \"=${version}\":" \ ${crate}/Cargo.toml git add ${crate}/Cargo.toml done diff --git a/pineappl/Cargo.toml b/pineappl/Cargo.toml index 75e1e3ca..bcf9be3c 100644 --- a/pineappl/Cargo.toml +++ b/pineappl/Cargo.toml @@ -26,8 +26,7 @@ git-version = "0.3.5" itertools = "0.10.1" lz4_flex = "0.11.6" ndarray = { features = ["serde"], version = "0.15.4" } -# TODO: opt out of default features in this crate to match line above -pineappl-v0 = { package = "pineappl", version = "0.8.2" } +pineappl_v0 = { path = "../pineappl_v0", version = "=1.3.3" } rayon = "1.5.1" rustc-hash = "1.1.0" serde = { features = ["derive"], version = "1.0.130" } diff --git a/pineappl/src/v0.rs b/pineappl/src/v0.rs index 5271cbd7..4b7b0c23 100644 --- a/pineappl/src/v0.rs +++ b/pineappl/src/v0.rs @@ -41,10 +41,26 @@ pub fn default_interps(flexible_scale: bool, convolutions: usize) -> Vec } pub fn read_uncompressed_v0(mut reader: impl BufRead) -> Result { + use pineappl_v0::convolutions::Convolution; + use pineappl_v0::pids::PidBasis as PidBasisV0; use pineappl_v0::subgrid::Subgrid as _; - let grid = GridV0::read(&mut reader).map_err(|err| Error::Other(err.into()))?; - let convolutions = read_convolutions_from_metadata(&grid); + let grid = GridV0::read_uncompressed(&mut reader).map_err(|err| Error::Other(err.into()))?; + let convolutions: Vec<_> = grid + .convolutions() + .into_iter() + .map(|old_conv| match old_conv { + Convolution::UnpolPDF(pid) => Some(Conv::new(ConvType::UnpolPDF, pid)), + Convolution::PolPDF(pid) => Some(Conv::new(ConvType::PolPDF, pid)), + Convolution::UnpolFF(pid) => Some(Conv::new(ConvType::UnpolFF, pid)), + Convolution::PolFF(pid) => Some(Conv::new(ConvType::PolFF, pid)), + Convolution::None => None, + }) + .collect(); + let pid_basis = match grid.pid_basis() { + PidBasisV0::Pdg => PidBasis::Pdg, + PidBasisV0::Evol => PidBasis::Evol, + }; let flexible_scale_grid = grid.subgrids().iter().any(|subgrid| { subgrid @@ -128,16 +144,7 @@ pub fn read_uncompressed_v0(mut reader: impl BufRead) -> Result { ) }) .collect(), - grid.key_values() - .and_then(|kv| kv.get("lumi_id_types")) - // TODO: use PidBasis::from_str - .map_or(PidBasis::Pdg, |lumi_id_types| { - match lumi_id_types.as_str() { - "pdg_mc_ids" => PidBasis::Pdg, - "evol" => PidBasis::Evol, - _ => panic!("unknown PID basis '{lumi_id_types}'"), - } - }), + pid_basis, convolutions.clone().into_iter().flatten().collect(), default_interps( flexible_scale_grid, @@ -261,74 +268,3 @@ pub fn read_uncompressed_v0(mut reader: impl BufRead) -> Result { Ok(result) } - -fn read_convolutions_from_metadata(grid: &GridV0) -> Vec> { - grid.key_values().map_or_else( - // if there isn't any metadata, we assume two unpolarized proton-PDFs are used - || vec![Some(Conv::new(ConvType::UnpolPDF, 2212)); 2], - |kv| { - // file format v0 only supports exactly two convolutions - (1..=2) - .map(|index| { - // if there are key-value pairs `convolution_particle_1` and - // `convolution_type_1` and the same with a higher index, we convert this - // metadata into `Convolution` - match ( - kv.get(&format!("convolution_particle_{index}")) - .map(|s| s.parse::()), - kv.get(&format!("convolution_type_{index}")) - .map(String::as_str), - ) { - (_, Some("None")) => None, - (Some(Ok(pid)), Some("UnpolPDF")) => Some(Conv::new(ConvType::UnpolPDF, pid)), - (Some(Ok(pid)), Some("PolPDF")) => Some(Conv::new(ConvType::PolPDF, pid)), - (Some(Ok(pid)), Some("UnpolFF")) => Some(Conv::new(ConvType::UnpolFF, pid)), - (Some(Ok(pid)), Some("PolFF")) => Some(Conv::new(ConvType::PolFF, pid)), - (None, None) => { - // if these key-value pairs are missing use the old metadata - match kv - .get(&format!("initial_state_{index}")) - .map(|s| s.parse::()) - { - Some(Ok(pid)) => { - // Addresses: https://github.com/NNPDF/pineappl/issues/334 - if grid.channels().is_empty() && pid == 2212 { - Some(Conv::new(ConvType::UnpolPDF, 2212)) - } else { - let condition = !grid.channels().iter().all(|entry| { - entry.entry().iter().all(|&(a, b, _)| - match index { - 1 => a, - 2 => b, - _ => unreachable!() - } == pid - ) - }); - - condition.then_some(Conv::new(ConvType::UnpolPDF, pid)) - } - } - None => Some(Conv::new(ConvType::UnpolPDF, 2212)), - Some(Err(err)) => panic!( - "metadata 'initial_state_{index}' could not be parsed: {err}" - ), - } - } - (None, Some(_)) => { - panic!("metadata 'convolution_type_{index}' is missing") - } - (Some(_), None) => { - panic!("metadata 'convolution_particle_{index}' is missing") - } - (Some(Ok(_)), Some(type_)) => { - panic!("metadata 'convolution_type_{index} = {type_}' is unknown") - } - (Some(Err(err)), Some(_)) => { - panic!("metadata 'convolution_particle_{index}' could not be parsed: {err}") - } - } - }) - .collect() - }, - ) -} diff --git a/pineappl_cli/tests/convolve.rs b/pineappl_cli/tests/convolve.rs index e20f894b..f9904fea 100644 --- a/pineappl_cli/tests/convolve.rs +++ b/pineappl_cli/tests/convolve.rs @@ -224,6 +224,95 @@ const NO_CHANNELS_GRID_STR: &str = 19 5 5 0.9 0.9 0.0000000e0 "; +const DYAA_SUBGRID_TEST: &str = "b x1 diff + [] [] +--+---+---+------------ + 0 0 0.1 5.4419693e-1 + 1 0.1 0.2 5.2673830e-1 + 2 0.2 0.3 5.5475101e-1 + 3 0.3 0.4 4.9474835e-1 + 4 0.4 0.5 4.6130464e-1 + 5 0.5 0.6 4.8110532e-1 + 6 0.6 0.7 4.7523978e-1 + 7 0.7 0.8 4.4445878e-1 + 8 0.8 0.9 4.2463353e-1 + 9 0.9 1 3.6203960e-1 +10 1 1.1 3.3418502e-1 +11 1.1 1.2 3.1297586e-1 +12 1.2 1.3 2.5858171e-1 +13 1.3 1.4 2.3601039e-1 +14 1.4 1.5 2.1117966e-1 +15 1.5 1.6 1.8150553e-1 +16 1.6 1.7 1.4598356e-1 +17 1.7 1.8 1.1183116e-1 +18 1.8 1.9 8.9304076e-2 +19 1.9 2 6.8197346e-2 +20 2 2.1 5.0815965e-2 +21 2.1 2.2 3.4465399e-2 +22 2.2 2.3 1.9973225e-2 +23 2.3 2.4 8.0989654e-3 +"; + +// TODO: why is this not the same as `DYAA_SUBGRID_V1_TEST`? Probably a bug in pineappl-0.4.1, +// which was used to generate the grid +const DYAA_SUBGRID_V1_TEST: &str = "b x1 diff + [] [] +--+---+---+------------ + 0 0 0.1 4.9421101e-1 + 1 0.1 0.2 5.3414733e-1 + 2 0.2 0.3 5.6047777e-1 + 3 0.3 0.4 4.9305119e-1 + 4 0.4 0.5 5.0444226e-1 + 5 0.5 0.6 4.8764547e-1 + 6 0.6 0.7 4.8822586e-1 + 7 0.7 0.8 4.3253468e-1 + 8 0.8 0.9 4.5631914e-1 + 9 0.9 1 4.4101667e-1 +10 1 1.1 3.7012152e-1 +11 1.1 1.2 3.2686684e-1 +12 1.2 1.3 2.8788858e-1 +13 1.3 1.4 2.5158100e-1 +14 1.4 1.5 1.9924476e-1 +15 1.5 1.6 1.6493339e-1 +16 1.6 1.7 1.5437771e-1 +17 1.7 1.8 1.2129587e-1 +18 1.8 1.9 9.4075605e-2 +19 1.9 2 6.2396273e-2 +20 2 2.1 4.9049203e-2 +21 2.1 2.2 3.5005153e-2 +22 2.2 2.3 2.0221966e-2 +23 2.3 2.4 5.9587519e-3 +"; + +const DYAA_IMPORT_ONLY_SUBGRID_V2_TEST: &str = "b x1 diff + [] [] +--+---+---+----------- + 0 0 0.1 2.6135143e2 + 1 0.1 0.2 2.6183087e2 + 2 0.2 0.3 2.6047594e2 + 3 0.3 0.4 2.6508580e2 + 4 0.4 0.5 2.5405348e2 + 5 0.5 0.6 2.4452454e2 + 6 0.6 0.7 2.4744331e2 + 7 0.7 0.8 2.4703713e2 + 8 0.8 0.9 2.3439723e2 + 9 0.9 1 2.6294172e2 +10 1 1.1 2.5651959e2 +11 1.1 1.2 1.9493319e2 +12 1.2 1.3 2.0817685e2 +13 1.3 1.4 1.8362573e2 +14 1.4 1.5 1.8185239e2 +15 1.5 1.6 1.5607187e2 +16 1.6 1.7 1.4492097e2 +17 1.7 1.8 9.8421742e1 +18 1.8 1.9 1.0014939e2 +19 1.9 2 7.8997062e1 +20 2 2.1 6.5452743e1 +21 2.1 2.2 5.2407224e1 +22 2.2 2.3 2.4198310e1 +23 2.3 2.4 8.6036652e0 +"; + #[test] fn help() { Command::cargo_bin("pineappl") @@ -463,7 +552,7 @@ fn three_convolutions() { } #[test] -fn no_channels_grid() { +fn issue_334() { Command::cargo_bin("pineappl") .unwrap() .args([ @@ -475,3 +564,59 @@ fn no_channels_grid() { .success() .stdout(NO_CHANNELS_GRID_STR); } + +#[test] +fn lagrange_subgrid_v1() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "convolve", + "../test-data/DYAA_0.4.1_LagrangeSubgridV1.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(DYAA_SUBGRID_V1_TEST); +} + +#[test] +fn lagrange_subgrid_v2() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "convolve", + "../test-data/DYAA_0.4.1_LagrangeSubgridV2.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(DYAA_SUBGRID_TEST); +} + +#[test] +fn lagrange_sparse_subgrid_v1() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "convolve", + "../test-data/DYAA_0.4.1_LagrangeSparseSubgridV1.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(DYAA_SUBGRID_TEST); +} + +#[test] +fn lagrange_subgrid_v1_mmv1() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "convolve", + "../test-data/DYAA_0.3.0_LagrangeSubgridV1.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(DYAA_SUBGRID_TEST); +} diff --git a/pineappl_v0/Cargo.toml b/pineappl_v0/Cargo.toml new file mode 100644 index 00000000..e7838b19 --- /dev/null +++ b/pineappl_v0/Cargo.toml @@ -0,0 +1,23 @@ +[package] +authors = ["Christopher Schwan "] +description = "PineAPPL is not an extension of APPLgrid" +name = "pineappl_v0" +readme = "README.md" + +categories.workspace = true +edition.workspace = true +keywords.workspace = true +license.workspace = true +repository.workspace = true +rust-version.workspace = true +version.workspace = true + +[lints] +workspace = true + +[dependencies] +bincode = "1.3.3" +enum_dispatch = "0.3.7" +ndarray = { features = ["serde"], version = "0.15.4" } +serde = { features = ["derive"], version = "1.0.130" } +thiserror = "1.0.30" diff --git a/pineappl_v0/README.md b/pineappl_v0/README.md new file mode 100644 index 00000000..f925b9c2 --- /dev/null +++ b/pineappl_v0/README.md @@ -0,0 +1,9 @@ +[![Rust](https://github.com/NNPDF/pineappl/workflows/Rust/badge.svg)](https://github.com/NNPDF/pineappl/actions?query=workflow%3ARust) +[![codecov](https://codecov.io/gh/NNPDF/pineappl/branch/master/graph/badge.svg)](https://codecov.io/gh/NNPDF/pineappl) +[![Documentation](https://docs.rs/pineappl/badge.svg)](https://docs.rs/pineappl) +[![crates.io](https://img.shields.io/crates/v/pineappl.svg)](https://crates.io/crates/pineappl) + +# PineAPPL + +PineAPPL is a library for recording and storing predictions for high-energy +physics observables independently of their parton distribution functions. diff --git a/pineappl_v0/src/bin.rs b/pineappl_v0/src/bin.rs new file mode 100644 index 00000000..6f764e9c --- /dev/null +++ b/pineappl_v0/src/bin.rs @@ -0,0 +1,160 @@ +//! Module that contains helpers for binning observables + +use super::convert::f64_from_usize; +use serde::Deserialize; +use std::f64; + +#[derive(Deserialize)] +enum Limits { + Equal { left: f64, right: f64, bins: usize }, + Unequal { limits: Vec }, +} + +/// Structure representing bin limits. +#[derive(Deserialize)] +pub struct BinLimits(Limits); + +/// Structure for remapping bin limits. +#[derive(Deserialize)] +pub struct BinRemapper { + normalizations: Vec, + limits: Vec<(f64, f64)>, +} + +/// Captures all information about the bins in a grid. +pub struct BinInfo<'a> { + limits: &'a BinLimits, + remapper: Option<&'a BinRemapper>, +} + +impl<'a> BinInfo<'a> { + /// Constructor. + #[must_use] + pub const fn new(limits: &'a BinLimits, remapper: Option<&'a BinRemapper>) -> Self { + Self { limits, remapper } + } + + /// Returns the number of bins. + #[must_use] + pub const fn bins(&self) -> usize { + self.limits.bins() + } + + /// Returns the number of dimensions. + #[must_use] + pub fn dimensions(&self) -> usize { + self.remapper.map_or(1, BinRemapper::dimensions) + } + + /// For each bin return a vector of `(left, right)` limits for each dimension. + #[must_use] + pub fn limits(&self) -> Vec> { + self.remapper.map_or_else( + || { + self.limits + .limits() + .windows(2) + .map(|window| vec![(window[0], window[1])]) + .collect() + }, + |remapper| { + remapper + .limits() + .to_vec() + .chunks_exact(self.dimensions()) + .map(<[(f64, f64)]>::to_vec) + .collect() + }, + ) + } + + /// Returns all normalization factors. + #[must_use] + pub fn normalizations(&self) -> Vec { + self.remapper.map_or_else( + || self.limits.bin_sizes(), + |remapper| remapper.normalizations().to_vec(), + ) + } +} + +impl BinRemapper { + /// Return the number of dimensions. + #[must_use] + pub const fn dimensions(&self) -> usize { + self.limits.len() / self.normalizations.len() + } + + /// Return tuples of left and right bin limits for all dimensions and all bins. + #[must_use] + pub fn limits(&self) -> &[(f64, f64)] { + &self.limits + } + + /// Return the normalization factors for all bins. + #[must_use] + pub fn normalizations(&self) -> &[f64] { + &self.normalizations + } +} + +impl BinLimits { + /// Returns the number of bins. + #[must_use] + pub const fn bins(&self) -> usize { + match &self.0 { + Limits::Equal { bins, .. } => *bins, + Limits::Unequal { limits } => limits.len() - 1, + } + } + + /// Returns the limits in a `Vec`. + /// + /// # Examples + /// + /// ```rust + /// use pineappl_v0::bin::BinLimits; + /// + /// // example with equally sized bins + /// let equal_bins = BinLimits::new(vec![0.25, 0.5, 0.75, 1.0]); + /// assert_eq!(equal_bins.limits(), vec![0.25, 0.5, 0.75, 1.0]); + /// + /// // example with unequally sized bins + /// let unequal_bins = BinLimits::new(vec![0.125, 0.25, 1.0, 1.5]); + /// assert_eq!(unequal_bins.limits(), vec![0.125, 0.25, 1.0, 1.5]); + /// ``` + #[must_use] + pub fn limits(&self) -> Vec { + match &self.0 { + Limits::Equal { left, right, bins } => (0..=*bins) + .map(|b| (*right - *left).mul_add(f64_from_usize(b) / f64_from_usize(*bins), *left)) + .collect(), + Limits::Unequal { limits } => limits.clone(), + } + } + + /// Returns the size for each bin. + /// + /// # Examples + /// + /// ```rust + /// use pineappl_v0::bin::BinLimits; + /// + /// // example with equally sized bins + /// let equal_bins = BinLimits::new(vec![0.25, 0.5, 0.75, 1.0]); + /// assert_eq!(equal_bins.bin_sizes(), vec![0.25, 0.25, 0.25]); + /// + /// // example with unequally sized bins + /// let unequal_bins = BinLimits::new(vec![0.125, 0.25, 1.0, 1.5]); + /// assert_eq!(unequal_bins.bin_sizes(), vec![0.125, 0.75, 0.5]); + /// ``` + #[must_use] + pub fn bin_sizes(&self) -> Vec { + match &self.0 { + Limits::Equal { left, right, bins } => { + vec![(*right - *left) / f64_from_usize(*bins); *bins] + } + Limits::Unequal { limits } => limits.windows(2).map(|x| x[1] - x[0]).collect(), + } + } +} diff --git a/pineappl_v0/src/boc.rs b/pineappl_v0/src/boc.rs new file mode 100644 index 00000000..7c811dbb --- /dev/null +++ b/pineappl_v0/src/boc.rs @@ -0,0 +1,45 @@ +//! Module containing structures for the 3 dimensions of a [`Grid`]: bins, [`Order`] and channels +//! (`boc`). + +use serde::Deserialize; + +/// Coupling powers for each grid. +#[derive(Deserialize)] +pub struct Order { + /// Exponent of the strong coupling. + pub alphas: u32, + /// Exponent of the electromagnetic coupling. + pub alpha: u32, + /// Exponent of the logarithm of the scale factor of the renomalization scale. + pub logxir: u32, + /// Exponent of the logarithm of the scale factor of the factorization scale. + pub logxif: u32, +} + +/// This structure represents a channel. Each channel consists of a tuple containing in the +/// following order, the particle ID of the first incoming parton, then the particle ID of the +/// second parton, and finally a numerical factor that will multiply the result for this specific +/// combination. +#[derive(Deserialize)] +pub struct Channel { + entry: Vec<(i32, i32, f64)>, +} + +impl Channel { + /// Returns a tuple representation of this entry. + /// + /// # Examples + /// + /// ```rust + /// use pineappl_v0::channel; + /// use pineappl_v0::boc::Channel; + /// + /// let entry = channel![4, 4, 1.0; 2, 2, 1.0]; + /// + /// assert_eq!(entry.entry(), [(2, 2, 1.0), (4, 4, 1.0)]); + /// ``` + #[must_use] + pub fn entry(&self) -> &[(i32, i32, f64)] { + &self.entry + } +} diff --git a/pineappl_v0/src/convert.rs b/pineappl_v0/src/convert.rs new file mode 100644 index 00000000..724e13c8 --- /dev/null +++ b/pineappl_v0/src/convert.rs @@ -0,0 +1,3 @@ +pub fn f64_from_usize(x: usize) -> f64 { + f64::from(u32::try_from(x).unwrap()) +} diff --git a/pineappl_v0/src/convolutions.rs b/pineappl_v0/src/convolutions.rs new file mode 100644 index 00000000..4fd82a34 --- /dev/null +++ b/pineappl_v0/src/convolutions.rs @@ -0,0 +1,19 @@ +//! Module for everything related to luminosity functions. + +/// Data type that indentifies different types of convolutions. +pub enum Convolution { + // TODO: eventually get rid of this value + /// No convolution. + None, + /// Unpolarized parton distribution function. The integer denotes the type of hadron with a PDG + /// MC ID. + UnpolPDF(i32), + /// Polarized parton distribution function. The integer denotes the type of hadron with a PDG + /// MC ID. + PolPDF(i32), + /// Unpolarized fragmentation function. The integer denotes the type of hadron with a PDG MC + /// ID. + UnpolFF(i32), + /// Polarized fragmentation function. The integer denotes the type of hadron with a PDG MC ID. + PolFF(i32), +} diff --git a/pineappl_v0/src/empty_subgrid.rs b/pineappl_v0/src/empty_subgrid.rs new file mode 100644 index 00000000..96d36b3c --- /dev/null +++ b/pineappl_v0/src/empty_subgrid.rs @@ -0,0 +1,32 @@ +//! TODO + +use super::subgrid::{Mu2, Subgrid, SubgridIndexedIter}; +use serde::Deserialize; +use std::borrow::Cow; +use std::iter; + +/// A subgrid type that is always empty. +#[derive(Deserialize)] +pub struct EmptySubgridV1; + +impl Subgrid for EmptySubgridV1 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + Cow::Borrowed(&[]) + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&[]) + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&[]) + } + + fn is_empty(&self) -> bool { + true + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + Box::new(iter::empty()) + } +} diff --git a/pineappl_v0/src/grid.rs b/pineappl_v0/src/grid.rs new file mode 100644 index 00000000..02c1808c --- /dev/null +++ b/pineappl_v0/src/grid.rs @@ -0,0 +1,228 @@ +//! Module containing all traits and supporting structures for grids. + +use super::bin::{BinInfo, BinLimits, BinRemapper}; +use super::boc::{Channel, Order}; +use super::convolutions::Convolution; +use super::pids::PidBasis; +use super::subgrid::{SubgridEnum, SubgridParams}; +use ndarray::{Array3, ArrayView3}; +use serde::Deserialize; +use std::collections::HashMap; +use std::io::BufRead; +use thiserror::Error; + +/// This structure represents a position (`x1`, `x2`, `q2`) in a `Subgrid` together with a +/// corresponding `weight`. The type `W` can either be a `f64` or `()`, which is used when multiple +/// weights should be signaled. +#[derive(Deserialize)] +pub struct Ntuple { + /// Momentum fraction of the first parton. + pub x1: f64, + /// Momentum fraction of the second parton. + pub x2: f64, + /// Squared scale. + pub q2: f64, + /// Weight of this entry. + pub weight: W, +} + +/// Error returned when merging two grids fails. +#[derive(Debug, Error)] +pub enum GridError { + /// Returned when failed to read a Grid. + #[error(transparent)] + ReadFailure(bincode::Error), +} + +#[derive(Deserialize)] +struct Mmv1; + +#[derive(Deserialize)] +struct Mmv2 { + remapper: Option, + key_value_db: HashMap, +} + +#[derive(Deserialize)] +struct Mmv3 { + remapper: Option, + key_value_db: HashMap, + _subgrid_template: SubgridEnum, +} + +// ALLOW: fixing the warning will break the file format +#[allow(clippy::large_enum_variant)] +#[derive(Deserialize)] +enum MoreMembers { + V1(Mmv1), + V2(Mmv2), + V3(Mmv3), +} + +/// Main data structure of `PineAPPL`. This structure contains a `Subgrid` for each `LumiEntry`, +/// bin, and coupling order it was created with. +#[derive(Deserialize)] +pub struct Grid { + subgrids: Array3, + channels: Vec, + bin_limits: BinLimits, + orders: Vec, + _subgrid_params: SubgridParams, + more_members: MoreMembers, +} + +impl Grid { + /// Return by which convention the particle IDs are encoded. + #[must_use] + pub fn pid_basis(&self) -> PidBasis { + if let Some(key_values) = self.key_values() + && let Some(lumi_id_types) = key_values.get("lumi_id_types") + { + match lumi_id_types.as_str() { + "pdg_mc_ids" => return PidBasis::Pdg, + "evol" => return PidBasis::Evol, + _ => unimplemented!("unknown particle ID convention {lumi_id_types}"), + } + } + + // if there's no basis explicitly set we're assuming to use PDG IDs + PidBasis::Pdg + } + + /// Construct a `Grid` by deserializing it from `reader`. + /// + /// # Errors + /// + /// If reading from the compressed or uncompressed stream fails an error is returned. + /// + /// # Panics + /// + /// Panics if the grid version is not `0`. + pub fn read_uncompressed(reader: impl BufRead) -> Result { + bincode::deserialize_from(reader).map_err(GridError::ReadFailure) + } + + /// Return the channels for this `Grid`. + #[must_use] + pub fn channels(&self) -> &[Channel] { + &self.channels + } + + /// Return a vector containing the type of convolutions performed with this grid. + /// + /// # Panics + /// + /// Panics if the metadata key--value pairs `convolution_particle_1` and `convolution_type_1`, + /// or `convolution_particle_2` and `convolution_type_2` are not correctly set. + #[must_use] + pub fn convolutions(&self) -> Vec { + self.key_values().map_or_else( + // if there isn't any metadata, we assume two unpolarized proton-PDFs are used + || vec![Convolution::UnpolPDF(2212), Convolution::UnpolPDF(2212)], + |kv| { + // the current file format only supports exactly two convolutions + (1..=2) + .map(|index| { + // if there are key-value pairs `convolution_particle_1` and + // `convolution_type_1` and the same with a higher index, we convert this + // metadata into `Convolution` + match ( + kv.get(&format!("convolution_particle_{index}")) + .map(|s| s.parse::()), + kv.get(&format!("convolution_type_{index}")) + .map(String::as_str), + ) { + (_, Some("None")) => Convolution::None, + (Some(Ok(pid)), Some("UnpolPDF")) => Convolution::UnpolPDF(pid), + (Some(Ok(pid)), Some("PolPDF")) => Convolution::PolPDF(pid), + (Some(Ok(pid)), Some("UnpolFF")) => Convolution::UnpolFF(pid), + (Some(Ok(pid)), Some("PolFF")) => Convolution::PolFF(pid), + (None, None) => { + // if these key-value pairs are missing use the old metadata + match kv + .get(&format!("initial_state_{index}")) + .map(|s| s.parse::()) + { + Some(Ok(pid)) => { + // Addresses: https://github.com/NNPDF/pineappl/issues/334 + if self.channels().is_empty() && pid == 2212 { + Convolution::UnpolPDF(2212) + } else { + let condition = !self.channels().iter().all(|entry| { + entry.entry().iter().all(|&(a, b, _)| + match index { + 1 => a, + 2 => b, + _ => unreachable!() + } == pid + ) + }); + + if condition { + Convolution::UnpolPDF(pid) + } else { + Convolution::None + } + } + } + None => Convolution::UnpolPDF(2212), + Some(Err(err)) => panic!("metadata 'initial_state_{index}' could not be parsed: {err}"), + } + } + (None, Some(_)) => { + panic!("metadata 'convolution_type_{index}' is missing") + } + (Some(_), None) => { + panic!("metadata 'convolution_particle_{index}' is missing") + } + (Some(Ok(_)), Some(type_)) => { + panic!("metadata 'convolution_type_{index} = {type_}' is unknown") + } + (Some(Err(err)), Some(_)) => panic!( + "metadata 'convolution_particle_{index}' could not be parsed: {err}" + ), + } + }) + .collect() + }, + ) + } + + /// Returns the subgrid parameters. + #[must_use] + pub fn orders(&self) -> &[Order] { + &self.orders + } + + /// Return all subgrids as an `ArrayView3`. + #[must_use] + pub fn subgrids(&self) -> ArrayView3<'_, SubgridEnum> { + self.subgrids.view() + } + + /// Return the currently set remapper, if there is any. + #[must_use] + pub const fn remapper(&self) -> Option<&BinRemapper> { + match &self.more_members { + MoreMembers::V1(_) => None, + MoreMembers::V2(mmv2) => mmv2.remapper.as_ref(), + MoreMembers::V3(mmv3) => mmv3.remapper.as_ref(), + } + } + + /// Returns all information about the bins in this grid. + #[must_use] + pub const fn bin_info(&self) -> BinInfo<'_> { + BinInfo::new(&self.bin_limits, self.remapper()) + } + + /// Returns a map with key-value pairs, if there are any stored in this grid. + #[must_use] + pub const fn key_values(&self) -> Option<&HashMap> { + match &self.more_members { + MoreMembers::V3(mmv3) => Some(&mmv3.key_value_db), + MoreMembers::V2(mmv2) => Some(&mmv2.key_value_db), + MoreMembers::V1(_) => None, + } + } +} diff --git a/pineappl_v0/src/import_only_subgrid.rs b/pineappl_v0/src/import_only_subgrid.rs new file mode 100644 index 00000000..b91e9b93 --- /dev/null +++ b/pineappl_v0/src/import_only_subgrid.rs @@ -0,0 +1,72 @@ +//! TODO + +use super::sparse_array3::SparseArray3; +use super::subgrid::{Mu2, Subgrid, SubgridIndexedIter}; +use serde::Deserialize; +use std::borrow::Cow; + +/// TODO +#[derive(Deserialize)] +pub struct ImportOnlySubgridV1 { + array: SparseArray3, + q2_grid: Vec, + x1_grid: Vec, + x2_grid: Vec, +} + +impl Subgrid for ImportOnlySubgridV1 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + self.q2_grid + .iter() + .copied() + .map(|q2| Mu2 { ren: q2, fac: q2 }) + .collect() + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&self.x1_grid) + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&self.x2_grid) + } + + fn is_empty(&self) -> bool { + self.array.is_empty() + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + Box::new(self.array.indexed_iter()) + } +} + +/// TODO +#[derive(Deserialize)] +pub struct ImportOnlySubgridV2 { + array: SparseArray3, + mu2_grid: Vec, + x1_grid: Vec, + x2_grid: Vec, +} + +impl Subgrid for ImportOnlySubgridV2 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + Cow::Borrowed(&self.mu2_grid) + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&self.x1_grid) + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&self.x2_grid) + } + + fn is_empty(&self) -> bool { + self.array.is_empty() + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + Box::new(self.array.indexed_iter()) + } +} diff --git a/pineappl_v0/src/lagrange_subgrid.rs b/pineappl_v0/src/lagrange_subgrid.rs new file mode 100644 index 00000000..23dce67c --- /dev/null +++ b/pineappl_v0/src/lagrange_subgrid.rs @@ -0,0 +1,297 @@ +//! Module containing the Lagrange-interpolation subgrid. + +use super::convert::f64_from_usize; +use super::sparse_array3::SparseArray3; +use super::subgrid::{Mu2, Subgrid, SubgridIndexedIter}; +use ndarray::Array3; +use serde::Deserialize; +use std::borrow::Cow; +use std::iter; + +fn weightfun(x: f64) -> f64 { + (x.sqrt() / (1.0 - 0.99 * x)).powi(3) +} + +fn fx(y: f64) -> f64 { + let mut yp = y; + + for _ in 0..100 { + let x = (-yp).exp(); + let delta = y - yp - 5.0 * (1.0 - x); + if (delta).abs() < 1e-12 { + return x; + } + let deriv = -1.0 - 5.0 * x; + yp -= delta / deriv; + } + + unreachable!(); +} + +fn fq2(tau: f64) -> f64 { + 0.0625 * tau.exp().exp() +} + +/// Subgrid which uses Lagrange-interpolation. +#[derive(Deserialize)] +pub struct LagrangeSubgridV1 { + grid: Option>, + ntau: usize, + ny: usize, + _yorder: usize, + _tauorder: usize, + itaumin: usize, + _itaumax: usize, + reweight: bool, + ymin: f64, + ymax: f64, + taumin: f64, + taumax: f64, +} + +impl LagrangeSubgridV1 { + fn deltay(&self) -> f64 { + (self.ymax - self.ymin) / f64_from_usize(self.ny - 1) + } + + fn deltatau(&self) -> f64 { + (self.taumax - self.taumin) / f64_from_usize(self.ntau - 1) + } + + fn gety(&self, iy: usize) -> f64 { + f64_from_usize(iy).mul_add(self.deltay(), self.ymin) + } + + fn gettau(&self, iy: usize) -> f64 { + f64_from_usize(iy).mul_add(self.deltatau(), self.taumin) + } +} + +impl Subgrid for LagrangeSubgridV1 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + (0..self.ntau) + .map(|itau| { + let q2 = fq2(self.gettau(itau)); + Mu2 { ren: q2, fac: q2 } + }) + .collect() + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + (0..self.ny).map(|iy| fx(self.gety(iy))).collect() + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + self.x1_grid() + } + + fn is_empty(&self) -> bool { + self.grid.is_none() + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + self.grid.as_ref().map_or_else( + || Box::new(iter::empty()) as Box>, + |grid| { + Box::new(grid.indexed_iter().filter(|&(_, &value)| value != 0.0).map( + |(tuple, &value)| { + ( + (self.itaumin + tuple.0, tuple.1, tuple.2), + value + * if self.reweight { + weightfun(fx(self.gety(tuple.1))) + * weightfun(fx(self.gety(tuple.2))) + } else { + 1.0 + }, + ) + }, + )) + }, + ) + } +} + +/// Subgrid which uses Lagrange-interpolation. +#[derive(Deserialize)] +pub struct LagrangeSubgridV2 { + grid: Option>, + ntau: usize, + ny1: usize, + ny2: usize, + _y1order: usize, + _y2order: usize, + _tauorder: usize, + itaumin: usize, + _itaumax: usize, + reweight1: bool, + reweight2: bool, + y1min: f64, + y1max: f64, + y2min: f64, + y2max: f64, + taumin: f64, + taumax: f64, + pub(crate) _static_q2: f64, +} + +impl LagrangeSubgridV2 { + fn deltay1(&self) -> f64 { + (self.y1max - self.y1min) / f64_from_usize(self.ny1 - 1) + } + + fn deltay2(&self) -> f64 { + (self.y1max - self.y2min) / f64_from_usize(self.ny2 - 1) + } + + fn deltatau(&self) -> f64 { + (self.taumax - self.taumin) / f64_from_usize(self.ntau - 1) + } + + fn gety1(&self, iy: usize) -> f64 { + if self.y1min == self.y1max { + debug_assert_eq!(iy, 0); + self.y1min + } else { + f64_from_usize(iy).mul_add(self.deltay1(), self.y1min) + } + } + + fn gety2(&self, iy: usize) -> f64 { + if self.y2min == self.y2max { + debug_assert_eq!(iy, 0); + self.y2min + } else { + f64_from_usize(iy).mul_add(self.deltay2(), self.y2min) + } + } + + fn gettau(&self, iy: usize) -> f64 { + if self.taumin == self.taumax { + debug_assert_eq!(iy, 0); + self.taumin + } else { + f64_from_usize(iy).mul_add(self.deltatau(), self.taumin) + } + } +} + +impl Subgrid for LagrangeSubgridV2 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + (0..self.ntau) + .map(|itau| { + let q2 = fq2(self.gettau(itau)); + Mu2 { ren: q2, fac: q2 } + }) + .collect() + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + (0..self.ny1).map(|iy| fx(self.gety1(iy))).collect() + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + (0..self.ny2).map(|iy| fx(self.gety2(iy))).collect() + } + + fn is_empty(&self) -> bool { + self.grid.is_none() + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + self.grid.as_ref().map_or_else( + || Box::new(iter::empty()) as Box>, + |grid| { + Box::new(grid.indexed_iter().filter(|&(_, &value)| value != 0.0).map( + |(tuple, &value)| { + ( + (self.itaumin + tuple.0, tuple.1, tuple.2), + value + * if self.reweight1 { + weightfun(fx(self.gety1(tuple.1))) + } else { + 1.0 + } + * if self.reweight2 { + weightfun(fx(self.gety2(tuple.2))) + } else { + 1.0 + }, + ) + }, + )) + }, + ) + } +} + +/// Subgrid which uses Lagrange-interpolation, but also stores its contents in a space-efficient +/// structure. +#[derive(Deserialize)] +pub struct LagrangeSparseSubgridV1 { + array: SparseArray3, + ntau: usize, + ny: usize, + _yorder: usize, + _tauorder: usize, + reweight: bool, + ymin: f64, + ymax: f64, + taumin: f64, + taumax: f64, +} + +impl LagrangeSparseSubgridV1 { + fn deltay(&self) -> f64 { + (self.ymax - self.ymin) / f64_from_usize(self.ny - 1) + } + + fn deltatau(&self) -> f64 { + (self.taumax - self.taumin) / f64_from_usize(self.ntau - 1) + } + + fn gety(&self, iy: usize) -> f64 { + f64_from_usize(iy).mul_add(self.deltay(), self.ymin) + } + + fn gettau(&self, iy: usize) -> f64 { + f64_from_usize(iy).mul_add(self.deltatau(), self.taumin) + } +} + +impl Subgrid for LagrangeSparseSubgridV1 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + (0..self.ntau) + .map(|itau| { + let q2 = fq2(self.gettau(itau)); + Mu2 { ren: q2, fac: q2 } + }) + .collect() + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + (0..self.ny).map(|iy| fx(self.gety(iy))).collect() + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + self.x1_grid() + } + + fn is_empty(&self) -> bool { + self.array.is_empty() + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + Box::new(self.array.indexed_iter().map(|(tuple, value)| { + ( + tuple, + value + * if self.reweight { + weightfun(fx(self.gety(tuple.1))) * weightfun(fx(self.gety(tuple.2))) + } else { + 1.0 + }, + ) + })) + } +} diff --git a/pineappl_v0/src/lib.rs b/pineappl_v0/src/lib.rs new file mode 100644 index 00000000..5d9db49a --- /dev/null +++ b/pineappl_v0/src/lib.rs @@ -0,0 +1,48 @@ +//! `PineAPPL` is not an extension of `APPLgrid`. +//! +//! # Overview +//! +//! The main type of this crate is [`Grid`], which represents the interpolation grids that +//! `PineAPPL` implements. Roughly speaking, a `Grid` is a three-dimensional array of [`Subgrid`] +//! objects together with metadata. The three dimensions are +//! 1. (perturbative) orders, represented by the type [`Order`] and accessible by +//! [`Grid::orders()`], +//! 2. bins, whose limits can be accessed by [`Grid::bin_info()`], and +//! 3. channels, whose definition is returned by [`Grid::channels()`]. +//! +//! `Subgrid` is a `trait` and objects that implement it are of the type [`SubgridEnum`]. The +//! latter is an `enum` of different types that are optimized to different scenarios: fast event +//! filling, small storage profile, etc. +//! +//! [`Grid`]: grid::Grid +//! [`Grid::bin_info()`]: grid::Grid::bin_info +//! [`Grid::channels()`]: grid::Grid::channels +//! [`Grid::orders()`]: grid::Grid::orders +//! [`Subgrid`]: subgrid::Subgrid +//! [`SubgridEnum`]: subgrid::SubgridEnum +//! [`Order`]: order::Order +//! +//! ## Metadata +//! +//! Metadata is a collection of key--value pairs, in which both keys and values are `String` +//! objects. In metadata anything a user whishes can be stored. However, there are [special keys], +//! which have meaning to `PineAPPL` and/or its CLI `pineappl`. This metadata enables the CLI to +//! automatically generate plots that are correctly labeled, for instance. For more applications +//! see also the [CLI tutorial]. +//! +//! [special keys]: https://nnpdf.github.io/pineappl/docs/metadata.html +//! [CLI tutorial]: https://nnpdf.github.io/pineappl/docs/cli-tutorial.html + +mod convert; + +pub mod bin; +pub mod boc; +pub mod convolutions; +pub mod empty_subgrid; +pub mod grid; +pub mod import_only_subgrid; +pub mod lagrange_subgrid; +pub mod ntuple_subgrid; +pub mod pids; +pub mod sparse_array3; +pub mod subgrid; diff --git a/pineappl_v0/src/ntuple_subgrid.rs b/pineappl_v0/src/ntuple_subgrid.rs new file mode 100644 index 00000000..56e8cf7d --- /dev/null +++ b/pineappl_v0/src/ntuple_subgrid.rs @@ -0,0 +1,34 @@ +//! Provides an implementation of the `Grid` trait with n-tuples. + +use super::grid::Ntuple; +use super::subgrid::{Mu2, Subgrid, SubgridIndexedIter}; +use serde::Deserialize; +use std::borrow::Cow; + +/// Structure holding a grid with an n-tuple as the storage method for weights. +#[derive(Deserialize)] +pub struct NtupleSubgridV1 { + ntuples: Vec>, +} + +impl Subgrid for NtupleSubgridV1 { + fn mu2_grid(&self) -> Cow<'_, [Mu2]> { + Cow::Borrowed(&[]) + } + + fn x1_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&[]) + } + + fn x2_grid(&self) -> Cow<'_, [f64]> { + Cow::Borrowed(&[]) + } + + fn is_empty(&self) -> bool { + self.ntuples.is_empty() + } + + fn indexed_iter(&self) -> SubgridIndexedIter<'_> { + panic!("NtupleSubgridV1 doesn't support the indexed_iter operation"); + } +} diff --git a/pineappl_v0/src/pids.rs b/pineappl_v0/src/pids.rs new file mode 100644 index 00000000..720edcca --- /dev/null +++ b/pineappl_v0/src/pids.rs @@ -0,0 +1,15 @@ +//! TODO + +/// Particle ID bases. In `PineAPPL` every particle is identified using a particle identifier +/// (PID), which is represented as an `i32`. The values of this `enum` specify how this value is +/// interpreted. +pub enum PidBasis { + /// This basis uses the [particle data group](https://pdg.lbl.gov/) (PDG) PIDs. For a complete + /// definition see the section 'Monte Carlo Particle Numbering Scheme' of the PDG Review, for + /// instance the [2023 review](https://pdg.lbl.gov/2023/mcdata/mc_particle_id_contents.html). + Pdg, + /// This basis specifies the evolution basis, which is the same as [`PidBasis::Pdg`], except + /// the following values have a special meaning: `100`, `103`, `108`, `115`, `124`, `135`, + /// `200`, `203`, `208`, `215`, `224`, `235`. + Evol, +} diff --git a/pineappl_v0/src/sparse_array3.rs b/pineappl_v0/src/sparse_array3.rs new file mode 100644 index 00000000..77f99a15 --- /dev/null +++ b/pineappl_v0/src/sparse_array3.rs @@ -0,0 +1,132 @@ +//! Module containing the `SparseArray3` struct. + +use serde::Deserialize; +use std::slice::Iter; + +/// Struct for a sparse three-dimensional array, which is optimized for the sparsity of +/// interpolation grids. +#[derive(Deserialize)] +pub struct SparseArray3 { + entries: Vec, + indices: Vec<(usize, usize)>, + start: usize, + dimensions: (usize, usize, usize), +} + +/// Immutable iterator over the elements of a `SparseArray3`. +pub struct IndexedIter<'a, T> { + entry_iter: Iter<'a, T>, + index_iter: Iter<'a, (usize, usize)>, + offset_a: Option<&'a (usize, usize)>, + offset_b: Option<&'a (usize, usize)>, + tuple: (usize, usize, usize), + dimensions: (usize, usize, usize), +} + +impl Iterator for IndexedIter<'_, T> { + type Item = ((usize, usize, usize), T); + + fn next(&mut self) -> Option { + if let Some(element) = self.entry_iter.next() { + let offset_a = self.offset_a.unwrap(); + let offset_b = self.offset_b.unwrap(); + + if self.dimensions.1 > self.dimensions.2 { + self.tuple.1 = self.tuple.1.max(offset_a.0); + + if self.tuple.1 >= (offset_b.1 - offset_a.1 + offset_a.0) { + loop { + self.offset_a = self.offset_b; + self.offset_b = self.index_iter.next(); + + let offset_a = self.offset_a.unwrap(); + let offset_b = self.offset_b?; + + self.tuple.2 += 1; + + if self.tuple.2 >= self.dimensions.2 { + self.tuple.0 += 1; + self.tuple.2 = 0; + } + + if (offset_b.1 - offset_a.1) != 0 { + self.tuple.1 = offset_a.0; + break; + } + } + } + + if *element == T::default() { + self.tuple.1 += 1; + self.next() + } else { + let result = Some((self.tuple, *element)); + self.tuple.1 += 1; + result + } + } else { + self.tuple.2 = self.tuple.2.max(offset_a.0); + + if self.tuple.2 >= (offset_b.1 - offset_a.1 + offset_a.0) { + loop { + self.offset_a = self.offset_b; + self.offset_b = self.index_iter.next(); + + let offset_a = self.offset_a.unwrap(); + let offset_b = self.offset_b?; + + self.tuple.1 += 1; + + if self.tuple.1 >= self.dimensions.1 { + self.tuple.0 += 1; + self.tuple.1 = 0; + } + + if (offset_b.1 - offset_a.1) != 0 { + self.tuple.2 = offset_a.0; + break; + } + } + } + + if *element == T::default() { + self.tuple.2 += 1; + self.next() + } else { + let result = Some((self.tuple, *element)); + self.tuple.2 += 1; + result + } + } + } else { + None + } + } +} + +impl SparseArray3 { + /// Returns `true` if the array contains no element. + #[must_use] + pub const fn is_empty(&self) -> bool { + self.entries.is_empty() + } + + /// Return an indexed `Iterator` over the non-zero elements of this array. The iterator element + /// type is `((usize, usize, usize), T)`. + #[must_use] + pub fn indexed_iter(&self) -> IndexedIter<'_, T> { + let mut result = IndexedIter { + entry_iter: self.entries.iter(), + index_iter: self.indices.iter(), + offset_a: None, + offset_b: None, + tuple: (self.start, 0, 0), + dimensions: self.dimensions, + }; + + result.offset_a = result.index_iter.next(); + result.offset_b = result.index_iter.next(); + + result + } +} diff --git a/pineappl_v0/src/subgrid.rs b/pineappl_v0/src/subgrid.rs new file mode 100644 index 00000000..a9819436 --- /dev/null +++ b/pineappl_v0/src/subgrid.rs @@ -0,0 +1,81 @@ +//! Module containing the trait `Subgrid` and supporting structs. + +use super::empty_subgrid::EmptySubgridV1; +use super::import_only_subgrid::{ImportOnlySubgridV1, ImportOnlySubgridV2}; +use super::lagrange_subgrid::{LagrangeSparseSubgridV1, LagrangeSubgridV1, LagrangeSubgridV2}; +use super::ntuple_subgrid::NtupleSubgridV1; +use enum_dispatch::enum_dispatch; +use serde::Deserialize; +use std::borrow::Cow; + +/// Enum which lists all possible `Subgrid` variants possible. +#[enum_dispatch(Subgrid)] +#[derive(Deserialize)] +pub enum SubgridEnum { + // WARNING: never change the order or content of this enum, only add to the end of it + /// Lagrange-interpolation subgrid. + LagrangeSubgridV1, + /// N-tuple subgrid. + NtupleSubgridV1, + /// Lagrange-interpolation subgrid. + LagrangeSparseSubgridV1, + /// Lagrange-interpolation subgrid with possibly different x1 and x2 bins. + LagrangeSubgridV2, + /// Import-only sparse subgrid with possibly different x1 and x2 bins. + ImportOnlySubgridV1, + /// Empty subgrid. + EmptySubgridV1, + /// Same as [`ImportOnlySubgridV1`], but with support for different renormalization and + /// factorization scales choices. + ImportOnlySubgridV2, +} + +/// Structure denoting renormalization and factorization scale values. +#[derive(Deserialize, Clone)] +pub struct Mu2 { + /// The (squared) renormalization scale value. + pub ren: f64, + /// The (squared) factorization scale value. + pub fac: f64, +} + +/// Trait each subgrid must implement. +#[enum_dispatch] +pub trait Subgrid { + /// Return a slice of [`Mu2`] values corresponding to the (squared) renormalization and + /// factorization values of the grid. If the subgrid does not use a grid, this method should + /// return an empty slice. + fn mu2_grid(&self) -> Cow<'_, [Mu2]>; + + /// Return a slice of values of `x1`. If the subgrid does not use a grid, this method should + /// return an empty slice. + fn x1_grid(&self) -> Cow<'_, [f64]>; + + /// Return a slice of values of `x2`. If the subgrid does not use a grid, this method should + /// return an empty slice. + fn x2_grid(&self) -> Cow<'_, [f64]>; + + /// Returns true if `fill` was never called for this grid. + fn is_empty(&self) -> bool; + + /// Return an iterator over all non-zero elements of the subgrid. + fn indexed_iter(&self) -> SubgridIndexedIter<'_>; +} + +/// Type to iterate over the non-zero contents of a subgrid. The tuple contains the indices of the +/// `mu2_grid`, the `x1_grid` and finally the `x2_grid`. +pub type SubgridIndexedIter<'a> = Box + 'a>; + +/// Subgrid creation parameters for subgrids that perform interpolation. +#[derive(Deserialize)] +pub struct SubgridParams { + _q2_bins: usize, + _q2_max: f64, + _q2_min: f64, + _q2_order: usize, + _reweight: bool, + _x_bins: usize, + _x_max: f64, + _x_min: f64, + _x_order: usize, +}