Merge pull request #394 from NNPDF/rust_sv

Recover SV, ev method, backwards method

.github/workflows/lha_bot_rust.yml

@@ -17,17 +17,27 @@ jobs:
   lhabench:
     name: LHA paper Benchmarks
     runs-on: ubuntu-latest
-    container:
-      image: ghcr.io/nnpdf/bench-evol:v2
-      credentials:
-        username: ${{ github.repository_owner }}
-        password: ${{ secrets.GITHUB_TOKEN }}
+    # container:
+    #   image: ghcr.io/nnpdf/bench-evol:v2
+    #   credentials:
+    #     username: ${{ github.repository_owner }}
+    #     password: ${{ secrets.GITHUB_TOKEN }}
+
+    strategy:
+      matrix:
+        python-version: ["3.12"]
+      fail-fast: false
 
     steps:
       - uses: actions/checkout@v2
+        # with:
+        #   # tags needed for dynamic versioning
+        #   fetch-depth: 0
+      - name: Set up Python ${{ matrix.python-version }} 🐍
+        id: setup-python
+        uses: actions/setup-python@v5
         with:
-          # tags needed for dynamic versioning
-          fetch-depth: 0
+          python-version: ${{ matrix.python-version }}
       - name: Install and configure Poetry
         uses: snok/install-poetry@v1
         with:
@@ -41,6 +51,7 @@ jobs:
       - name: Run benchmark
         run: |
           ./rustify.sh
+          poe compile
           poe lha -m "nnlo and sv"
           # TODO wait for polarized to reactivate
           # poe lha -m "ffns_pol and sv"

src/eko/evolution_operator/__init__.py.patch

