[sshin23]

Adds bench/ordering/ — the GPU-hybrid nested-dissection ordering prototype developed alongside #126 (not part of the #126 squash) — with two performance changes on top of the version used in that work, and the end-to-end MadNLP measurements that motivated them.

Contents

• MetisExact.jl — bitwise-faithful Julia translation of METIS 5.1.0's METIS_NodeND path (all-default options); validated perm- and iperm-identical to libmetis on lap2d/lap3d and the case78484 KKT dumps up to n = 2.09M.
• GpuND.jl — KernelAbstractions device kernels for the coarsening (hashed-handshake HEM matching incl. 2-hop, device-resident contraction), stage-validated against the exact port.
• Hybrid.jl — nodend_hybrid(backend, A): device coarsening for blocks ≥ GPU_MIN_NVTXS, exact host initial separator / FM refinement / MMD below, task-parallel recursion.
• GpuRefine.jl, GpuBatch.jl — device refinement rounds (off by default: fill cost) and level-batched coarsening, kept as measured evidence: batched device coarsening of the sub-50k blocks is 4.2x slower than 16 host threads on the case78484 population, which is why those blocks deliberately stay on the CPU.
• gpu_nodend.jl loader and the nodend_smoke.jl gate.

Performance changes in this PR

1. Packed single-sort contraction (GpuND.gpu_contract): the sortperm over (crow<<32 | ccol) keys was ~90% of the device coarsening cost (72–77 ms of 85 ms at the 7.3M-edge top level). When the coarse ids fit 20 bits and the edge count 24 bits (guarded; the original sortperm path remains as fallback for larger graphs), one radix sort! of (crow<<44 | ccol<<24 | edge index) replaces sortperm + two gathers, recovering the weight through the embedded edge index. Output is element-identical; the coarsening chain over the 674k graph is 1.8x faster on a GV100.
2. Task-parallel sub-threshold recursion (Hybrid.mlnd!): previously each top subtree ran one serial MlevelNestedDissection! task, so the largest subtree set the wall. The recursion now spawns above SPAWN_MIN (own CtrlT copy, task-local RNG, disjoint order ranges — the same invariants the existing top-level spawns rely on). Ordering 6.3 s → 4.3 s at fill 1.002x of exact METIS.
3. Opt-in parallel top (Hybrid.PARALLEL_TOP + TASK_SETUP): overlaps the host work (initial separator, FM, split, downloads) of the ≥ threshold sibling bisections while the device coarsening serializes per device. Off by default. On a multi-GPU host the TASK_SETUP hook must pin every task to one device — unpinned tasks default to device 0 while the caller may be on another device, and that cross-device mix was observed to corrupt the ordering (caught by the Int32 overflow guard in the symbolic maps). The module header documents the recipe.

Measurements (78k ACOPF condensed KKT, n = 674,562, Quadro GV100, 16 Julia threads)

configuration	MadNLP init_time
algo4 (host METIS, wrapper default)	8.62 s
GPU ND + packed contraction	7.42 s (−14%)
GPU ND + packed contraction + parallel top	5.64 s (−35%)

All three runs: 101 IPM iterations, identical objective, linsolve flat (fill identical: 1.48e7 nnz_L vs METIS 1.476e7). Interleaved same-load comparison: the hybrid ordering runs at 0.63x of host METIS per repetition. No runtime data is cached; a fresh analysis runs every time.

Requirement: the recursion gains need a multi-threaded Julia session (the measurements used JULIA_NUM_THREADS=16). Single-threaded, only the contraction gain remains and the Julia host recursion is slower than libmetis C — single-threaded sessions should keep algo4.

🤖 Generated with Claude Code

[claude]

nit (non-blocking): this flips the module-global MetisExact.RNG_TASK[] and never restores it, and MetisExact.NodeND does not force it back to false. In a session with JULIA_NUM_THREADS > 1, any nodend_exact call made after a nodend_hybrid call silently runs on the Julia task RNG, so the "bitwise-identical to libmetis" oracle is no longer exact (the smoke script only avoids this because it calls the exact port first). Fix: either save/restore the previous value in a try … finally here, or have NodeND set RNG_TASK[] = false for its own run (it has no spawned tasks).

[claude]

nit (non-blocking): Main.GpuRefine (also lines 125 and 127) ties this module to being included at top level of Main; the sibling modules are reached with using ..MetisExact / using ..GpuND, and gpu_nodend.jl already includes GpuRefine.jl before Hybrid.jl, so using ..GpuRefine at the top of the module would work the same way and keep the loader order the only requirement.

[claude]

nit (non-blocking): with the default grid = 60 the graph has n = 3600 < Hybrid.GPU_MIN_NVTXS[] = 50_000, so node_bisection! takes the MlevelNodeBisectionL2! host branch everywhere and the default smoke run never executes GpuND.gpu_match / gpu_contract (nor the packed-sort path this PR adds); it only gates the task-parallel host recursion. Line 49 already prints "below threshold" for the refine path, but the coarsening path is silently skipped too. Suggest a default grid ≥ 224 (n ≥ 50_176), or temporarily lowering Hybrid.GPU_MIN_NVTXS[] (and REFINE_HOST_BELOW[]) inside the smoke so the device kernels are exercised on KA.CPU() at any grid size.

[claude]

