Warp-aggregate the attractive-force atomics in t-SNE - #8565
Conversation
Both t-SNE attractive-force kernels apply one atomicAdd per non-zero to attr_forces[row], and RepulsionKernel applies one per body to the single Z_norm scalar. from_knn_symmetrize_matrix lays the symmetrized COO out row by row, so a warp normally spans one or two rows: all 32 lanes hit the same address and the adds serialize in L2. Sum each run of equal row indices inside the warp first and let the run's last lane issue one atomic per row, and fold Z_norm into one atomic per warp. The reduction is driven by lane position rather than by comparing row indices pairwise, so a row that appears in two separate runs within a warp is still counted once, and lanes past the end of the COO stay in the warp with a row index no edge can match rather than returning early. Measured on an RTX 5090 (sm_120, CUDA 13.2, driver 580.126.09) with the parameters cuml.TSNE configures by default, baseline and candidate runs interleaved, median of 3 rounds of 5: Barnes-Hut n= 20,000 604 ms -> 444 ms (1.36x) Barnes-Hut n= 50,000 760 ms -> 601 ms (1.27x) Barnes-Hut n=200,000 2105 ms -> 1551 ms (1.36x) FFT n= 20,000 392 ms -> 193 ms (2.03x) FFT n= 50,000 424 ms -> 268 ms (1.58x) FFT n=200,000 1699 ms -> 1185 ms (1.43x) attractive_kernel_bh alone goes from 225 us to 38 us per launch, and compute_Pij_x_Qij_kernel was 61% of unseeded FFT t-SNE GPU time before the change. Both kernels were already non-deterministic (atomicAdd ordering, atomicCAS tree construction), so no reproducibility contract changes. The seeded FFT path uses compute_Pij_x_Qij_deterministic_rows, which walks rows without atomics; it is untouched and stays bit-identical. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HMTsqoAxrKxVQR6RfcUZhz
|
No actionable comments were generated in the recent review. 🎉 ℹ️ Recent review info⚙️ Run configurationConfiguration used: Path: .coderabbit.yaml Review profile: CHILL Plan: Enterprise Run ID: 📒 Files selected for processing (2)
Included review availability: Your plan provides up to 12 included reviews per hour; 11 remain after this review. 📝 SummarySummary by CodeRabbit
WalkthroughThe t-SNE Barnes-Hut and FFT CUDA kernels now use warp-level reductions to aggregate normalization and force updates. The changes add CUDA utilities, preserve inactive lanes for warp operations, and reduce atomic updates for contiguous rows and warp participants. Changest-SNE CUDA warp aggregation
Estimated code review effort: 4 (Complex) | ~45 minutes Merge Risk: ⚪ Minimal · up to The CUDA optimization reduces atomic contention while preserving reported t-SNE output quality and deterministic seeded behavior. No merge-blocking risk remains. Suggested reviewers: 🚥 Pre-merge checks | ✅ 5✅ Passed checks (5 passed)
✨ Finishing Touches 💡 1🛠️ Fix failing CI checks 💡
🧪 Generate unit tests (beta)
Comment |
Problem
Both t-SNE attractive-force kernels apply one
atomicAddper non-zero toattr_forces[row]:BH::attractive_kernel_bh(Barnes-Hut)FFT::compute_Pij_x_Qij_kernel(FFT, unseeded path)raft::sparse::linalg::from_knn_symmetrize_matrixwrites the symmetrized COOinto per-row spans, so it comes out grouped by row. A warp therefore normally
spans one or two rows, all 32 lanes target the same address, and the adds
serialize in L2.
BH::RepulsionKernelhas the same problem in a sharper form:one
atomicAddper body onto the singleZ_normscalar, i.e. ~n serializedadds to one address per launch.
Profiling
TSNE_fiton an RTX 5090 (n=50,000, 1000 iterations, the parameterscuml.TSNEconfigures by default) putcompute_Pij_x_Qij_kernelat 61% ofall GPU time on the FFT path and
attractive_kernel_bhat 26% on theBarnes-Hut path.
Change
Sum each run of equal row indices inside the warp first and let the run's last
lane issue one atomic per row; fold
Z_norminto one atomic per warp.Two details worth reviewing:
(
run_start, derived from a ballot of run heads), not by comparing rowindices pairwise. A naive
other == itest would double-count a row thathappens to appear in two separate runs within one warp; this does not.
carrying a row index no real edge can match, so every warp-wide primitive
sees a full mask.
Results
RTX 5090 (sm_120), CUDA 13.2, driver 580.126.09,
-O3,CMAKE_CUDA_ARCHITECTURES=120-real.Baseline and candidate runs interleaved in the same process sequence, median of
3 rounds of 5 timed iterations; run-to-run spread was under 1%. Workload is 20
Gaussian clusters,
random_stateunset, with then_neighbors/learning-rate schedule
cuml.TSNEderives fromn(its defaultlearning_rate_method="adaptive"), not the rawTSNEParamsdefaults.At the kernel level,
attractive_kernel_bhgoes from 225 µs to 38 µs perlaunch (5.9x).
The mechanism is not architecture-specific — same-address
atomicAddisserialized on every supported architecture, and the intrinsics used
(
__shfl_up_sync,__ballot_sync,__clz) are available on all of them — butI only have sm_120 hardware, so the numbers above are the only measured ones.
The magnitude will differ elsewhere: this GPU's 96 MB L2 holds most of the COO
working set, so the post-change kernel runs largely out of L2 here.
Correctness
{Barnes-Hut, FFT} × 2 seeds: maximum delta 0.008, inside the baseline's own
seed-to-seed spread. No non-finite outputs.
(
atomicAddordering; Barnes-Hut also builds its tree withatomicCAS), sorandom_statenever made these paths reproducible. The seeded FFT path usescompute_Pij_x_Qij_deterministic_rows, which walks rows without atomics — itis untouched and its output is bit-identical before and after. The warp
reduction is a tree sum, so it is if anything slightly better conditioned
than the previous serialized adds.
compute-sanitizer:memcheckandsynccheckreport 0 errors for bothalgorithms.
racecheckreports exactly the same hazard counts as thebaseline build (Barnes-Hut 9, FFT 16 — all pre-existing warnings, 0 errors),
so no new hazards.
p, across Barnes-Hut,FFT and exact — finite KL divergence, no failures. Partial final warps
(NNZ not a multiple of 32) are exercised by every one of these.
pre-commit runis clean on both files.Not addressed here
FFT::compute_interpolated_indices(~9% of unseeded FFT GPU time) has the samecontention pattern, but its colliding keys are strided rather than contiguous,
so it needs
__match_any_syncrather than this run-based reduction. Left for afollow-up.
The ~18-line reduction idiom is duplicated between the two kernels, which live
in different headers. I kept it inline rather than adding a shared device helper
because that is the exact code I benchmarked and validated; happy to factor it
into
cpp/src/tsne/utils.cuh(which both already include) if reviewers prefer.🤖 Generated with Claude Code
https://claude.ai/code/session_01HMTsqoAxrKxVQR6RfcUZhz