@@ -586,7 +586,7 @@ index bd1b19d6..f543f7bc 100644
      def initialize_op_members(self):
          """Init all operators with the identity or zeros."""
          eye = OpMember(
-@@ -857,10 +318,13 @@ class Operator(sv.ScaleVariationModeMixin):
+@@ -857,10 +318,14 @@ class Operator(sv.ScaleVariationModeMixin):
              else:
                  self.op_members[n] = zero.copy()
 
@@ -599,12 +599,13 @@ index bd1b19d6..f543f7bc 100644
 +        cfg.as1 = self.as_list[1]
 +        cfg.as0 = self.as_list[0]
 +        cfg.py = ekors.ffi.cast("void *", cb_quad_ker_qcd.address)
++        cfg.method_num = self.ev_method
 +
 +    def run_op_integration(self, log_grid):
          """Run the integration for each grid point.
 
          Parameters
-@@ -875,18 +333,54 @@ class Operator(sv.ScaleVariationModeMixin):
+@@ -875,18 +339,53 @@ class Operator(sv.ScaleVariationModeMixin):
          """
          column = []
          k, logx = log_grid
@@ -619,11 +620,10 @@ index bd1b19d6..f543f7bc 100644
 +        cfg.nf = self.nf
 +        cfg.is_log = self.int_disp.log
 +        cfg.logx = logx
-+        cfg.L = np.log(self.xif2)
-+        cfg.method_num = 1
++        cfg.Lsv = np.log(self.xif2)
 +        cfg.ev_op_iterations = self.config["ev_op_iterations"]
 +        cfg.ev_op_max_order_qcd = self.config["ev_op_max_order"][0]
-+        cfg.sv_mode_num = 1
++        cfg.sv_mode_num = self.sv_mode
 +        cfg.is_threshold = self.is_threshold
 +        self.update_cfg(cfg)
 +

src/eko/evolution_operator/operator_matrix_element.py.patch

@@ -143,7 +143,7 @@ index 098b52db..b56e27d1 100644
  class OperatorMatrixElement(Operator):
      r"""Internal representation of a single |OME|.
 
-@@ -290,41 +186,14 @@ class OperatorMatrixElement(Operator):
+@@ -290,41 +186,15 @@ class OperatorMatrixElement(Operator):
              )
          return labels
 
@@ -190,6 +190,7 @@ index 098b52db..b56e27d1 100644
 +        cfg.as1 = self.a_s
 +        cfg.as0 = 0.0
 +        cfg.Lsv = np.log(self.xif2)
++        cfg.method_num = self.backward_method
 
      @property
      def a_s(self):

src/eko/evolution_operator/quad_ker.py

@@ -90,15 +90,15 @@ def cb_quad_ker_qcd(
     areas_raw,
     areas_x,
     areas_y,
-    L,
-    _method_num,
+    _L,
+    ev_method,
     as1,
     as0,
     ev_op_iterations,
     ev_op_max_order_qcd,
-    _sv_mode_num,
+    sv_mode,
     is_threshold,
-    _Lsv,
+    Lsv,
 ):
     """C Callback inside integration kernel."""
     # recover complex variables
@@ -107,9 +107,6 @@ def cb_quad_ker_qcd(
     # combute basis functions
     areas = nb.carray(areas_raw, (areas_x, areas_y))
     pj = interpolation.evaluate_grid(n, is_log, logx, areas)
-    # TODO recover parameters
-    method = "iterate-exact"
-    sv_mode = sv.Modes.exponentiated
     order = (order_qcd, 0)
     ev_op_max_order = (ev_op_max_order_qcd, 0)
     if is_singlet:
@@ -119,12 +116,12 @@ def cb_quad_ker_qcd(
         gamma_singlet = re_gamma_singlet + im_gamma_singlet * 1j
         if sv_mode == sv.Modes.exponentiated:
             gamma_singlet = sv_exponentiated.gamma_variation(
-                gamma_singlet, order, nf, L
+                gamma_singlet, order, nf, Lsv
             )
         # construct eko
         ker = s.dispatcher(
             order,
-            method,
+            ev_method,
             gamma_singlet,
             as1,
             as0,
@@ -135,7 +132,7 @@ def cb_quad_ker_qcd(
         # scale var expanded is applied on the kernel
         if sv_mode == sv.Modes.expanded and not is_threshold:
             ker = np.ascontiguousarray(
-                sv_expanded.singlet_variation(gamma_singlet, as1, order, nf, L, dim=2)
+                sv_expanded.singlet_variation(gamma_singlet, as1, order, nf, Lsv, dim=2)
             ) @ np.ascontiguousarray(ker)
         ker = select_singlet_element(ker, mode0, mode1)
     else:
@@ -144,19 +141,19 @@ def cb_quad_ker_qcd(
         im_gamma_ns = nb.carray(im_gamma_raw, order_qcd)
         gamma_ns = re_gamma_ns + im_gamma_ns * 1j
         if sv_mode == sv.Modes.exponentiated:
-            gamma_ns = sv_exponentiated.gamma_variation(gamma_ns, order, nf, L)
+            gamma_ns = sv_exponentiated.gamma_variation(gamma_ns, order, nf, Lsv)
         # construct eko
         ker = ns.dispatcher(
             order,
-            method,
+            ev_method,
             gamma_ns,
             as1,
             as0,
             nf,
             ev_op_iterations,
         )
         if sv_mode == sv.Modes.expanded and not is_threshold:
-            ker = sv_expanded.non_singlet_variation(gamma_ns, as1, order, nf, L) * ker
+            ker = sv_expanded.non_singlet_variation(gamma_ns, as1, order, nf, Lsv) * ker
     # recombine everything
     res = ker * pj * jac
     return np.real(res)
@@ -238,12 +235,12 @@ def cb_quad_ker_ome(
     areas_x,
     areas_y,
     L,
-    _method_num,
+    backward_method,
     as1,
     _as0,
     _ev_op_iterations,
     _ev_op_max_order_qcd,
-    _sv_mode_num,
+    sv_mode,
     _is_threshold,
     Lsv,
 ):
@@ -254,8 +251,6 @@ def cb_quad_ker_ome(
     # compute basis functions
     areas = nb.carray(areas_raw, (areas_x, areas_y))
     pj = interpolation.evaluate_grid(n, is_log, logx, areas)
-    # TODO recover parameters
-    sv_mode = sv.Modes.exponentiated
     order = (order_qcd, 0)
     if is_singlet:
         indices = {21: 0, 100: 1, 90: 2}
@@ -274,9 +269,6 @@ def cb_quad_ker_ome(
     if sv_mode == sv.Modes.exponentiated:
         A = sv.exponentiated.gamma_variation(A, order, nf, Lsv)
 
-    # TODO recover InversionMethod
-    backward_method = ""
-
     # build the expansion in alpha_s depending on the strategy
     ker = build_ome(A, order, as1, backward_method)