diff --git a/pineappl/src/fk_table.rs b/pineappl/src/fk_table.rs index 0ec4e64e..9300b42d 100644 --- a/pineappl/src/fk_table.rs +++ b/pineappl/src/fk_table.rs @@ -1,10 +1,12 @@ //! Provides the [`FkTable`] type. -use super::boc::{Kinematics, Order}; +use super::boc::{Channel, Kinematics, Order}; use super::convolutions::ConvolutionCache; +use super::empty_subgrid::EmptySubgridV1; use super::grid::Grid; +use super::pids::OptRules; use super::subgrid::{self, Subgrid}; -use ndarray::ArrayD; +use ndarray::{s, ArrayD}; use std::fmt::{self, Display, Formatter}; use std::iter; use std::str::FromStr; @@ -240,70 +242,40 @@ impl FkTable { ) } - /// Optimizes the storage of FK tables based of assumptions of the PDFs at the FK table's - /// scale. - /// - /// # Panics - /// - /// TODO + /// Optimize the size of this FK-table by throwing away heavy quark flavors assumed to be zero + /// at the FK-table's scales and calling [`Grid::optimize`]. pub fn optimize(&mut self, assumptions: FkAssumptions) { - let mut add = Vec::new(); + let OptRules(sum, delete) = self.grid.pid_basis().opt_rules(assumptions); - match assumptions { - FkAssumptions::Nf6Ind => { - // nothing to do here - } - FkAssumptions::Nf6Sym => { - add.push((235, 200)); - } - FkAssumptions::Nf5Ind => { - add.extend_from_slice(&[(235, 200), (135, 100)]); - } - FkAssumptions::Nf5Sym => { - add.extend_from_slice(&[(235, 200), (135, 100), (224, 200)]); - } - FkAssumptions::Nf4Ind => { - add.extend_from_slice(&[(235, 200), (135, 100), (224, 200), (124, 100)]); - } - FkAssumptions::Nf4Sym => { - add.extend_from_slice(&[ - (235, 200), - (135, 100), - (224, 200), - (124, 100), - (215, 200), - ]); - } - FkAssumptions::Nf3Ind => { - add.extend_from_slice(&[ - (235, 200), - (135, 100), - (224, 200), - (124, 100), - (215, 200), - (115, 100), - ]); - } - FkAssumptions::Nf3Sym => { - add.extend_from_slice(&[ - (235, 200), - (135, 100), - (224, 200), - (124, 100), - (215, 200), - (115, 100), - (208, 200), - ]); + for idx in 0..self.grid.channels().len() { + let &[(ref pids, factor)] = self.grid.channels()[idx].entry() else { + // every FK-table must have a trivial channel definition + unreachable!() + }; + let mut pids = pids.clone(); + + for pid in &mut pids { + if delete.iter().any(|&delete| *pid == delete) { + for subgrid in self.grid.subgrids_mut().slice_mut(s![.., .., idx]) { + *subgrid = EmptySubgridV1.into(); + } + } else if let Some(replace) = sum + .iter() + .find_map(|&(search, replace)| (*pid == search).then_some(replace)) + { + *pid = replace; + } } + + self.grid.channels_mut()[idx] = Channel::new(vec![(pids, factor)]); } - self.grid.rewrite_channels(&add, &[]); + self.grid.optimize(); // store the assumption so that we can check it later on self.grid .metadata_mut() .insert("fk_assumptions".to_owned(), assumptions.to_string()); - self.grid.optimize(); } } diff --git a/pineappl/src/grid.rs b/pineappl/src/grid.rs index 9b064d55..55f9c769 100644 --- a/pineappl/src/grid.rs +++ b/pineappl/src/grid.rs @@ -1505,53 +1505,6 @@ impl Grid { } } - pub(crate) fn rewrite_channels(&mut self, add: &[(i32, i32)], del: &[i32]) { - // TODO: generalize this method to n convolutions - assert_eq!(self.convolutions().len(), 2); - - self.channels = self - .channels() - .iter() - .map(|entry| { - Channel::new( - entry - .entry() - .iter() - .map(|(pids, f)| { - ( - vec![ - // if `a` is to be added to another pid replace it with this pid - add.iter().fold(pids[0], |id, &(source, target)| { - if id == source { - target - } else { - id - } - }), - // if `b` is to be added to another pid replace it with this pid - add.iter().fold(pids[1], |id, &(source, target)| { - if id == source { - target - } else { - id - } - }), - ], - // if any of the pids `a` or `b` are to b deleted set the factor to - // zero - if del.iter().any(|&id| id == pids[0] || id == pids[1]) { - 0.0 - } else { - *f - }, - ) - }) - .collect(), - ) - }) - .collect(); - } - /// Splits the grid such that each channel contains only a single tuple of PIDs. pub fn split_channels(&mut self) { let indices: Vec<_> = self diff --git a/pineappl/src/pids.rs b/pineappl/src/pids.rs index e7453857..329cc6fd 100644 --- a/pineappl/src/pids.rs +++ b/pineappl/src/pids.rs @@ -1,6 +1,7 @@ //! TODO use super::boc::Channel; +use super::fk_table::FkAssumptions; use float_cmp::approx_eq; use serde::{Deserialize, Serialize}; use std::str::FromStr; @@ -114,8 +115,66 @@ impl PidBasis { (&Self::Evol, Self::Pdg) => channel.translate(&evol_to_pdg_mc_ids), } } + + /// TODO + #[must_use] + pub fn opt_rules(&self, assumptions: FkAssumptions) -> OptRules { + match (*self, assumptions) { + (Self::Evol | Self::Pdg, FkAssumptions::Nf6Ind) => OptRules(Vec::new(), Vec::new()), + (Self::Evol, FkAssumptions::Nf6Sym) => OptRules(vec![(235, 200)], Vec::new()), + (Self::Evol, FkAssumptions::Nf5Ind) => { + OptRules(vec![(235, 200), (135, 100)], Vec::new()) + } + (Self::Evol, FkAssumptions::Nf5Sym) => { + OptRules(vec![(235, 200), (135, 100), (224, 200)], Vec::new()) + } + (Self::Evol, FkAssumptions::Nf4Ind) => OptRules( + vec![(235, 200), (135, 100), (224, 200), (124, 100)], + Vec::new(), + ), + (Self::Evol, FkAssumptions::Nf4Sym) => OptRules( + vec![(235, 200), (135, 100), (224, 200), (124, 100), (215, 200)], + Vec::new(), + ), + (Self::Evol, FkAssumptions::Nf3Ind) => OptRules( + vec![ + (235, 200), + (135, 100), + (224, 200), + (124, 100), + (215, 200), + (115, 100), + ], + Vec::new(), + ), + (Self::Evol, FkAssumptions::Nf3Sym) => OptRules( + vec![ + (235, 200), + (135, 100), + (224, 200), + (124, 100), + (215, 200), + (115, 100), + (208, 200), + ], + Vec::new(), + ), + (Self::Pdg, FkAssumptions::Nf6Sym) => OptRules(vec![(-6, 6)], Vec::new()), + (Self::Pdg, FkAssumptions::Nf5Ind) => OptRules(Vec::new(), vec![-6, 6]), + (Self::Pdg, FkAssumptions::Nf5Sym) => OptRules(vec![(-5, 5)], vec![-6, 6]), + (Self::Pdg, FkAssumptions::Nf4Ind) => OptRules(Vec::new(), vec![-6, 6, -5, 5]), + (Self::Pdg, FkAssumptions::Nf4Sym) => OptRules(vec![(-4, 4)], vec![-6, 6, -5, 5]), + (Self::Pdg, FkAssumptions::Nf3Ind) => OptRules(Vec::new(), vec![-6, 6, -5, 5, -4, 4]), + (Self::Pdg, FkAssumptions::Nf3Sym) => { + OptRules(vec![(-3, 3)], vec![-6, 6, -5, 5, -4, 4]) + } + } + } } +/// Return type of [`PidBasis::optimization_rules`]. +pub struct OptRules(pub Vec<(i32, i32)>, pub Vec); + /// Error returned by [`PidBasis::from_str`] when passed with an unknown argument. #[derive(Debug, Error)] #[error("unknown PID basis: {basis}")]