From 2a4f79c87d13b31010e7125b7d72fcf9a8c11cb0 Mon Sep 17 00:00:00 2001 From: t7phy Date: Fri, 4 Oct 2024 15:47:05 +0200 Subject: [PATCH] Save convolution type in `ConvolutionCache` --- pineappl/src/convolutions.rs | 12 +++---- pineappl/tests/drell_yan_lo.rs | 13 +++++-- pineappl_capi/src/lib.rs | 16 +++++++-- pineappl_cli/src/helpers.rs | 66 +++++++++++++++++++++++++++------- pineappl_cli/tests/import.rs | 6 +++- 5 files changed, 87 insertions(+), 26 deletions(-) diff --git a/pineappl/src/convolutions.rs b/pineappl/src/convolutions.rs index 2f88034d..f4a7f5b6 100644 --- a/pineappl/src/convolutions.rs +++ b/pineappl/src/convolutions.rs @@ -20,14 +20,14 @@ pub struct ConvolutionCache<'a> { imur2: Vec, imuf2: Vec, ix: Vec>, - pdg: Vec, + pdg: Vec, perm: Vec>, } impl<'a> ConvolutionCache<'a> { /// TODO pub fn new( - pdg: Vec, + pdg: Vec, xfx: Vec<&'a mut dyn FnMut(i32, f64, f64) -> f64>, alphas: &'a mut dyn FnMut(f64) -> f64, ) -> Self { @@ -58,16 +58,16 @@ impl<'a> ConvolutionCache<'a> { .iter() .enumerate() .map(|(max_idx, conv)| { - conv.pid().map(|pid| { + conv.pid().map(|_| { self.pdg .iter() .take(max_idx + 1) .enumerate() .rev() - .find_map(|(idx, &pdg)| { - if pid == pdg { + .find_map(|(idx, pdg)| { + if conv == pdg { Some((idx, false)) - } else if pid == pids::charge_conjugate_pdg_pid(pdg) { + } else if *conv == pdg.charge_conjugate() { Some((idx, true)) } else { None diff --git a/pineappl/tests/drell_yan_lo.rs b/pineappl/tests/drell_yan_lo.rs index 1dd30879..3b411ebd 100644 --- a/pineappl/tests/drell_yan_lo.rs +++ b/pineappl/tests/drell_yan_lo.rs @@ -351,7 +351,11 @@ fn perform_grid_tests( grid.scale_by_order(10.0, 1.0, 10.0, 10.0, 4.0); // TEST 5: `convolve` - let mut convolution_cache = ConvolutionCache::new(vec![2212], vec![&mut xfx], &mut alphas); + let mut convolution_cache = ConvolutionCache::new( + vec![Convolution::UnpolPDF(2212)], + vec![&mut xfx], + &mut alphas, + ); let bins = grid.convolve(&mut convolution_cache, &[], &[], &[], &[(1.0, 1.0, 1.0)]); for (result, reference) in bins.iter().zip(reference.iter()) { @@ -362,8 +366,11 @@ fn perform_grid_tests( let mut xfx1 = |id, x, q2| pdf.xfx_q2(id, x, q2); let mut xfx2 = |id, x, q2| pdf.xfx_q2(id, x, q2); let mut alphas2 = |_| 0.0; - let mut convolution_cache2 = - ConvolutionCache::new(vec![2212, 2212], vec![&mut xfx1, &mut xfx2], &mut alphas2); + let mut convolution_cache2 = ConvolutionCache::new( + vec![Convolution::UnpolPDF(2212), Convolution::UnpolPDF(2212)], + vec![&mut xfx1, &mut xfx2], + &mut alphas2, + ); let bins2 = grid.convolve(&mut convolution_cache2, &[], &[], &[], &[(1.0, 1.0, 1.0)]); for (result, reference) in bins2.iter().zip(reference.iter()) { diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 290c7d04..7343e8b6 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -455,7 +455,11 @@ pub unsafe extern "C" fn pineappl_grid_convolve_with_one( unsafe { slice::from_raw_parts(channel_mask, grid.channels().len()) }.to_vec() }; let results = unsafe { slice::from_raw_parts_mut(results, grid.bin_info().bins()) }; - let mut convolution_cache = ConvolutionCache::new(vec![pdg_id], vec![&mut xfx], &mut als); + let mut convolution_cache = ConvolutionCache::new( + vec![Convolution::UnpolPDF(pdg_id)], + vec![&mut xfx], + &mut als, + ); results.copy_from_slice(&grid.convolve( &mut convolution_cache, @@ -517,8 +521,14 @@ pub unsafe extern "C" fn pineappl_grid_convolve_with_two( unsafe { slice::from_raw_parts(channel_mask, grid.channels().len()) }.to_vec() }; let results = unsafe { slice::from_raw_parts_mut(results, grid.bin_info().bins()) }; - let mut convolution_cache = - ConvolutionCache::new(vec![pdg_id1, pdg_id2], vec![&mut xfx1, &mut xfx2], &mut als); + let mut convolution_cache = ConvolutionCache::new( + vec![ + Convolution::UnpolPDF(pdg_id1), + Convolution::UnpolPDF(pdg_id2), + ], + vec![&mut xfx1, &mut xfx2], + &mut als, + ); results.copy_from_slice(&grid.convolve( &mut convolution_cache, diff --git a/pineappl_cli/src/helpers.rs b/pineappl_cli/src/helpers.rs index 30f6acf4..9dbcf19e 100644 --- a/pineappl_cli/src/helpers.rs +++ b/pineappl_cli/src/helpers.rs @@ -2,7 +2,7 @@ use super::GlobalConfiguration; use anyhow::{anyhow, ensure, Context, Error, Result}; use lhapdf::{Pdf, PdfSet}; use ndarray::{Array3, Ix3}; -use pineappl::convolutions::ConvolutionCache; +use pineappl::convolutions::{Convolution, ConvolutionCache}; use pineappl::grid::Grid; use prettytable::format::{FormatBuilder, LinePosition, LineSeparator}; use prettytable::Table; @@ -279,19 +279,39 @@ pub fn convolve_scales( .iter() .map(|fun| move |q2| fun.alphas_q2(q2)) .collect(); - let pdg_ids: Vec<_> = conv_funs + let convolutions: Vec<_> = conv_funs .iter() - .map(|fun| { - // if the field 'Particle' is missing we assume it's a proton PDF - fun.set() + .zip(grid.convolutions()) + .map(|(fun, convolution)| { + 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() + .unwrap(); + + match fun.set().entry("SetType").unwrap_or_default().as_str() { + "fragfn" => Convolution::UnpolFF(pid), + "" => { + // if we can not figure out the type of the convolution from the PDF set, we + // assume it from the grid convolution at the same index + match convolution { + // if the grid convolution is None, we assume it's a proton PDF + Convolution::None | Convolution::UnpolPDF(_) => Convolution::UnpolPDF(pid), + Convolution::PolPDF(_) => Convolution::PolPDF(pid), + Convolution::UnpolFF(_) => Convolution::UnpolFF(pid), + Convolution::PolFF(_) => Convolution::PolFF(pid), + } + } + // TODO: convince the LHAPDF maintainers to make SetType necessary for polarized + // PDFs and all FFs + _ => unimplemented!(), + } }) .collect(); - let mut cache = ConvolutionCache::new(pdg_ids, xfx, &mut alphas_funs[cfg.use_alphas_from]); + let mut cache = ConvolutionCache::new(convolutions, xfx, &mut alphas_funs[cfg.use_alphas_from]); let mut results = grid.convolve(&mut cache, &orders, bins, channels, scales); match mode { @@ -430,19 +450,39 @@ pub fn convolve_subgrid( .iter() .map(|fun| move |q2| fun.alphas_q2(q2)) .collect(); - let pdg_ids: Vec<_> = conv_funs + let convolutions: Vec<_> = conv_funs .iter() - .map(|fun| { - // if the field 'Particle' is missing we assume it's a proton PDF - fun.set() + .zip(grid.convolutions()) + .map(|(fun, convolution)| { + 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() + .unwrap(); + + match fun.set().entry("SetType").unwrap_or_default().as_str() { + "fragfn" => Convolution::UnpolFF(pid), + "" => { + // if we can not figure out the type of the convolution from the PDF set, we + // assume it from the grid convolution at the same index + match convolution { + // if the grid convolution is None, we assume it's a proton PDF + Convolution::None | Convolution::UnpolPDF(_) => Convolution::UnpolPDF(pid), + Convolution::PolPDF(_) => Convolution::PolPDF(pid), + Convolution::UnpolFF(_) => Convolution::UnpolFF(pid), + Convolution::PolFF(_) => Convolution::PolFF(pid), + } + } + // TODO: convince the LHAPDF maintainers to make SetType necessary for polarized + // PDFs and all FFs + _ => unimplemented!(), + } }) .collect(); - let mut cache = ConvolutionCache::new(pdg_ids, xfx, &mut alphas_funs[cfg.use_alphas_from]); + let mut cache = ConvolutionCache::new(convolutions, xfx, &mut alphas_funs[cfg.use_alphas_from]); let subgrid = grid.convolve_subgrid(&mut cache, order, bin, lumi, (1.0, 1.0, 1.0)); subgrid diff --git a/pineappl_cli/tests/import.rs b/pineappl_cli/tests/import.rs index 6501c9c0..871b951b 100644 --- a/pineappl_cli/tests/import.rs +++ b/pineappl_cli/tests/import.rs @@ -893,7 +893,11 @@ fn import_hadronic_fktable() { let pdf = Pdf::with_setname_and_member("NNPDF31_nlo_as_0118_luxqed", 0).unwrap(); let mut xfx = |id, x, q2| pdf.xfx_q2(id, x, q2); let mut alphas = |_| 0.0; - let mut convolution_cache = ConvolutionCache::new(vec![2212], vec![&mut xfx], &mut alphas); + let mut convolution_cache = ConvolutionCache::new( + vec![Convolution::UnpolPDF(2212)], + vec![&mut xfx], + &mut alphas, + ); let results = grid.convolve(&mut convolution_cache, &[], &[], &[], &[(1.0, 1.0, 1.0)]); let mut fk_table = FkTable::try_from(grid).unwrap();