diff --git a/.github/workflows/rust.yml b/.github/workflows/rust.yml index 3ab65727..ed6c64fa 100644 --- a/.github/workflows/rust.yml +++ b/.github/workflows/rust.yml @@ -36,6 +36,7 @@ jobs: # `-C link-dead-code` is needed to prevent 'warning: XX functions have mismatched data' warnings RUSTFLAGS: '-Cinstrument-coverage -Clink-dead-code' run: | + export NEOPDF_DATA_PATH=/usr/local/share/LHAPDF # we need stderr, but we can't run test twice because it'll regenerate/modify the binaries which interferes with `llvm-cov` cargo test --features=applgrid,evolve,fastnlo,fktable --no-fail-fast 2> >(tee stderr 1>&2) # from https://stackoverflow.com/a/51141872/812178 diff --git a/CHANGELOG.md b/CHANGELOG.md index 12004603..a0006352 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,7 +9,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added -- Added a `pineappl_grid_convolution_pids` to the C API +- added a `pineappl_grid_convolution_pids` to the C API +- added support for the [NeoPDF](https://qcdlab.github.io/neopdf/) library as + an alternative to LHAPDF in the CLI; use `--backend=neopdf` to switch to it ## [1.4.2] - 31/05/2026 diff --git a/Cargo.lock b/Cargo.lock index dd682171..36de53bc 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -3,10 +3,10 @@ version = 4 [[package]] -name = "adler" -version = "1.0.2" +name = "adler2" +version = "2.0.1" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "f26201604c87b1e01bd3d98f8d5d9a8fcbb815e8cedb41ffccbeb4bf593a35fe" +checksum = "320119579fcad9c21884f5c4861d16174d0e06250625266f50fe6898340abefa" [[package]] name = "aho-corasick" @@ -245,7 +245,7 @@ checksum = "af491d569909a7e4dee0ad7db7f5341fef5c614d5b8ec8cf765732aba3cff681" dependencies = [ "serde", "termcolor", - "unicode-width", + "unicode-width 0.2.2", ] [[package]] @@ -254,6 +254,17 @@ version = "1.0.5" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "1d07550c9036bf2ae0c684c4297d503f838287c83c53686d05370d0e139ae570" +[[package]] +name = "console" +version = "0.16.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d64e8af5551369d19cf50138de61f1c42074ab970f74e99be916646777f8fc87" +dependencies = [ + "encode_unicode", + "libc", + "windows-sys 0.61.2", +] + [[package]] name = "cpufeatures" version = "0.2.17" @@ -438,6 +449,12 @@ dependencies = [ "winapi", ] +[[package]] +name = "dyn-clone" +version = "1.0.20" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "d0881ea181b1df73ff77ffaaf9c7544ecc11e82fba9b5f27b262a3c73a332555" + [[package]] name = "either" version = "1.15.0" @@ -478,7 +495,7 @@ source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "39cab71617ae0d63f51a36d69f866391735b51691dbda63cf6f96d042b63efeb" dependencies = [ "libc", - "windows-sys 0.52.0", + "windows-sys 0.61.2", ] [[package]] @@ -505,9 +522,9 @@ checksum = "5baebc0774151f905a1a2cc41989300b1e6fbb29aff0ceffa1064fdd3088d582" [[package]] name = "flate2" -version = "1.0.28" +version = "1.1.9" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "46303f565772937ffe1d394a4fac6f411c6013172fadde9dcdb1e147a086940e" +checksum = "843fba2746e448b37e26a819579957415c8cef339bf08564fe8b7ddbd959573c" dependencies = [ "crc32fast", "miniz_oxide", @@ -555,6 +572,17 @@ dependencies = [ "wasi", ] +[[package]] +name = "getrandom" +version = "0.4.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "300e883d756b2e4ec94e02791f39b04b522276138852cfc41d9fb7e904106099" +dependencies = [ + "cfg-if", + "libc", + "r-efi", +] + [[package]] name = "git-version" version = "0.3.9" @@ -660,6 +688,17 @@ dependencies = [ "hashbrown", ] +[[package]] +name = "indicatif" +version = "0.18.4" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "25470f23803092da7d239834776d653104d551bc4d7eacaf31e6837854b8e9eb" +dependencies = [ + "console", + "portable-atomic", + "unit-prefix", +] + [[package]] name = "is-terminal" version = "0.4.17" @@ -668,14 +707,14 @@ checksum = "3640c1c38b8e4e43584d8df18be5fc6b0aa314ce6ebf51b53313d4306cca8e46" dependencies = [ "hermit-abi", "libc", - "windows-sys 0.52.0", + "windows-sys 0.61.2", ] [[package]] name = "itertools" -version = "0.14.0" +version = "0.13.0" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "2b192c782037fadd9cfa75548310488aabdbf3d2da73885b31bd0abd03351285" +checksum = "413ee7dfc52ee1a4949ceeb7dbc8a33f2d6c088194d9f922fb8318faf1f01186" dependencies = [ "either", ] @@ -775,11 +814,12 @@ checksum = "f8ca58f447f06ed17d5fc4043ce1b10dd205e060fb3ce5b979b8ed8e59ff3f79" [[package]] name = "miniz_oxide" -version = "0.7.1" +version = "0.8.9" source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "e7810e0be55b428ada41041c41f32c9f1a42817901b4ccf45fa3d4b6561e74c7" +checksum = "1fa76a2c86f704bdb222d66965fb3d63269ce38518b83cb0575fca855ebb6316" dependencies = [ - "adler", + "adler2", + "simd-adler32", ] [[package]] @@ -811,6 +851,68 @@ dependencies = [ "zip", ] +[[package]] +name = "neopdf" +version = "0.4.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "1f4081bcb68db0532fcbb63cff45a7abb90aa30b3cd3a8cf3c8085918e47d729" +dependencies = [ + "bincode", + "flate2", + "git-version", + "indicatif", + "itertools", + "lz4_flex", + "ndarray", + "neopdf_legacy", + "ninterp", + "rayon", + "regex", + "serde", + "serde_yaml", + "tar", + "tempfile", + "thiserror", + "ureq", +] + +[[package]] +name = "neopdf_legacy" +version = "0.4.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "6d176650a696abb440cc1bd576934433f96a989ac4b4a414673f3568ff544626" +dependencies = [ + "bincode", + "flate2", + "git-version", + "indicatif", + "itertools", + "lz4_flex", + "ndarray", + "ninterp", + "rayon", + "regex", + "serde", + "serde_yaml", + "tar", + "tempfile", + "thiserror", + "ureq", +] + +[[package]] +name = "ninterp" +version = "0.8.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "67588d36f7a821f32edf2ad9699804b65f65b835634ccc267f1c2bcdaca1acf4" +dependencies = [ + "dyn-clone", + "itertools", + "ndarray", + "num-traits", + "thiserror", +] + [[package]] name = "num-bigint" version = "0.4.6" @@ -957,7 +1059,6 @@ dependencies = [ "cc", "cxx", "cxx-build", - "managed-lhapdf", "pkg-config", ] @@ -990,6 +1091,7 @@ dependencies = [ "managed-lhapdf", "ndarray", "ndarray-npy", + "neopdf", "pineappl", "pineappl_applgrid", "pineappl_fastnlo", @@ -1091,7 +1193,7 @@ dependencies = [ "is-terminal", "lazy_static", "term", - "unicode-width", + "unicode-width 0.1.11", ] [[package]] @@ -1183,6 +1285,12 @@ dependencies = [ "proc-macro2", ] +[[package]] +name = "r-efi" +version = "6.0.0" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f8dcc9c7d52a811697d2151c701e0d08956f92b0e24136cf4cf27b57a6a0d9bf" + [[package]] name = "rand" version = "0.10.1" @@ -1239,11 +1347,23 @@ version = "0.4.6" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "ba009ff324d1fc1b900bd1fdb31564febe58a8ccc8a6fdbb93b543d33b13ca43" dependencies = [ - "getrandom", + "getrandom 0.2.17", "libredox", "thiserror", ] +[[package]] +name = "regex" +version = "1.12.3" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "e10754a14b9137dd7b1e3e5b0493cc9171fdd105e0ab477f51b72e7f3ac0e276" +dependencies = [ + "aho-corasick", + "memchr", + "regex-automata", + "regex-syntax", +] + [[package]] name = "regex-automata" version = "0.4.14" @@ -1269,7 +1389,7 @@ checksum = "a4689e6c2294d81e88dc6261c768b63bc4fcdb852be6d1352498b114f61383b7" dependencies = [ "cc", "cfg-if", - "getrandom", + "getrandom 0.2.17", "libc", "untrusted", "windows-sys 0.52.0", @@ -1297,7 +1417,7 @@ dependencies = [ "errno", "libc", "linux-raw-sys", - "windows-sys 0.52.0", + "windows-sys 0.61.2", ] [[package]] @@ -1431,6 +1551,12 @@ version = "1.3.0" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "0fda2ff0d084019ba4d7c6f371c95d8fd75ce3524c3cb8fb653a3023f6323e64" +[[package]] +name = "simd-adler32" +version = "0.3.9" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "703d5c7ef118737c72f1af64ad2f6f8c5e1921f818cdcb97b8fe6fc69bf66214" + [[package]] name = "smallvec" version = "1.15.1" @@ -1483,9 +1609,10 @@ source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "32497e9a4c7b38532efcdebeef879707aa9f794296a4f0244f6f69e9bc8574bd" dependencies = [ "fastrand", + "getrandom 0.4.3", "once_cell", "rustix", - "windows-sys 0.52.0", + "windows-sys 0.61.2", ] [[package]] @@ -1598,6 +1725,18 @@ version = "0.1.11" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "e51733f11c9c4f72aa0c160008246859e340b00807569a0da0e7a1079b27ba85" +[[package]] +name = "unicode-width" +version = "0.2.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "b4ac048d71ede7ee76d585517add45da530660ef4390e49b098733c6e897f254" + +[[package]] +name = "unit-prefix" +version = "0.5.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "81e544489bf3d8ef66c953931f56617f423cd4b5494be343d9b9d3dda037b9a3" + [[package]] name = "unsafe-libyaml" version = "0.2.11" @@ -1713,7 +1852,7 @@ version = "0.1.11" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "c2a7b1c03c876122aa43f3020e6c3c3ee5c05081c9a00739faf7503aeba10d22" dependencies = [ - "windows-sys 0.52.0", + "windows-sys 0.61.2", ] [[package]] @@ -1722,6 +1861,12 @@ version = "0.4.0" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "712e227841d057c1ee1cd2fb22fa7e5a5461ae8e48fa2ca79ec42cfc1931183f" +[[package]] +name = "windows-link" +version = "0.2.1" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "f0805222e57f7521d6a62e36fa9163bc891acd422f971defe97d64e70d0a4fe5" + [[package]] name = "windows-sys" version = "0.48.0" @@ -1740,6 +1885,15 @@ dependencies = [ "windows-targets 0.52.6", ] +[[package]] +name = "windows-sys" +version = "0.61.2" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "ae137229bcbd6cdf0f7b80a31df61766145077ddf49416a728b02cb3921ff3fc" +dependencies = [ + "windows-link", +] + [[package]] name = "windows-targets" version = "0.48.5" diff --git a/Cargo.toml b/Cargo.toml index 1098c4dc..9bb18305 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -35,11 +35,12 @@ enum_dispatch = "0.3.7" flate2 = "1.0.22" float-cmp = { default-features = false, version = "0.10.0" } git-version = "0.3.5" -itertools = "0.14.0" +itertools = "0.13.0" lhapdf = { package = "managed-lhapdf", version = "0.4.1" } lz4_flex = "0.13.0" ndarray = { version = "0.17.2" } ndarray-npy = { default-features = false, features = ["npz"], version = "0.10.0" } +neopdf = { version = "0.4.0" } num-complex = "0.4.4" numpy = "0.28.0" pkg-config = "0.3.26" diff --git a/README.md b/README.md index 886318c2..b94cf5d3 100644 --- a/README.md +++ b/README.md @@ -46,7 +46,8 @@ If you use PineAPPL, please cite By using PineAPPL, you're probably also using - [APPLgrid], -- [fastNLO] and/or +- [fastNLO], +- [NeoPDF] and/or - [LHAPDF]. If that is the case, please cite these accordingly. @@ -54,6 +55,7 @@ If that is the case, please cite these accordingly. [APPLgrid]: https://applgrid.hepforge.org [fastNLO]: https://fastnlo.hepforge.org [LHAPDF]: https://lhapdf.hepforge.org +[NeoPDF]: https://qcdlab.github.io/neopdf/ [zenodo DOI]: https://zenodo.org/badge/latestdoi/248306479 [version 0]: https://inspirehep.net/literature/1814432 [version 1]: https://inspirehep.net/literature/3169212 diff --git a/pineappl_applgrid/Cargo.toml b/pineappl_applgrid/Cargo.toml index c7a0ebf2..cc26b521 100644 --- a/pineappl_applgrid/Cargo.toml +++ b/pineappl_applgrid/Cargo.toml @@ -17,7 +17,6 @@ workspace = true [dependencies] cxx.workspace = true -lhapdf.workspace = true [build-dependencies] cc.workspace = true @@ -25,4 +24,4 @@ cxx-build.workspace = true pkg-config.workspace = true [features] -static = ["lhapdf/static"] +static = [] diff --git a/pineappl_applgrid/src/lib.rs b/pineappl_applgrid/src/lib.rs index d238bb94..bf805053 100644 --- a/pineappl_applgrid/src/lib.rs +++ b/pineappl_applgrid/src/lib.rs @@ -194,7 +194,6 @@ pub mod ffi { } } -use lhapdf::Pdf; use std::mem; use std::pin::Pin; use std::ptr; @@ -206,14 +205,22 @@ static MUTEX: OnceLock> = OnceLock::new(); /// TODO. pub fn grid_convolve_with_one( grid: Pin<&mut ffi::grid>, - pdf: &mut Pdf, + xfx_q2: &mut impl FnMut(i32, f64, f64) -> f64, + alphas_q2: &mut impl FnMut(f64) -> f64, nloops: i32, rscale: f64, fscale: f64, escale: f64, ) -> Vec { - let xfx = |x: &f64, q: &f64, results: *mut f64, pdf: *mut ffi::c_void| { - let pdf = unsafe { &mut *pdf.cast::() }; + struct FnMuts<'a> { + xfx_q2: &'a mut dyn FnMut(i32, f64, f64) -> f64, + alphas_q2: &'a mut dyn FnMut(f64) -> f64, + } + + let mut fn_muts = FnMuts { xfx_q2, alphas_q2 }; + + let xfx = |x: &f64, q: &f64, results: *mut f64, ptr: *mut ffi::c_void| { + let fn_muts = unsafe { &mut *ptr.cast::>() }; let results = unsafe { slice::from_raw_parts_mut(results, 14) }; for (pid, result) in [-6, -5, -4, -3, -2, -1, 21, 1, 2, 3, 4, 5, 6, 22] .into_iter() @@ -221,13 +228,13 @@ pub fn grid_convolve_with_one( { // some grids have x-nodes at slightly larger values than `1.0`; in that cases these // are numerical problems which we 'fix' by evaluating at exactly `1.0` instead - *result = pdf.xfx_q2(pid, x.min(1.0), *q * *q); + *result = (fn_muts.xfx_q2)(pid, x.min(1.0), *q * *q); } }; - let alphas = |q: &f64, pdf: *mut ffi::c_void| -> f64 { - let pdf = unsafe { &mut *pdf.cast::() }; - pdf.alphas_q2(*q * *q) + let alphas = |q: &f64, ptr: *mut ffi::c_void| -> f64 { + let fn_muts = unsafe { &mut *ptr.cast::>() }; + (fn_muts.alphas_q2)(*q * *q) }; let lock = MUTEX @@ -241,7 +248,7 @@ pub fn grid_convolve_with_one( grid, xfx, alphas, - ptr::from_mut(pdf).cast::(), + ptr::from_mut(&mut fn_muts).cast::(), nloops, rscale, fscale, diff --git a/pineappl_cli/Cargo.toml b/pineappl_cli/Cargo.toml index 65e8041e..854cca95 100644 --- a/pineappl_cli/Cargo.toml +++ b/pineappl_cli/Cargo.toml @@ -30,6 +30,7 @@ lhapdf.workspace = true lz4_flex = { optional = true, workspace = true } ndarray.workspace = true ndarray-npy = { optional = true, workspace = true } +neopdf.workspace = true pineappl.workspace = true pineappl_applgrid = { optional = true, workspace = true } pineappl_fastnlo = { optional = true, workspace = true } diff --git a/pineappl_cli/src/analyze.rs b/pineappl_cli/src/analyze.rs index cc5f8738..4fa01737 100644 --- a/pineappl_cli/src/analyze.rs +++ b/pineappl_cli/src/analyze.rs @@ -64,7 +64,7 @@ pub struct CkfOpts { impl Subcommand for CkfOpts { fn run(&self, cfg: &GlobalConfiguration) -> Result { let grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; let orders_den = if self.orders_den.is_empty() { grid.orders() diff --git a/pineappl_cli/src/channels.rs b/pineappl_cli/src/channels.rs index 0ba04375..a92722fa 100644 --- a/pineappl_cli/src/channels.rs +++ b/pineappl_cli/src/channels.rs @@ -64,7 +64,7 @@ pub struct Opts { impl Subcommand for Opts { fn run(&self, cfg: &GlobalConfiguration) -> Result { let grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; let mut channels: Vec<_> = self.channels.iter().cloned().flatten().collect(); channels.sort_unstable(); diff --git a/pineappl_cli/src/convolve.rs b/pineappl_cli/src/convolve.rs index f6ceae11..d33d166f 100644 --- a/pineappl_cli/src/convolve.rs +++ b/pineappl_cli/src/convolve.rs @@ -58,7 +58,7 @@ pub struct Opts { impl Subcommand for Opts { fn run(&self, cfg: &GlobalConfiguration) -> Result { let grid = helpers::read_grid(&self.input)?; - let mut conv_funs_0 = helpers::create_conv_funs(&self.conv_funs[0])?; + let mut conv_funs_0 = helpers::create_conv_funs(&self.conv_funs[0], cfg.pdf_backend)?; let bins: Vec<_> = self.bins.iter().cloned().flatten().collect(); let results = helpers::convolve_scales( @@ -91,7 +91,7 @@ impl Subcommand for Opts { .iter() .flat_map(|conv_funs| { let conv_types = &conv_funs.conv_types; - let mut conv_funs = helpers::create_conv_funs(conv_funs).unwrap(); + let mut conv_funs = helpers::create_conv_funs(conv_funs, cfg.pdf_backend).unwrap(); helpers::convolve( &grid, &mut conv_funs, diff --git a/pineappl_cli/src/diff.rs b/pineappl_cli/src/diff.rs index 151549be..0b8e9172 100644 --- a/pineappl_cli/src/diff.rs +++ b/pineappl_cli/src/diff.rs @@ -120,7 +120,7 @@ impl Subcommand for Opts { bail!("channels differ"); } - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; let mut table = helpers::create_table(); let mut title = Row::empty(); diff --git a/pineappl_cli/src/evolve/mod.rs b/pineappl_cli/src/evolve/mod.rs index dc4413c9..437c31bc 100644 --- a/pineappl_cli/src/evolve/mod.rs +++ b/pineappl_cli/src/evolve/mod.rs @@ -2,10 +2,10 @@ mod eko; use super::helpers::{self, ConvFuns, ConvoluteMode}; +use super::pdf_backend::PdfBackend; use super::{GlobalConfiguration, Subcommand}; use anyhow::{Result, anyhow}; use clap::{Parser, ValueHint}; -use lhapdf::Pdf; use pineappl::fk_table::FkTable; use pineappl::grid::Grid; use std::path::{Path, PathBuf}; @@ -59,7 +59,7 @@ impl Subcommand for Opts { use prettytable::row; let grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; let results = helpers::convolve_scales( &grid, &mut conv_funs, @@ -133,7 +133,7 @@ impl Subcommand for Opts { fn evolve_grid( grid: &Grid, ekos: &[&Path], - use_alphas_from: &Pdf, + use_alphas_from: &dyn PdfBackend, orders: &[(u8, u8)], xir: f64, xif: f64, @@ -167,7 +167,7 @@ fn evolve_grid( fn evolve_grid( _: &Grid, _: &[&Path], - _: &Pdf, + _: &dyn PdfBackend, _: &[(u8, u8)], _: f64, _: f64, diff --git a/pineappl_cli/src/export/applgrid.rs b/pineappl_cli/src/export/applgrid.rs index 910cf2cf..30efceba 100644 --- a/pineappl_cli/src/export/applgrid.rs +++ b/pineappl_cli/src/export/applgrid.rs @@ -2,7 +2,6 @@ use anyhow::{Result, bail}; use cxx::{UniquePtr, let_cxx_string}; use float_cmp::approx_eq; use itertools::izip; -use lhapdf::Pdf; use ndarray::{Axis, s}; use pineappl::boc::{Channel, Kinematics, Order}; use pineappl::grid::Grid; @@ -13,7 +12,6 @@ use pineappl_applgrid::ffi::{self, grid}; use std::f64::consts::TAU; use std::iter; use std::path::Path; -use std::pin::Pin; fn reconstruct_subgrid_params(grid: &Grid, order: usize, bin: usize) -> Result> { if grid @@ -392,13 +390,3 @@ pub(super) fn convert_into_applgrid( Ok((applgrid, order_mask)) } - -// TODO: deduplicate this function from import -pub(super) fn convolve_applgrid(grid: Pin<&mut grid>, conv_funs: &mut [Pdf]) -> Vec { - let nloops = grid.nloops(); - - // TODO: add support for convolving an APPLgrid with two functions - assert_eq!(conv_funs.len(), 1); - - pineappl_applgrid::grid_convolve_with_one(grid, &mut conv_funs[0], nloops, 1.0, 1.0, 1.0) -} diff --git a/pineappl_cli/src/export/mod.rs b/pineappl_cli/src/export/mod.rs index 32232eaf..9917566b 100644 --- a/pineappl_cli/src/export/mod.rs +++ b/pineappl_cli/src/export/mod.rs @@ -2,11 +2,11 @@ mod applgrid; use super::helpers::{self, ConvFuns, ConvoluteMode}; +use super::pdf_backend::ConvFunBackend; use super::{GlobalConfiguration, Subcommand}; use anyhow::{Result, anyhow}; use clap::builder::{PossibleValuesParser, TypedValueParser as _}; use clap::{Parser, ValueHint}; -use lhapdf::Pdf; use pineappl::boc::Order; use pineappl::grid::Grid; use std::path::{Path, PathBuf}; @@ -50,13 +50,13 @@ impl Subcommand for Opts { use prettytable::{cell, row}; let mut grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; // TODO: figure out `member` from `self.pdfset` let (grid_type, results, scale_variations, order_mask) = convert_into_grid( &self.output, &mut grid, - &mut conv_funs, + &conv_funs, self.scales, self.discard_non_matching_values, )?; @@ -190,15 +190,30 @@ impl Subcommand for Opts { fn convert_into_applgrid( output: &Path, grid: &mut Grid, - conv_funs: &mut [Pdf], + conv_funs: &[ConvFunBackend], _: usize, discard_non_matching_values: bool, ) -> Result<(&'static str, Vec, usize, Vec)> { + use super::pdf_backend::PdfBackend as _; + // TODO: check also scale-varied results let (mut applgrid, order_mask) = applgrid::convert_into_applgrid(grid, output, discard_non_matching_values)?; - let results = applgrid::convolve_applgrid(applgrid.pin_mut(), conv_funs); + let nloops = applgrid.nloops(); + + // TODO: add support for convolving an APPLgrid with two functions + assert_eq!(conv_funs.len(), 1); + + let results = pineappl_applgrid::grid_convolve_with_one( + applgrid.pin_mut(), + &mut |pid, x, q2| conv_funs[0].xfx_q2(pid, x, q2), + &mut |q2| conv_funs[0].alphas_q2(q2), + nloops, + 1.0, + 1.0, + 1.0, + ); Ok(("APPLgrid", results, 1, order_mask)) } @@ -207,7 +222,7 @@ fn convert_into_applgrid( fn convert_into_applgrid( _: &Path, _: &mut Grid, - _: &mut [Pdf], + _: &[ConvFunBackend], _: usize, _: bool, ) -> Result<(&'static str, Vec, usize, Vec)> { @@ -219,7 +234,7 @@ fn convert_into_applgrid( fn convert_into_grid( output: &Path, grid: &mut Grid, - conv_funs: &mut [Pdf], + conv_funs: &[ConvFunBackend], scales: usize, discard_non_matching_values: bool, ) -> Result<(&'static str, Vec, usize, Vec)> { diff --git a/pineappl_cli/src/helpers.rs b/pineappl_cli/src/helpers.rs index 227fcc7c..8e2ef8b9 100644 --- a/pineappl_cli/src/helpers.rs +++ b/pineappl_cli/src/helpers.rs @@ -1,7 +1,10 @@ use super::GlobalConfiguration; +use super::pdf_backend::{ + self, Backend, ConvFunBackend, ConvFunSetBackend, ForcePositive, PdfBackend as _, + PdfSetBackend as _, +}; use anyhow::{Context as _, Error, Result, anyhow, bail}; use itertools::Itertools as _; -use lhapdf::{Pdf, PdfSet}; use pineappl::boc::{ScaleFuncForm, Scales}; use pineappl::convolutions::{Conv, ConvType, ConvolutionCache}; use pineappl::grid::Grid; @@ -136,50 +139,33 @@ pub enum ConvoluteMode { Normal, } -pub fn create_conv_funs(funs: &ConvFuns) -> Result> { - Ok(funs - .lhapdf_names +/// Creates convolution functions using the specified backend. +pub fn create_conv_funs(funs: &ConvFuns, backend: Backend) -> Result> { + funs.lhapdf_names .iter() .zip(&funs.members) - .map(|(lhapdf_name, member)| { - lhapdf_name.parse().map_or_else( - |_| { - let member = member.unwrap_or(0); - // UNWRAP: we don't support sets with more members than `i32` - Pdf::with_setname_and_member(lhapdf_name, member.try_into().unwrap()) - }, - Pdf::with_lhaid, - ) + .map(|(name, member)| { + let member = member.unwrap_or(0); + pdf_backend::create_pdf(name, member, backend) }) - .collect::>()?) + .collect() } +/// Creates convolution functions for a PDF set using the specified backend. pub fn create_conv_funs_for_set( funs: &ConvFuns, index_of_set: usize, -) -> Result<(PdfSet, Vec>)> { + backend: Backend, +) -> Result<(ConvFunSetBackend, Vec>)> { let setname = &funs.lhapdf_names[index_of_set]; - let set = setname.parse().map_or_else( - |_| Ok::<_, Error>(PdfSet::new(setname)?), - |lhaid| { - Ok(PdfSet::new( - &lhapdf::lookup_pdf(lhaid) - .map(|(set, _)| set) - .ok_or_else(|| { - anyhow!("no convolution function for LHAID = `{lhaid}` found") - })?, - )?) - }, - )?; + let set = pdf_backend::create_pdf_set(setname, backend)?; - let conv_funs = set - .mk_pdfs()? + let set_members = set.mk_pdfs()?; + let conv_funs = set_members .into_iter() - .map(|conv_fun| { - // TODO: do not create objects that are getting overwritten in any case - let mut conv_funs = create_conv_funs(funs)?; - conv_funs[index_of_set] = conv_fun; - + .map(|member_pdf| { + let mut conv_funs = create_conv_funs(funs, backend)?; + conv_funs[index_of_set] = member_pdf; Ok::<_, Error>(conv_funs) }) .collect::>()?; @@ -249,9 +235,10 @@ pub fn labels_and_units(grid: &Grid, integrated: bool) -> (Vec<(String, &str)>, ) } +/// Performs convolution with scale variations using the backend abstraction. pub fn convolve_scales( grid: &Grid, - conv_funs: &mut [Pdf], + conv_funs: &mut [ConvFunBackend], conv_types: &[ConvType], orders: &[(u8, u8)], bins: &[usize], @@ -273,7 +260,7 @@ pub fn convolve_scales( if cfg.force_positive { for fun in conv_funs.iter_mut() { - fun.set_force_positive(1); + fun.set_force_positive(ForcePositive::ClipNegative); } } @@ -289,11 +276,12 @@ pub fn convolve_scales( .iter_mut() .map(|fun| (fun.x_min(), fun.x_max())) .collect(); + let mut funs: Vec<_> = conv_funs .iter() - .zip(x_min_max) - .map(|(fun, (x_min, x_max))| { - move |id, x, q2| { + .zip(&x_min_max) + .map(|(fun, &(x_min, x_max))| { + move |id: i32, x: f64, q2: f64| { if !cfg.allow_extrapolation && (x < x_min || x > x_max) { 0.0 } else { @@ -302,26 +290,19 @@ pub fn convolve_scales( } }) .collect(); + let xfx: Vec<_> = funs .iter_mut() .map(|fun| fun as &mut dyn FnMut(i32, f64, f64) -> f64) .collect(); - let mut alphas_funs: Vec<_> = conv_funs - .iter() - .map(|fun| move |q2| fun.alphas_q2(q2)) - .collect(); + + let mut alphas_funs: Vec<_> = conv_funs.iter().map(|fun| |q2| fun.alphas_q2(q2)).collect(); + let convolutions: Vec<_> = conv_funs .iter() .zip(conv_types) .map(|(fun, &conv_type)| { - let pid = fun - .set() - .entry("Particle") - // if the field 'Particle' is missing we assume it's a proton PDF - .map_or(Ok(2212), |string| string.parse::()) - // UNWRAP: if this fails, there's a non-integer string in the LHAPDF info file - .unwrap(); - + let pid = fun.particle_id(); Conv::new(conv_type, pid) }) .collect(); @@ -384,9 +365,10 @@ pub fn scales_vector(grid: &Grid, scales: usize) -> &[(f64, f64, f64)] { } } +/// Performs convolution using the backend abstraction. pub fn convolve( grid: &Grid, - conv_funs: &mut [Pdf], + conv_funs: &mut [ConvFunBackend], conv_types: &[ConvType], orders: &[(u8, u8)], bins: &[usize], diff --git a/pineappl_cli/src/import/applgrid.rs b/pineappl_cli/src/import/applgrid.rs index 5652bd37..efa4bce6 100644 --- a/pineappl_cli/src/import/applgrid.rs +++ b/pineappl_cli/src/import/applgrid.rs @@ -1,6 +1,5 @@ use anyhow::Result; use float_cmp::assert_approx_eq; -use lhapdf::Pdf; use pineappl::boc::{BinsWithFillLimits, Channel, Kinematics, Order, ScaleFuncForm, Scales}; use pineappl::convolutions::{Conv, ConvType}; use pineappl::grid::Grid; @@ -299,12 +298,3 @@ pub(super) fn convert_applgrid(grid: Pin<&mut grid>, alpha: u8) -> Result Ok(grid0) } - -pub(super) fn convolve_applgrid(grid: Pin<&mut grid>, conv_funs: &mut [Pdf]) -> Vec { - let nloops = grid.nloops(); - - // TODO: add support for convolving an APPLgrid with two functions - assert_eq!(conv_funs.len(), 1); - - pineappl_applgrid::grid_convolve_with_one(grid, &mut conv_funs[0], nloops, 1.0, 1.0, 1.0) -} diff --git a/pineappl_cli/src/import/mod.rs b/pineappl_cli/src/import/mod.rs index 01df27ad..f124f7d7 100644 --- a/pineappl_cli/src/import/mod.rs +++ b/pineappl_cli/src/import/mod.rs @@ -6,11 +6,11 @@ mod fastnlo; mod fktable; use super::helpers::{self, ConvFuns, ConvoluteMode}; +use super::pdf_backend::ConvFunBackend; use super::{GlobalConfiguration, Subcommand}; use anyhow::{Result, anyhow}; use clap::builder::{PossibleValuesParser, TypedValueParser as _}; use clap::{Parser, ValueHint}; -use lhapdf::Pdf; use pineappl::grid::Grid; use std::path::{Path, PathBuf}; use std::process::ExitCode; @@ -63,13 +63,13 @@ impl Subcommand for Opts { fn run(&self, cfg: &GlobalConfiguration) -> Result { use prettytable::{cell, row}; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; // TODO: figure out `member` from `self.pdfset` let (grid_type, mut grid, reference_results, scale_variations) = convert_grid( &self.input, self.alpha, - &mut conv_funs, + &conv_funs, &self.conv_funs, 0, self.scales, @@ -172,16 +172,30 @@ impl Subcommand for Opts { fn convert_applgrid( input: &Path, alpha: u8, - conv_funs: &mut [Pdf], + conv_funs: &[ConvFunBackend], _: usize, ) -> Result<(&'static str, Grid, Vec, usize)> { + use super::pdf_backend::PdfBackend as _; use pineappl_applgrid::ffi; // TODO: check AMCATNLO scale variations let mut grid = ffi::make_grid(input.to_str().unwrap())?; let pgrid = applgrid::convert_applgrid(grid.pin_mut(), alpha)?; - let results = applgrid::convolve_applgrid(grid.pin_mut(), conv_funs); + let nloops = grid.nloops(); + + // TODO: add support for convolving an APPLgrid with two functions + assert_eq!(conv_funs.len(), 1); + + let results = pineappl_applgrid::grid_convolve_with_one( + grid.pin_mut(), + &mut |pid, x, q2| conv_funs[0].xfx_q2(pid, x, q2), + &mut |q2| conv_funs[0].alphas_q2(q2), + nloops, + 1.0, + 1.0, + 1.0, + ); Ok(("APPLgrid", pgrid, results, 1)) } @@ -190,7 +204,7 @@ fn convert_applgrid( fn convert_applgrid( _: &Path, _: u8, - _: &mut [Pdf], + _: &[ConvFunBackend], _: usize, ) -> Result<(&'static str, Grid, Vec, usize)> { Err(anyhow!( @@ -299,7 +313,7 @@ fn convert_fktable(_: &Path) -> Result<(&'static str, Grid, Vec, usize)> { fn convert_grid( input: &Path, alpha: u8, - conv_funs: &mut [Pdf], + conv_funs: &[ConvFunBackend], fun_names: &ConvFuns, member: usize, scales: usize, diff --git a/pineappl_cli/src/lib.rs b/pineappl_cli/src/lib.rs index 50abbb09..98c50987 100644 --- a/pineappl_cli/src/lib.rs +++ b/pineappl_cli/src/lib.rs @@ -11,6 +11,7 @@ mod helpers; mod import; mod merge; mod orders; +pub mod pdf_backend; mod plot; mod pull; mod read; @@ -39,6 +40,9 @@ pub struct GlobalConfiguration { /// Choose the PDF/FF set for the strong coupling. #[arg(default_value = "0", long, value_name = "IDX")] pub use_alphas_from: usize, + /// Select the PDF interpolation backend: 'lhapdf' or 'neopdf'. + #[arg(default_value = "lhapdf", long, value_name = "BACKEND")] + pub pdf_backend: pdf_backend::Backend, } /// TODO. diff --git a/pineappl_cli/src/orders.rs b/pineappl_cli/src/orders.rs index 7000d738..b16ab589 100644 --- a/pineappl_cli/src/orders.rs +++ b/pineappl_cli/src/orders.rs @@ -41,7 +41,7 @@ pub struct Opts { impl Subcommand for Opts { fn run(&self, cfg: &GlobalConfiguration) -> Result { let grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; let mut orders: Vec<_> = grid .orders() diff --git a/pineappl_cli/src/pdf_backend.rs b/pineappl_cli/src/pdf_backend.rs new file mode 100644 index 00000000..b0d5f158 --- /dev/null +++ b/pineappl_cli/src/pdf_backend.rs @@ -0,0 +1,399 @@ +//! PDF backend abstraction layer. +//! +//! This module provides a unified interface for different PDF interpolation backends, +//! currently supporting `LHAPDF` and `NeoPDF`. It allows runtime selection of the backend +//! and provides type-safe access to PDF metadata. + +use anyhow::{Context as _, Result, anyhow}; +use enum_dispatch::enum_dispatch; +// use std::fmt; +use lhapdf::PdfSet; +use neopdf::gridpdf::ForcePositive as NeoPDFForcePositive; +use neopdf::pdf::PDF; +use std::str::FromStr; + +pub use neopdf::uncertainty::CL_1_SIGMA; + +/// Method for handling negative PDF values. +#[derive(Clone, Copy, Debug)] +pub enum ForcePositive { + /// No clipping - return values as-is. + None, + /// Clip negative values to zero. + ClipNegative, +} + +/// Available PDF backends. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq)] +pub enum Backend { + /// `LHAPDF` backend (C++ library with Rust bindings). + #[default] + Lhapdf, + /// `NeoPDF` backend (pure Rust implementation). + Neopdf, +} + +impl FromStr for Backend { + type Err = anyhow::Error; + + fn from_str(s: &str) -> Result { + match s.to_lowercase().as_str() { + "lhapdf" => Ok(Self::Lhapdf), + "neopdf" => Ok(Self::Neopdf), + _ => Err(anyhow!( + "unknown PDF backend '{s}'; must be 'lhapdf' or 'neopdf'" + )), + } + } +} + +/// Unified PDF backend interface. +/// +/// This trait provides a common interface for PDF interpolation backends, +/// allowing the CLI to work with different implementations transparently. +#[enum_dispatch] +pub trait PdfBackend: Send { + /// Evaluates the strong coupling constant `alpha_s(Q^2)`. + fn alphas_q2(&self, q2: f64) -> f64; + + // /// Returns the error type of the PDF set (e.g., "replicas", "hessian"). + // fn error_type(&self) -> String; + + // /// Returns whether this is a polarized PDF. + // fn is_polarized(&self) -> bool; + + /// Returns the hadron particle ID (PDG code). + fn particle_id(&self) -> i32; + + /// Sets the method for handling negative PDF values. + fn set_force_positive(&mut self, method: ForcePositive); + + // /// Returns the PDF set type (space-like or time-like). + // fn set_type(&self) -> SetType; + + /// Returns the maximum valid x value. + fn x_max(&mut self) -> f64; + + /// Returns the minimum valid x value. + fn x_min(&mut self) -> f64; + + /// Evaluates xf(x, Q^2) for a given flavor. + /// + /// # Arguments + /// * `id` - The Monte Carlo PDG flavor ID. + /// * `x` - The momentum fraction. + /// * `q2` - The squared energy scale. + /// + /// # Returns + /// The PDF value xf(x, Q^2). + /// + /// # TODO + /// Extend to support `NeoPDF` multi-parameters interpolation. + fn xfx_q2(&self, id: i32, x: f64, q2: f64) -> f64; +} + +/// TODO. +#[expect( + clippy::large_enum_variant, + reason = "LhapdfPdf is artificially small due to indirection" +)] +#[enum_dispatch(PdfBackend)] +pub enum ConvFunBackend { + /// TODO. + LhapdfPdf, + /// TODO. + NeopdfPdf, +} + +/// Unified PDF set interface for uncertainty calculations. +#[enum_dispatch] +pub trait PdfSetBackend { + /// Creates all PDF members in the set. + /// + /// # Errors + /// + /// TODO. + fn mk_pdfs(&self) -> Result>; + + /// Calculates the uncertainty for a set of values. + /// + /// # Errors + /// + /// TODO. + fn uncertainty(&self, values: &[f64], cl: f64, alternative: bool) -> Result; +} + +/// TODO. +#[enum_dispatch(PdfSetBackend)] +pub enum ConvFunSetBackend { + /// TODO. + LhapdfSet, + /// TODO. + NeopdfSet, +} + +/// Uncertainty information from a PDF set. +#[derive(Clone, Debug)] +pub struct Uncertainty { + /// Central value. + pub central: f64, + /// Negative error (absolute value). + pub errminus: f64, + /// Positive error (absolute value). + pub errplus: f64, +} + +/// LHAPDF backend wrapper. +pub struct LhapdfPdf { + pdf: lhapdf::Pdf, +} + +impl LhapdfPdf { + /// Creates a new LHAPDF PDF by LHAID. + /// + /// # Errors + /// + /// TODO. + pub fn with_lhaid(lhaid: i32) -> Result { + Ok(Self { + pdf: lhapdf::Pdf::with_lhaid(lhaid) + .with_context(|| format!("failed to load LHAPDF with LHAID {lhaid}"))?, + }) + } + + /// Creates a new LHAPDF PDF by set name and member index. + /// + /// # Errors + /// + /// TODO. + pub fn with_setname_and_member(setname: &str, member: i32) -> Result { + Ok(Self { + pdf: lhapdf::Pdf::with_setname_and_member(setname, member).with_context(|| { + format!("failed to load LHAPDF set '{setname}' member {member}") + })?, + }) + } +} + +impl PdfBackend for LhapdfPdf { + fn alphas_q2(&self, q2: f64) -> f64 { + self.pdf.alphas_q2(q2) + } + + fn particle_id(&self) -> i32 { + self.pdf + .set() + .entry("Particle") + .map_or(Ok(2212), |s| s.parse()) + .unwrap_or(2212) + } + + fn set_force_positive(&mut self, method: ForcePositive) { + match method { + ForcePositive::None => self.pdf.set_force_positive(0), + ForcePositive::ClipNegative => self.pdf.set_force_positive(1), + } + } + + fn x_max(&mut self) -> f64 { + self.pdf.x_max() + } + + fn x_min(&mut self) -> f64 { + self.pdf.x_min() + } + + fn xfx_q2(&self, id: i32, x: f64, q2: f64) -> f64 { + self.pdf.xfx_q2(id, x, q2) + } +} + +/// LHAPDF set wrapper. +pub struct LhapdfSet { + set: lhapdf::PdfSet, +} + +impl LhapdfSet { + /// Creates a new LHAPDF set by name. + /// + /// # Errors + /// + /// TODO. + pub fn new(setname: &str) -> Result { + Ok(Self { + set: PdfSet::new(setname) + .with_context(|| format!("failed to load LHAPDF set '{setname}'"))?, + }) + } +} + +impl PdfSetBackend for LhapdfSet { + fn mk_pdfs(&self) -> Result> { + let pdfs = self.set.mk_pdfs()?; + Ok(pdfs + .into_iter() + .map(|pdf| LhapdfPdf { pdf }.into()) + .collect()) + } + + fn uncertainty(&self, values: &[f64], cl: f64, alternative: bool) -> Result { + let unc = self.set.uncertainty(values, cl, alternative)?; + Ok(Uncertainty { + central: unc.central, + errminus: unc.errminus, + errplus: unc.errplus, + }) + } +} + +/// `NeoPDF` backend wrapper. +pub struct NeopdfPdf { + pdf: PDF, +} + +impl NeopdfPdf { + /// Loads a NeoPDF member by set name and member index. + #[must_use] + pub fn load(pdf_name: &str, member: usize) -> Self { + Self { + pdf: PDF::load(pdf_name, member), + } + } + + /// Loads a NeoPDF member by its LHAPDF ID (LHAID). + #[must_use] + pub fn load_by_lhaid(lhaid: u32) -> Self { + Self { + pdf: PDF::load_by_lhaid(lhaid), + } + } +} + +impl PdfBackend for NeopdfPdf { + fn alphas_q2(&self, q2: f64) -> f64 { + self.pdf.alphas_q2(q2) + } + + fn particle_id(&self) -> i32 { + self.pdf.metadata().hadron_pid + } + + fn set_force_positive(&mut self, method: ForcePositive) { + match method { + ForcePositive::None => { + self.pdf.set_force_positive(NeoPDFForcePositive::NoClipping); + } + ForcePositive::ClipNegative => { + self.pdf + .set_force_positive(NeoPDFForcePositive::ClipNegative); + } + } + } + + fn x_max(&mut self) -> f64 { + self.pdf.metadata().x_max + } + + fn x_min(&mut self) -> f64 { + self.pdf.metadata().x_min + } + + fn xfx_q2(&self, id: i32, x: f64, q2: f64) -> f64 { + self.pdf.xfxq2(id, &[x, q2]) + } +} + +/// `NeoPDF` set wrapper. +pub struct NeopdfSet { + pdf_name: String, + // num_members: usize, + error_type: String, +} + +impl NeopdfSet { + /// Creates a new NeoPDF set by name. + #[must_use] + pub fn new(pdf_name: &str) -> Self { + // Load member 0 to get metadata + let pdf0 = PDF::load(pdf_name, 0); + let metadata = pdf0.metadata(); + + Self { + pdf_name: pdf_name.to_owned(), + // num_members: metadata.num_members as usize, + error_type: metadata.error_type.clone(), + } + } +} + +impl PdfSetBackend for NeopdfSet { + // fn num_members(&self) -> usize { + // self.num_members + // } + + fn mk_pdfs(&self) -> Result> { + let pdfs = PDF::load_pdfs(&self.pdf_name); + Ok(pdfs + .into_iter() + .map(|pdf| NeopdfPdf { pdf }.into()) + .collect()) + } + + fn uncertainty(&self, values: &[f64], cl: f64, alternative: bool) -> Result { + let unc = + neopdf::uncertainty::uncertainty(values, &self.error_type, CL_1_SIGMA, cl, alternative) + .map_err(|e| anyhow::anyhow!(e))?; + + Ok(Uncertainty { + central: unc.central, + errminus: unc.errminus, + errplus: unc.errplus, + }) + } +} + +/// Creates a single PDF from a set name and member index. +/// +/// # Errors +/// +/// TODO. +pub fn create_pdf(name: &str, member: usize, backend: Backend) -> Result { + match backend { + Backend::Lhapdf => { + if let Ok(lhaid) = name.parse::() { + Ok(LhapdfPdf::with_lhaid(lhaid)?.into()) + } else { + Ok(LhapdfPdf::with_setname_and_member(name, member.try_into().unwrap())?.into()) + } + } + Backend::Neopdf => { + if let Ok(lhaid) = name.parse::() { + Ok(NeopdfPdf::load_by_lhaid(lhaid).into()) + } else { + Ok(NeopdfPdf::load(name, member).into()) + } + } + } +} + +/// Creates a PDF set for uncertainty calculations. +/// +/// # Errors +/// +/// TODO. +pub fn create_pdf_set(name: &str, backend: Backend) -> Result { + match backend { + Backend::Lhapdf => { + // Try parsing as LHAID first + let setname = if let Ok(lhaid) = name.parse::() { + lhapdf::lookup_pdf(lhaid) + .map(|(set, _)| set) + .ok_or_else(|| anyhow!("no convolution function for LHAID = `{lhaid}` found"))? + } else { + name.to_owned() + }; + Ok(LhapdfSet::new(&setname)?.into()) + } + Backend::Neopdf => Ok(NeopdfSet::new(name).into()), + } +} diff --git a/pineappl_cli/src/plot.rs b/pineappl_cli/src/plot.rs index d661a016..825c9e7a 100644 --- a/pineappl_cli/src/plot.rs +++ b/pineappl_cli/src/plot.rs @@ -1,4 +1,5 @@ use super::helpers::{self, ConvFuns, ConvoluteMode}; +use super::pdf_backend::PdfSetBackend as _; use super::{GlobalConfiguration, Subcommand}; use anyhow::Result; use clap::builder::{PossibleValuesParser, TypedValueParser as _}; @@ -63,7 +64,7 @@ impl Subcommand for Opts { }; let grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs[0])?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs[0], cfg.pdf_backend)?; let slices = grid.bwfl().slices(); let mut data_string = String::new(); @@ -159,7 +160,7 @@ impl Subcommand for Opts { .map(|conv_funs| { if self.no_conv_fun_unc { let conv_types = &conv_funs.conv_types; - let mut conv_funs = helpers::create_conv_funs(conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(conv_funs, cfg.pdf_backend)?; let results = helpers::convolve( &grid, @@ -178,6 +179,7 @@ impl Subcommand for Opts { let (set, funs) = helpers::create_conv_funs_for_set( conv_funs, self.conv_fun_uncert_from, + cfg.pdf_backend, )?; let pdf_results: Vec<_> = funs diff --git a/pineappl_cli/src/pull.rs b/pineappl_cli/src/pull.rs index 0b40586d..dd739213 100644 --- a/pineappl_cli/src/pull.rs +++ b/pineappl_cli/src/pull.rs @@ -1,8 +1,8 @@ use super::helpers::{self, ConvFuns, ConvoluteMode}; +use super::pdf_backend::{ConvFunBackend, PdfSetBackend}; use super::{GlobalConfiguration, Subcommand}; use anyhow::{Error, Result}; use clap::{Parser, ValueHint}; -use lhapdf::{Pdf, PdfSet}; use prettytable::{Row, cell}; use rayon::{ThreadPoolBuilder, prelude::*}; use std::num::NonZeroUsize; @@ -53,9 +53,9 @@ impl Subcommand for Opts { let grid = helpers::read_grid(&self.input)?; let (set1, mut conv_funs1) = - helpers::create_conv_funs_for_set(&self.conv_funs1, self.pull_from)?; + helpers::create_conv_funs_for_set(&self.conv_funs1, self.pull_from, cfg.pdf_backend)?; let (set2, mut conv_funs2) = - helpers::create_conv_funs_for_set(&self.conv_funs2, self.pull_from)?; + helpers::create_conv_funs_for_set(&self.conv_funs2, self.pull_from, cfg.pdf_backend)?; ThreadPoolBuilder::new() .num_threads(self.threads) @@ -141,73 +141,75 @@ impl Subcommand for Opts { (diff / unc1.hypot(unc2), unc1, unc2) }; - let channel_results = - |conv_funs: &ConvFuns, pdfset: &mut [Vec], set: &PdfSet| -> Vec { - if let Some(member) = conv_funs.members[self.pull_from] { - (0..grid.channels().len()) - .map(|channel| { - let mut channel_mask = vec![false; grid.channels().len()]; - channel_mask[channel] = true; - match helpers::convolve( - &grid, - &mut pdfset[member], - &conv_funs.conv_types, - &self.orders, - &[bin], - &channel_mask, - 1, - ConvoluteMode::Normal, - cfg, - ) - .as_slice() - { - [value] => *value, - _ => unreachable!(), - } - }) - .collect() - } else { - let results: Vec<_> = pdfset - .iter_mut() - .flat_map(|fun| { - (0..grid.channels().len()) - .map(|channel| { - let mut channel_mask = vec![false; grid.channels().len()]; - channel_mask[channel] = true; - match helpers::convolve( - &grid, - fun, - &conv_funs.conv_types, - &self.orders, - &[bin], - &channel_mask, - 1, - ConvoluteMode::Normal, - cfg, - ) - .as_slice() - { - [value] => *value, - _ => unreachable!(), - } - }) - .collect::>() - }) - .collect(); - - (0..grid.channels().len()) - .map(|channel| { - let central: Vec<_> = results - .iter() - .skip(channel) - .step_by(grid.channels().len()) - .copied() - .collect(); - set.uncertainty(¢ral, self.cl, false).unwrap().central - }) - .collect() - } - }; + let channel_results = |conv_funs: &ConvFuns, + pdfset: &mut [Vec], + set: &dyn PdfSetBackend| + -> Vec { + if let Some(member) = conv_funs.members[self.pull_from] { + (0..grid.channels().len()) + .map(|channel| { + let mut channel_mask = vec![false; grid.channels().len()]; + channel_mask[channel] = true; + match helpers::convolve( + &grid, + &mut pdfset[member], + &conv_funs.conv_types, + &self.orders, + &[bin], + &channel_mask, + 1, + ConvoluteMode::Normal, + cfg, + ) + .as_slice() + { + [value] => *value, + _ => unreachable!(), + } + }) + .collect() + } else { + let results: Vec<_> = pdfset + .iter_mut() + .flat_map(|fun| { + (0..grid.channels().len()) + .map(|channel| { + let mut channel_mask = vec![false; grid.channels().len()]; + channel_mask[channel] = true; + match helpers::convolve( + &grid, + fun, + &conv_funs.conv_types, + &self.orders, + &[bin], + &channel_mask, + 1, + ConvoluteMode::Normal, + cfg, + ) + .as_slice() + { + [value] => *value, + _ => unreachable!(), + } + }) + .collect::>() + }) + .collect(); + + (0..grid.channels().len()) + .map(|channel| { + let central: Vec<_> = results + .iter() + .skip(channel) + .step_by(grid.channels().len()) + .copied() + .collect(); + set.uncertainty(¢ral, self.cl, false).unwrap().central + }) + .collect() + } + }; let mut pull_tuples = if self.limit == 0 { Vec::new() diff --git a/pineappl_cli/src/uncert.rs b/pineappl_cli/src/uncert.rs index f116c496..063218e6 100644 --- a/pineappl_cli/src/uncert.rs +++ b/pineappl_cli/src/uncert.rs @@ -1,4 +1,5 @@ use super::helpers::{self, ConvFuns, ConvoluteMode}; +use super::pdf_backend::PdfSetBackend as _; use super::{GlobalConfiguration, Subcommand}; use anyhow::{Error, Result}; use clap::builder::{PossibleValuesParser, TypedValueParser as _}; @@ -66,7 +67,7 @@ pub struct Opts { #[command(flatten)] group: Group, /// Confidence level in per cent, for convolution function uncertainties. - #[arg(default_value_t = lhapdf::CL_1_SIGMA, long)] + #[arg(default_value_t = super::pdf_backend::CL_1_SIGMA, long)] cl: f64, /// Show integrated numbers (without bin widths) instead of differential ones. #[arg(long, short)] @@ -94,7 +95,7 @@ pub struct Opts { impl Subcommand for Opts { fn run(&self, cfg: &GlobalConfiguration) -> Result { let grid = helpers::read_grid(&self.input)?; - let mut conv_funs = helpers::create_conv_funs(&self.conv_funs)?; + let mut conv_funs = helpers::create_conv_funs(&self.conv_funs, cfg.pdf_backend)?; let limits = helpers::convolve_limits( &grid, @@ -116,7 +117,8 @@ impl Subcommand for Opts { .conv_fun .iter() .map(|&index| { - let (set, funs) = helpers::create_conv_funs_for_set(&self.conv_funs, index)?; + let (set, funs) = + helpers::create_conv_funs_for_set(&self.conv_funs, index, cfg.pdf_backend)?; let results: Vec<_> = funs .into_par_iter() .map(|mut funs| { @@ -145,7 +147,7 @@ impl Subcommand for Opts { .map(|pdf| results[pdf][bin]) .collect::>() }) - .map(|values| Ok(set.uncertainty(&values, self.cl, false)?)) + .map(|values| set.uncertainty(&values, self.cl, false)) .collect::>() }) .collect::>()?; diff --git a/pineappl_cli/tests/convolve.rs b/pineappl_cli/tests/convolve.rs index 95c66fd9..fcca14ba 100644 --- a/pineappl_cli/tests/convolve.rs +++ b/pineappl_cli/tests/convolve.rs @@ -536,6 +536,117 @@ fn issue_334() { .stdout(NO_CHANNELS_GRID_STR); } +#[test] +fn neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "convolve", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(str::contains(DEFAULT_STR)); +} + +#[test] +fn lhaid_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "convolve", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "324900", + ]) + .assert() + .success() + .stdout(DEFAULT_STR); +} + +#[test] +fn neopdf_multiple_pdfs() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "convolve", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + "324900=NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(DEFAULT_MULTIPLE_PDFS_STR); +} + +#[test] +fn integrated_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "convolve", + "--integrated", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(INTEGRATED_STR); +} + +#[test] +fn three_pdfs_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "convolve", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed/0", + "NNPDF31_nlo_as_0118_luxqed/1", + "NNPDF31_nlo_as_0118_luxqed/2", + ]) + .assert() + .success() + .stdout(THREE_PDFS_STR); +} + +#[test] +fn bins_13567_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "convolve", + "--bins=1,3,5-7", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(BINS_13567_STR); +} + +#[test] +fn force_positive_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--force-positive", + "--pdf-backend=neopdf", + "convolve", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(FORCE_POSITIVE_STR); +} + #[test] fn lagrange_subgrid_v1() { Command::cargo_bin("pineappl") diff --git a/pineappl_cli/tests/export.rs b/pineappl_cli/tests/export.rs index d1bc766e..5a176efc 100644 --- a/pineappl_cli/tests/export.rs +++ b/pineappl_cli/tests/export.rs @@ -1,10 +1,11 @@ #![expect(missing_docs, reason = "non-public items will not be documented")] use assert_cmd::Command; -use predicates::str as pred_str; #[cfg(feature = "applgrid")] use assert_fs::NamedTempFile; +#[cfg(feature = "applgrid")] +use predicates::str as pred_str; const HELP_STR: &str = "Converts PineAPPL grids to APPLgrid files diff --git a/pineappl_cli/tests/import.rs b/pineappl_cli/tests/import.rs index 55282537..d9920ad9 100644 --- a/pineappl_cli/tests/import.rs +++ b/pineappl_cli/tests/import.rs @@ -1,10 +1,11 @@ #![expect(missing_docs, reason = "non-public items will not be documented")] use assert_cmd::Command; -use predicates::str as pred_str; #[cfg(any(feature = "applgrid", feature = "fastnlo", feature = "fktable"))] use assert_fs::NamedTempFile; +#[cfg(any(feature = "applgrid", feature = "fastnlo"))] +use predicates::str as pred_str; #[cfg(feature = "fastnlo")] const HELP_STR: &str = "Converts APPLgrid/fastNLO/FastKernel files to PineAPPL grids diff --git a/pineappl_cli/tests/main.rs b/pineappl_cli/tests/main.rs index 2663622e..2ddc5869 100644 --- a/pineappl_cli/tests/main.rs +++ b/pineappl_cli/tests/main.rs @@ -29,6 +29,7 @@ Options: --force-positive Forces negative PDF values to zero --allow-extrapolation Allow extrapolation of PDFs outside their region of validity --use-alphas-from Choose the PDF/FF set for the strong coupling [default: 0] + --pdf-backend Select the PDF interpolation backend: 'lhapdf' or 'neopdf' [default: lhapdf] -h, --help Print help -V, --version Print version "; diff --git a/pineappl_cli/tests/uncert.rs b/pineappl_cli/tests/uncert.rs index dfcd0067..30ac3c39 100644 --- a/pineappl_cli/tests/uncert.rs +++ b/pineappl_cli/tests/uncert.rs @@ -1,6 +1,7 @@ #![expect(missing_docs, reason = "non-public items will not be documented")] use assert_cmd::Command; +use predicates::str; use std::num::NonZeroUsize; use std::thread; @@ -204,6 +205,28 @@ const SCALE_ENV_9_STR: &str = "b etal dsig/detal 9pt-svar (env) 7 4 4.5 2.7517266e1 -5.36 5.22 "; +const NEOPDF_CONV_FUN_STR: &str = "-+----+----+-----------+-----------+---------+--------- +0 2 2.25 7.5459110e2 7.5461655e2 -1.14 1.14 +1 2.25 2.5 6.9028342e2 6.9027941e2 -1.16 1.16 +2 2.5 2.75 6.0025198e2 6.0022595e2 -1.18 1.18 +3 2.75 3 4.8552235e2 4.8548211e2 -1.22 1.22 +4 3 3.25 3.6195456e2 3.6191001e2 -1.27 1.27 +5 3.25 3.5 2.4586691e2 2.4582640e2 -1.35 1.35 +6 3.5 4 1.1586851e2 1.1584074e2 -1.51 1.51 +7 4 4.5 2.7517266e1 2.7504644e1 -2.77 2.77 +"; + +const NEOPDF_CONV_FUN_CL_90_STR: &str = "-+----+----+-----------+-----------+---------+--------- +0 2 2.25 7.5459110e2 7.5461655e2 -1.87 1.87 +1 2.25 2.5 6.9028342e2 6.9027941e2 -1.90 1.90 +2 2.5 2.75 6.0025198e2 6.0022595e2 -1.95 1.95 +3 2.75 3 4.8552235e2 4.8548211e2 -2.00 2.00 +4 3 3.25 3.6195456e2 3.6191001e2 -2.08 2.08 +5 3.25 3.5 2.4586691e2 2.4582640e2 -2.22 2.22 +6 3.5 4 1.1586851e2 1.1584074e2 -2.48 2.48 +7 4 4.5 2.7517266e1 2.7504644e1 -4.56 4.56 +"; + #[test] fn help() { Command::cargo_bin("pineappl") @@ -437,3 +460,88 @@ fn scale_env_9() { .success() .stdout(SCALE_ENV_9_STR); } + +#[test] +fn conv_fun_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "uncert", + "--conv-fun", + "--threads=1", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(str::contains(NEOPDF_CONV_FUN_STR)); +} + +#[test] +fn conv_fun_cl_90_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "uncert", + "--conv-fun", + "--cl=90", + "--threads=1", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(str::contains(NEOPDF_CONV_FUN_CL_90_STR)); +} + +#[test] +fn scale_env_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "uncert", + "--scale-env", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(str::contains(SCALE_ENV_STR)); +} + +#[test] +fn conv_fun_orders_a2_as1a2_neopdf_beckend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "uncert", + "--conv-fun=0", + "--orders=a2,as1a2", + "--threads=1", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(ORDERS_A2_AS1A2_STR); +} + +#[test] +fn scale_abs_9_neopdf_backend() { + Command::cargo_bin("pineappl") + .unwrap() + .args([ + "--pdf-backend=neopdf", + "uncert", + "--scale-abs=9", + "../test-data/LHCB_WP_7TEV_opt.pineappl.lz4", + "NNPDF31_nlo_as_0118_luxqed", + ]) + .assert() + .success() + .stdout(SCALE_ABS_9_STR); +}