diff --git a/pineappl/src/grid.rs b/pineappl/src/grid.rs
index d754d521..3784e1d5 100644
--- a/pineappl/src/grid.rs
+++ b/pineappl/src/grid.rs
@@ -1429,6 +1429,50 @@ impl Grid {
             "[{elapsed_precise}] {bar:50.cyan/blue} {pos:>7}/{len:7} - ETA: {eta_precise} {msg}",
         ));
 
+        // which (tgt_pid, src_pid) tuples are non-zero in general?
+        let non_zero_pid_indices: Vec<_> = (0..operator.dim().0)
+            .cartesian_product(0..operator.dim().1)
+            .filter(|&(tgt_pid_idx, src_pid_idx)| {
+                operator
+                    .slice(s![tgt_pid_idx, src_pid_idx, .., .., ..])
+                    .iter()
+                    .any(|&value| value != 0.0)
+            })
+            .collect();
+
+        // the x1 and x2 grid values are the same across all non-zero subgrids
+        let first_non_zero_grid = self
+            .subgrids
+            .iter()
+            .filter(|subgrid| !subgrid.is_empty())
+            .next()
+            .unwrap_or_else(|| unreachable!());
+        // source x1/x2 grid might differ and be differently sorted than the operator
+        let x1_grid: Vec<_> = first_non_zero_grid
+            .x1_grid()
+            .iter()
+            .map(|x| {
+                eko_info
+                    .grid_axes
+                    .x_grid
+                    .iter()
+                    .position(|xi| approx_eq!(f64, *xi, *x, ulps = 64))
+                    .unwrap_or_else(|| unreachable!())
+            })
+            .collect();
+        let x2_grid: Vec<_> = first_non_zero_grid
+            .x2_grid()
+            .iter()
+            .map(|x| {
+                eko_info
+                    .grid_axes
+                    .x_grid
+                    .iter()
+                    .position(|xi| approx_eq!(f64, *xi, *x, ulps = 64))
+                    .unwrap_or_else(|| unreachable!())
+            })
+            .collect();
+
         // iterate over all bins, which are mapped one-to-one from the target to the source grid
         for bin in 0..self.bin_info().bins() {
             // iterate over the source grid luminosities
@@ -1477,32 +1521,6 @@ impl Grid {
 
                     let src_subgrid = &self.subgrids[[order, bin, src_lumi]];
 
-                    // source x1/x2 grid might differ and be differently sorted than the operator
-                    let x1_grid: Vec<_> = src_subgrid
-                        .x1_grid()
-                        .iter()
-                        .map(|x| {
-                            eko_info
-                                .grid_axes
-                                .x_grid
-                                .iter()
-                                .position(|xi| approx_eq!(f64, *xi, *x, ulps = 64))
-                                .unwrap_or_else(|| unreachable!())
-                        })
-                        .collect();
-                    let x2_grid: Vec<_> = src_subgrid
-                        .x2_grid()
-                        .iter()
-                        .map(|x| {
-                            eko_info
-                                .grid_axes
-                                .x_grid
-                                .iter()
-                                .position(|xi| approx_eq!(f64, *xi, *x, ulps = 64))
-                                .unwrap_or_else(|| unreachable!())
-                        })
-                        .collect();
-
                     for ((iq2, ix1, ix2), &value) in src_subgrid.iter() {
                         let scale = src_subgrid.mu2_grid()[iq2].fac;
                         let src_iq2 = src_array_q2_grid
@@ -1595,6 +1613,20 @@ impl Grid {
                             0
                         };
 
+                        // if `op1` and `op2` below are zero there's no work to do
+                        // TODO: ideally we change the for loops instead of vetoing here
+                        if (has_pdf1
+                            && !non_zero_pid_indices
+                                .iter()
+                                .any(|&tuple| tuple == (tgt_pid1_idx, src_pid1_idx)))
+                            || (has_pdf2
+                                && !non_zero_pid_indices
+                                    .iter()
+                                    .any(|&tuple| tuple == (tgt_pid2_idx, src_pid2_idx)))
+                        {
+                            continue;
+                        }
+
                         // create target subgrid
                         let mut tgt_array =
                             SparseArray3::new(1, tgt_x1_grid.len(), tgt_x2_grid.len());