From f5fee88f3ad964d4e42952bff190108dd7776959 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Mon, 7 Apr 2025 12:03:13 +0200 Subject: [PATCH 01/15] Add functions to get and rotate PID basis in C-API --- examples/cpp/advanced-convolution.cpp | 9 +++++++ pineappl_capi/src/lib.rs | 36 +++++++++++++++++++++++++++ 2 files changed, 45 insertions(+) diff --git a/examples/cpp/advanced-convolution.cpp b/examples/cpp/advanced-convolution.cpp index e16f9e30..95fdfb4a 100644 --- a/examples/cpp/advanced-convolution.cpp +++ b/examples/cpp/advanced-convolution.cpp @@ -123,5 +123,14 @@ int main(int argc, char* argv[]) { << normalizations.at(bin) << '\n'; } + // modify the particle ID representation of the Grid + char* ref_pid_repr_basis = pineappl_get_pid_basis(grid); + assert(std::string(ref_pid_repr_basis) == "PINEAPPL_PID_BASIS_EVOL"); + + pineappl_rotate_pid_basis(grid, PINEAPPL_PID_BASIS_PDG); + + char* mod_pid_repr_basis = pineappl_get_pid_basis(grid); + assert(std::string(mod_pid_repr_basis) == "PINEAPPL_PID_BASIS_PDG"); + pineappl_grid_delete(grid); } diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 20b792ef..263f3cf6 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -895,6 +895,42 @@ pub unsafe extern "C" fn pineappl_grid_split_lumi(grid: *mut Grid) { grid.split_channels(); } +/// Change the particle ID representation of a given Grid. +/// +/// # Safety +/// +/// If `grid` does not point to a valid `Grid` object (for example when `grid` is the null pointer) +/// or if `pid_basis` doe not refer to a correct basis, then this function is not safe to call. +#[no_mangle] +pub unsafe extern "C" fn pineappl_rotate_pid_basis(grid: *mut Grid, pid_basis: PidBasis) { + let grid = unsafe { &mut *grid }; + + grid.rotate_pid_basis(pid_basis); +} + +/// Get the particle ID representation of a Grid. +/// +/// # Safety +/// +/// If `grid` does not point to a valid `Grid` object, for example when `grid` is the `NULL` +/// pointer, this function is not safe to call. The parameter `key` must be non-`NULL` and a valid +/// C string. +/// +/// # Panics +/// +/// This function might panic if `pid_basis` contains interior null bytes. This occurs if a new +/// variant is added to the `PidBasis` enum in the future and the match statement isn't updated. +#[no_mangle] +pub unsafe extern "C" fn pineappl_get_pid_basis(grid: *mut Grid) -> *mut c_char { + let grid = unsafe { &mut *grid }; + let pid_basis = match grid.pid_basis() { + PidBasis::Pdg => "PINEAPPL_PID_BASIS_PDG", + PidBasis::Evol => "PINEAPPL_PID_BASIS_EVOL", + }; + + CString::new(pid_basis).unwrap().into_raw() +} + /// Optimizes the grid representation for space efficiency. /// /// # Safety From 9ded95b07e0bdc0b7ce7f2e0c1300a88fbcca825 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Tue, 8 Apr 2025 12:44:29 +0200 Subject: [PATCH 02/15] Fix memory leak in `advanced-convolution.cpp` --- examples/cpp/advanced-convolution.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/examples/cpp/advanced-convolution.cpp b/examples/cpp/advanced-convolution.cpp index 95fdfb4a..3ddb796c 100644 --- a/examples/cpp/advanced-convolution.cpp +++ b/examples/cpp/advanced-convolution.cpp @@ -126,11 +126,13 @@ int main(int argc, char* argv[]) { // modify the particle ID representation of the Grid char* ref_pid_repr_basis = pineappl_get_pid_basis(grid); assert(std::string(ref_pid_repr_basis) == "PINEAPPL_PID_BASIS_EVOL"); + pineappl_string_delete(ref_pid_repr_basis); pineappl_rotate_pid_basis(grid, PINEAPPL_PID_BASIS_PDG); char* mod_pid_repr_basis = pineappl_get_pid_basis(grid); assert(std::string(mod_pid_repr_basis) == "PINEAPPL_PID_BASIS_PDG"); + pineappl_string_delete(mod_pid_repr_basis); pineappl_grid_delete(grid); } From e03260cf3eeb89b051315cd34b90f7c55af32522 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Tue, 8 Apr 2025 21:32:22 +0200 Subject: [PATCH 03/15] Address review comments --- examples/cpp/advanced-convolution.cpp | 12 +++++------- pineappl_capi/src/lib.rs | 24 +++++++----------------- 2 files changed, 12 insertions(+), 24 deletions(-) diff --git a/examples/cpp/advanced-convolution.cpp b/examples/cpp/advanced-convolution.cpp index 3ddb796c..16f67eed 100644 --- a/examples/cpp/advanced-convolution.cpp +++ b/examples/cpp/advanced-convolution.cpp @@ -124,15 +124,13 @@ int main(int argc, char* argv[]) { } // modify the particle ID representation of the Grid - char* ref_pid_repr_basis = pineappl_get_pid_basis(grid); - assert(std::string(ref_pid_repr_basis) == "PINEAPPL_PID_BASIS_EVOL"); - pineappl_string_delete(ref_pid_repr_basis); + pineappl_pid_basis ref_pid_repr_basis = pineappl_grid_pid_basis(grid); + assert(ref_pid_repr_basis == PINEAPPL_PID_BASIS_EVOL); - pineappl_rotate_pid_basis(grid, PINEAPPL_PID_BASIS_PDG); + pineappl_grid_rotate_pid_basis(grid, PINEAPPL_PID_BASIS_PDG); - char* mod_pid_repr_basis = pineappl_get_pid_basis(grid); - assert(std::string(mod_pid_repr_basis) == "PINEAPPL_PID_BASIS_PDG"); - pineappl_string_delete(mod_pid_repr_basis); + pineappl_pid_basis mod_pid_repr_basis = pineappl_grid_pid_basis(grid); + assert(mod_pid_repr_basis == PINEAPPL_PID_BASIS_PDG); pineappl_grid_delete(grid); } diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 263f3cf6..4baaa1d1 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -895,40 +895,30 @@ pub unsafe extern "C" fn pineappl_grid_split_lumi(grid: *mut Grid) { grid.split_channels(); } -/// Change the particle ID representation of a given Grid. +/// Change the particle ID basis of a given Grid. /// /// # Safety /// /// If `grid` does not point to a valid `Grid` object (for example when `grid` is the null pointer) -/// or if `pid_basis` doe not refer to a correct basis, then this function is not safe to call. +/// or if `pid_basis` does not refer to a correct basis, then this function is not safe to call. #[no_mangle] -pub unsafe extern "C" fn pineappl_rotate_pid_basis(grid: *mut Grid, pid_basis: PidBasis) { +pub unsafe extern "C" fn pineappl_grid_rotate_pid_basis(grid: *mut Grid, pid_basis: PidBasis) { let grid = unsafe { &mut *grid }; grid.rotate_pid_basis(pid_basis); } -/// Get the particle ID representation of a Grid. +/// Get the particle ID basis of a Grid. /// /// # Safety /// /// If `grid` does not point to a valid `Grid` object, for example when `grid` is the `NULL` -/// pointer, this function is not safe to call. The parameter `key` must be non-`NULL` and a valid -/// C string. -/// -/// # Panics -/// -/// This function might panic if `pid_basis` contains interior null bytes. This occurs if a new -/// variant is added to the `PidBasis` enum in the future and the match statement isn't updated. +/// pointer, this function is not safe to call. #[no_mangle] -pub unsafe extern "C" fn pineappl_get_pid_basis(grid: *mut Grid) -> *mut c_char { +pub unsafe extern "C" fn pineappl_grid_pid_basis(grid: *mut Grid) -> PidBasis { let grid = unsafe { &mut *grid }; - let pid_basis = match grid.pid_basis() { - PidBasis::Pdg => "PINEAPPL_PID_BASIS_PDG", - PidBasis::Evol => "PINEAPPL_PID_BASIS_EVOL", - }; - CString::new(pid_basis).unwrap().into_raw() + *grid.pid_basis() } /// Optimizes the grid representation for space efficiency. From d3c7b40f7d884835a4c16af33457a7989ebb275e Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Thu, 10 Apr 2025 21:49:55 +0200 Subject: [PATCH 04/15] Proposal to expose the subgrid contents of a Grid --- examples/cpp/Makefile | 4 ++ examples/cpp/get-subgrids | Bin 0 -> 16912 bytes examples/cpp/output | 26 ++++++++++++ pineappl/src/packed_array.rs | 7 +++- pineappl_capi/src/lib.rs | 76 +++++++++++++++++++++++++++++++++++ 5 files changed, 112 insertions(+), 1 deletion(-) create mode 100755 examples/cpp/get-subgrids diff --git a/examples/cpp/Makefile b/examples/cpp/Makefile index 7f44b4f6..aef4689b 100644 --- a/examples/cpp/Makefile +++ b/examples/cpp/Makefile @@ -13,6 +13,7 @@ PROGRAMS = \ advanced-filling \ convolve-grid-deprecated \ convolve-grid \ + get-subgrids \ deprecated \ display-channels-deprecated \ display-channels \ @@ -47,6 +48,9 @@ convolve-grid: convolve-grid.cpp deprecated: deprecated.cpp $(CXX) $(CXXFLAGS) $< $(LHAPDF_DEPS) $(PINEAPPL_DEPS) -o $@ +get-subgrids: get-subgrids.cpp + $(CXX) $(CXXFLAGS) $< $(PINEAPPL_DEPS) -o $@ + display-channels-deprecated: display-channels-deprecated.cpp $(CXX) $(CXXFLAGS) $< $(PINEAPPL_DEPS) -o $@ diff --git a/examples/cpp/get-subgrids b/examples/cpp/get-subgrids new file mode 100755 index 0000000000000000000000000000000000000000..32dacc697f6dfd8bfc8b4b25d0499a50f38434f1 GIT binary patch literal 16912 zcmeHO4RBo5b-t?~YmJRoc5D!)!3(TG4rr~FY#|vOc4f&tnU$?ra&V#Md9}NewybuQ z-6wnP5~3iR@kS_!kPr_rjhji@)R_j_rkNJfR@h+SHf{Ln6lc& zT-2B{*;P_e$%TA(m8r*^>vQcn&3aM5e?bZ$=2klzo5^}vEhg9H^BUD}jXR^_qTwRD zgUarpvSaR4TRMi~Od)qsQM)GeE~in38LMHt)LgiXrxff1D{nd?b+2@V=|vOjXD{--CNO1<^R zuV43r{^J`v55H_1Asy-m*^mwqVRSH2_?O`z9dmqYWSKCcr3V(SH%oK>-vD^=993T_ zg2&O|RdeY7YZ3i>pg))YUl-BeUPS*2bZ{;|-z|baRm8u)i2vJ)_&HMqUr_|FEP}Tc z!S5&%=XZ+Wa3h>J3Ta)Deto}4yGb~hYdk*<+*#-xGF^{z;Z10Fkys|Cyifqupgl=? zXK$At3Maw?k)#n$^meU}#$(~$KwmU0bbVkb9@CRXAYtgbkg9Vk4o70)!0>QXA4o(( z@D>dlVUZ8^MPhm|J`yu#@#ILK)K3luhAryBcx-1pI$}iPv7{ai$3AQjNF)NgXJP3i z5e|feetVD6*r4miU?RSYYzK@%J)B6y6WzB4QM;)p-l})$kytXEFgk}iI=1x0d(p9i ziqLFDDB8W#bnA)4lX_nu8P>7nBm#-ua1$|5sSl)5`fxarjK>1eh_PGW*&r=E&AJ{> z8VGENTMYIX9#3#EkkE}pAYvptgVd7YX@cKKkhME|_1*@(qo-N#zO|=GkHTAQ*N}|T zxIx#0sg%d#X*5mXKN1`0oWpTXfV@+8g9yP0cWfE@gg!hCa1x3nhvUg`hdCOnaVia= zZ;kDW#6tR(aB^gbrVD0D@Gd<#c$eNEh#+nns-T)bLxD&PNH7)9`y(_`d%|KUJd_L@ z0$RDr5)2GSX!-UKnlwVe>#ifQnZ%)7X&?x$h$5wRP$qcaylq z-`NHc@M?EcW1iUP$rD?xl3CUp-K%q@rrBbndkrof+&OT!D92HOyGuF#Y*^VVL>2DY z6=J@uw~31&lcr7HK}@S3;##Epg%d~jHM0-<3sSBUrubw|CDJe{~x#Ds!55TwBWTC{R3W8pyiMH%QWRROZ1$;_hgAEoMPc6 zlW@BQ$8hGT#)4b#|FsreKHn+}rv*2k3(Sf-3(jLoH4PS=pIHcRvEUaGpvdbr3lu*l zT8=cYStNX64#j-5;OY@w)^4-l4vYQ{3%d@9fwvr2yr+%*+OCaP9$O`Zwr|QPn|WRvd&Yi5j>62Ee*!$S;wSj; zSnfrM<~6?K3{7HSjg zUVw4tZPaBpL%(a^DXz6y~(5TD2PmPTCx?1`j@k%$E zRdv0tZpR`);7WjwMc#{PC}7m<(vDWt2^>HU4dy?ie6wYCu zAcD-prtJf;1$xRv|A)BxOq)BXJr&pEY-o*Zt*O^UdNeK3z?Td3WL%3aaouSuWG?)M zCfJN)?+RRS+TUWX5IOvHCc#mM!T z-P+i*woC&dtx=cVvG*Goz^;AgMkYN!@9@l=dZwQzrN8xk$6nf$YU!WJQFbgJm8}O4 zYmUcl+St?f^x@9dpN^dN2hUr3bl!3AB>cIixjSzG@~6-HSDw_KK40F|`YXrY z@4`Pi?|HG`ac1h&CGvXG#tz#$$ENLC>tV+jU1$A1hr4P|neofc@eXwAG)-8FM;q5L zYqR4%%v#^C`Xgk_MTH-ihW%ABSDk*j>u=cuX5Rmv%_io#A$^X%M) zYa#PK46w#H_U^>BMJM(UfNTVUjJLz{3Qj=QLDPcP1L}3=25LkuMgJm8c}r0<7XOrI z>W=NcKk#k$-RaYh{3Bg1nLbj@jW11=?BAi4&h$+7FTu_F^t0Io2r|2d{_T+AKAYVR z&gY-*H=fr!(=XhqrGKG~y>hm@x5G2#Ifgaub5NUE+O$%Lw`UhI;O_#_CMw?~)_UGp zii`acD&tl)f99l`KSy{TO{_a6_17FGw(^3OKC3-_`er(*n(e6e!Ubaq9K4FvnjF|? zPGV|eT-cuC8(sG|pFxjDnzm|V>mH$XT}!`eV0~M6Cw$H<-47@eM6gFI?}wso`;pmr zkMz%$>3&K6vIEeCm)yNJcXtI?qS6SIMxZnTr4cBNKxqU@Bk(&O0eaVrpHM_{z!?dp zi0uqSN5Xh>TyG&Ef`M2pZaB?X{!Y9x){C&IlSw-dO ze37lJJZBf*oA>v&({}N7dnGE4*$IO3#g|rgi(6svN0p=0lLY0ia)IQ zpDO+%#g8if?~4CY@pFpLQ%O9RD1MFN&5CnelU^>F9JzOEops25Ug`9#c6;0&XJbR- z>V`E9H;7_Yjg9#(@+6uj72*y0IVO)Nm#9os3$KdjK}+pi`sTb}l*?ad^;Ia?kWo}P z-YpWZ&COfZ--H6ST}VLybKzeEUV}`sgWh}#dQ;#8zP|w93$;a*k>)h@@#;d-L->D{ zDqGnPAS}i>oK)jL?>dP78pmAm{AUsTEbv8Tmx)O)3KpW@Y_9}l3@$P((ji1b+zQL;d3KTI4$-INMo7e}57D(IWVd2rt9xs^&Y#`D5TT4u$9I zxgvgCGETOh<-35-m1jKxoaRH}d4ELW%fxZttPd)EOW{-AeEj>s=NkW)i{KZE;8$QF zoXbyh5xfhyGd~n5MT+Q;7QtUCg1=V8|GPzSF%UNDRRVd^9ULAOlwa%)#-qpr*Rc-i z9h$yji?6Fg*V?xTz5Vw$`MNsS3;h=V<~E;S-@IW%PY3d>eQo{@y&#pmFn|3+3F&$` z6fgqFhu@6c<`5FwN0MP`)|%Vwm8r)_A*O6$GxzzUlG9CfGgDcmMdx|Y^-wada+6Y0+S?Uej(wKEgZ^ALVVY-awON%#e6 z;5nagNh$Ma6<-saZyq2{0{Il{ApI(?>Qhek$AbC0d zUe2-qL}Q!(F+yISFt)j%^S=xE*dL(?cXIcT5$FRq5+)zyQY40NTZV-jzm0_51F;eJ za3YRRHjLdCsBa{K&&VPn1>qBddLu9(r21eWIVjwr-7z#_awB2Z?8Jv3_);ZL&{2oq zQvzh5O2bh@xaI738{rhn^7?Zp;<%hX?ikL|`ba4~AgN5}e{Yhj1zw9yWCJ2$#O} zp_`{o`rtB}3JgVpFp0xE{<&q7rh80AfwggHD2(rKek(-#4cr>>EM`(I)@dTd&kF2$Kf`=wE-Ge{ zVYMH1nwre%w_+<5wtozjR{JNMrZ)4pt;Q(x`D6Gc6s-2&xzem*{#jO3oOG@9Duikr z>aw5x?z`U9Dr}Es=5IjFVX^1;EC)WOJ?*dP6m$E$59RO4o$AEU^+Uw%@_Cqt4pE=j zp1-G>qKz%6&Sbg&Y|nf_fjz&c$te4SIZLy^cFdO)*z>+-e}k#d*D1SUJC>Khp0Hf} zUbE&P+f(N(Y(FaQ7m literal 0 HcmV?d00001 diff --git a/examples/cpp/output b/examples/cpp/output index e7262e30..d7f24048 100644 --- a/examples/cpp/output +++ b/examples/cpp/output @@ -146,6 +146,32 @@ idx left right dsig/dx dx 21 2.1 2.2 3.507294e-02 0.1 22 2.2 2.3 1.976542e-02 0.1 23 2.3 2.4 5.590552e-03 0.1 + bin sg idx sg value + --- ------ ------------ + 0 41020 -4.93616e-07 + 1 41020 -4.5499e-08 + 2 40971 -4.10008e-08 + 3 40971 -2.78597e-07 + 4 40971 -2.01924e-07 + 5 40922 -3.39322e-10 + 6 40922 -8.43358e-08 + 7 40922 -2.99247e-07 + 8 40922 -1.5117e-07 + 9 40872 -9.95013e-11 + 10 40873 -1.30578e-07 + 11 40873 -2.55869e-07 + 12 40823 -1.0823e-09 + 13 40823 -4.22686e-09 + 14 40824 -1.68149e-07 + 15 40774 -7.58702e-10 + 16 40774 -2.41887e-08 + 17 40774 -8.83139e-09 + 18 40725 -1.14932e-09 + 19 40725 -3.03857e-08 + 20 40725 -2.66414e-08 + 21 40675 -3.60243e-10 + 22 40676 -1.22785e-08 + 23 40626 -7.04493e-11 idx p-p c#0 l#0 p-p~ c#0 l# p-d c#0 l#0 p-d dx --- ------------ ----------- ------------- ------------ ------ 0 5.263109e-01 5.263109e-01 5.263109e-01 5.263109e-01 0.1 diff --git a/pineappl/src/packed_array.rs b/pineappl/src/packed_array.rs index 0c5d9587..f93c544d 100644 --- a/pineappl/src/packed_array.rs +++ b/pineappl/src/packed_array.rs @@ -171,7 +171,12 @@ impl From> for PackedArray { } /// Converts a `multi_index` into a flat index. -fn ravel_multi_index(multi_index: &[usize], shape: &[usize]) -> usize { +/// +/// # Panics +/// +/// TODO +#[must_use] +pub fn ravel_multi_index(multi_index: &[usize], shape: &[usize]) -> usize { assert_eq!(multi_index.len(), shape.len()); multi_index diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 4baaa1d1..8aa0cf75 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -60,7 +60,9 @@ use pineappl::boc::{Bin, BinsWithFillLimits, Channel, Kinematics, Order, ScaleFu use pineappl::convolutions::{Conv, ConvType, ConvolutionCache}; use pineappl::grid::{Grid, GridOptFlags}; use pineappl::interpolation::{Interp, InterpMeth, Map, ReweightMeth}; +use pineappl::packed_array::ravel_multi_index; use pineappl::pids::PidBasis; +use pineappl::subgrid::Subgrid; use std::collections::HashMap; use std::ffi::{CStr, CString}; use std::fs::File; @@ -1893,3 +1895,77 @@ pub unsafe extern "C" fn pineappl_grid_convolve( mu_scales, )); } + +/// Get the number of convolutions for a given Grid. +/// +/// # Safety +/// +/// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, +/// this function is not safe to call. +#[no_mangle] +pub unsafe extern "C" fn pineappl_convolutions_len(grid: *mut Grid) -> usize { + let grid = unsafe { &mut *grid }; + grid.convolutions().len() +} + +/// Get the shape of a subgrid for a given bin, channel, and order. +/// +/// # Safety +/// +/// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, +/// this function is not safe to call. Additionally, the pointer that specifies the shape of the +/// subgrid must be an array whose size must be `(n+1)` where `n` represents the number of +/// convolutions as given by `pineappl_convolutions_len`. +#[no_mangle] +pub unsafe extern "C" fn pineappl_subgrid_shape( + grid: *mut Grid, + bin: usize, + order: usize, + channel: usize, + shape: *mut usize, +) { + let grid = unsafe { &mut *grid }; + let mut subgrid = grid.subgrids()[[order, bin, channel]].clone(); + + let subgrid_shape = if subgrid.is_empty() { + let subgrid_dim = grid.convolutions().len() + 1; + &vec![0; subgrid_dim] + } else { + subgrid.shape() + }; + + let shape = unsafe { slice::from_raw_parts_mut(shape, grid.convolutions().len() + 1) }; + shape.copy_from_slice(subgrid_shape); +} + +/// Get the subgrid for a given bin, channel, and order +/// +/// # Safety +/// +/// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, +/// this function is not safe to call. Additionally, the pointer that specifies the size of the subgrid +/// when flattened must be an array; its size must be computed by multiplying the shape dimension as +/// given by `pineappl_subgrid_shape`. +#[no_mangle] +pub unsafe extern "C" fn pineappl_subgrid_array( + grid: *mut Grid, + bin: usize, + order: usize, + channel: usize, + subgrid_array: *mut f64, +) { + let grid = unsafe { &mut *grid }; + let mut subgrid = grid.subgrids()[[order, bin, channel]].clone(); + + let flattened_shape = subgrid.shape().iter().copied().product(); + let mut flattened_subgrid_array = vec![0.0; flattened_shape]; + + for (index, value) in subgrid.indexed_iter() { + let ravel_index = ravel_multi_index(index.as_slice(), subgrid.clone().shape()); + flattened_subgrid_array[ravel_index] = value; + } + + let subgrid_array = unsafe { slice::from_raw_parts_mut(subgrid_array, flattened_shape) }; + + subgrid_array.copy_from_slice(&flattened_subgrid_array); +} From 9ed7d3cb48c36145293b94911f6a94ef41337cc5 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Thu, 10 Apr 2025 22:00:49 +0200 Subject: [PATCH 05/15] Add the correct `C++` file --- examples/cpp/get-subgrids | Bin 16912 -> 0 bytes examples/cpp/get-subgrids.cpp | 54 ++++++++++++++++++++++++++++++++++ 2 files changed, 54 insertions(+) delete mode 100755 examples/cpp/get-subgrids create mode 100644 examples/cpp/get-subgrids.cpp diff --git a/examples/cpp/get-subgrids b/examples/cpp/get-subgrids deleted file mode 100755 index 32dacc697f6dfd8bfc8b4b25d0499a50f38434f1..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 16912 zcmeHO4RBo5b-t?~YmJRoc5D!)!3(TG4rr~FY#|vOc4f&tnU$?ra&V#Md9}NewybuQ z-6wnP5~3iR@kS_!kPr_rjhji@)R_j_rkNJfR@h+SHf{Ln6lc& zT-2B{*;P_e$%TA(m8r*^>vQcn&3aM5e?bZ$=2klzo5^}vEhg9H^BUD}jXR^_qTwRD zgUarpvSaR4TRMi~Od)qsQM)GeE~in38LMHt)LgiXrxff1D{nd?b+2@V=|vOjXD{--CNO1<^R zuV43r{^J`v55H_1Asy-m*^mwqVRSH2_?O`z9dmqYWSKCcr3V(SH%oK>-vD^=993T_ zg2&O|RdeY7YZ3i>pg))YUl-BeUPS*2bZ{;|-z|baRm8u)i2vJ)_&HMqUr_|FEP}Tc z!S5&%=XZ+Wa3h>J3Ta)Deto}4yGb~hYdk*<+*#-xGF^{z;Z10Fkys|Cyifqupgl=? zXK$At3Maw?k)#n$^meU}#$(~$KwmU0bbVkb9@CRXAYtgbkg9Vk4o70)!0>QXA4o(( z@D>dlVUZ8^MPhm|J`yu#@#ILK)K3luhAryBcx-1pI$}iPv7{ai$3AQjNF)NgXJP3i z5e|feetVD6*r4miU?RSYYzK@%J)B6y6WzB4QM;)p-l})$kytXEFgk}iI=1x0d(p9i ziqLFDDB8W#bnA)4lX_nu8P>7nBm#-ua1$|5sSl)5`fxarjK>1eh_PGW*&r=E&AJ{> z8VGENTMYIX9#3#EkkE}pAYvptgVd7YX@cKKkhME|_1*@(qo-N#zO|=GkHTAQ*N}|T zxIx#0sg%d#X*5mXKN1`0oWpTXfV@+8g9yP0cWfE@gg!hCa1x3nhvUg`hdCOnaVia= zZ;kDW#6tR(aB^gbrVD0D@Gd<#c$eNEh#+nns-T)bLxD&PNH7)9`y(_`d%|KUJd_L@ z0$RDr5)2GSX!-UKnlwVe>#ifQnZ%)7X&?x$h$5wRP$qcaylq z-`NHc@M?EcW1iUP$rD?xl3CUp-K%q@rrBbndkrof+&OT!D92HOyGuF#Y*^VVL>2DY z6=J@uw~31&lcr7HK}@S3;##Epg%d~jHM0-<3sSBUrubw|CDJe{~x#Ds!55TwBWTC{R3W8pyiMH%QWRROZ1$;_hgAEoMPc6 zlW@BQ$8hGT#)4b#|FsreKHn+}rv*2k3(Sf-3(jLoH4PS=pIHcRvEUaGpvdbr3lu*l zT8=cYStNX64#j-5;OY@w)^4-l4vYQ{3%d@9fwvr2yr+%*+OCaP9$O`Zwr|QPn|WRvd&Yi5j>62Ee*!$S;wSj; zSnfrM<~6?K3{7HSjg zUVw4tZPaBpL%(a^DXz6y~(5TD2PmPTCx?1`j@k%$E zRdv0tZpR`);7WjwMc#{PC}7m<(vDWt2^>HU4dy?ie6wYCu zAcD-prtJf;1$xRv|A)BxOq)BXJr&pEY-o*Zt*O^UdNeK3z?Td3WL%3aaouSuWG?)M zCfJN)?+RRS+TUWX5IOvHCc#mM!T z-P+i*woC&dtx=cVvG*Goz^;AgMkYN!@9@l=dZwQzrN8xk$6nf$YU!WJQFbgJm8}O4 zYmUcl+St?f^x@9dpN^dN2hUr3bl!3AB>cIixjSzG@~6-HSDw_KK40F|`YXrY z@4`Pi?|HG`ac1h&CGvXG#tz#$$ENLC>tV+jU1$A1hr4P|neofc@eXwAG)-8FM;q5L zYqR4%%v#^C`Xgk_MTH-ihW%ABSDk*j>u=cuX5Rmv%_io#A$^X%M) zYa#PK46w#H_U^>BMJM(UfNTVUjJLz{3Qj=QLDPcP1L}3=25LkuMgJm8c}r0<7XOrI z>W=NcKk#k$-RaYh{3Bg1nLbj@jW11=?BAi4&h$+7FTu_F^t0Io2r|2d{_T+AKAYVR z&gY-*H=fr!(=XhqrGKG~y>hm@x5G2#Ifgaub5NUE+O$%Lw`UhI;O_#_CMw?~)_UGp zii`acD&tl)f99l`KSy{TO{_a6_17FGw(^3OKC3-_`er(*n(e6e!Ubaq9K4FvnjF|? zPGV|eT-cuC8(sG|pFxjDnzm|V>mH$XT}!`eV0~M6Cw$H<-47@eM6gFI?}wso`;pmr zkMz%$>3&K6vIEeCm)yNJcXtI?qS6SIMxZnTr4cBNKxqU@Bk(&O0eaVrpHM_{z!?dp zi0uqSN5Xh>TyG&Ef`M2pZaB?X{!Y9x){C&IlSw-dO ze37lJJZBf*oA>v&({}N7dnGE4*$IO3#g|rgi(6svN0p=0lLY0ia)IQ zpDO+%#g8if?~4CY@pFpLQ%O9RD1MFN&5CnelU^>F9JzOEops25Ug`9#c6;0&XJbR- z>V`E9H;7_Yjg9#(@+6uj72*y0IVO)Nm#9os3$KdjK}+pi`sTb}l*?ad^;Ia?kWo}P z-YpWZ&COfZ--H6ST}VLybKzeEUV}`sgWh}#dQ;#8zP|w93$;a*k>)h@@#;d-L->D{ zDqGnPAS}i>oK)jL?>dP78pmAm{AUsTEbv8Tmx)O)3KpW@Y_9}l3@$P((ji1b+zQL;d3KTI4$-INMo7e}57D(IWVd2rt9xs^&Y#`D5TT4u$9I zxgvgCGETOh<-35-m1jKxoaRH}d4ELW%fxZttPd)EOW{-AeEj>s=NkW)i{KZE;8$QF zoXbyh5xfhyGd~n5MT+Q;7QtUCg1=V8|GPzSF%UNDRRVd^9ULAOlwa%)#-qpr*Rc-i z9h$yji?6Fg*V?xTz5Vw$`MNsS3;h=V<~E;S-@IW%PY3d>eQo{@y&#pmFn|3+3F&$` z6fgqFhu@6c<`5FwN0MP`)|%Vwm8r)_A*O6$GxzzUlG9CfGgDcmMdx|Y^-wada+6Y0+S?Uej(wKEgZ^ALVVY-awON%#e6 z;5nagNh$Ma6<-saZyq2{0{Il{ApI(?>Qhek$AbC0d zUe2-qL}Q!(F+yISFt)j%^S=xE*dL(?cXIcT5$FRq5+)zyQY40NTZV-jzm0_51F;eJ za3YRRHjLdCsBa{K&&VPn1>qBddLu9(r21eWIVjwr-7z#_awB2Z?8Jv3_);ZL&{2oq zQvzh5O2bh@xaI738{rhn^7?Zp;<%hX?ikL|`ba4~AgN5}e{Yhj1zw9yWCJ2$#O} zp_`{o`rtB}3JgVpFp0xE{<&q7rh80AfwggHD2(rKek(-#4cr>>EM`(I)@dTd&kF2$Kf`=wE-Ge{ zVYMH1nwre%w_+<5wtozjR{JNMrZ)4pt;Q(x`D6Gc6s-2&xzem*{#jO3oOG@9Duikr z>aw5x?z`U9Dr}Es=5IjFVX^1;EC)WOJ?*dP6m$E$59RO4o$AEU^+Uw%@_Cqt4pE=j zp1-G>qKz%6&Sbg&Y|nf_fjz&c$te4SIZLy^cFdO)*z>+-e}k#d*D1SUJC>Khp0Hf} zUbE&P+f(N(Y(FaQ7m diff --git a/examples/cpp/get-subgrids.cpp b/examples/cpp/get-subgrids.cpp new file mode 100644 index 00000000..f0255943 --- /dev/null +++ b/examples/cpp/get-subgrids.cpp @@ -0,0 +1,54 @@ +#include + +#include +#include +#include +#include +#include +#include +#include + +int main() { + std::string filename = "drell-yan-rap-ll.pineappl.lz4"; + + // read the grid from a file + auto* grid = pineappl_grid_read(filename.c_str()); + + // Determine the number of bins and the index of order and channel + std::size_t n_bins = pineappl_grid_bin_count(grid); + std::size_t order = 0; + std::size_t channel = 0; + + // Get the dimension of the subgrids + std::size_t subgrid_dim = pineappl_convolutions_len(grid) + 1; + + std::cout << std::right << std::setw(10) << "bin" << std::setw(10) << "sg idx" + << std::setw(16) << "sg value" << "\n"; + std::cout << std::right << std::setw(10) << "---" << std::setw(10) << "------" + << std::setw(16) << "------------" << "\n"; + + for (std::size_t b = 0; b < n_bins; ++b) { + std::vector subgrid_shape(subgrid_dim); + pineappl_subgrid_shape(grid, b, order, channel, subgrid_shape.data()); + + // Check if the subgrid is not empty + if (subgrid_shape[0] != 0) { + std::size_t flat_shape = std::accumulate(subgrid_shape.begin(), + subgrid_shape.end(), 1, std::multiplies()); + std::vector subgrid_array(flat_shape); + + pineappl_subgrid_array(grid, b, order, channel, subgrid_array.data()); + + // TODO: Unravel of the index & check against some reference values + for (std::size_t value = 0; value < subgrid_array.size(); ++value) { + if (subgrid_array[value] != 0) { + std::cout << std::right << std::setw(10) << b << std::setw(10) + << value << std::setw(16) << subgrid_array[value] << "\n"; + break; + } + } + } + } + + pineappl_grid_delete(grid); +} From bd68162dd7c88b9961ccbbf9641804cf2a0a2764 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Thu, 10 Apr 2025 23:13:13 +0200 Subject: [PATCH 06/15] Address comments from Code Review --- examples/cpp/get-subgrids.cpp | 6 +++--- pineappl_capi/src/lib.rs | 34 ++++++++++++++++++++++------------ 2 files changed, 25 insertions(+), 15 deletions(-) diff --git a/examples/cpp/get-subgrids.cpp b/examples/cpp/get-subgrids.cpp index f0255943..a69c76d1 100644 --- a/examples/cpp/get-subgrids.cpp +++ b/examples/cpp/get-subgrids.cpp @@ -20,7 +20,7 @@ int main() { std::size_t channel = 0; // Get the dimension of the subgrids - std::size_t subgrid_dim = pineappl_convolutions_len(grid) + 1; + std::size_t subgrid_dim = pineappl_grid_kinematics_len(grid); std::cout << std::right << std::setw(10) << "bin" << std::setw(10) << "sg idx" << std::setw(16) << "sg value" << "\n"; @@ -29,7 +29,7 @@ int main() { for (std::size_t b = 0; b < n_bins; ++b) { std::vector subgrid_shape(subgrid_dim); - pineappl_subgrid_shape(grid, b, order, channel, subgrid_shape.data()); + pineappl_grid_subgrid_shape(grid, b, order, channel, subgrid_shape.data()); // Check if the subgrid is not empty if (subgrid_shape[0] != 0) { @@ -37,7 +37,7 @@ int main() { subgrid_shape.end(), 1, std::multiplies()); std::vector subgrid_array(flat_shape); - pineappl_subgrid_array(grid, b, order, channel, subgrid_array.data()); + pineappl_grid_subgrid_array(grid, b, order, channel, subgrid_array.data()); // TODO: Unravel of the index & check against some reference values for (std::size_t value = 0; value < subgrid_array.size(); ++value) { diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 8aa0cf75..f950ea78 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -1903,21 +1903,32 @@ pub unsafe extern "C" fn pineappl_grid_convolve( /// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, /// this function is not safe to call. #[no_mangle] -pub unsafe extern "C" fn pineappl_convolutions_len(grid: *mut Grid) -> usize { +pub unsafe extern "C" fn pineappl_grid_convolutions_len(grid: *mut Grid) -> usize { let grid = unsafe { &mut *grid }; grid.convolutions().len() } +/// Get the number of different kinematics for a given Grid. +/// +/// # Safety +/// +/// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, +/// this function is not safe to call. +#[no_mangle] +pub unsafe extern "C" fn pineappl_grid_kinematics_len(grid: *mut Grid) -> usize { + let grid = unsafe { &mut *grid }; + grid.kinematics().len() +} + /// Get the shape of a subgrid for a given bin, channel, and order. /// /// # Safety /// /// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, /// this function is not safe to call. Additionally, the pointer that specifies the shape of the -/// subgrid must be an array whose size must be `(n+1)` where `n` represents the number of -/// convolutions as given by `pineappl_convolutions_len`. +/// subgrid has to be an array whose size must be as given by `pineappl_grid_kinematics_len`. #[no_mangle] -pub unsafe extern "C" fn pineappl_subgrid_shape( +pub unsafe extern "C" fn pineappl_grid_subgrid_shape( grid: *mut Grid, bin: usize, order: usize, @@ -1928,13 +1939,13 @@ pub unsafe extern "C" fn pineappl_subgrid_shape( let mut subgrid = grid.subgrids()[[order, bin, channel]].clone(); let subgrid_shape = if subgrid.is_empty() { - let subgrid_dim = grid.convolutions().len() + 1; + let subgrid_dim = grid.kinematics().len(); &vec![0; subgrid_dim] } else { subgrid.shape() }; - let shape = unsafe { slice::from_raw_parts_mut(shape, grid.convolutions().len() + 1) }; + let shape = unsafe { slice::from_raw_parts_mut(shape, grid.kinematics().len()) }; shape.copy_from_slice(subgrid_shape); } @@ -1945,9 +1956,9 @@ pub unsafe extern "C" fn pineappl_subgrid_shape( /// If `grid` does not point to a valid `Grid` object, for example when `grid` is the null pointer, /// this function is not safe to call. Additionally, the pointer that specifies the size of the subgrid /// when flattened must be an array; its size must be computed by multiplying the shape dimension as -/// given by `pineappl_subgrid_shape`. +/// given by `pineappl_grid_subgrid_shape`. #[no_mangle] -pub unsafe extern "C" fn pineappl_subgrid_array( +pub unsafe extern "C" fn pineappl_grid_subgrid_array( grid: *mut Grid, bin: usize, order: usize, @@ -1957,15 +1968,14 @@ pub unsafe extern "C" fn pineappl_subgrid_array( let grid = unsafe { &mut *grid }; let mut subgrid = grid.subgrids()[[order, bin, channel]].clone(); - let flattened_shape = subgrid.shape().iter().copied().product(); - let mut flattened_subgrid_array = vec![0.0; flattened_shape]; + let subgrid_array = + unsafe { slice::from_raw_parts_mut(subgrid_array, subgrid.shape().iter().product()) }; + let mut flattened_subgrid_array = vec![0.0; subgrid_array.len()]; for (index, value) in subgrid.indexed_iter() { let ravel_index = ravel_multi_index(index.as_slice(), subgrid.clone().shape()); flattened_subgrid_array[ravel_index] = value; } - let subgrid_array = unsafe { slice::from_raw_parts_mut(subgrid_array, flattened_shape) }; - subgrid_array.copy_from_slice(&flattened_subgrid_array); } From 770b76309084d064f6f876f3a6119fff2be0445c Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 11 Apr 2025 11:10:44 +0200 Subject: [PATCH 07/15] Add `unravel_index` and tests to example --- examples/cpp/get-subgrids.cpp | 53 +++++++++++++++++++++++++++++++---- examples/cpp/output | 52 +++++++++++++++++----------------- 2 files changed, 73 insertions(+), 32 deletions(-) diff --git a/examples/cpp/get-subgrids.cpp b/examples/cpp/get-subgrids.cpp index a69c76d1..b37f0db4 100644 --- a/examples/cpp/get-subgrids.cpp +++ b/examples/cpp/get-subgrids.cpp @@ -1,13 +1,44 @@ #include +#include #include #include #include #include +#include #include #include #include + +std::vector unravel_index(std::size_t flat_index, const std::vector& shape) { + std::size_t ndim = shape.size(); + std::vector coords(ndim); + + for (int i = ndim - 1; i >= 0; --i) { + coords[i] = flat_index % shape[i]; + flat_index /= shape[i]; + } + + return coords; +} + +std::string coords_to_string(const std::vector& coords) { + std::ostringstream osstream; + + osstream << "("; + for (std::size_t i = 0; i < coords.size(); ++i) { + osstream << coords[i]; + if (i != coords.size() - 1) { + osstream << ", "; + } + } + osstream << ")"; + + return osstream.str(); +} + + int main() { std::string filename = "drell-yan-rap-ll.pineappl.lz4"; @@ -23,9 +54,11 @@ int main() { std::size_t subgrid_dim = pineappl_grid_kinematics_len(grid); std::cout << std::right << std::setw(10) << "bin" << std::setw(10) << "sg idx" - << std::setw(16) << "sg value" << "\n"; + << std::setw(6 * subgrid_dim) << "sg coordinates" << std::setw(16) + << "sg value" << "\n"; std::cout << std::right << std::setw(10) << "---" << std::setw(10) << "------" - << std::setw(16) << "------------" << "\n"; + << std::setw(6 * subgrid_dim) << "--------------" << std::setw(16) + << "------------" << "\n"; for (std::size_t b = 0; b < n_bins; ++b) { std::vector subgrid_shape(subgrid_dim); @@ -39,11 +72,19 @@ int main() { pineappl_grid_subgrid_array(grid, b, order, channel, subgrid_array.data()); - // TODO: Unravel of the index & check against some reference values - for (std::size_t value = 0; value < subgrid_array.size(); ++value) { - if (subgrid_array[value] != 0) { + for (std::size_t index = 0; index < subgrid_array.size(); ++index) { + if (subgrid_array[index] != 0) { + // Unravel the index to recover the standard coordinates + std::vector coords = unravel_index(index, subgrid_shape); + std::cout << std::right << std::setw(10) << b << std::setw(10) - << value << std::setw(16) << subgrid_array[value] << "\n"; + << index << std::setw(6 * subgrid_dim) << coords_to_string(coords) + << std::setw(16) << subgrid_array[index] << "\n"; + + // Compare to some reference value + if (b==0 && index==41020) { + assert(subgrid_array[index] == -4.936156925096015e-07); + } break; } } diff --git a/examples/cpp/output b/examples/cpp/output index d7f24048..d5ad2d4d 100644 --- a/examples/cpp/output +++ b/examples/cpp/output @@ -146,32 +146,32 @@ idx left right dsig/dx dx 21 2.1 2.2 3.507294e-02 0.1 22 2.2 2.3 1.976542e-02 0.1 23 2.3 2.4 5.590552e-03 0.1 - bin sg idx sg value - --- ------ ------------ - 0 41020 -4.93616e-07 - 1 41020 -4.5499e-08 - 2 40971 -4.10008e-08 - 3 40971 -2.78597e-07 - 4 40971 -2.01924e-07 - 5 40922 -3.39322e-10 - 6 40922 -8.43358e-08 - 7 40922 -2.99247e-07 - 8 40922 -1.5117e-07 - 9 40872 -9.95013e-11 - 10 40873 -1.30578e-07 - 11 40873 -2.55869e-07 - 12 40823 -1.0823e-09 - 13 40823 -4.22686e-09 - 14 40824 -1.68149e-07 - 15 40774 -7.58702e-10 - 16 40774 -2.41887e-08 - 17 40774 -8.83139e-09 - 18 40725 -1.14932e-09 - 19 40725 -3.03857e-08 - 20 40725 -2.66414e-08 - 21 40675 -3.60243e-10 - 22 40676 -1.22785e-08 - 23 40626 -7.04493e-11 + bin sg idx sg coordinates sg value + --- ------ -------------- ------------ + 0 41020 (16, 20, 20) -4.93616e-07 + 1 41020 (16, 20, 20) -4.5499e-08 + 2 40971 (16, 19, 21) -4.10008e-08 + 3 40971 (16, 19, 21) -2.78597e-07 + 4 40971 (16, 19, 21) -2.01924e-07 + 5 40922 (16, 18, 22) -3.39322e-10 + 6 40922 (16, 18, 22) -8.43358e-08 + 7 40922 (16, 18, 22) -2.99247e-07 + 8 40922 (16, 18, 22) -1.5117e-07 + 9 40872 (16, 17, 22) -9.95013e-11 + 10 40873 (16, 17, 23) -1.30578e-07 + 11 40873 (16, 17, 23) -2.55869e-07 + 12 40823 (16, 16, 23) -1.0823e-09 + 13 40823 (16, 16, 23) -4.22686e-09 + 14 40824 (16, 16, 24) -1.68149e-07 + 15 40774 (16, 15, 24) -7.58702e-10 + 16 40774 (16, 15, 24) -2.41887e-08 + 17 40774 (16, 15, 24) -8.83139e-09 + 18 40725 (16, 14, 25) -1.14932e-09 + 19 40725 (16, 14, 25) -3.03857e-08 + 20 40725 (16, 14, 25) -2.66414e-08 + 21 40675 (16, 13, 25) -3.60243e-10 + 22 40676 (16, 13, 26) -1.22785e-08 + 23 40626 (16, 12, 26) -7.04493e-11 idx p-p c#0 l#0 p-p~ c#0 l# p-d c#0 l#0 p-d dx --- ------------ ----------- ------------- ------------ ------ 0 5.263109e-01 5.263109e-01 5.263109e-01 5.263109e-01 0.1 From 1d194a382688734410e2f8a14af95365b62566cd Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 11 Apr 2025 11:31:28 +0200 Subject: [PATCH 08/15] Fix `pineappl_channels_entry` --- examples/cpp/display-channels.cpp | 7 +++++-- pineappl_capi/src/lib.rs | 3 ++- 2 files changed, 7 insertions(+), 3 deletions(-) diff --git a/examples/cpp/display-channels.cpp b/examples/cpp/display-channels.cpp index c22da961..d3fd3f2e 100644 --- a/examples/cpp/display-channels.cpp +++ b/examples/cpp/display-channels.cpp @@ -22,6 +22,9 @@ int main(int argc, char* argv[]) { // read the grid from a file auto* grid = pineappl_grid_read(filename.c_str()); + // How many convolutions are there? + auto n_conv = pineappl_grid_convolutions_len(grid); + // extract all channels auto* channels = pineappl_grid_channels(grid); @@ -36,11 +39,11 @@ int main(int argc, char* argv[]) { auto combinations = pineappl_channels_combinations(channels, channel); std::vector factors(combinations); - std::vector pids(2 * combinations); + std::vector pids(n_conv * combinations); // read out the channel with index given by `channel`, writing the particle identifiers into // `pids` and the corresponding factors into `factors` - pineappl_channels_entry(channels, channel, pids.data(), factors.data()); + pineappl_channels_entry(channels, channel, n_conv, pids.data(), factors.data()); for (std::size_t combination = 0; combination != combinations; ++combination) { auto factor = factors.at(combination); diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index f950ea78..853e51d4 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -1773,12 +1773,13 @@ pub unsafe extern "C" fn pineappl_grid_split_channels(grid: *mut Grid) { pub unsafe extern "C" fn pineappl_channels_entry( channels: *const Channels, entry: usize, + nb_convolutions: usize, pdg_ids: *mut i32, factors: *mut f64, ) { let channels = unsafe { &*channels }; let entry = channels.0[entry].entry(); - let pdg_ids = unsafe { slice::from_raw_parts_mut(pdg_ids, 3 * entry.len()) }; + let pdg_ids = unsafe { slice::from_raw_parts_mut(pdg_ids, nb_convolutions * entry.len()) }; let factors = unsafe { slice::from_raw_parts_mut(factors, entry.len()) }; entry From 233b9a4430d5de4ac85132808b73895af5766023 Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 11 Apr 2025 11:32:26 +0200 Subject: [PATCH 09/15] Add empty newline for consistency --- pineappl_capi/src/lib.rs | 2 ++ 1 file changed, 2 insertions(+) diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 853e51d4..e532524e 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -1906,6 +1906,7 @@ pub unsafe extern "C" fn pineappl_grid_convolve( #[no_mangle] pub unsafe extern "C" fn pineappl_grid_convolutions_len(grid: *mut Grid) -> usize { let grid = unsafe { &mut *grid }; + grid.convolutions().len() } @@ -1918,6 +1919,7 @@ pub unsafe extern "C" fn pineappl_grid_convolutions_len(grid: *mut Grid) -> usiz #[no_mangle] pub unsafe extern "C" fn pineappl_grid_kinematics_len(grid: *mut Grid) -> usize { let grid = unsafe { &mut *grid }; + grid.kinematics().len() } From 4d2f0c9ba2186abeef3405501bbf3f8349dfd4fd Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 11 Apr 2025 15:43:06 +0200 Subject: [PATCH 10/15] Hit missing part in the coverage related to empty grid --- examples/cpp/advanced-filling.cpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/examples/cpp/advanced-filling.cpp b/examples/cpp/advanced-filling.cpp index 01a1fbfe..d85272e6 100644 --- a/examples/cpp/advanced-filling.cpp +++ b/examples/cpp/advanced-filling.cpp @@ -158,6 +158,13 @@ int main() { } //-----------------------------------------------------------------------// + // Check that the Grid contains an empty subgrid at (b, o, c)=(0, 0, 0) + auto subgrid_dim = pineappl_grid_kinematics_len(grid); + std::vector subgrid_shape(subgrid_dim); + std::vector zero_vector(subgrid_dim, 0); + pineappl_grid_subgrid_shape(grid, 0, 0, 0, subgrid_shape.data()); + assert(subgrid_shape == zero_vector); + pineappl_grid_write(grid, "advanced-filling.pineappl.lz4"); // release memory From 3a9f23b3a22a213b4f614e05386909f65051e06b Mon Sep 17 00:00:00 2001 From: Christopher Schwan Date: Fri, 11 Apr 2025 17:24:10 +0200 Subject: [PATCH 11/15] Remove `mut` from `Subgrid::shape` --- pineappl/src/subgrid.rs | 8 ++++---- pineappl_capi/src/lib.rs | 6 +++--- 2 files changed, 7 insertions(+), 7 deletions(-) diff --git a/pineappl/src/subgrid.rs b/pineappl/src/subgrid.rs index c115f1f9..799bb51a 100644 --- a/pineappl/src/subgrid.rs +++ b/pineappl/src/subgrid.rs @@ -33,7 +33,7 @@ impl Subgrid for EmptySubgridV1 { Vec::new() } - fn shape(&mut self) -> &[usize] { + fn shape(&self) -> &[usize] { panic!("EmptySubgridV1 doesn't have a shape"); } @@ -166,7 +166,7 @@ impl Subgrid for ImportSubgridV1 { Box::new(self.array.indexed_iter()) } - fn shape(&mut self) -> &[usize] { + fn shape(&self) -> &[usize] { self.array.shape() } @@ -275,7 +275,7 @@ impl Subgrid for InterpSubgridV1 { self.array.is_empty() } - fn shape(&mut self) -> &[usize] { + fn shape(&self) -> &[usize] { self.array.shape() } @@ -471,7 +471,7 @@ pub trait Subgrid { fn scale(&mut self, factor: f64); /// Return the shape of the subgrid - fn shape(&mut self) -> &[usize]; + fn shape(&self) -> &[usize]; /// Assume that the convolution functions for indices `a` and `b` for this grid are the same /// and use this to optimize the size of the grid. diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index e532524e..de3c275a 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -1939,7 +1939,7 @@ pub unsafe extern "C" fn pineappl_grid_subgrid_shape( shape: *mut usize, ) { let grid = unsafe { &mut *grid }; - let mut subgrid = grid.subgrids()[[order, bin, channel]].clone(); + let subgrid = &grid.subgrids()[[order, bin, channel]]; let subgrid_shape = if subgrid.is_empty() { let subgrid_dim = grid.kinematics().len(); @@ -1969,14 +1969,14 @@ pub unsafe extern "C" fn pineappl_grid_subgrid_array( subgrid_array: *mut f64, ) { let grid = unsafe { &mut *grid }; - let mut subgrid = grid.subgrids()[[order, bin, channel]].clone(); + let subgrid = &grid.subgrids()[[order, bin, channel]]; let subgrid_array = unsafe { slice::from_raw_parts_mut(subgrid_array, subgrid.shape().iter().product()) }; let mut flattened_subgrid_array = vec![0.0; subgrid_array.len()]; for (index, value) in subgrid.indexed_iter() { - let ravel_index = ravel_multi_index(index.as_slice(), subgrid.clone().shape()); + let ravel_index = ravel_multi_index(index.as_slice(), subgrid.shape()); flattened_subgrid_array[ravel_index] = value; } From a02beae383409569e78c4684d1939ae9af69883a Mon Sep 17 00:00:00 2001 From: Christopher Schwan Date: Fri, 11 Apr 2025 17:35:48 +0200 Subject: [PATCH 12/15] Use less `mut` and copy directly --- pineappl_capi/src/lib.rs | 31 ++++++++++++++++--------------- 1 file changed, 16 insertions(+), 15 deletions(-) diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index de3c275a..48a8824e 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -1932,24 +1932,24 @@ pub unsafe extern "C" fn pineappl_grid_kinematics_len(grid: *mut Grid) -> usize /// subgrid has to be an array whose size must be as given by `pineappl_grid_kinematics_len`. #[no_mangle] pub unsafe extern "C" fn pineappl_grid_subgrid_shape( - grid: *mut Grid, + grid: *const Grid, bin: usize, order: usize, channel: usize, shape: *mut usize, ) { - let grid = unsafe { &mut *grid }; + let grid = unsafe { &*grid }; let subgrid = &grid.subgrids()[[order, bin, channel]]; - let subgrid_shape = if subgrid.is_empty() { + // avoid calling `Subgrid::shape()` for empty grids, which may panic let subgrid_dim = grid.kinematics().len(); &vec![0; subgrid_dim] } else { subgrid.shape() }; - let shape = unsafe { slice::from_raw_parts_mut(shape, grid.kinematics().len()) }; - shape.copy_from_slice(subgrid_shape); + + shape.copy_from_slice(&subgrid_shape); } /// Get the subgrid for a given bin, channel, and order @@ -1962,23 +1962,24 @@ pub unsafe extern "C" fn pineappl_grid_subgrid_shape( /// given by `pineappl_grid_subgrid_shape`. #[no_mangle] pub unsafe extern "C" fn pineappl_grid_subgrid_array( - grid: *mut Grid, + grid: *const Grid, bin: usize, order: usize, channel: usize, subgrid_array: *mut f64, ) { - let grid = unsafe { &mut *grid }; + let grid = unsafe { &*grid }; let subgrid = &grid.subgrids()[[order, bin, channel]]; - let subgrid_array = - unsafe { slice::from_raw_parts_mut(subgrid_array, subgrid.shape().iter().product()) }; - let mut flattened_subgrid_array = vec![0.0; subgrid_array.len()]; + // avoid calling `Subgrid::shape()` for empty grids, which may panic + if !subgrid.is_empty() { + let shape = subgrid.shape(); + let subgrid_array = + unsafe { slice::from_raw_parts_mut(subgrid_array, shape.iter().product()) }; - for (index, value) in subgrid.indexed_iter() { - let ravel_index = ravel_multi_index(index.as_slice(), subgrid.shape()); - flattened_subgrid_array[ravel_index] = value; + for (index, value) in subgrid.indexed_iter() { + let ravel_index = ravel_multi_index(index.as_slice(), &shape); + subgrid_array[ravel_index] = value; + } } - - subgrid_array.copy_from_slice(&flattened_subgrid_array); } From 77979408bbb8bc5636973e51a70b705ee790cf3c Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 11 Apr 2025 18:03:13 +0200 Subject: [PATCH 13/15] Check that calling `pineappl_grid_subgrid_array` on an empty subgrid does not panic --- examples/cpp/advanced-filling.cpp | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/examples/cpp/advanced-filling.cpp b/examples/cpp/advanced-filling.cpp index d85272e6..caf5801e 100644 --- a/examples/cpp/advanced-filling.cpp +++ b/examples/cpp/advanced-filling.cpp @@ -6,6 +6,7 @@ #include #include #include +#include int main() { // Construct the channel object based on the nb of convolutions @@ -165,6 +166,10 @@ int main() { pineappl_grid_subgrid_shape(grid, 0, 0, 0, subgrid_shape.data()); assert(subgrid_shape == zero_vector); + // Check that a call to an empty subgrid does not panic + std::vector subgrid_array(0); + pineappl_grid_subgrid_array(grid, 0, 0, 0, subgrid_array.data()); + pineappl_grid_write(grid, "advanced-filling.pineappl.lz4"); // release memory From 9c08d4f8c80f9ee797f6168adafef663b5ca70bc Mon Sep 17 00:00:00 2001 From: Radonirinaunimi Date: Fri, 11 Apr 2025 20:12:11 +0200 Subject: [PATCH 14/15] Remove no longer needed include --- examples/cpp/advanced-filling.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/examples/cpp/advanced-filling.cpp b/examples/cpp/advanced-filling.cpp index caf5801e..8c7a4a72 100644 --- a/examples/cpp/advanced-filling.cpp +++ b/examples/cpp/advanced-filling.cpp @@ -6,7 +6,6 @@ #include #include #include -#include int main() { // Construct the channel object based on the nb of convolutions From 9d6f89761f17805ccc47bc42a9358e5cd8b4c27b Mon Sep 17 00:00:00 2001 From: Christopher Schwan Date: Sat, 12 Apr 2025 09:20:08 +0200 Subject: [PATCH 15/15] Undo interface change to `pineappl_channels_entry` --- examples/cpp/display-channels.cpp | 2 +- pineappl_capi/src/lib.rs | 18 ++++++++++-------- 2 files changed, 11 insertions(+), 9 deletions(-) diff --git a/examples/cpp/display-channels.cpp b/examples/cpp/display-channels.cpp index d3fd3f2e..636fd328 100644 --- a/examples/cpp/display-channels.cpp +++ b/examples/cpp/display-channels.cpp @@ -43,7 +43,7 @@ int main(int argc, char* argv[]) { // read out the channel with index given by `channel`, writing the particle identifiers into // `pids` and the corresponding factors into `factors` - pineappl_channels_entry(channels, channel, n_conv, pids.data(), factors.data()); + pineappl_channels_entry(channels, channel, pids.data(), factors.data()); for (std::size_t combination = 0; combination != combinations; ++combination) { auto factor = factors.at(combination); diff --git a/pineappl_capi/src/lib.rs b/pineappl_capi/src/lib.rs index 48a8824e..c28ee961 100644 --- a/pineappl_capi/src/lib.rs +++ b/pineappl_capi/src/lib.rs @@ -1760,26 +1760,28 @@ pub unsafe extern "C" fn pineappl_grid_split_channels(grid: *mut Grid) { grid.split_channels(); } -/// Similar to `pineappl_lumi_entry` but for luminosity channels that involve 3 or more partons, ie. -/// in the case of multiple convolutions. +/// Read out the channel with index `entry` of the given `channels`. The PDG ids and factors will +/// be copied into `pdg_ids` and `factors`. /// /// # Safety /// -/// The parameter `channels` must point to a valid `Channels` object created by `pineappl_channels_new` or -/// `pineappl_grid_channels`. The parameter `factors` must point to an array as long as the size -/// returned by `pineappl_channels_combinations` and `pdg_ids` must point to an array that is twice as -/// long. +/// The parameter `channels` must point to a valid [`Channels`] object created by +/// [`pineappl_channels_new`] or [`pineappl_grid_channels`]. The parameter `factors` must point to +/// an array as long as the result of [`pineappl_channels_combinations`] and `pdg_ids` must +/// point to an array as long as the result multiplied with the number of convolutions. #[no_mangle] pub unsafe extern "C" fn pineappl_channels_entry( channels: *const Channels, entry: usize, - nb_convolutions: usize, pdg_ids: *mut i32, factors: *mut f64, ) { let channels = unsafe { &*channels }; let entry = channels.0[entry].entry(); - let pdg_ids = unsafe { slice::from_raw_parts_mut(pdg_ids, nb_convolutions * entry.len()) }; + // if the channel has no entries we assume no convolutions, which is OK we don't copy anything + // in this case + let convolutions = entry.get(0).map_or(0, |x| x.0.len()); + let pdg_ids = unsafe { slice::from_raw_parts_mut(pdg_ids, convolutions * entry.len()) }; let factors = unsafe { slice::from_raw_parts_mut(factors, entry.len()) }; entry