Merge pull request #61 from N3PDF/vectorized_example

Vectorized example C-API

.gitignore

@@ -134,3 +134,5 @@ benchmarks/*.png
 # cffi stuff
 *.pc
 cpdfflow.c
+# example executables
+example

capi/README.md

@@ -7,11 +7,9 @@ This repository contains a C library to access PDFFlow from programming language
 
 Make sure you have installed the `pdfflow` module in Python, then in order to install proceed with the usual cmake steps:
 ```bash
-mkdir build
-cd build
-cmake .. -DCMAKE_INSTALL_PREFIX=<your install prefix>
-make
-make install
+cmake -S . -B build -DCMAKE_INSTALL_PREFIX=<your install prefix>
+cmake --build build
+cmake --install build
 ```
 
 ## Usage
@@ -23,4 +21,11 @@ pkg-config pdfflow --cflags
 pkg-config pdfflow --libs
 ```
 
+If you installed to a non-standard location, you need to set up the `PKG_CONFIG_PATH` and `LD_LIBRARY_PATH`, e.g.:
+```bash
+export PKG_CONFIG_PATH=${VIRTUAL_ENV}/lib/pkgconfig/:${PKG_CONFIG_PATH}:
+export LD_LIBRARY_PATH=${VIRTUAL_ENV}/lib/:${LD_LIBRARY_PATH}:
+```
+
+
 Sample programs using this library are provided in the `capi/examples/` directory.

capi/examples/example.c

@@ -4,13 +4,21 @@
 
 int main() {
     // load pdf
-    mkpdf("NNPDF31_nlo_as_0118/0", "/usr/share/lhapdf/LHAPDF/");
+    mkpdf("NNPDF40_nnlo_as_01180/0", "/usr/share/lhapdf/LHAPDF/");
 
     // test xfxq2 and alphasq2
-    const double x = 0.1, q2 = 1.65;
-    for (int fl=-5; fl <=5; fl++)
-        printf("flv=%d - xfx = %f\n", fl, xfxq2(fl, x, q2));
-    printf("alphas(q2=%f) = %f\n", q2, alphasq2(q2));
+    int pid[] = {-5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5};
+    double xs[] = {0.1};
+    double q2s[] = {1.65};
+    double *xf_vectorized = xfxq2(pid, 11, xs, 1, q2s, 1);
+    for (int fl = 0; fl < 11; fl++)
+        printf("flv=%d - xfx = %f\n", fl-5, xf_vectorized[fl]);
+
+    double q2_vectorized[] = {1.65, 10.65};
+    double* as_vectorized = alphasq2(q2_vectorized, 2);
+
+    printf("alphas(q2=%f) = %f\n", q2_vectorized[0], as_vectorized[0]);
+    printf("alphas(q2=%f) = %f\n", q2_vectorized[1], as_vectorized[1]);
 
     return 0;
 }

capi/examples/fortran/example.f90

@@ -6,27 +6,35 @@ program example
 
    integer, parameter :: dp = kind(1.d0)
 
-   integer :: pid
-   real(dp) :: alpha_s, q2, x, xfx
+   integer :: pid(0:10), i
+   real(dp) :: alpha_s, q2(0:1), x(0:1), xfx(0:11)
+   real(dp) :: q2_vectorized(0:2), as_vectorized(0:2)
 
-   character(kind=c_char, len=21) :: ss = "NNPDF31_nlo_as_0118/0"
+   character(kind=c_char, len=23) :: ss = "NNPDF40_nnlo_as_01180/0"
    character(kind=c_char, len=24) :: pp = "/usr/share/lhapdf/LHAPDF/"
 
    call mkPDF(ss, pp)
 
-   q2 = 1.65
-   x = 0.1
-   call alphasq2(q2, alpha_s)
+   do i = 0, 10
+      pid(i) = i - 5
+   enddo
+   x(0) = 0.1
+   q2(0) = 1.65
 
-   write(*, '(A, F7.5)')"The value of alpha_s is: ", alpha_s
+   call xfxq2(pid, 11, x, 1, q2, 1, xfx)
 
-   write(*, fmt=100)"Value of x*f(x) at q2=", q2, " x=",x
-   do pid = -5, 5
-      call xfxq2(pid, x, q2, xfx)
-      write(*, fmt=200)"Flavour: ", pid, " value: ", xfx
+   do i = 0, 10
+      write(*, fmt=200)"Flavour: ", i - 5, " value: ", xfx(i)
    enddo
 
-100 format (A, F6.2, A, F4.2)
+   q2_vectorized(0) = 1.65
+   q2_vectorized(1) = 10.65
+   call alphasq2(q2_vectorized, 2, as_vectorized);
+
+   write(*, fmt=100)"alphas(q2=",q2_vectorized(0),") = ", as_vectorized(0)
+   write(*, fmt=100)"alphas(q2=",q2_vectorized(1),") = ", as_vectorized(1)
+
+100 format (A, F10.7, A, F10.7)
 200 format ("    ", A, I0, A, F10.7)
 
 end program

capi/examples/fortran/pdfflow_f_interface.c

@@ -3,22 +3,26 @@
 #include "pdfflow/pdfflow.h"
 
 // Some notes about the Fortran interface (to be moved to the docs)
-// 
+//
 // There are two caveats to take into account when interfacing with Fortran:
 // 1 - In fortran all arguments are passed by reference, therefore this file
 // does little more than dereference the arguments to pass it through the cffi interface
 // 2 - Fortran assumes all functions to have a `_` at the end. This is not true in C
-// it is possible to change this behaviour with the  -fno-underscoring gcc flag, but since 
+// it is possible to change this behaviour with the  -fno-underscoring gcc flag, but since
 // it is not standard we prefer to provide an interface between the "cnaming" and the "fnaming"
 
 void mkpdf_(const char *fname, const char *dirname) {
    mkpdf(fname, dirname);
 }
 
-void alphasq2_(const double *q2, double *alphas) {
-   *alphas = alphasq2(*q2);
+void alphasq2_(double *q2, const int *n, double *alphas) {
+   double *as = alphasq2(q2, *n);
+   for (int i = 0; i < *n; i++)
+      alphas[i] = as[i];
 }
 
-void xfxq2_(const int *pid, const double *x, const double *q2, double *f) {
-   *f = xfxq2(*pid, *x, *q2);
+void xfxq2_(int *pid, const int *n, double *x, const int *m, double *q2, const int *o, double *f) {
+   double *xf = xfxq2(pid, *n, x, *m, q2, *o);
+   for (int i = 0; i < *n * *m * *o; i++)
+      f[i] = xf[i];
 }

capi/src/pdfflow/pdfflow.h

@@ -4,6 +4,6 @@
 
 extern void mkpdf(const char *fname, const char * dirname);
 
-extern double xfxq2(int pid, double x, double q2);
+extern double *xfxq2(int *pid, int n, double *x, int m, double *q2, int o);
 
-extern double alphasq2(double q2);
+extern double *alphasq2(double *q2, int n);

capi/src/wrapper.py

@@ -1,6 +1,7 @@
 # This file is part of
 from cpdfflow import ffi
 from pdfflow import pflow
+import numpy as np
 
 pdf = None
 
@@ -14,14 +15,20 @@ def mkpdf(fname, dirname):
 
 
 @ffi.def_extern()
-def xfxq2(pid, x, q2):
-    """Returns the xfxQ2 value for a given PID at specific x and q2."""
+def xfxq2(pid, n, x, m, q2, o):
+    """Returns the xfxQ2 value for arrays of PID, x and q2 values."""
     global pdf
-    return pdf.xfxQ2([pid], [x], [q2])
+    pid_numpy = np.frombuffer(ffi.buffer(pid, n*ffi.sizeof('int')), dtype='int32')
+    x_numpy = np.frombuffer(ffi.buffer(x, m*ffi.sizeof('double')), dtype='double')
+    q2_numpy = np.frombuffer(ffi.buffer(q2, o*ffi.sizeof('double')), dtype='double')
+    ret = pdf.xfxQ2(pid_numpy, x_numpy, q2_numpy).numpy()
+    return ffi.cast("double*", ffi.from_buffer(ret))
 
 
 @ffi.def_extern()
-def alphasq2(q2):
-    """Returns the alpha strong coupling at specific q2 value."""
+def alphasq2(q2, n):
+    """Returns the alpha strong coupling at for an array of q2 values."""
     global pdf
-    return pdf.alphasQ2([q2])
+    q2_numpy = np.frombuffer(ffi.buffer(q2, n*ffi.sizeof('double')), dtype='double')
+    ret = pdf.alphasQ2(q2_numpy).numpy()
+    return ffi.cast("double*", ffi.from_buffer(ret))

setup.cfg

@@ -30,7 +30,7 @@ install_requires =
 where = src
 
 [options.extras_require]
-capi = ciffi
+capi = cffi
 docs =
     sphinx_rtd_theme
     recommonmark

src/pdfflow/pflow.py

@@ -108,6 +108,7 @@ def mkPDF(fname, dirname=None):
     ----------
         fname: str
             PDF name and member in the format '<set_name>/<set member number>'
+            If the set number is not given, assume member 0
         dirname: str
             LHAPDF datadir, if None will try to guess from LHAPDF
 
@@ -116,8 +117,12 @@ def mkPDF(fname, dirname=None):
         PDF: pdfflow.PDF
             instantiated member of the PDF class
     """
-    fname_sp, member = fname.split("/")
-    member = int(member)
+    try:
+        fname_sp, member = fname.split("/")
+        member = int(member)
+    except ValueError:
+        fname_sp = fname
+        member = 0
     return mkPDFs(fname_sp, [member], dirname=dirname)