From 23e807fb483caa1d9476cfcc757502a4756358f6 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 12 Sep 2025 16:40:26 +0200 Subject: [PATCH 1/9] Start adding drop-in interface for pure LHAPDF compatibility --- neopdf_capi/src/include/NeoPDF.hpp | 44 +++++++ neopdf_capi/src/lib.rs | 179 +++++++++++++++++++++++++++++ neopdf_capi/tests/Makefile | 7 +- neopdf_fapi/neopdf.f90 | 44 +++++++ 4 files changed, 272 insertions(+), 2 deletions(-) diff --git a/neopdf_capi/src/include/NeoPDF.hpp b/neopdf_capi/src/include/NeoPDF.hpp index 66d0c79b..e3893850 100644 --- a/neopdf_capi/src/include/NeoPDF.hpp +++ b/neopdf_capi/src/include/NeoPDF.hpp @@ -454,3 +454,47 @@ class GridWriter { }; } // namespace neopdf + +// Everything below this line is added for LHAPDF compatibility. + +namespace LHAPDF { + + class PDF { + protected: + std::unique_ptr _neopdf; + + PDF() = default; + explicit PDF(std::unique_ptr&& pdf) : _neopdf(std::move(pdf)) {} + + public: + virtual ~PDF() = default; + + virtual double xfxQ2(int id, double x, double q2) const { + return _neopdf->xfxQ2(id, x, q2); + } + + virtual double alphasQ2(double q2) const { + return _neopdf->alphasQ2(q2); + } + + double xMin() const { return _neopdf->x_min(); } + double xMax() const { return _neopdf->x_max(); } + double q2Min() const { return _neopdf->q2_min(); } + double q2Max() const { return _neopdf->q2_max(); } + }; + + class GridPDF : public PDF { + public: + explicit GridPDF(std::unique_ptr&& pdf) : PDF(std::move(pdf)) {} + }; + + inline PDF* mkPDF(const std::string& name, int member = 0) { + auto neopdf_ptr = std::make_unique(name, member); + return new GridPDF(std::move(neopdf_ptr)); + } + + inline void setVerbosity(int /*verbosity*/) { + // neopdf doesn't have a verbosity setting, so this is a no-op. + } + +} // namespace LHAPDF diff --git a/neopdf_capi/src/lib.rs b/neopdf_capi/src/lib.rs index d4b590a9..b4703ca9 100644 --- a/neopdf_capi/src/lib.rs +++ b/neopdf_capi/src/lib.rs @@ -1073,3 +1073,182 @@ pub unsafe extern "C" fn neopdf_grid_compress( Err(_) => NeopdfResult::ErrorMemoryError, } } + +// LHAPDF C-API drop-in compatibility layer. +/// +/// This global state stores the loaded PDF set and the currently selected +/// member index. It is not thread-safe and mirrors the behavior of legacy +/// LHAPDF interfaces. +struct LhapdfState { + pdf_set: Option>, + member: usize, +} + +static mut LHAPDF_STATE: LhapdfState = LhapdfState { + pdf_set: None, + member: 0, +}; + +/// Sets LHAPDF runtime parameters from a string (no-op). +/// +/// This function is provided for compatibility with the LHAPDF C API and does +/// nothing in this implementation. +/// +/// # Safety +/// +/// The pointer `line` must be a valid, null-terminated C string if non-null. +#[no_mangle] +pub const unsafe extern "C" fn setlhaparm(_line: *const c_char) {} + +/// Fortran name-mangled variant of `setlhaparm` (no-op). +/// +/// This function is provided for compatibility with the LHAPDF Fortran API and +/// does nothing in this implementation. +/// +/// # Safety +/// +/// The pointer `line` must be valid for reads of `len` bytes if non-null. +#[no_mangle] +pub const unsafe extern "C" fn setlhaparm_(_line: *const c_char, _len: isize) {} + +/// Initializes a PDF set by its name/path and loads all members. +/// +/// The loaded set is stored in a global state used by the other LHAPDF-compatible +/// functions, with the current member index reset to 0. +/// +/// # Panics +/// +/// Panics if `name` is not valid UTF-8 or is not a valid, null-terminated C string. +/// +/// # Safety +/// +/// `name` must be a valid, null-terminated C string pointing to a UTF-8 path or set name. +#[no_mangle] +pub unsafe extern "C" fn initpdfsetbyname(name: *const c_char) { + unsafe { + let c_str = CStr::from_ptr(name); + let pdf_name = c_str.to_str().expect("Invalid UTF-8 string"); + let pdfs = PDF::load_pdfs(pdf_name); + LHAPDF_STATE.pdf_set = Some(pdfs); + LHAPDF_STATE.member = 0; + } +} + +/// Fortran name-mangled variant of `initpdfsetbyname`. +/// +/// Reads a fixed-length Fortran character buffer, trims trailing spaces, and loads +/// the corresponding PDF set into the global state. The current member index is +/// reset to 0. +/// +/// # Panics +/// +/// TODO +/// +/// # Safety +/// +/// `name` must be valid for reads of `len` bytes. The buffer may not be +/// null-terminated; trailing spaces are trimmed. +#[no_mangle] +#[allow(clippy::cast_sign_loss)] +pub unsafe extern "C" fn initpdfsetbyname_(name: *const c_char, len: c_int) { + unsafe { + let name_slice = slice::from_raw_parts(name.cast::(), len as usize); + let pdf_name = std::str::from_utf8(name_slice).unwrap().trim_end(); + let pdfs = PDF::load_pdfs(pdf_name); + LHAPDF_STATE.pdf_set = Some(pdfs); + LHAPDF_STATE.member = 0; + } +} + +/// Selects the active member of the currently loaded PDF set. +/// +/// # Safety +/// +/// This function does not perform bounds checking here; subsequent calls that +/// use the member will return early if the index is out of range. +#[no_mangle] +#[allow(clippy::cast_sign_loss)] +pub unsafe extern "C" fn initpdf(member: c_int) { + unsafe { LHAPDF_STATE.member = member as usize }; +} + +/// Fortran name-mangled variant of `initpdf`. +/// +/// # Safety +/// +/// `member` must be a valid pointer to an integer. +#[no_mangle] +#[allow(clippy::cast_sign_loss)] +pub unsafe extern "C" fn initpdf_(member: *const c_int) { + unsafe { LHAPDF_STATE.member = *member as usize }; +} + +/// Evaluates parton distribution functions at given `(x, q)` for the active member. +/// +/// # Safety +/// +/// - `f` must point to writable memory for at least 13 `c_double` values. +/// - Requires that a PDF set has been initialized via `initpdfsetbyname` or its +/// Fortran variant. If no set is loaded or the member index is out of range, +/// the function returns without writing. +#[no_mangle] +pub unsafe extern "C" fn evolvepdf(x: c_double, q: c_double, f: *mut c_double) { + unsafe { + if let Some(pdfs) = &LHAPDF_STATE.pdf_set { + if LHAPDF_STATE.member >= pdfs.len() { + return; + } + let pdf = &pdfs[LHAPDF_STATE.member]; + let q2 = q * q; + + let (available_pids, _) = pdf.pids().clone().into_raw_vec_and_offset(); + + let pids_to_calc = [-6, -5, -4, -3, -2, -1, 21, 1, 2, 3, 4, 5, 6]; + let out_slice = slice::from_raw_parts_mut(f, 13); + + for i in 0..13 { + let pid = pids_to_calc[i]; + if available_pids.contains(&pid) { + out_slice[i] = pdf.xfxq2(pid, &[x, q2]); + } else { + out_slice[i] = 0.0; + } + } + } + } +} + +/// Fortran name-mangled variant of `evolvepdf`. +/// +/// # Safety +/// +/// - `x`, `q` must be valid pointers to `c_double` values. +/// - `f` must point to writable memory for at least 13 `c_double` values. +/// - A PDF set must have been initialized and the member index must be valid or +/// the function will return without writing. +#[no_mangle] +pub unsafe extern "C" fn evolvepdf_(x: *const c_double, q: *const c_double, f: *mut c_double) { + unsafe { + if let Some(pdfs) = &LHAPDF_STATE.pdf_set { + if LHAPDF_STATE.member >= pdfs.len() { + return; + } + let pdf = &pdfs[LHAPDF_STATE.member]; + let q2 = (*q) * (*q); + + let (available_pids, _) = pdf.pids().clone().into_raw_vec_and_offset(); + + let pids_to_calc = [-6, -5, -4, -3, -2, -1, 21, 1, 2, 3, 4, 5, 6]; + let out_slice = slice::from_raw_parts_mut(f, 13); + + for i in 0..13 { + let pid = pids_to_calc[i]; + if available_pids.contains(&pid) { + out_slice[i] = pdf.xfxq2(pid, &[*x, q2]); + } else { + out_slice[i] = 0.0; + } + } + } + } +} diff --git a/neopdf_capi/tests/Makefile b/neopdf_capi/tests/Makefile index 4024f217..f7f0c39e 100644 --- a/neopdf_capi/tests/Makefile +++ b/neopdf_capi/tests/Makefile @@ -5,12 +5,12 @@ NEOPDF_DEPS != pkg-config --cflags --libs neopdf_capi LHAPDF_DEPS != pkg-config --cflags --libs lhapdf MATH_LIBS = -lm -PROGRAMS = check-capi check-oop check-writer check-writer-oop check-xapi check-xwriter +PROGRAMS = check-capi check-oop check-writer check-writer-oop check-xapi check-xwriter check-lhapdf-compat all: $(PROGRAMS) test-examples: $(PROGRAMS) - set -e && for i in $(PROGRAMS); do ./$${i} > output; diff -u $${i}.output output; done; rm output + set -e && for i in $(PROGRAMS); do ./${i} > output; diff -u ${i}.output output; done; rm output check-capi: check-capi.cpp $(CXX) $(CXXFLAGS) $< $(LHAPDF_DEPS) $(NEOPDF_DEPS) -o $@ @@ -30,6 +30,9 @@ check-xapi: check-capi.c check-xwriter: check-writer.c $(CC) $(CXXFLAGS) $< $(NEOPDF_DEPS) $(MATH_LIBS) -o $@ +check-lhapdf-compat: check-lhapdf-compat.cpp + $(CXX) $(CXXFLAGS) $< $(NEOPDF_DEPS) -o $@ + .PHONY: clean clean: diff --git a/neopdf_fapi/neopdf.f90 b/neopdf_fapi/neopdf.f90 index 02d2e786..0e2c4853 100644 --- a/neopdf_fapi/neopdf.f90 +++ b/neopdf_fapi/neopdf.f90 @@ -300,6 +300,29 @@ function c_neopdf_grid_compress(collection, metadata, output_path) bind(c, name= character (c_char) :: output_path(*) integer (c_int) :: c_neopdf_grid_compress end function + + subroutine c_setlhaparm(line, len) bind(c, name="setlhaparm_") + use iso_c_binding + character(kind=c_char) :: line(*) + integer(c_int), value :: len + end subroutine + + subroutine c_initpdfsetbyname(name, len) bind(c, name="initpdfsetbyname_") + use iso_c_binding + character(kind=c_char) :: name(*) + integer(c_int), value :: len + end subroutine + + subroutine c_initpdf(member) bind(c, name="initpdf_") + use iso_c_binding + integer(c_int) :: member + end subroutine + + subroutine c_evolvepdf(x, q, f) bind(c, name="evolvepdf_") + use iso_c_binding + real(c_double) :: x, q + real(c_double) :: f(*) + end subroutine end interface contains @@ -586,4 +609,25 @@ function neopdf_grid_compress(collection, metadata, output_path) result(res) res = c_neopdf_grid_compress(collection%ptr, metadata, output_path // c_null_char) end function + subroutine setlhaparm(line) + character(len=*), intent(in) :: line + call c_setlhaparm(line, int(len(trim(line)), c_int)) + end subroutine + + subroutine initpdfsetbyname(name) + character(len=*), intent(in) :: name + call c_initpdfsetbyname(name, int(len(trim(name)), c_int)) + end subroutine + + subroutine initpdf(member) + integer, intent(in) :: member + call c_initpdf(member) + end subroutine + + subroutine evolvepdf(x, q, f) + real(dp), intent(in) :: x, q + real(dp), intent(out) :: f(-6:6) + call c_evolvepdf(x, q, f) + end subroutine + end module From 4911498ffb365c2e82dadf379b37c557fcd94040 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 12 Sep 2025 16:47:04 +0200 Subject: [PATCH 2/9] Add tests --- neopdf_capi/tests/Makefile | 4 +- .../tests/check-lhapdf-compatibility.cpp | 79 +++++++++++++++++++ 2 files changed, 81 insertions(+), 2 deletions(-) create mode 100644 neopdf_capi/tests/check-lhapdf-compatibility.cpp diff --git a/neopdf_capi/tests/Makefile b/neopdf_capi/tests/Makefile index f7f0c39e..c4a74c71 100644 --- a/neopdf_capi/tests/Makefile +++ b/neopdf_capi/tests/Makefile @@ -5,7 +5,7 @@ NEOPDF_DEPS != pkg-config --cflags --libs neopdf_capi LHAPDF_DEPS != pkg-config --cflags --libs lhapdf MATH_LIBS = -lm -PROGRAMS = check-capi check-oop check-writer check-writer-oop check-xapi check-xwriter check-lhapdf-compat +PROGRAMS = check-capi check-oop check-writer check-writer-oop check-xapi check-xwriter check-lhapdf-compatibility all: $(PROGRAMS) @@ -30,7 +30,7 @@ check-xapi: check-capi.c check-xwriter: check-writer.c $(CC) $(CXXFLAGS) $< $(NEOPDF_DEPS) $(MATH_LIBS) -o $@ -check-lhapdf-compat: check-lhapdf-compat.cpp +check-lhapdf-compatibility: check-lhapdf-compatibility.cpp $(CXX) $(CXXFLAGS) $< $(NEOPDF_DEPS) -o $@ .PHONY: clean diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.cpp b/neopdf_capi/tests/check-lhapdf-compatibility.cpp new file mode 100644 index 00000000..f147cc7f --- /dev/null +++ b/neopdf_capi/tests/check-lhapdf-compatibility.cpp @@ -0,0 +1,79 @@ +#include +#include +#include +#include +#include +#include + +using namespace LHAPDF; + +const double TOLERANCE = 1e-12; + +void test_lhapdf_compatibility_oop() { + std::cout << "=== Test LHAPDF C++ OOP Compatibility Layer ===\n"; + + std::string pdfname = "NNPDF40_nnlo_as_01180"; + + // Use the compatibility layer + setVerbosity(0); + PDF* neo_lha_pdf = mkPDF(pdfname, 0); + + // Use the native NeoPDF C API for comparison + NeoPDFWrapper* neo_pdf = neopdf_pdf_load(pdfname.c_str(), 0); + + double x = 1e-4; + double q2 = 100.0; + int pid = 21; + + double result = neo_lha_pdf->xfxQ2(pid, x, q2); + double expected = neopdf_pdf_xfxq2(neo_pdf, pid, x, q2); + double reldif = std::abs(result - expected) / expected; + + std::cout << "pid=" << pid << ", x=" << x << ", Q2=" << q2 << std::endl; + std::cout << "Native NeoPDF: " << expected << std::endl; + std::cout << "NeoPDF (LHAPDF compat): " << result << std::endl; + std::cout << "Relative difference: " << reldif << std::endl; + + assert(reldif < TOLERANCE); + + delete neo_lha_pdf; + neopdf_pdf_free(neo_pdf); + + std::cout << "LHAPDF C++ OOP Compatibility test passed." << std::endl; +} + +void test_lhapdf_compatibility_c() { + printf("=== Test LHAPDF C Compatibility Layer ===\n"); + + const char* pdfname = "NNPDF40_nnlo_as_01180"; + + // Use the compatibility layer + initpdfsetbyname(pdfname); + initpdf(0); + double x = 0.1; + double q = 10.0; + double f_compat[13]; + evolvepdf(x, q, f_compat); + + // Use the native NeoPDF API for comparison + NeoPDFWrapper* neo_pdf = neopdf_pdf_load(pdfname, 0); + double expected_g = neopdf_pdf_xfxq2(neo_pdf, 21, x, q*q); + neopdf_pdf_free(neo_pdf); + + printf("Comparing gluon PDF at x=0.1, Q=10.0\n"); + printf("Native NeoPDF: %e\n", expected_g); + printf("NeoPDF (LHAPDF compat): %e\n", f_compat[6]); + + double reldiff = fabs(f_compat[6] - expected_g) / expected_g; + printf("Relative difference: %e\n", reldiff); + + assert(reldiff < TOLERANCE); + + printf("LHAPDF C Compatibility test passed.\n"); +} + +int main() { + test_lhapdf_compatibility_oop(); + test_lhapdf_compatibility_c(); + return 0; +} From a0e9d8a5a275f40b7b132aa10e528b0af8b41897 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 12 Sep 2025 20:44:20 +0200 Subject: [PATCH 3/9] Fix error in `NeoPDF.hpp` --- neopdf_capi/src/include/NeoPDF.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/neopdf_capi/src/include/NeoPDF.hpp b/neopdf_capi/src/include/NeoPDF.hpp index e3893850..00dde89a 100644 --- a/neopdf_capi/src/include/NeoPDF.hpp +++ b/neopdf_capi/src/include/NeoPDF.hpp @@ -489,7 +489,7 @@ namespace LHAPDF { }; inline PDF* mkPDF(const std::string& name, int member = 0) { - auto neopdf_ptr = std::make_unique(name, member); + std::unique_ptr neopdf_ptr(new neopdf::NeoPDF(name, member)); return new GridPDF(std::move(neopdf_ptr)); } From f112e5d0f95253d60cc5911af535926da8ae3ce5 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Sat, 13 Sep 2025 17:45:04 +0200 Subject: [PATCH 4/9] Fix errors & clean up a bit --- neopdf_capi/src/include/NeoPDF.hpp | 25 ++++++++----------- neopdf_capi/tests/Makefile | 2 +- .../tests/check-lhapdf-compatibility.cpp | 3 ++- 3 files changed, 14 insertions(+), 16 deletions(-) diff --git a/neopdf_capi/src/include/NeoPDF.hpp b/neopdf_capi/src/include/NeoPDF.hpp index 00dde89a..a3de2bc9 100644 --- a/neopdf_capi/src/include/NeoPDF.hpp +++ b/neopdf_capi/src/include/NeoPDF.hpp @@ -455,11 +455,10 @@ class GridWriter { } // namespace neopdf -// Everything below this line is added for LHAPDF compatibility. +/** @brief LHAPDF compatibility for no-code migration. */ +namespace NEOLHAPDF { -namespace LHAPDF { - - class PDF { +class PDF { protected: std::unique_ptr _neopdf; @@ -481,20 +480,18 @@ namespace LHAPDF { double xMax() const { return _neopdf->x_max(); } double q2Min() const { return _neopdf->q2_min(); } double q2Max() const { return _neopdf->q2_max(); } - }; +}; - class GridPDF : public PDF { +class GridPDF : public PDF { public: explicit GridPDF(std::unique_ptr&& pdf) : PDF(std::move(pdf)) {} - }; +}; - inline PDF* mkPDF(const std::string& name, int member = 0) { - std::unique_ptr neopdf_ptr(new neopdf::NeoPDF(name, member)); - return new GridPDF(std::move(neopdf_ptr)); - } +inline PDF* mkPDF(const std::string& name, int member = 0) { + std::unique_ptr neopdf_ptr(new neopdf::NeoPDF(name, member)); + return new GridPDF(std::move(neopdf_ptr)); +} - inline void setVerbosity(int /*verbosity*/) { - // neopdf doesn't have a verbosity setting, so this is a no-op. - } +inline void setVerbosity(int /*verbosity*/) { } } // namespace LHAPDF diff --git a/neopdf_capi/tests/Makefile b/neopdf_capi/tests/Makefile index c4a74c71..46eb724f 100644 --- a/neopdf_capi/tests/Makefile +++ b/neopdf_capi/tests/Makefile @@ -10,7 +10,7 @@ PROGRAMS = check-capi check-oop check-writer check-writer-oop check-xapi check-x all: $(PROGRAMS) test-examples: $(PROGRAMS) - set -e && for i in $(PROGRAMS); do ./${i} > output; diff -u ${i}.output output; done; rm output + set -e && for i in $(PROGRAMS); do ./$${i} > output; diff -u $${i}.output output; done; rm output check-capi: check-capi.cpp $(CXX) $(CXXFLAGS) $< $(LHAPDF_DEPS) $(NEOPDF_DEPS) -o $@ diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.cpp b/neopdf_capi/tests/check-lhapdf-compatibility.cpp index f147cc7f..494cf126 100644 --- a/neopdf_capi/tests/check-lhapdf-compatibility.cpp +++ b/neopdf_capi/tests/check-lhapdf-compatibility.cpp @@ -5,7 +5,7 @@ #include #include -using namespace LHAPDF; +using namespace NEOLHAPDF; const double TOLERANCE = 1e-12; @@ -75,5 +75,6 @@ void test_lhapdf_compatibility_c() { int main() { test_lhapdf_compatibility_oop(); test_lhapdf_compatibility_c(); + return 0; } From 70ff6a9ef2ef41e69d84ba9b7f3326d684fdf8ad Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Sat, 13 Sep 2025 17:51:38 +0200 Subject: [PATCH 5/9] Add missing file --- neopdf_capi/tests/check-lhapdf-compatibility.output | 12 ++++++++++++ 1 file changed, 12 insertions(+) create mode 100644 neopdf_capi/tests/check-lhapdf-compatibility.output diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.output b/neopdf_capi/tests/check-lhapdf-compatibility.output new file mode 100644 index 00000000..0137760a --- /dev/null +++ b/neopdf_capi/tests/check-lhapdf-compatibility.output @@ -0,0 +1,12 @@ +=== Test LHAPDF C++ OOP Compatibility Layer === +pid=21, x=0.0001, Q2=100 +Native NeoPDF: 30.8985 +NeoPDF (LHAPDF compat): 30.8985 +Relative difference: 0 +LHAPDF C++ OOP Compatibility test passed. +=== Test LHAPDF C Compatibility Layer === +Comparing gluon PDF at x=0.1, Q=10.0 +Native NeoPDF: 1.132265e+00 +NeoPDF (LHAPDF compat): 1.132265e+00 +Relative difference: 0.000000e+00 +LHAPDF C Compatibility test passed. From 7c207799d0b0076393fe6e9b24fb32d5ff301a3e Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Sat, 13 Sep 2025 21:51:35 +0200 Subject: [PATCH 6/9] Change default in PIDs computed with `evolvepdf` (include photons) --- neopdf_capi/src/lib.rs | 18 ++++++++---------- .../tests/check-lhapdf-compatibility.cpp | 10 +++++----- .../tests/check-lhapdf-compatibility.output | 4 ++-- 3 files changed, 15 insertions(+), 17 deletions(-) diff --git a/neopdf_capi/src/lib.rs b/neopdf_capi/src/lib.rs index b4703ca9..290f46ae 100644 --- a/neopdf_capi/src/lib.rs +++ b/neopdf_capi/src/lib.rs @@ -10,6 +10,8 @@ use neopdf::parser::SubgridData; use neopdf::pdf::PDF; use neopdf::writer::GridArrayCollection; +const DEFAULT_PIDS: [i32; 14] = [21, -6, -5, -4, -3, -2, -1, 1, 2, 3, 4, 5, 6, 22]; + /// Result codes for `NeoPDF` operations #[repr(C)] #[derive(Debug, Clone, Copy, PartialEq, Eq)] @@ -1202,12 +1204,10 @@ pub unsafe extern "C" fn evolvepdf(x: c_double, q: c_double, f: *mut c_double) { let q2 = q * q; let (available_pids, _) = pdf.pids().clone().into_raw_vec_and_offset(); + let out_slice = slice::from_raw_parts_mut(f, 14); - let pids_to_calc = [-6, -5, -4, -3, -2, -1, 21, 1, 2, 3, 4, 5, 6]; - let out_slice = slice::from_raw_parts_mut(f, 13); - - for i in 0..13 { - let pid = pids_to_calc[i]; + for i in 0..14 { + let pid = DEFAULT_PIDS[i]; if available_pids.contains(&pid) { out_slice[i] = pdf.xfxq2(pid, &[x, q2]); } else { @@ -1237,12 +1237,10 @@ pub unsafe extern "C" fn evolvepdf_(x: *const c_double, q: *const c_double, f: * let q2 = (*q) * (*q); let (available_pids, _) = pdf.pids().clone().into_raw_vec_and_offset(); + let out_slice = slice::from_raw_parts_mut(f, 14); - let pids_to_calc = [-6, -5, -4, -3, -2, -1, 21, 1, 2, 3, 4, 5, 6]; - let out_slice = slice::from_raw_parts_mut(f, 13); - - for i in 0..13 { - let pid = pids_to_calc[i]; + for i in 0..14 { + let pid = DEFAULT_PIDS[i]; if available_pids.contains(&pid) { out_slice[i] = pdf.xfxq2(pid, &[*x, q2]); } else { diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.cpp b/neopdf_capi/tests/check-lhapdf-compatibility.cpp index 494cf126..9b99836f 100644 --- a/neopdf_capi/tests/check-lhapdf-compatibility.cpp +++ b/neopdf_capi/tests/check-lhapdf-compatibility.cpp @@ -52,19 +52,19 @@ void test_lhapdf_compatibility_c() { initpdf(0); double x = 0.1; double q = 10.0; - double f_compat[13]; - evolvepdf(x, q, f_compat); + double xfxs[14]; + evolvepdf(x, q, xfxs); // Use the native NeoPDF API for comparison NeoPDFWrapper* neo_pdf = neopdf_pdf_load(pdfname, 0); - double expected_g = neopdf_pdf_xfxq2(neo_pdf, 21, x, q*q); + double expected_g = neopdf_pdf_xfxq2(neo_pdf, -1, x, q*q); neopdf_pdf_free(neo_pdf); printf("Comparing gluon PDF at x=0.1, Q=10.0\n"); printf("Native NeoPDF: %e\n", expected_g); - printf("NeoPDF (LHAPDF compat): %e\n", f_compat[6]); + printf("NeoPDF (LHAPDF compat): %e\n", xfxs[6]); - double reldiff = fabs(f_compat[6] - expected_g) / expected_g; + double reldiff = fabs(xfxs[6] - expected_g) / expected_g; printf("Relative difference: %e\n", reldiff); assert(reldiff < TOLERANCE); diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.output b/neopdf_capi/tests/check-lhapdf-compatibility.output index 0137760a..5eeadc38 100644 --- a/neopdf_capi/tests/check-lhapdf-compatibility.output +++ b/neopdf_capi/tests/check-lhapdf-compatibility.output @@ -6,7 +6,7 @@ Relative difference: 0 LHAPDF C++ OOP Compatibility test passed. === Test LHAPDF C Compatibility Layer === Comparing gluon PDF at x=0.1, Q=10.0 -Native NeoPDF: 1.132265e+00 -NeoPDF (LHAPDF compat): 1.132265e+00 +Native NeoPDF: 1.420505e-01 +NeoPDF (LHAPDF compat): 1.420505e-01 Relative difference: 0.000000e+00 LHAPDF C Compatibility test passed. From a20e1741fa6371c7436eaf7aedcbcc91a18326ed Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Sat, 13 Sep 2025 21:54:26 +0200 Subject: [PATCH 7/9] Update changelogs --- CHANGELOG.md | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 38ac9f48..7f4a56f6 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,7 +9,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added -- Add an interface to the Wolfram Language to allow Rust APIs to be called in +- Added proper LHAPDF drop-in compatibility layer for no-code migration. +- Added an interface to the Wolfram Language to allow Rust APIs to be called in Mathematica. - Added the logic to determine the Euclidean distance of a point to the closest subgrid in order to allow extrapolation. From ecc9a96bbc44d3ff241fc4530a961f6be138fd76 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Sat, 13 Sep 2025 22:23:34 +0200 Subject: [PATCH 8/9] Address warnings --- neopdf_capi/src/lib.rs | 41 +++++++++++++++++++++++------------------ 1 file changed, 23 insertions(+), 18 deletions(-) diff --git a/neopdf_capi/src/lib.rs b/neopdf_capi/src/lib.rs index 290f46ae..ac3d6ff8 100644 --- a/neopdf_capi/src/lib.rs +++ b/neopdf_capi/src/lib.rs @@ -1196,23 +1196,26 @@ pub unsafe extern "C" fn initpdf_(member: *const c_int) { #[no_mangle] pub unsafe extern "C" fn evolvepdf(x: c_double, q: c_double, f: *mut c_double) { unsafe { - if let Some(pdfs) = &LHAPDF_STATE.pdf_set { - if LHAPDF_STATE.member >= pdfs.len() { + let state_ptr = &raw const LHAPDF_STATE; + let pdf_set_ptr = &raw const (*state_ptr).pdf_set; + + if let Some(pdfs) = &*pdf_set_ptr { + let member = (*state_ptr).member; + if member >= pdfs.len() { return; } - let pdf = &pdfs[LHAPDF_STATE.member]; + let pdf = &pdfs[member]; let q2 = q * q; let (available_pids, _) = pdf.pids().clone().into_raw_vec_and_offset(); let out_slice = slice::from_raw_parts_mut(f, 14); - for i in 0..14 { - let pid = DEFAULT_PIDS[i]; - if available_pids.contains(&pid) { - out_slice[i] = pdf.xfxq2(pid, &[x, q2]); + for (out, &pid) in out_slice.iter_mut().zip(DEFAULT_PIDS.iter()) { + *out = if available_pids.contains(&pid) { + pdf.xfxq2(pid, &[x, q2]) } else { - out_slice[i] = 0.0; - } + 0.0 + }; } } } @@ -1229,23 +1232,25 @@ pub unsafe extern "C" fn evolvepdf(x: c_double, q: c_double, f: *mut c_double) { #[no_mangle] pub unsafe extern "C" fn evolvepdf_(x: *const c_double, q: *const c_double, f: *mut c_double) { unsafe { - if let Some(pdfs) = &LHAPDF_STATE.pdf_set { - if LHAPDF_STATE.member >= pdfs.len() { + let pdf_set_ptr = &raw const LHAPDF_STATE.pdf_set; + + if let Some(pdfs) = &*pdf_set_ptr { + let member = (&raw const LHAPDF_STATE.member).read(); + if member >= pdfs.len() { return; } - let pdf = &pdfs[LHAPDF_STATE.member]; + let pdf = &pdfs[member]; let q2 = (*q) * (*q); let (available_pids, _) = pdf.pids().clone().into_raw_vec_and_offset(); let out_slice = slice::from_raw_parts_mut(f, 14); - for i in 0..14 { - let pid = DEFAULT_PIDS[i]; - if available_pids.contains(&pid) { - out_slice[i] = pdf.xfxq2(pid, &[*x, q2]); + for (out, &pid) in out_slice.iter_mut().zip(DEFAULT_PIDS.iter()) { + *out = if available_pids.contains(&pid) { + pdf.xfxq2(pid, &[*x, q2]) } else { - out_slice[i] = 0.0; - } + 0.0 + }; } } } From 46cb2bb23e5d019bc290d779e3dfb45d7bf87486 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Sat, 13 Sep 2025 22:49:40 +0200 Subject: [PATCH 9/9] Add `alphaspdf` LHAPDF compatibility --- neopdf_capi/src/lib.rs | 54 +++++++++++++++++++ .../tests/check-lhapdf-compatibility.cpp | 23 +++++++- .../tests/check-lhapdf-compatibility.output | 8 +++ 3 files changed, 83 insertions(+), 2 deletions(-) diff --git a/neopdf_capi/src/lib.rs b/neopdf_capi/src/lib.rs index ac3d6ff8..a50d73ea 100644 --- a/neopdf_capi/src/lib.rs +++ b/neopdf_capi/src/lib.rs @@ -1255,3 +1255,57 @@ pub unsafe extern "C" fn evolvepdf_(x: *const c_double, q: *const c_double, f: * } } } + +/// Evaluates the strong coupling `alpha_s` at scale `q` for the active member. +/// +/// # Safety +/// +/// - Requires that a PDF set has been initialized via `initpdfsetbyname` or its +/// Fortran variant. If no set is loaded or the member index is out of range, +/// the function returns 0.0. +#[no_mangle] +pub unsafe extern "C" fn alphaspdf(q: c_double) -> c_double { + unsafe { + let state_ptr = &raw const LHAPDF_STATE; + let pdf_set_ptr = &raw const (*state_ptr).pdf_set; + + if let Some(pdfs) = &*pdf_set_ptr { + let member = (*state_ptr).member; + if member >= pdfs.len() { + return 0.0; + } + let pdf = &pdfs[member]; + let q2 = q * q; + pdf.alphas_q2(q2) + } else { + 0.0 + } + } +} + +/// Fortran name-mangled variant of `alphaspdf`. +/// +/// # Safety +/// +/// - `q` must be a valid pointer to a `c_double` value. +/// - A PDF set must have been initialized and the member index must be valid or +/// the function will return 0.0. +#[no_mangle] +pub unsafe extern "C" fn alphaspdf_(q: *const c_double) -> c_double { + unsafe { + let state_ptr = &raw const LHAPDF_STATE; + let pdf_set_ptr = &raw const (*state_ptr).pdf_set; + + if let Some(pdfs) = &*pdf_set_ptr { + let member = (*state_ptr).member; + if member >= pdfs.len() { + return 0.0; + } + let pdf = &pdfs[member]; + let q2 = (*q) * (*q); + pdf.alphas_q2(q2) + } else { + 0.0 + } + } +} diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.cpp b/neopdf_capi/tests/check-lhapdf-compatibility.cpp index 9b99836f..496c3c55 100644 --- a/neopdf_capi/tests/check-lhapdf-compatibility.cpp +++ b/neopdf_capi/tests/check-lhapdf-compatibility.cpp @@ -34,7 +34,17 @@ void test_lhapdf_compatibility_oop() { std::cout << "NeoPDF (LHAPDF compat): " << result << std::endl; std::cout << "Relative difference: " << reldif << std::endl; + double as_result = neo_lha_pdf->alphasQ2(q2); + double as_expected = neopdf_pdf_alphas_q2(neo_pdf, q2); + double as_reldif = std::abs(as_result - as_expected) / as_expected; + + std::cout << "Q2=" << q2 << std::endl; + std::cout << "Native NeoPDF: " << as_expected << std::endl; + std::cout << "NeoPDF (LHAPDF compat): " << as_result << std::endl; + std::cout << "Relative difference: " << as_reldif << std::endl; + assert(reldif < TOLERANCE); + assert(as_reldif < TOLERANCE); delete neo_lha_pdf; neopdf_pdf_free(neo_pdf); @@ -54,20 +64,29 @@ void test_lhapdf_compatibility_c() { double q = 10.0; double xfxs[14]; evolvepdf(x, q, xfxs); + double a_s = alphaspdf(q); // Use the native NeoPDF API for comparison NeoPDFWrapper* neo_pdf = neopdf_pdf_load(pdfname, 0); double expected_g = neopdf_pdf_xfxq2(neo_pdf, -1, x, q*q); neopdf_pdf_free(neo_pdf); + double as_expected = neopdf_pdf_alphas_q2(neo_pdf, q*q); + + double reldiff = fabs(xfxs[6] - expected_g) / expected_g; + double as_reldiff = fabs(a_s - as_expected) / as_expected; printf("Comparing gluon PDF at x=0.1, Q=10.0\n"); printf("Native NeoPDF: %e\n", expected_g); printf("NeoPDF (LHAPDF compat): %e\n", xfxs[6]); - - double reldiff = fabs(xfxs[6] - expected_g) / expected_g; printf("Relative difference: %e\n", reldiff); + printf("Comparing AlphaS at Q=10.0\n"); + printf("Native NeoPDF: %e\n", as_expected); + printf("NeoPDF (LHAPDF compat): %e\n", a_s); + printf("Relative difference: %e\n", as_reldiff); + assert(reldiff < TOLERANCE); + assert(as_reldiff < TOLERANCE); printf("LHAPDF C Compatibility test passed.\n"); } diff --git a/neopdf_capi/tests/check-lhapdf-compatibility.output b/neopdf_capi/tests/check-lhapdf-compatibility.output index 5eeadc38..7b8c4a53 100644 --- a/neopdf_capi/tests/check-lhapdf-compatibility.output +++ b/neopdf_capi/tests/check-lhapdf-compatibility.output @@ -3,10 +3,18 @@ pid=21, x=0.0001, Q2=100 Native NeoPDF: 30.8985 NeoPDF (LHAPDF compat): 30.8985 Relative difference: 0 +Q2=100 +Native NeoPDF: 0.178123 +NeoPDF (LHAPDF compat): 0.178123 +Relative difference: 0 LHAPDF C++ OOP Compatibility test passed. === Test LHAPDF C Compatibility Layer === Comparing gluon PDF at x=0.1, Q=10.0 Native NeoPDF: 1.420505e-01 NeoPDF (LHAPDF compat): 1.420505e-01 Relative difference: 0.000000e+00 +Comparing AlphaS at Q=10.0 +Native NeoPDF: 1.781227e-01 +NeoPDF (LHAPDF compat): 1.781227e-01 +Relative difference: 0.000000e+00 LHAPDF C Compatibility test passed.