Skip to content

Bugfix and refactor Grid::evolve #354

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
cschwan merged 17 commits into from
Jul 4, 2025
Merged

Bugfix and refactor Grid::evolve #354

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
17 commits
Select commit Hold shift + click to select a range
1403d79
Vary fragmentation scale with the correct factor
cschwan Jun 26, 2025
a22caf6
Remove `op_fac1` and `op_frg1` only used in error messages
cschwan Jun 26, 2025
e40fdfe
Shift some calls closer together
cschwan Jun 26, 2025
b00e368
Put grid scales into an array
cschwan Jun 26, 2025
3bd08d7
Put operator scales into an array
cschwan Jun 26, 2025
f47e31f
Put initial scales into an array
cschwan Jun 26, 2025
6069e1e
Put yet another pair of scales into an array
cschwan Jun 26, 2025
540fdf2
Rename `used_op_scales1` to `op_scales1`
cschwan Jun 26, 2025
4934248
Merge two code paths
cschwan Jun 26, 2025
45d9261
Merge two more code paths
cschwan Jun 26, 2025
b2bbab8
Prepare for merging yet another code path
cschwan Jun 26, 2025
a842aa2
Fix a useless check - this didn't do what it was supposed to
cschwan Jun 26, 2025
5bb428f
Merge two more code paths
cschwan Jun 26, 2025
227124c
Merge more code
cschwan Jun 26, 2025
d8edb0d
Fix clippy warning
cschwan Jun 26, 2025
fdc3df5
Fix failing test
cschwan Jun 26, 2025
b2eda63
Fix another failing test
cschwan Jun 27, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
214 changes: 71 additions & 143 deletions pineappl/src/grid.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1107,34 +1107,24 @@ impl Grid {
// the EKOs' indices matching this Grid's convolutions
let mut eko_map = Vec::new();

// initial factorization scale
let mut fac0 = None;
// factorization slices we use
let mut used_op_fac1 = Vec::new();
// factorization slices we encounter, but possibly don't use
let mut op_fac1 = Vec::new();

// initial fragmentation scale
let mut frg0 = None;
// fragmentation slices we use
let mut used_op_frg1 = Vec::new();
// fragmentation slices we encounter, but possibly don't use
let mut op_frg1 = Vec::new();

// factorization scales needed by the grid
let grid_fac1: Vec<_> = self
.evolve_info(order_mask)
.fac1
.into_iter()
.map(|fac| xi.1 * xi.1 * fac)
.collect();
// fragmentation scales needed by the grid
let grid_frg1: Vec<_> = self
.evolve_info(order_mask)
.frg1
.into_iter()
.map(|frg| xi.1 * xi.1 * frg)
.collect();
// initial factorization and fragmentation scales
let mut scales0 = [None, None];
// factorization and fragmentation slices we use
let mut op_scales1 = [Vec::new(), Vec::new()];
// names for the scales in the same ordering
let names = ["fac", "frg"];

let EvolveInfo {
fac1: grid_fac1,
frg1: grid_frg1,
..
} = self.evolve_info(order_mask);
let grid_scales1: [Vec<_>; 2] = [
// factorization scales needed by the grid
grid_fac1.into_iter().map(|fac| xi.1 * xi.1 * fac).collect(),
// fragmentation scales needed by the grid
grid_frg1.into_iter().map(|frg| xi.2 * xi.2 * frg).collect(),
];

for slice in zip_n(slices) {
let (infos, operators): (Vec<_>, Vec<_>) = slice
Expand All @@ -1150,9 +1140,9 @@ impl Grid {
));
}

let mut fac1 = None;
let mut frg1 = None;
let mut scales1 = [None, None];

// check that all operators are compatible with each other
for (info, operator) in infos.iter().zip(&operators) {
let dim_op_info = (
info.pids1.len(),
Expand All @@ -1168,49 +1158,24 @@ impl Grid {
)));
}

