I have completed my review. Here is my assessment.


Adversarial Review: Ridge GCV → cuSOLVER SPD solve (Phase 1)

VERDICT: fix-then-ship

The actual Phase‑1 deliverable — spd_solve, the cuSOLVER handle, the gating helper, the 3 call‑site swaps, the CMake link/compile‑def edits, and the equivalence test — is correct and faithful to the approved design. I could not find a correctness, ABI, or fold‑isolation bug in any of it.

There is one real problem: an unrelated, out‑of‑scope rewrite of the hot gemm() dispatch path is bundled into the same diff. It isn’t in the approved design, it touches a far broader blast radius than the isolated solve, and it’s not covered by the equivalence test. Split it out (or justify+test it separately) before merging; the Ridge work itself is shippable.


Findings

1. [SCOPE] gemm() was rewritten to use static thread_local ReusableDeviceBuffer scratch — not part of this design. cpp/src/core/cuda_dispatch.cpp:1559-1572 (and the dA/dB/dC deref changes at :1597-1603). The approved design (§3) scopes Phase 1 to spd_solve + the cuSOLVER handle + CMake + the 3 dual call‑site swaps + the test. This gemm change is a separate performance optimization (replacing per‑call DevicePtr with grow‑only per‑thread scratch to avoid cudaMalloc/cudaFree device‑sync stalls). It is plausibly correctdA/dB are fully overwritten by copy_h2d each call, dC is only read when beta != 0 (same contract as the original uninitialized DevicePtr), grow‑only ensure() is sound, and thread_local gives per‑thread isolation — but:

  • It modifies the primitive used by every CUDA GEMM in libn4m, including the K = Xc·Xcᵀ build for the dual design and all PLS moment builds — much wider than the isolated spd_solve.

  • It is not exercised by the equivalence test (which only calls spd_solve directly), and the design’s optional moment‑sweep parity sanity (§5) wasn’t added.

  • §4.2 of the design explicitly cautioned against introducing reused/persistent device buffers in Phase 1.

  • Fix: extract this into its own PR with its own benchmark + a run_moment_sweep parity check (N4M_CUDA_RIDGE_DISABLE=1 vs not, scores agree ≤1e‑8). If it must ship together, call it out explicitly and add the parity gate.

