REVIEW 2 major objections 6 minor 42 references
A GPU graph engine makes BAO-scale 3- and 4-point galaxy clustering, including parity-odd 4-point signals, practical on a single consumer card.
Reviewed by Pith at T0; open to challenge. T0 means a machine referee read the full paper against a public rubric. the ladder, T0–T4 →
T0 review · grok-4.5
2026-07-11 02:09 UTC pith:76HYUZ3K
load-bearing objection Solid engineering v2 of GRAMSCI: merge-walk + OpenACC GPU with out-of-core tiling delivers measured 2.6–9× speedups and usable connected/parity 4pCF on DESI-scale graphs, with careful validation. the 2 major comments →
Fast Graph-based Higher-Order Clustering Statistics on the GPU
The pith
A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.
Core claim
The complete GRAMSCI query engine for the 3-point, 4-point, and parity-decomposed 4-point correlation functions can be ported to a single consumer GPU via OpenACC, delivering measured speedups of 2.6 imes (3pCF) to 9 imes (4pCF) over a 64-thread CPU node on identical binned estimators, while an automatic window-tiling scheme processes graphs larger than device memory (demonstrated on a 9 imes10^9-edge DESI BAO graph on a 24 GB card at ~20% overhead) and reproduces CPU counts to floating-point precision.
What carries the argument
The sorted adjacency graph plus merge-walk (and three-pointer walk for tetrahedra): hub and neighbor lists are already ordered by node index, so their intersection is a single linear two- or three-pointer pass of cost O(m) rather than repeated binary searches; the same CSR layout is mapped to OpenACC gangs (one hub per block) with vector lanes walking a shared list, enabling the GPU kernels and out-of-core tiling.
Load-bearing premise
That the GPU mapping rules (one block per hub, shared-list vector lanes, gang-slot scratch, and atomic slot-strided accumulation) stay free of race or ordering bugs on the heavy high-multiplicity tails of real clustered, FKP-weighted survey graphs where four-point cost is dominated by rare dense hubs.
What would settle it
Run the GPU and CPU 3pCF/4pCF/parity kernels on the same production DESI-scale FKP-weighted graph (or a controlled high-multiplicity synthetic with known analytic counts) and check whether every binned configuration still agrees to floating-point tolerance after out-of-core tiling; any systematic discrepancy falsifies the claim of exact, race-free acceleration.
If this is right
- Connected 4-point measurements across the full BAO range become feasible for DESI-scale catalogues, including the thousands of mocks needed for covariance estimation.
- A configuration-space parity-odd 4-point estimator with deterministic treatment of zero-volume tetrahedra is available as an independent cross-check of harmonic-basis parity searches.
- Joint 2+3+4-point BAO analyses that bin anisotropic multipoles can extract extra signal-to-noise from the same galaxy sample.
- Graphs exceeding device memory (billions of edges) can be measured on a single consumer card with only modest PCIe overhead, removing the need for multi-node CPU clusters for many survey analyses.
Where Pith is reading between the lines
- Because the tile decomposition partitions work by hub window with no inter-tile communication, a multi-GPU extension that assigns one hub window per device is a natural next step for even larger surveys.
- The same CSR-plus-merge-walk structure is geometry-agnostic and could be reused for higher-order clustering of any discrete point set (halos, voids, radio sources) once the graph is built.
- If the connected 4-point residual remains consistent with mocks while the parity-odd channel is non-zero, the code supplies a clean separation between late-time non-Gaussianity and a true parity-violating signal.
Editorial analysis
A structured set of objections, weighed in public.
Referee Report
Summary. This paper presents GRAMSCI Version 2, a substantial update to the 2019 graph-database algorithm for configuration-space N-point correlation functions. The main contributions are: (i) a merge-walk sorted-set intersection that reduces the 3pCF/4pCF inner loops from O(m log m) to O(m); (ii) a parity-decomposed 4pCF with deterministic handling of zero-volume tetrahedra; (iii) internal estimation of the disconnected 4pCF on the same graph to return the connected 4pCF; (iv) a NumPy Python interface; and (v) an OpenACC GPU port of the full query engine with out-of-core window tiling. Measured speedups on identical estimators are 2.6× (3pCF) to 9× (4pCF) versus a 64-thread CPU node, with a 9×10^9-edge BAO-scale DESI 3pCF run on a 24 GB card at ~20% overhead. Correctness is checked by exact CPU–GPU count agreement (max relative difference 7×10^{-7}), analytic chiral-tetrahedron injection (0.3% recovery), a Quijote parity-null (χ²/dof ≈ 1.10), and DESI DR1 LRG 3pCF/connected-4pCF consistency with the EZmock ensemble.
Significance. If the reported performance and correctness hold, the work is a practically useful methods contribution for large-scale structure. Configuration-space connected and parity-odd 4pCFs at BAO scales have been expensive; the combination of merge-walk, internal disconnected subtraction, and a single-consumer-GPU out-of-core engine makes DESI-scale covariance estimation and parity cross-checks more routine. Strengths that should be credited explicitly include: public code, exact binned-count agreement between CPU and GPU (including chunked production-scale DESI), analytic injection tests for the parity-odd channel, a field-level Quijote null, and end-to-end application to DESI DR1 LRG versus EZmocks. The paper deliberately avoids cross-code timing claims and reports wall-clock measurements on fixed hardware and identical estimators, which is the right standard for an engineering result.
major comments (2)
- §5.2 notes that survey 4pCF cost is dominated by the high-multiplicity tail (⟨m³⟩ ≫ ⟨m⟩³) and that NGC at r_max=80 h^{-1} Mpc takes 2.6 h, an order of magnitude beyond the uniform-catalogue scaling of §4.3. Table 2 and Figure 2 report GPU speedups primarily on the quasi-uniform test catalogue (and on DESI 3pCF). For the central 4pCF claim, a short production-scale GPU timing (or at least a clear statement of whether the 9× figure still holds under FKP-weighted clustered tails) would make the headline speedup more representative of the science use case the paper itself emphasizes.
- §3.3 and Eq. (4)–(5): the disconnected subtraction uses an internally computed ξ on the same graph and bins, which is the right design. The manuscript should state more explicitly whether the same Landy–Szalay-style signed-weight estimator and RR normalization are used for every pair partition, and whether any edge-effect or bin-center approximation enters when evaluating ξ(r_ab) for the six edges of a binned tetrahedron. A one-sentence clarification would remove residual ambiguity about self-consistency of ζ_conn.
minor comments (6)
- Figure 2 caption and §4.3: the fixed ~0.7 s GPU offset is important for interpreting the crossover scales; consider marking the transfer/launch floor explicitly on the plot or in the legend.
- §3.2: the |V|<10^{-9} threshold for degenerate tetrahedra is free but well-motivated; a brief note on sensitivity (or that the even channel is unaffected) would help readers who re-implement the estimator.
- §3.5.3 / Table 2: the out-of-core DESI 3pCF (W=5, ~20% overhead) is a strong result; stating the window-size selection rule or environment-variable override more prominently would aid reproducibility.
- Listing 1 and §3.4: the Python interface is welcome; a one-line note on whether temporary-directory I/O is the only path or whether an in-memory/ctypes path is planned would set expectations.
- §4.7: expected A100/H100 scaling is framed as bandwidth-based expectation, not measurement — keep that distinction clear in the abstract/conclusions if those numbers are quoted elsewhere.
- Minor typography: title and running headers show “F ast” (space in Fast); arXiv abstract has “Sabiu et.al” (missing space/period consistency).
Circularity Check
No significant circularity: engineering measurements and self-consistent estimators, not definitional or fitted predictions.
full rationale
This is a methods/engineering paper whose central claims are wall-clock speedups of identical estimators (merge-walk vs binary search; OpenACC GPU vs 64-thread CPU), out-of-core tiling overhead, and exact reproduction of binned counts. The merge-walk is a standard two/three-pointer intersection on already-sorted adjacency lists (Eq. 1 invariant); its O(m) cost is algorithmic, not definitional of the measured times. The connected 4pCF (Eq. 5) subtracts the Wick disconnected term built from an internally computed 2pCF on the same graph and bins—an intentional self-consistency choice, not a circular prediction. Parity even/odd channels are defined by the sign of the signed volume V (Eq. 3) with a deterministic zero-weight rule for degenerates; they are validated by analytic chiral-tetrahedron injection and a Quijote field-level null, not by construction. DESI applications compare measured 3pCF/connected-4pCF to an external EZmock ensemble; no cosmological parameter is fitted and re-predicted. Self-citation of Sabiu et al. 2019 is ordinary lineage for the base graph construction and is not load-bearing for the new kernels, speedups, or science results. No uniqueness theorem, ansatz smuggling, or renaming of a known empirical pattern appears. The derivation chain is therefore self-contained against external catalogues and analytic tests.
Axiom & Free-Parameter Ledger
free parameters (3)
- degenerate tetrahedron volume threshold |V|<10^{-9}
- out-of-core window count W / window size
- direction-pixel quantization for parity
axioms (5)
- standard math Sorted adjacency lists by node index remain a valid invariant after kd-tree neighbor search, enabling two- and three-pointer merge-walks.
- domain assumption Landy–Szalay / Szapudi–Szalay-style signed data-minus-random weights yield unbiased NpCF estimators when data weights sum to +1 and randoms to −1.
- domain assumption Disconnected 4pCF is exactly the three Wick pairings of ξ products and is parity-even, so subtracting it isolates the connected (and leaves the odd channel untouched).
- ad hoc to paper OpenACC gang/vector mapping with one hub per gang and vector lanes on a shared adjacency list produces correct, race-free counts when combined with atomic slot-strided accumulation.
- domain assumption EZmock and Quijote catalogues are adequate external references for consistency and parity-null tests of the DESI LRG measurement pipeline.
read the original abstract
We present a significant update to GRAMSCI (GRAph Made Statistics for Cosmological Information; Sabiu et.al 2019), an algorithm for the fast computation of the general $N$-point spatial correlation function of any discrete point set embedded in $\mathbb{R}^n$. Utilizing the concepts of kd-trees and graph databases, we count all possible $N$-tuples in binned configurations within a given length scale. In this {\em Version 2 update} we describe several additions to the original code. We replace the binary-search inner loop, which cost $O(m\log m)$ per hub--spoke pair, where $m$ is the mean neighbor count, with a merge-walk algorithm that reduces the inner loop to $O(m)$. We implement a parity-decomposed 4pCF that separates the signal into even and odd channels, enabling direct tests of parity violation in the galaxy distribution. We estimate the disconnected 4pCF internally on the same graph to return the connected 4pCF. We provide a Python interface so the Fortran engine can be called directly from NumPy arrays. Finally, and principally, we present a GPU port of the full query engine (OpenACC): the 3pCF, 4pCF, and parity-decomposed 4pCF kernels run on a single consumer GPU with measured speedups of $2.6\times$ (3pCF) to $9\times$ (4pCF) over a 64-thread CPU node, and an out-of-core tiling scheme allows graphs far exceeding device memory. We measure a $9\times10^9$-edge BAO-scale 3pCF on a 24\,GB card with ${\sim}20\%$ overhead. We validate the code against its CPU reference, against analytic injection tests, and demonstrate BAO-scale applications on the DESI DR1 LRG sample compared against the EZmock ensemble.
Figures
Reference graph
Works this paper leans on
-
[1]
G., Aguilar , J., Ahlen , S., et al
Adame , A. G., Aguilar , J., Ahlen , S., et al. 2025, title DESI 2024 VI: cosmological constraints from the measurements of baryon acoustic oscillations , , 2025, 021, 10.1088/1475-7516/2025/02/021
-
[2]
Alam , S., Ata , M., Bailey , S., et al. 2017, title The clustering of galaxies in the completed SDSS-III Baryon Oscillation Spectroscopic Survey: cosmological analysis of the DR12 galaxy sample , , 470, 2617, 10.1093/mnras/stx721
-
[3]
Anderson , L., Aubourg , \'E ., Bailey , S., et al. 2014, title The clustering of galaxies in the SDSS-III Baryon Oscillation Spectroscopic Survey: baryon acoustic oscillations in the Data Releases 10 and 11 Galaxy samples , , 441, 24, 10.1093/mnras/stu523
-
[4]
Astropy Collaboration , Robitaille , T. P., Tollerud , E. J., et al. 2013, title Astropy: A community Python package for astronomy , , 558, A33, 10.1051/0004-6361/201322068
-
[5]
Astropy Collaboration , Price-Whelan , A. M., Sip o cz , B. M., et al. 2018, title The Astropy Project: Building an Open-science Project and Status of the v2.0 Core Package , , 156, 123, 10.3847/1538-3881/aabc4f
-
[6]
Astropy Collaboration , Price-Whelan , A. M., Lim , P. L., et al. 2022, title The Astropy Project: Sustaining and Growing a Community-oriented Open-source Project and the Latest Major Release (v5.0) of the Core Package , , 935, 167, 10.3847/1538-4357/ac7c74
work page internal anchor Pith review doi:10.3847/1538-4357/ac7c74 2022
-
[7]
Cahn , R. N., Slepian , Z., & Hou , J. 2023, title Test for Cosmological Parity Violation Using the 3D Distribution of Galaxies , , 130, 201002, 10.1103/PhysRevLett.130.201002
-
[8]
Chuang , C.-H., Kitaura , F.-S., Prada , F., Zhao , C., & Yepes , G. 2015, title EZmocks: extending the Zel'dovich approximation to generate mock galaxy catalogues with accurate clustering statistics , , 446, 2621, 10.1093/mnras/stu2301
-
[9]
Davis , M., & Peebles , P. J. E. 1983, title A survey of galaxy redshifts. V. The two-point position and velocity correlations. , , 267, 465, 10.1086/160884
-
[10]
DESI Collaboration , Abdul Karim , M., Adame , A. G., et al. 2026, title Data Release 1 of the Dark Energy Spectroscopic Instrument , , 171, 285, 10.3847/1538-3881/ae4c43
-
[11]
Eisenstein , D. J., Zehavi , I., Hogg , D. W., et al. 2005, title Detection of the Baryon Acoustic Peak in the Large-Scale Correlation Function of SDSS Luminous Red Galaxies , , 633, 560, 10.1086/466512
-
[12]
Frieman , J. A., & Gaztanaga , E. 1994, title The Three-Point Function as a Probe of Models for Large-Scale Structure , , 425, 392, 10.1086/173995
-
[13]
Fry , J. N., & Gaztanaga , E. 1993, title Biasing and Hierarchical Statistics in Large-Scale Structure , , 413, 447, 10.1086/173015
work page internal anchor Pith review doi:10.1086/173015 1993
-
[14]
2009, title Clustering of luminous red galaxies - III
Gazta \ n aga , E., Cabr \'e , A., Castander , F., Crocce , M., & Fosalba , P. 2009, title Clustering of luminous red galaxies - III. Baryon acoustic peak in the three-point correlation , , 399, 801, 10.1111/j.1365-2966.2009.15313.x
-
[15]
Gazta \ n aga , E., & Scoccimarro , R. 2005, title The three-point function in large-scale structure: redshift distortions and galaxy bias , , 361, 824, 10.1111/j.1365-2966.2005.09234.x
-
[16]
Gil-Mar \' n , H., Percival , W. J., Verde , L., et al. 2017, title The clustering of galaxies in the SDSS-III Baryon Oscillation Spectroscopic Survey: RSD measurement from the power spectrum and bispectrum of the DR12 BOSS galaxies , , 465, 1757, 10.1093/mnras/stw2679
-
[17]
Groth , E. J., & Peebles , P. J. E. 1977, title Statistical analysis of catalogs of extragalactic objects. VII. Two- and three-point correlation functions for the high-resolution Shane-Wirtanen catalog of galaxies. , , 217, 385, 10.1086/155588
-
[18]
Guo , H., Zheng , Z., Jing , Y. P., et al. 2015, title Modelling the redshift-space three-point correlation function in SDSS-III. , , 449, L95, 10.1093/mnrasl/slv020
-
[19]
Hamilton , A. J. S. 1993, title Toward Better Ways to Measure the Galaxy Correlation Function , , 417, 19, 10.1086/173288
-
[20]
Harris , C. R., Millman , K. J., van der Walt , S. J., et al. 2020, title Array programming with NumPy , , 585, 357, 10.1038/s41586-020-2649-2
-
[21]
Hou , J., Slepian , Z., & Cahn , R. N. 2023, title Measurement of parity-odd modes in the large-scale 4-point correlation function of Sloan Digital Sky Survey Baryon Oscillation Spectroscopic Survey twelfth data release CMASS and LOWZ galaxies , , 522, 5701, 10.1093/mnras/stad1062
-
[22]
Hunter , J. D. 2007, title Matplotlib: A 2D Graphics Environment , Computing in Science and Engineering, 9, 90, 10.1109/MCSE.2007.55
-
[23]
2004, title The skewness of the aperture mass statistic, MNRAS, 352, 338
Jarvis , M., Bernstein , G., & Jain , B. 2004, title The skewness of the aperture mass statistic, MNRAS, 352, 338
work page 2004
-
[24]
Kennel , M. B. 2004, title KDTREE 2: Fortran 95 and C++ software to efficiently search for near neighbors in a multi-dimensional Euclidean space , arXiv e-prints, physics/0408067, 10.48550/arXiv.physics/0408067
work page internal anchor Pith review Pith/arXiv arXiv doi:10.48550/arxiv.physics/0408067 2004
-
[25]
Landy , S. D., & Szalay , A. S. 1993, title Bias and Variance of Angular Correlation Functions , , 412, 64, 10.1086/172900
-
[26]
Peebles , P. J. E., & Hauser , M. G. 1974, title Statistical Analysis of Catalogs of Extragalactic Objects. III. The Shane-Wirtanen and Zwicky Catalogs , , 28, 19, 10.1086/190308
-
[27]
Philcox , O. H. E. 2022, title Probing parity violation with the four-point correlation function of BOSS galaxies , , 106, 063501, 10.1103/PhysRevD.106.063501
-
[28]
Philcox , O. H. E., Slepian , Z., Hou , J., et al. 2022, title ENCORE: an O (N _ g ^ 2 ) estimator for galaxy N-point correlation functions , , 509, 2457, 10.1093/mnras/stab3025
-
[29]
G., Hoyle , B., Kim , J., & Li , X.-D
Sabiu , C. G., Hoyle , B., Kim , J., & Li , X.-D. 2019, title Graph Database Solution for Higher-order Spatial Statistics in the Era of Big Data , , 242, 29, 10.3847/1538-4365/ab22b5
-
[30]
2006, title Cosmology and the bispectrum , , 74, 023522, 10.1103/PhysRevD.74.023522
Sefusatti , E., Crocce , M., Pueblas , S., & Scoccimarro , R. 2006, title Cosmology and the bispectrum , , 74, 023522, 10.1103/PhysRevD.74.023522
-
[31]
Sefusatti , E., & Komatsu , E. 2007, title Bispectrum of galaxies from high-redshift galaxy surveys: Primordial non-Gaussianity and nonlinear galaxy bias , , 76, 083004, 10.1103/PhysRevD.76.083004
-
[32]
Sinha , M., & Garrison , L. H. 2020, title CORRFUNC - a suite of blazing fast correlation functions on the CPU , , 491, 3022, 10.1093/mnras/stz3157
-
[33]
Slepian , Z., & Eisenstein , D. J. 2015, title Computing the three-point correlation function of galaxies in O(N\^2) time , , 454, 4142, 10.1093/mnras/stv2119
-
[34]
Slepian , Z., & Eisenstein , D. J. 2016, title Accelerating the two-point and three-point galaxy correlation functions using Fourier transforms , , 455, L31, 10.1093/mnrasl/slv133
-
[36]
Slepian , Z., Eisenstein , D. J., Brownstein , J. R., et al. 2017 b , title Detection of baryon acoustic oscillation features in the large-scale three-point correlation function of SDSS BOSS DR12 CMASS galaxies , , 469, 1738, 10.1093/mnras/stx488
-
[37]
Slepian , Z., Krolewski , A., Greco , A., et al. 2025, title Measurement of Parity-Violating Modes of the Dark Energy Spectroscopic Instrument (DESI) Year 1 Luminous Red Galaxies' 4-Point Correlation Function , arXiv e-prints, arXiv:2508.09133, 10.48550/arXiv.2508.09133
work page internal anchor Pith review Pith/arXiv arXiv doi:10.48550/arxiv.2508.09133 2025
-
[38]
Szapudi , I., & Szalay , A. S. 1998, title A New Class of Estimators for the N-Point Correlations , , 494, L41, 10.1086/311146
-
[39]
Totsuji , H., & Kihara , T. 1969, title The Correlation Function for the Distribution of Galaxies , , 21, 221, 10.1093/pasj/21.3.221
-
[40]
Verde , L., Heavens , A. F., Percival , W. J., et al. 2002, title The 2dF Galaxy Redshift Survey: the bias of galaxies and the density of the Universe , , 335, 432, 10.1046/j.1365-8711.2002.05620.x
-
[41]
2020, title The Quijote Simulations , , 250, 2, 10.3847/1538-4365/ab9d82
Villaescusa-Navarro , F., Hahn , C., Massara , E., et al. 2020, title The Quijote Simulations , , 250, 2, 10.3847/1538-4365/ab9d82
-
[42]
SciPy 1.0: Fundamental Algorithms for Scientific Computing in Python
Virtanen , P., Gommers , R., Oliphant , T. E., et al. 2020, title SciPy 1.0: fundamental algorithms for scientific computing in Python , Nature Medicine, 17, 261, 10.1038/s41592-019-0686-2
-
[43]
Zehavi , I., Zheng , Z., Weinberg , D. H., et al. 2011, title Galaxy Clustering in the Completed SDSS Redshift Survey: The Dependence on Color and Luminosity , , 736, 59, 10.1088/0004-637X/736/1/59
discussion (0)
Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.