From 366ace1f5592f2afb9481530c56533c3d4d2ac06 Mon Sep 17 00:00:00 2001 From: t7phy Date: Mon, 28 Oct 2024 15:40:02 +0100 Subject: [PATCH] Generalize static node value detection --- pineappl/src/empty_subgrid.rs | 8 ++-- pineappl/src/grid.rs | 13 ++----- pineappl/src/import_subgrid.rs | 36 +++-------------- pineappl/src/interp_subgrid.rs | 71 +++++++++++++++++++++------------- pineappl/src/interpolation.rs | 25 +++++++++--- pineappl/src/subgrid.rs | 4 +- 6 files changed, 80 insertions(+), 77 deletions(-) diff --git a/pineappl/src/empty_subgrid.rs b/pineappl/src/empty_subgrid.rs index 1adf05e8..de33d72f 100644 --- a/pineappl/src/empty_subgrid.rs +++ b/pineappl/src/empty_subgrid.rs @@ -1,7 +1,7 @@ //! TODO use super::interpolation::Interp; -use super::subgrid::{Mu2, Stats, Subgrid, SubgridEnum, SubgridIndexedIter}; +use super::subgrid::{Stats, Subgrid, SubgridEnum, SubgridIndexedIter}; use serde::{Deserialize, Serialize}; use std::iter; @@ -47,9 +47,7 @@ impl Subgrid for EmptySubgridV1 { } } - fn static_scale(&self) -> Option { - None - } + fn optimize_static_nodes(&mut self) {} } #[cfg(test)] @@ -66,6 +64,7 @@ mod tests { subgrid.merge(&EmptySubgridV1.into(), None); subgrid.scale(2.0); subgrid.symmetrize(1, 2); + subgrid.optimize_static_nodes(); assert_eq!( subgrid.stats(), Stats { @@ -76,7 +75,6 @@ mod tests { bytes_per_value: 0, } ); - assert_eq!(subgrid.static_scale(), None); } #[test] diff --git a/pineappl/src/grid.rs b/pineappl/src/grid.rs index b9a5a924..6d73cb1a 100644 --- a/pineappl/src/grid.rs +++ b/pineappl/src/grid.rs @@ -850,7 +850,7 @@ impl Grid { } } - fn optimize_subgrid_type(&mut self, static_scale_detection: bool) { + fn optimize_subgrid_type(&mut self, optimize_static_nodes: bool) { for subgrid in &mut self.subgrids { match subgrid { // replace empty subgrids of any type with `EmptySubgridV1` @@ -858,15 +858,10 @@ impl Grid { *subgrid = EmptySubgridV1.into(); } _ => { - // TODO: this requires a `pub(crate)` in `InterpSubgridV1`; we should - // replace this with a method - if !static_scale_detection { - if let SubgridEnum::InterpSubgridV1(subgrid) = subgrid { - // disable static-scale detection - subgrid.static_q2 = -1.0; - } + if optimize_static_nodes { + subgrid.optimize_static_nodes(); } - + // TODO: check if we should remove this *subgrid = ImportSubgridV1::from(&*subgrid).into(); } } diff --git a/pineappl/src/import_subgrid.rs b/pineappl/src/import_subgrid.rs index d181118a..314ec1cb 100644 --- a/pineappl/src/import_subgrid.rs +++ b/pineappl/src/import_subgrid.rs @@ -3,8 +3,8 @@ use super::evolution::EVOLVE_INFO_TOL_ULPS; use super::interpolation::Interp; use super::packed_array::PackedArray; -use super::subgrid::{Mu2, Stats, Subgrid, SubgridEnum, SubgridIndexedIter}; -use float_cmp::{approx_eq, assert_approx_eq}; +use super::subgrid::{Stats, Subgrid, SubgridEnum, SubgridIndexedIter}; +use float_cmp::approx_eq; use itertools::izip; use serde::{Deserialize, Serialize}; use std::mem; @@ -130,17 +130,7 @@ impl Subgrid for ImportSubgridV1 { } } - fn static_scale(&self) -> Option { - if let &[static_scale] = self.node_values()[0].as_slice() { - Some(Mu2 { - ren: static_scale, - fac: static_scale, - frg: -1.0, - }) - } else { - None - } - } + fn optimize_static_nodes(&mut self) {} } impl From<&SubgridEnum> for ImportSubgridV1 { @@ -161,32 +151,18 @@ impl From<&SubgridEnum> for ImportSubgridV1 { }, ); - let mut new_node_values: Vec<_> = subgrid + let new_node_values: Vec<_> = subgrid .node_values() .iter() .zip(&ranges) .map(|(values, range)| values[range.clone()].to_vec()) .collect(); - let static_scale = if let Some(Mu2 { ren, fac, frg }) = subgrid.static_scale() { - assert_approx_eq!(f64, ren, fac, ulps = 4); - assert_approx_eq!(f64, frg, -1.0, ulps = 4); - new_node_values[0] = vec![fac]; - true - } else { - false - }; let mut array = PackedArray::new(new_node_values.iter().map(Vec::len).collect()); for (mut indices, value) in subgrid.indexed_iter() { - for (idx, (index, range)) in indices.iter_mut().zip(&ranges).enumerate() { - // TODO: generalize static scale detection - *index = if static_scale && idx == 0 { - // if there's a static scale we want every value to be added to same grid point - 0 - } else { - *index - range.start - }; + for (index, range) in indices.iter_mut().zip(&ranges) { + *index = *index - range.start; } array[indices.as_slice()] += value; diff --git a/pineappl/src/interp_subgrid.rs b/pineappl/src/interp_subgrid.rs index d8b1b6a0..da749470 100644 --- a/pineappl/src/interp_subgrid.rs +++ b/pineappl/src/interp_subgrid.rs @@ -2,7 +2,7 @@ use super::interpolation::{self, Interp}; use super::packed_array::PackedArray; -use super::subgrid::{Mu2, Stats, Subgrid, SubgridEnum, SubgridIndexedIter}; +use super::subgrid::{Stats, Subgrid, SubgridEnum, SubgridIndexedIter}; use float_cmp::approx_eq; use serde::{Deserialize, Serialize}; use std::mem; @@ -12,7 +12,7 @@ use std::mem; pub struct InterpSubgridV1 { array: PackedArray, interps: Vec, - pub(crate) static_q2: f64, + static_nodes: Vec>, } impl InterpSubgridV1 { @@ -22,7 +22,7 @@ impl InterpSubgridV1 { Self { array: PackedArray::new(interps.iter().map(Interp::nodes).collect()), interps: interps.to_vec(), - static_q2: 0.0, + static_nodes: vec![Some(-1.0); interps.len()], } } } @@ -32,14 +32,14 @@ impl Subgrid for InterpSubgridV1 { debug_assert_eq!(interps.len(), ntuple.len()); if interpolation::interpolate(interps, ntuple, weight, &mut self.array) { - // TODO: make this more general - let q2 = ntuple[0]; - if self.static_q2 == 0.0 { - self.static_q2 = q2; - } else if !approx_eq!(f64, self.static_q2, -1.0, ulps = 4) - && !approx_eq!(f64, self.static_q2, q2, ulps = 4) - { - self.static_q2 = -1.0; + for (value, previous_node) in ntuple.iter().zip(&mut self.static_nodes) { + if let Some(previous_value) = previous_node { + if *previous_value < 0.0 { + *previous_value = *value; + } else if !approx_eq!(f64, *previous_value, *value, ulps = 4) { + *previous_node = None; + } + } } } } @@ -110,12 +110,40 @@ impl Subgrid for InterpSubgridV1 { } } - fn static_scale(&self) -> Option { - (self.static_q2 > 0.0).then_some(Mu2 { - ren: self.static_q2, - fac: self.static_q2, - frg: -1.0, - }) + fn optimize_static_nodes(&mut self) { + let mut new_array = PackedArray::new( + self.array + .shape() + .iter() + .zip(&self.static_nodes) + .map(|(&dim, static_node)| if static_node.is_some() { 1 } else { dim }) + .collect(), + ); + + for (mut index, value) in self.array.indexed_iter() { + for (idx, static_node) in index.iter_mut().zip(&self.static_nodes) { + if static_node.is_some() { + *idx = 0; + } + } + new_array[index.as_slice()] += value; + } + + self.array = new_array; + + for (static_node, interp) in self.static_nodes.iter_mut().zip(&mut self.interps) { + if let &mut Some(value) = static_node { + *interp = Interp::new( + value, + value, + 1, + 0, + interp.reweight_meth(), + interp.map(), + interp.interp_meth(), + ); + } + } } } @@ -175,14 +203,6 @@ mod tests { assert!(!subgrid.is_empty()); assert_eq!(subgrid.indexed_iter().count(), 4 * 4 * 4); - assert_eq!( - subgrid.static_scale(), - Some(Mu2 { - ren: 1000.0, - fac: 1000.0, - frg: -1.0 - }) - ); assert_eq!( subgrid.stats(), Stats { @@ -198,7 +218,6 @@ mod tests { assert!(!subgrid.is_empty()); assert_eq!(subgrid.indexed_iter().count(), 2 * 4 * 4 * 4); - assert_eq!(subgrid.static_scale(), None); assert_eq!( subgrid.stats(), Stats { diff --git a/pineappl/src/interpolation.rs b/pineappl/src/interpolation.rs index 367747e7..13ea92de 100644 --- a/pineappl/src/interpolation.rs +++ b/pineappl/src/interpolation.rs @@ -207,9 +207,13 @@ impl Interp { /// TODO #[must_use] pub fn node_values(&self) -> Vec { - (0..self.nodes) - .map(|node| self.map_y_to_x(self.gety(node))) - .collect() + if self.nodes == 1 { + vec![self.map_y_to_x(self.min)] + } else { + (0..self.nodes) + .map(|node| self.map_y_to_x(self.gety(node))) + .collect() + } } fn map_y_to_x(&self, y: f64) -> f64 { @@ -249,6 +253,18 @@ impl Interp { pub const fn map(&self) -> Map { self.map } + + /// TODO + #[must_use] + pub const fn interp_meth(&self) -> InterpMeth { + self.interp_meth + } + + /// TODO + #[must_use] + pub const fn reweight_meth(&self) -> ReweightMeth { + self.reweight + } } /// TODO @@ -682,7 +698,6 @@ mod tests { assert_eq!(node_values.len(), 1); - // TODO: the return value is not the one expected (90^2), because `deltay` is zero - assert!(node_values[0].is_nan()); + assert_approx_eq!(f64, node_values[0], 90.0 * 90.0, ulps = 16); } } diff --git a/pineappl/src/subgrid.rs b/pineappl/src/subgrid.rs index 027bbd8d..fc2165cb 100644 --- a/pineappl/src/subgrid.rs +++ b/pineappl/src/subgrid.rs @@ -85,8 +85,8 @@ pub trait Subgrid { /// Return statistics for this subgrid. fn stats(&self) -> Stats; - /// Return the static (single) scale, if this subgrid has one. - fn static_scale(&self) -> Option; + /// TODO + fn optimize_static_nodes(&mut self); } /// Type to iterate over the non-zero contents of a subgrid. The tuple contains the indices of the