2. [STYLE] static thread_local lifetime footnote (only if finding #1 stays). cuda_dispatch.cpp:1564-1566. The scratch destructors run cudaFree at thread/process exit, where the CUDA runtime may already be torn down (cudaErrorCudartUnloading). The return is ignored, so it’s benign, and it mirrors the existing CublasState dtor pattern — but it’s a new lifetime that the per‑call DevicePtr didn’t have, and the buffers stay resident at high‑water mark for the process lifetime. Acceptable as a documented tradeoff; noting for completeness.

3. [OK] spd_solve correctness (cuda_dispatch.cpp:1156-1260). Verified against §3.1 and the §4 traps:

  • A uploaded as‑is, no transpose — correct (A symmetric ⇒ row‑major image == column‑major image; uplo immaterial, documented at :1207).

  • B packed column‑major Bcm[target*n+i]=B[i*q+target] (:1200-1204) and X unpacked inverse X[i*q+target]=Xcm[target*n+i] (:1242-1245) — both correct.

  • devInfo checked after both potrf (:1219-1228) and potrs (:1234-1238).

  • Return codes exactly per spec: potrf info>0 1 (not‑PD), info<0 2, potrs info!=0 2, and check_cusolver throws → caught by the bad_alloc/exception/... envelope → 2 (:1248-1259).

  • Workspace sized via cusolverDnDpotrf_bufferSize (:1210-1213), no assumed size.

  • Overflow/validity guards are thorough: null/zero, n,q INT_MAX, mul_overflows(n,n), mul_overflows(n,q), and byte‑size overflow (:1170-1189).

  • No leaks (RAII DevicePtr for dA/dB/dInfo/dWork, std::vector staging). Default‑stream ordering is correct — copy_d2h of info is a synchronous cudaMemcpy on the default stream, so it waits for potrf/potrs (no explicit sync needed).

4. [OK] Scope of the swaps. Exactly the 3 dual sites route to solve_dual_spd: predict_ridge_from_dual_design (sweep.cpp:529), ridge_dual_cross_heldout_sse (:577), fit_ridge_from_dual_design (:623). The k×k PᵀW inverse (sweep.cpp:721), the primal moment solve (:913), and the PLS‑prefix inverse (:1403) all remain on host solve_square_qr. ridge.cpp is untouched (confirmed via diff‑stat). No Phase‑2 eigh code was added.

5. [OK] Fallback + gating (sweep.cpp:435-469). solve_dual_spd falls through to host QR on s==1 and s==2 (the if (s==0) return; leaves both other codes falling out of the #if block to solve_square_qr). ridge_cuda_dual_enabled: n<256 false, N4M_CUDA_RIDGE_DISABLE honored, cuda_runtime_available() checked, #else false. The λ=0 generic case is handled correctly: potrf reports not‑PD → 1 → host QR, which itself rejects the singular R with N4M_ERR_NUMERICAL_FAILURE — identical to the pre‑change behavior. The [[maybe_unused]] on ridge_cuda_dual_enabled correctly prevents an unused‑function -Werror failure in the non‑CUDA build.

6. [OK] Fold‑clean. spd_solve uses only per‑call DevicePtrs and the passed‑in A/B; no persistent or cross‑fold device buffer. K/Y_work are fold‑local. No fold‑to‑fold state in the solve path. (Caveat: finding #1’s gemm scratch is reused, but it’s fully overwritten each call so it carries no numerical state — not a leakage vector, but it’s the kind of reuse §4.2 flagged.)

7. [OK] ABI. spd_solve is declared only in the internal cuda_dispatch.hpp, has no N4M_API, and the project compiles with -fvisibility=hidden (CMakeLists.txt:77-79, n4m_targets.cmake:111-112) plus the n4m_linux.map version script — so it stays out of the dynamic symbol table, exactly like the existing gemm/pls1_* internals. No cpp/abi/, cpp/include/, bindings/, or catalog/ changes (confirmed via diff‑stat). cuSOLVER handle lifecycle is correct: created after cuBLAS, available gated on both succeeding, cuBLAS rolled back (cublasDestroy_v2 + handle=nullptr) if cusolverDnCreate fails, and destroyed before cuBLAS in the dtor (cuda_dispatch.cpp:47-76).

8. [OK] The test (test_internal_linalg.cpp:304-371). It’s a real test: builds a deterministic SPD A = (GᵀG)/n + I (eigenvalues ≥ 1, well‑conditioned), solves via spd_solve and a host QR reference (qr_square_multi, byte‑for‑byte the same algorithm as the production solve_square_qr), and asserts max abs diff ≤ 1e-8. It self‑skips with a printed message when no GPU is present. The not‑PD path is covered with a diag(1,…,1,0) matrix (a clean zero pivot → potrf info=n>0), asserting spd_solve(...) == 1. The #if defined(N4M_USE_CUDA) guards (helper, test, and invocation) plus the cpp/tests/CMakeLists.txt:94-95 generator‑expression compile‑def correctly make the TU see the macro that n4m_c_static keeps PRIVATE. Note: it tests n=64 < 256 but calls spd_solve directly, bypassing the n<256 gate — correct, so the device path is actually exercised.

9. [STYLE] N4M_CUDA_RIDGE_DISABLE semantics (sweep.cpp:439-442). Any non‑empty value other than a leading '0' disables — so N4M_CUDA_RIDGE_DISABLE=false would also disable CUDA, diverging from the truthy_env idiom used elsewhere (cuda_dispatch.cpp:386). This matches the design’s §3.4 spec verbatim and is defensible for a kill‑switch, so it’s not a bug — just a minor inconsistency worth a one‑line comment if you care about env‑var uniformity.


Bottom line: the Ridge spd_solve change is clean, correct, ABI‑safe, fold‑isolated, and well‑tested. Pull the gemm() thread_local‑scratch optimization out of this diff (finding #1) and the Phase‑1 work ships as designed. The orchestrator’s GPU run should confirm the ≤1e-8 equivalence and the not‑PD ==1 assertions, which my static read indicates will pass.