You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
CostMatrix uses alignas(64) storage and pads rows to BLOCK_SIZE = 16, so vector loads stay within cache lines. A transposed copy gives sequential column access.
Bit counting uses hardware POPCNT (TreeTools::count_bits, plus inline asm in hmi.cpp and spr_lookup.cpp).
LAPJV vectorisation was analysed in .AGENTS/memory/LAP_profiling.md. The reduction-transfer loop doesn't auto-vectorise because int_fast64_t has no vector type and the loop nest blocks it. The Dijkstra phase is gather-bound through col_list.
SIMD was ruled out for reduce_trees (pointer-chasing, branchy code).
Outstanding:
Try the scanned[] mask rewrite of the LAPJV Dijkstra loop so it auto-vectorises, and benchmark on a release build.
Try explicit or #pragma omp simd vectorisation of the split-overlap inner loops (a_row[bin] & b_row[bin] with popcount) for trees with many bins.
Any intrinsics must stay portable (CRAN: no -march=native). This probably means runtime dispatch or relying on compiler auto-vectorisation.
Both outstanding items tested; both discounted. Closing.
Setup: main @ 62897d4, R 4.3.3, GCC 13, Xeon with AVX2 (no AVX-512 VPOPCNTDQ), OMP_NUM_THREADS=1. Release flags (-O2, as on CRAN) unless stated. Median of 7 runs per workload, 2–3 interleaved rounds.
1. scanned[] mask in the LAPJV Dijkstra phase: slower
I tried two variants. Both give results identical to main on CID/MSD/PID/Nye and on 400 random LAPJV() matrices.
V1: a single pass over 0..dim-1 that skips scanned columns.
V2: a branchless, vectorisable update pass (d/pred via select), then a separate pass to collect new d == min columns.
workload
main
V1
V2
V2, -O3 -mavx2 (main, same flags)
CID, 60×50 tips
0.15 s
0.20
0.29
0.22 (0.14)
MSD, 60×50 tips
0.07
0.09
0.11
0.08 (0.06)
CID, 20×200 tips
0.47
0.70
1.06
0.81 (0.42)
MSD, 20×200 tips
0.15
0.19
0.30
0.20 (0.14)
CID, 5×1000 tips
2.3
3.0
5.7
4.9 (2.0)
LAPJV, 50 × 200² uniform
0.05
0.05
0.05
0.05 (0.05)
On uniform random matrices it is neutral, which matches the earlier note in LAP_profiling.md. On tree cost matrices it is 1.3–2.9× slower. Tree cost matrices produce long Dijkstra searches, and col_list shrinks the unscanned set with each step, so scanning all dim columns costs far more than vectorisation saves. This holds even with AVX2 enabled. Also, cost is int64, and the SSE2 baseline has no packed 64-bit compare, so under CRAN flags the loop would not vectorise anyway.
2. Vectorising the split-overlap popcount loop: ≤ 2 % ceiling
Callgrind on CID + PID for 1000-tip trees attributes about 75 % of instructions to lap() and about 0.9 % to the count_bits(a & b) lines. I microbenchmarked the exact loop shape (all split pairs × n_bins), then compared it with the end-to-end time for one pair of trees:
tips
CID per pair
scalar popcnt loop
4-accumulator unroll
AVX2 pshufb
1000
0.31 s
4.6 ms (1.5 %)
3.7 ms
2.9 ms
4000
19.6 s
840 ms (4.3 %)
775 ms
525 ms
The best case is AVX2 with runtime dispatch. It saves about 0.5 % at 1000 tips and about 1.6 % at 4000 tips, and less for PID, which does more work per cell. #pragma omp simd cannot help under baseline x86-64 flags because SSE2 has no vector popcount. The unroll alone saves at most 0.3–0.5 %. Neither justifies target attributes and dispatch code in a CRAN package.
Conclusion: the remaining cost is in LAP's comparison count, not in memory access or popcount throughput. No code change.
Activity
[-]Performance gain opportunity?[/-][+]Use SIMD to improve performance (?)[/+]ms609 commentedon Oct 12, 2026
Status: partly explored, no explicit SIMD yet.
Done:
CostMatrixusesalignas(64)storage and pads rows toBLOCK_SIZE = 16, so vector loads stay within cache lines. A transposed copy gives sequential column access.TreeTools::count_bits, plus inline asm inhmi.cppandspr_lookup.cpp)..AGENTS/memory/LAP_profiling.md. The reduction-transfer loop doesn't auto-vectorise becauseint_fast64_thas no vector type and the loop nest blocks it. The Dijkstra phase is gather-bound throughcol_list.reduce_trees(pointer-chasing, branchy code).Outstanding:
scanned[]mask rewrite of the LAPJV Dijkstra loop so it auto-vectorises, and benchmark on a release build.#pragma omp simdvectorisation of the split-overlap inner loops (a_row[bin] & b_row[bin]with popcount) for trees with many bins.-march=native). This probably means runtime dispatch or relying on compiler auto-vectorisation.Generated by Claude Code
ms609 commentedon Oct 12, 2026
Both outstanding items tested; both discounted. Closing.
Setup: main @ 62897d4, R 4.3.3, GCC 13, Xeon with AVX2 (no AVX-512 VPOPCNTDQ),
OMP_NUM_THREADS=1. Release flags (-O2, as on CRAN) unless stated. Median of 7 runs per workload, 2–3 interleaved rounds.1.
scanned[]mask in the LAPJV Dijkstra phase: slowerI tried two variants. Both give results identical to main on CID/MSD/PID/Nye and on 400 random
LAPJV()matrices.0..dim-1that skips scanned columns.d/predvia select), then a separate pass to collect newd == mincolumns.-O3 -mavx2(main, same flags)LAPJV, 50 × 200² uniformOn uniform random matrices it is neutral, which matches the earlier note in
LAP_profiling.md. On tree cost matrices it is 1.3–2.9× slower. Tree cost matrices produce long Dijkstra searches, andcol_listshrinks the unscanned set with each step, so scanning alldimcolumns costs far more than vectorisation saves. This holds even with AVX2 enabled. Also,costisint64, and the SSE2 baseline has no packed 64-bit compare, so under CRAN flags the loop would not vectorise anyway.2. Vectorising the split-overlap popcount loop: ≤ 2 % ceiling
Callgrind on CID + PID for 1000-tip trees attributes about 75 % of instructions to
lap()and about 0.9 % to thecount_bits(a & b)lines. I microbenchmarked the exact loop shape (all split pairs ×n_bins), then compared it with the end-to-end time for one pair of trees:popcntlooppshufbThe best case is AVX2 with runtime dispatch. It saves about 0.5 % at 1000 tips and about 1.6 % at 4000 tips, and less for PID, which does more work per cell.
#pragma omp simdcannot help under baseline x86-64 flags because SSE2 has no vector popcount. The unroll alone saves at most 0.3–0.5 %. Neither justifiestargetattributes and dispatch code in a CRAN package.Conclusion: the remaining cost is in LAP's comparison count, not in memory access or popcount throughput. No code change.
Generated by Claude Code