nit (non-blocking): view(segid, m:m) (and line 321 in the fallback) is a BoundsError when m == 0. The coarsening loop in Hybrid.coarsen_gpu! cannot reach this on later levels (the nedges > nvtxs ÷ 2 test breaks first), but the root call can: a graph with ≥ GPU_MIN_NVTXS vertices and no off-diagonal entries (nedges == 0) goes straight to gpu_contract, where METIS would take the Match_RM/RandomBisection route. A m == 0 && return (cn, …, empty arrays, cmap) guard (or routing nedges == 0 to the host port in node_bisection!) closes it.

[claude]

VERDICT: APPROVE (head cb21b8d)

Scope. This is not a TNN task PR: it adds the GPU-hybrid NodeND prototype under bench/ordering/ only (7 new files, 3675 lines). src/, ext/, test/, PLAN.md, TASKS.md, Project.toml and Manifest.toml are untouched, so the task-section/Report checks do not apply and no public name, parameter string or phase string is affected. bench/Project.toml already carries Metis; the new code reaches KernelAbstractions and Atomix through SparseDirectSolver.*, and Printf/Libdl/Random are stdlibs, so no new dependency is introduced.

What I checked.

• GpuND.gpu_contract packed path: the (crow << 44) | (ccol << 24) | edge layout fits under the cn < 2^20, m < 2^24 guard; the typemax(UInt64) self-edge sentinel cannot collide with a real key (a non-self edge never has both ids equal to 2^20-1); segment boundaries compare >> 24, the weight is recovered through the embedded edge index, and the rowcount[crow + 2] / cumsum offsetting yields a cn+1 row pointer with cxadj[1] = 0. Representative scan, cmap, cvwgt are consistent with the host CreateCoarseGraph contract; the only atomics are integer adds (deterministic) in a bench prototype, not in library assembly/extend-add.
• gpu_match: mutual-proposal resolve is symmetric and atomic-free; the maxvwgt cap matches METIS SHEM semantics (heavy vertices never propose and end as singletons); the 2-hop stage mirrors Match_2HopAny gating (>10% unmatched, degree classes 2/3/any) and, like METIS, does not apply the weight cap.
• Hybrid.coarsen_gpu! / node_bisection! follow CoarsenGraph / MlevelNodeBisectionL1 control flow (maxvwgt = 1.5·tvwgt/CoarsenTo, the COARSEN_FRACTION / nedges > nvtxs ÷ 2 break, niparts ÷ 2 on early stop); finer/coarser/cmap links are set so the host Refine2WayNode! walks the device-built chain. Replacing L2 by an L1-shaped bisection above the threshold is a documented quality trade (1.003–1.015 fill).
• Task-parallel recursion (mlnd!, dissect!): spawned tasks get their own CtrlT copy (all 15 fields, pijbm copied), a reseeded task-local RNG keyed on (nvtxs, lastvtx), and disjoint order ranges; MetisExact module state touched from tasks is read-only (RNG_TASK, DBG, GENMMD_PTR) or locked (prof!), and genmmd is reentrant. The spawn sequence is a deterministic function of the graph, so the ordering is reproducible for a fixed thread count. The multi-GPU pinning hazard of PARALLEL_TOP is documented at the knob.
• nodend_hybrid mirrors METIS_NodeND (compression, nseps bump, uncompress, 0→1-based return); nodend_smoke.jl uses the existing SymmetricPattern(rowptr, colval, n, structure; view, index) and evaluate_ordering API correctly.
• Int32 device casts go through Int32.(…), which throws InexactError on overflow rather than wrapping.

Not verified here. The runner lets me neither execute Julia nor read the workflow logs, so I could not run nodend_smoke.jl or confirm the perm-identity and timing claims of the description; CI is pending on this head (fix shows cancelled, which is the pipeline job, not a test). The claims are consistent with the code.

Non-blocking findings (inline).

1. Hybrid.jl:298 leaves MetisExact.RNG_TASK[] set, so nodend_exact after a hybrid call in a threaded session is no longer bitwise libmetis.
2. Hybrid.jl:123 uses Main.GpuRefine instead of using ..GpuRefine.
3. nodend_smoke.jl default grid (n = 3600) is below GPU_MIN_NVTXS, so the default run never executes the GpuND kernels.
4. GpuND.gpu_contract view(segid, m:m) errors for m == 0 (a ≥ 50k-vertex graph without off-diagonal entries).
5. bench/README.md file table does not list bench/ordering/.

No blocking findings.

[github-actions]

CI-FIX round 1 of 3: https://github.com/exanauts/SparseDirectSolver.jl/actions/runs/38343843966 failed on cb21b8d; Claude is diagnosing it.

[claude]

CI-FIX round 1 — no code change.

Root cause: the self-hosted CUDA runner, not this PR. In run 38343843966, only test-gpu (… cuda) failed. test-github-cpuonly and test-gpu (… amdgpu) passed. Every CUDA error is Out of GPU memory (CUDA error: out of memory (code 2, ERROR_OUT_OF_MEMORY)), thrown from cuDevicePrimaryCtxRetain → CuContext → create_context, along with many Failed to load compiled kernel image messages. The device had no memory left even to create a CUDA context in the test workers, so the GPU was already full, most likely from another job or a leftover process on the shared cuda machine. The failures cover unrelated files (test_api, test_schur, test_numeric_*, test_symbolic_device, ported cuDSS tests…), which fits an environment problem.

This PR (cb21b8d) adds only files under bench/ordering/. Nothing in src/, ext/ or test/ includes them, and main at c5e2d86 passed the same CUDA job a few hours earlier.

Fix: none in the code. Please free the GPU on the cuda runner (check nvidia-smi for stale processes) and re-run the failed job: gh run rerun 38343843966 --failed.

🤖 Generated with Claude Code