if info.conv_type.is_pdf() {
if let Some(fac0) = fac0 {
// check that this EKO slice is compatible with all previous slices
if !approx_eq!(f64, fac0, info.fac0, ulps = 8) {
let idx = usize::from(!info.conv_type.is_pdf());

// check that both initial- and process-level scales agree among all operators of
// the same convolution type
for (scale, op_scale) in [
(&mut scales0[idx], info.fac0),
(&mut scales1[idx], info.fac1),
] {
if let &mut Some(scale) = scale {
// check that the initial scale of all EKOs in this slice agree with each other
if !approx_eq!(f64, scale, op_scale, ulps = 8) {
return Err(Error::General(format!(
"EKO slice's fac0 = '{}' is incompatible with previous slices' fac0 = '{fac0}'",
info.fac0
"EKO slice's {0}{idx} = '{op_scale}' is incompatible with previous slices' {0}{idx} = '{scale}'",
names[idx],
)));
}
} else {
fac0 = Some(info.fac0);
}

if let Some(fac1) = fac1 {
// we assume that all EKO slices always share the same factorization and/or
// fragmentation scale at process level. If this isn't the case, for
// instance when the fragmentation scale is functionally different from the
// factorization scale, this implementation isn't general enough and has to
// be changed
if !subgrid::node_value_eq(info.fac1, fac1) {
unimplemented!();
}
} else {
fac1 = Some(info.fac1);
}
} else {
if let Some(frg0) = frg0 {
if !approx_eq!(f64, frg0, info.fac0, ulps = 8) {
return Err(Error::General(format!(
"EKO slice's frg0 = '{}' is incompatible with previous slices' frg0 = '{frg0}'",
info.fac0
)));
}
} else {
frg0 = Some(info.fac0);
}

if let Some(frg1) = frg1 {
if !subgrid::node_value_eq(info.fac1, frg1) {
unimplemented!();
}
} else {
frg1 = Some(info.fac1);
*scale = Some(op_scale);
}
}
}
Expand All @@ -1235,58 +1200,40 @@ impl Grid {
.collect::<Result<_>>()?;
}

if let Some(fac1) = fac1 {
op_fac1.push(fac1);

// it's possible that due to small numerical differences we get two slices which
// are almost the same. We have to skip those in order not to evolve the 'same'
// slice twice
if used_op_fac1
.iter()
.any(|&fac| subgrid::node_value_eq(fac, fac1))
{
continue;
}

// skip slices that the grid doesn't use
if !grid_fac1
.iter()
.any(|&fac| subgrid::node_value_eq(fac, fac1))
{
continue;
}
}

if let Some(frg1) = frg1 {
op_frg1.push(frg1);
for ((&scale1, op_scales1), grid_scales1) in
scales1.iter().zip(&mut op_scales1).zip(&grid_scales1)
{
if let Some(scale1) = scale1 {
// it's possible that due to small numerical differences we get two slices which
// are almost the same. We have to skip those in order not to evolve the 'same'
// slice twice
if op_scales1
.iter()
.any(|&s| subgrid::node_value_eq(s, scale1))
{
continue;
}

// it's possible that due to small numerical differences we get two slices which
// are almost the same. We have to skip those in order not to evolve the 'same'
// slice twice
if used_op_frg1
.iter()
.any(|&frg| subgrid::node_value_eq(frg, frg1))
{
continue;
}
// skip slices that the grid doesn't use
if !grid_scales1
.iter()
.any(|&s| subgrid::node_value_eq(s, scale1))
{
continue;
}

// skip slices that the grid doesn't use
if !grid_frg1
.iter()
.any(|&frg| subgrid::node_value_eq(frg, frg1))
{
continue;
op_scales1.push(scale1);
}
}

let operators: Vec<_> = eko_map.iter().map(|&idx| operators[idx].view()).collect();
let infos: Vec<_> = eko_map.iter().map(|&idx| infos[idx].clone()).collect();

let (fac, frg, scale_values) = match (fac0, frg0) {
(None, None) => unreachable!(),
(Some(fac0), None) => (ScaleFuncForm::Scale(0), ScaleFuncForm::NoScale, vec![fac0]),
(None, Some(frg0)) => (ScaleFuncForm::NoScale, ScaleFuncForm::Scale(0), vec![frg0]),
(Some(fac0), Some(frg0)) => {
let (fac, frg, scale_values) = match scales0 {
[None, None] => unreachable!(),
[Some(fac0), None] => (ScaleFuncForm::Scale(0), ScaleFuncForm::NoScale, vec![fac0]),
[None, Some(frg0)] => (ScaleFuncForm::NoScale, ScaleFuncForm::Scale(0), vec![frg0]),
[Some(fac0), Some(frg0)] => {
if approx_eq!(f64, fac0, frg0, ulps = 8) {
(ScaleFuncForm::Scale(0), ScaleFuncForm::Scale(0), vec![fac0])
} else {
Expand Down Expand Up @@ -1344,41 +1291,22 @@ impl Grid {
} else {
result = Some(evolved_slice);
}

if let Some(fac1) = fac1 {
used_op_fac1.push(fac1);
}

if let Some(frg1) = frg1 {
used_op_frg1.push(frg1);
}
}

op_fac1.sort_by(f64::total_cmp);
op_frg1.sort_by(f64::total_cmp);

// make sure we've evolved all slices
if let Some(fac1) = grid_fac1.into_iter().find(|&grid_fac1| {
!used_op_fac1
.iter()
.any(|&eko_fac1| subgrid::node_value_eq(grid_fac1, eko_fac1))
}) {
return Err(Error::General(format!(
"no operator for fac1 = {fac1} found in {op_fac1:?}"
)));
}

// make sure we've evolved all slices
if let Some(frg1) = grid_frg1.into_iter().find(|&grid_frg1| {
!used_op_frg1
// make sure we've evolved all slices in the grid
for ((grid_scales1, op_scales1), name) in grid_scales1.iter().zip(&op_scales1).zip(&names) {
if let Some(scale1) = grid_scales1
.iter()
.any(|&eko_frg1| subgrid::node_value_eq(grid_frg1, eko_frg1))
}) {
return Err(Error::General(format!(
"no operator for frg1 = {frg1} found in {op_frg1:?}"
)));
.find(|&&a| !op_scales1.iter().any(|&b| subgrid::node_value_eq(a, b)))
{
return Err(Error::General(format!(
"no operator for {name}1 = '{scale1}' found"
)));
}
}

// if we didn't actually perform an evolution an empty grid isn't guaranteed to be an
// FK-table
let result =
result.ok_or_else(|| Error::General("no evolution was performed".to_owned()))?;

Expand Down
2 changes: 1 addition & 1 deletion pineappl_cli/tests/evolve.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -103,7 +103,7 @@ const LHCB_WP_7TEV_V2_XIF_2_STR: &str = "b Grid FkTable
";

const LHCB_WP_7TEV_V2_XIF_2_ERROR_STR: &str =
"Error: no operator for fac1 = 25825.775616000003 found in [6456.443904000001]
"Error: no operator for fac1 = '25825.775616000003' found
";

const LHCB_WP_7TEV_OPTIMIZED_STR: &str = "b etal dsig/detal
Expand Down