Support momentum-flavour manipulation of grids?

[felixhekhorn]

This issue is similar to #346, but not the same.

Physics case

We are interested in dijet cross sections in photo-production at HERA, which have been measures here and the corresponding theory is available here to which all equations refer henceforth. Photo-production cross sections can be split into a "direct" contribution (called "point" in the paper) and a "resolved" contribution (called "hadr" in the paper), see Eq. 2.4 there. Combining Eq. 2.2 and Eq. 2.6 we can write the resolved part conceptually as

$$\sigma_{ep} = f_{\gamma/e} \otimes f_{i/\gamma} \otimes \hat \sigma_{ij} \otimes f_{j/p}$$

where $f_{a/b}$ is the probability to find a in b. While $f_{j/p}$ is well known 🙃 , $f_{i/\gamma}$ is not and we are hunting for it.

The way the underlying program works (if I understand correctly) is to first pre-compute

$$f_{i/\gamma/e} = f_{\gamma/e} \otimes f_{i/\gamma} $$

and then convolute that with the partonic matrix elements.

Question

can I use pineappl to directly produce a grid, which takes $f_{i/\gamma}$ and $f_{j/p}$ and yields $\sigma_{ep}$?

For that one would need to support the manipulation of flavour+momentum, in the following sense:

1. one would grid the underlying partonic matrix elements (which just collide two QCD partons - as in other colliders), call it $\Sigma_{kj}$ where $k$ and $j$ are multi-indices denoting both momentum and flavour
2. one would grid the additional convolution, call it $F_{ik}$ (with also $i$ multi-index) (F for photon flux)
3. one would glue them together $F\Sigma_{ij} = F_{ik} \Sigma_{kj}$

In this specific case $F_{ik}$ only mixes momentum fractions (it's an integral), but is diagonal in flavour.

Algorithmically, $F$ behaves similar to an EKO, i.e. it has two indices, but in contrast to EKO it does not collapse a convolution slot, but changes its type in a sense (specifically it's parent particle).

Alternatives

If we are not going for the full-fledged solution straight away: could I do at least some subgrid manipulation? i.e. we compute $\Sigma_{kj}$, then I turn the subgrids somehow into ImportOnlySubgridVX (or something) such that they are plain arrays of weights (and not Lagrange weight or stuff), then do a np.einsum on those and write them back?

As a poor man solution, I can of course

1. for grids: do the additional multiplication/convolution outside pineappl and directly feed the grid the correct thing
2. for FK tables: add the convolution to the EKO as I know that EKOs are nothing but numpy arrays

[Radonirinaunimi]

At some point, we discussed about the ability to collapse a convolution along a given direction and returns a grid with lesser convolutions, but no one has had the time to look into it yet. For instance, if one has a grid which requires two convolutions, a PDF and a FF, one can convolve the FF part and get a grid with the PDF dependence (which can be convolved with a PDF set in the standard way to get the cross-sections).

I believe this is a generalization of what you need, right?

[felixhekhorn]

At some point, we discussed about the ability to collapse a convolution along a given direction and returns a grid with lesser convolutions,

that's #346

I believe this is a generalization of what you need, right?

mmm, not really - it is related, but different. The above would remove an (multi-)index, this would transform it. In a sense #346 is a vector operation, but this is a matrix operation.

Let me ask more explicitly about the short term solution:

• Is there an operation (e.g. optimize), which will cast any subgrid to a dense numpy grid with plain weights, regardless of how the grid was generated (e.g. in MC generators)?
• can I then iterate on each subgrid and generate a new grid replacing the subgrids one by one with a modified version?

[cschwan]

• Is there an operation (e.g. optimize), which will cast any subgrid to a dense numpy grid with plain weights, regardless of how the grid was generated (e.g. in MC generators)?

In Rust you can use Grid::subgrids to get a reference to all subgrids, and then you can iterate over each non-zero entry of a subgrid with Subgrid::indexed_iter.

• can I then iterate on each subgrid and generate a new grid replacing the subgrids one by one with a modified version?

Yes, using Grid::set_subgrid, but that will only make sense if the convolution and all other properties of the grid don't change.

[Radonirinaunimi]

* Is there an operation (e.g. `optimize`), which will cast any subgrid to a dense numpy grid with plain weights, regardless of how the grid was generated (e.g. in MC generators)?

Yes, this can done with all of the APIs except Fortran. With Python in particular, you can use grid.subgrid() which will return a SubgridEnum object from which you can call .to_array() to return a numpy array.

* can I then iterate on each subgrid and generate a new grid replacing the subgrids one by one with a modified version?

Yes, you could look at:

pineappl/pineappl_py/tests/test_grid.py

Line 177 in e442048

def test_set_subgrid(self, fake_grids):

But as @cschwan said, it makes sense if the properties don't change (which if needed could be faked).