Add sparse Lanczos SVD solver - #3034
Conversation
aamijar
added
non-breaking
|
Hi @Intron7, I haven't looked at this PR closely yet, but one thing to think about is that we should try and reuse parts of the existing lanczos eigensolver as much as possible. They should have some things in common right? |
|
/ok to test eabeaa4 |
|
@aamijar the algorithm is different even is the name is similar. The two paths share the Lanczos name but the kernels are different algorithms: the eigensolver builds a symmetric tridiagonal via a one-vector recurrence and ritz-solves with syevd; the SVD builds a bidiagonal via Golub-Kahan with two coupled bases (A @ v and Aᵀ @ u) and ritz-solves with gesvdj. Restart is also different, the SVD path locks converged singular triplets and restarts on the unconverged subspace. |
|
/ok to test b80ca30 |
|
/ok to test 93e9a19 |
|
@coderabbitai review |
✅ Action performedReview finished.
|
|
No actionable comments were generated in the recent review. 🎉 ℹ️ Recent review info⚙️ Run configurationConfiguration used: Path: .coderabbit.yaml Review profile: CHILL Plan: Pro Plus Run ID: 📒 Files selected for processing (8)
🚧 Files skipped from review as they are similar to previous changes (5)
📝 WalkthroughSummary by CodeRabbit
WalkthroughAdds a sparse Lanczos SVD solver across C++, Python, benchmarks, tests, and documentation. It adds shared configuration types, runtime wrappers, Python solver selection, validation tooling, and build integration for Lanczos and randomized sparse SVD paths. ChangesSparse Lanczos SVD solver
Estimated code review effort: 4 (Complex) | ~75 minutes Suggested reviewers: 🚥 Pre-merge checks | ✅ 5✅ Passed checks (5 passed)
✨ Finishing Touches🧪 Generate unit tests (beta)
Thanks for using CodeRabbit! It's free for OSS, and your support helps us grow. If you like it, consider giving us a shout-out. Comment |
![coderabbitai[bot]](https://avatars.githubusercontent.com/in/347564?s=60&v=4){kind=small}
There was a problem hiding this comment.
Actionable comments posted: 2
🧹 Nitpick comments (6)
python/pylibraft/pylibraft/tests/test_sparse.py (1)
228-405: 📐 Maintainability & Code Quality | 🔵 Trivial | 💤 Low valuePrefer standalone test functions over class methods.
The new Lanczos tests are added as methods on
class TestSvds. The repo guideline asks for standalonetest_foo_bar()functions rather thanclass TestFootest classes. These new tests follow the file's existing class convention, so this is non-blocking, but consider extracting them (or the whole suite) into module-level functions in a follow-up.As per coding guidelines: "Write tests as standalone
test_foo_bar()functions rather thanclass TestFootest classes."🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the rest with a brief reason, keep changes minimal, and validate. In `@python/pylibraft/pylibraft/tests/test_sparse.py` around lines 228 - 405, The new Lanczos tests are implemented as methods on TestSvds, but the test guideline prefers standalone test_foo_bar() functions. Move these cases out of the class and make them module-level tests, preserving the same assertions and parameters in functions like test_lanczos_diagonal_spectrum, test_lanczos_return_singular_vectors_modes, test_lanczos_rotated_clustered_spectrum, test_lanczos_rank_deficient_and_edge_shapes, test_lanczos_non_convergence_raises, and test_invalid_solver.Source: Coding guidelines
cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh (2)
251-254: 🎯 Functional Correctness | 🔵 Trivial | ⚡ Quick winBreakdown epsilon is not type-aware.
eps = 1e-9is fixed regardless ofValueTypeT, butfloat's machine epsilon (~1.19e-7) is larger than this threshold, so breakdown detection is comparatively looser/less reliable for thefloatinstantiation than fordouble.🔢 Scale epsilon by machine precision
- ValueTypeT eps = ValueTypeT(1e-9); + ValueTypeT eps = std::sqrt(std::numeric_limits<ValueTypeT>::epsilon());🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the rest with a brief reason, keep changes minimal, and validate. In `@cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh` around lines 251 - 254, The breakdown threshold in lanczos_svds.cuh is hardcoded in the Lanczos SVD setup, so it is not scaled to ValueTypeT and is too small for float. Update the epsilon initialization near the op.rows()/op.cols() setup to derive the threshold from machine precision for the current ValueTypeT (for example via std::numeric_limits<ValueTypeT>::epsilon() with an appropriate scale) so the breakdown check in the Lanczos routine behaves consistently for both float and double.
40-50: 🚀 Performance & Scalability | 🔵 Trivial | 🏗️ Heavy liftRedundant host-device sync in the hot Lanczos step loop.
cublasnrm2in the default host pointer mode already synchronizes internally to populateresult; the followingresource::sync_stream(handle, stream)on Line 48 is a second blocking barrier. Combined with the fact thatnormalize_or_randomize(and thusvector_norm) is invoked twice per Lanczos step (Lines 296-306, 320-331), this yields up to2 * active_ncvexplicit host syncs per restart iteration, which can materially hurt performance on large sparse matrices in this hot path.⚡ Remove the redundant explicit sync
template <typename ValueTypeT> ValueTypeT vector_norm(raft::resources const& handle, ValueTypeT const* x, int n) { common::nvtx::range<common::nvtx::domain::raft> scope("lanczos_svds::vector_norm"); auto cublas_handle = resource::get_cublas_handle(handle); auto stream = resource::get_cuda_stream(handle); ValueTypeT result{}; RAFT_CUBLAS_TRY(raft::linalg::detail::cublasnrm2(cublas_handle, n, x, 1, &result, stream)); - resource::sync_stream(handle, stream); return result; }If the extra sync was added deliberately to guard against non-default cuBLAS stream/pointer-mode configurations, please confirm with a comment; otherwise this is safe to drop.
Also applies to: 192-215, 275-332
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the rest with a brief reason, keep changes minimal, and validate. In `@cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh` around lines 40 - 50, The `vector_norm` helper in `lanczos_svds.cuh` is adding a redundant blocking `resource::sync_stream(handle, stream)` after `raft::linalg::detail::cublasnrm2`, which already completes before `result` is read in the default host pointer mode. Remove that extra sync from `vector_norm` and keep the rest of the Lanczos flow (`normalize_or_randomize` call sites) unchanged unless there is a documented need for a non-default cuBLAS mode; if so, add a brief comment near `vector_norm` explaining the requirement.Source: Coding guidelines
cpp/bench/prims/sparse/svds.cu (1)
312-371: 📐 Maintainability & Code Quality | 🔵 Trivial | 💤 Low valueDuplicate benchmark classes differ by a single config line.
lanczos_svds_benchandlanczos_mgs2_svds_benchare otherwise identical. Consider parameterizing a single class with a constructor/bool flag instead of duplicating the fullrun_benchmarkbody.♻️ Proposed consolidation
template <typename value_t> -class lanczos_svds_bench : public svds_bench_base<value_t> { +class lanczos_svds_bench : public svds_bench_base<value_t> { public: - using svds_bench_base<value_t>::svds_bench_base; + explicit lanczos_svds_bench(svds_input const& p, bool use_mgs2 = false) + : svds_bench_base<value_t>(p), use_mgs2_(use_mgs2) + { + } void run_benchmark(::benchmark::State& state) override { auto csr_matrix = this->csr_view(); raft::sparse::solver::sparse_lanczos_svd_config<value_t> config; config.n_components = this->params.k; config.ncv = this->params.ncv; config.tolerance = value_t(1e-4); config.max_iterations = this->params.max_iterations; config.seed = 1234; + config.use_mgs2_orthogonalization = use_mgs2_; this->loop_on_state( state, [&]() { raft::sparse::solver::sparse_lanczos_svd(this->handle, config, csr_matrix, this->singular_values.view(), this->U.view(), this->Vt.view()); }, false); } + + private: + bool use_mgs2_; }; - -template <typename value_t> -class lanczos_mgs2_svds_bench : public svds_bench_base<value_t> { - public: - using svds_bench_base<value_t>::svds_bench_base; - ... -};🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the rest with a brief reason, keep changes minimal, and validate. In `@cpp/bench/prims/sparse/svds.cu` around lines 312 - 371, Consolidate the duplicated benchmark logic in lanczos_svds_bench and lanczos_mgs2_svds_bench by extracting the shared run_benchmark flow into one implementation and parameterizing the orthogonalization choice with a constructor flag or similar field. Keep the common sparse_lanczos_svd_config setup in one place, only toggling config.use_mgs2_orthogonalization for the MGS2 variant, and preserve the existing loop_on_state and sparse_lanczos_svd call sites.cpp/tests/sparse/solver/lanczos_svds.cu (2)
88-125: 📐 Maintainability & Code Quality | 🔵 Trivial | ⚡ Quick winMissing coverage for
U/Vt=std::nullopt.
sparse_lanczos_svd_configdocumentsstd::nulloptas a valid way to skip storingU/Vt(persolver_types.hpp), but no test in this file exercises that path. Given RAFT's testing guideline to maintain high coverage of public APIs, consider adding a case that passesstd::nulloptfor one or both outputs.As per coding guidelines, "It's important for RAFT to maintain a high test coverage of the public APIs in order to minimize the potential for downstream projects to encounter unexpected build or runtime behavior."
Also applies to: 205-241, 315-330
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the rest with a brief reason, keep changes minimal, and validate. In `@cpp/tests/sparse/solver/lanczos_svds.cu` around lines 88 - 125, Add a test case in lanczos_svds.cu that exercises sparse_lanczos_svd with std::nullopt for U and/or Vt, since sparse_lanczos_svd_config and the solver API support skipping these outputs. Update the Run helper or add a new helper around sparse_lanczos_svd to call it with nullopt, then verify the singular values still match expectations without materializing the omitted outputs. Use the existing sparse_lanczos_svd, sparse_lanczos_svd_config, and Run symbols to keep the new coverage aligned with the current test structure.Source: Coding guidelines
25-544: 📐 Maintainability & Code Quality | 🔵 Trivial | ⚡ Quick winConsider extracting shared test harness.
The CSR buffer setup → config construction → solve → host copy-back → sync pattern is repeated near-identically across
LanczosSvdsTest::Run,LanczosClusteredSpectrumTest::Run,run_diagonal_hardening_case, andrun_rotated_block_hardening_case. Extracting a common helper (e.g., taking indptr/indices/values + config and returning host S/U/Vt) would reduce duplication and ease future maintenance of this test suite.🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the rest with a brief reason, keep changes minimal, and validate. In `@cpp/tests/sparse/solver/lanczos_svds.cu` around lines 25 - 544, The test file repeats the same CSR setup, sparse_lanczos_svd config construction, solve, host copy-back, and stream sync flow in LanczosSvdsTest::Run, LanczosClusteredSpectrumTest::Run, run_diagonal_hardening_case, and run_rotated_block_hardening_case. Extract that sequence into a shared helper that accepts the CSR buffers and sparse_lanczos_svd_config and returns the host S/U/Vt results, then update each of those call sites to use it while keeping only the case-specific assertions and matrix construction in place.
🤖 Prompt for all review comments with AI agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.
Inline comments:
In `@cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh`:
- Around line 124-162: The cuBLAS pointer mode in mgs2_orthogonalize is only
restored on the success path, so an exception from any RAFT_CUBLAS_TRY or
RAFT_CUDA_TRY can leave the shared cublas handle in device mode. Update
mgs2_orthogonalize to use an exception-safe scoped guard, or reuse an existing
RAFT canonical pointer-mode RAII helper if one exists, so the handle is always
returned to CUBLAS_POINTER_MODE_HOST on every exit path. Keep the existing
cublasSetPointerMode setup and ensure the restore happens automatically even if
the loop aborts early.
In `@cpp/include/raft/sparse/solver/lanczos_svds.cuh`:
- Around line 1-16: Move the Lanczos SVDS implementation templates out of the
public header by splitting the existing wrapper in `lanczos_svds.cuh` into the
usual inline/explicit-instantiation pair. Keep only the public wrapper
declarations in the header, move the template definitions from
`detail/lanczos_svds.cuh` into an `-inl.cuh`/`-ext.cuh` layout, and add
`RAFT_COMPILED` / `RAFT_EXPLICIT_INSTANTIATE_ONLY` guards so downstream
translation units do not instantiate the implementation directly. Make sure the
wrapper overloads remain the only public entry points and that the float/double
explicit instantiations used by the `.cu` files are reused.
---
Nitpick comments:
In `@cpp/bench/prims/sparse/svds.cu`:
- Around line 312-371: Consolidate the duplicated benchmark logic in
lanczos_svds_bench and lanczos_mgs2_svds_bench by extracting the shared
run_benchmark flow into one implementation and parameterizing the
orthogonalization choice with a constructor flag or similar field. Keep the
common sparse_lanczos_svd_config setup in one place, only toggling
config.use_mgs2_orthogonalization for the MGS2 variant, and preserve the
existing loop_on_state and sparse_lanczos_svd call sites.
In `@cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh`:
- Around line 251-254: The breakdown threshold in lanczos_svds.cuh is hardcoded
in the Lanczos SVD setup, so it is not scaled to ValueTypeT and is too small for
float. Update the epsilon initialization near the op.rows()/op.cols() setup to
derive the threshold from machine precision for the current ValueTypeT (for
example via std::numeric_limits<ValueTypeT>::epsilon() with an appropriate
scale) so the breakdown check in the Lanczos routine behaves consistently for
both float and double.
- Around line 40-50: The `vector_norm` helper in `lanczos_svds.cuh` is adding a
redundant blocking `resource::sync_stream(handle, stream)` after
`raft::linalg::detail::cublasnrm2`, which already completes before `result` is
read in the default host pointer mode. Remove that extra sync from `vector_norm`
and keep the rest of the Lanczos flow (`normalize_or_randomize` call sites)
unchanged unless there is a documented need for a non-default cuBLAS mode; if
so, add a brief comment near `vector_norm` explaining the requirement.
In `@cpp/tests/sparse/solver/lanczos_svds.cu`:
- Around line 88-125: Add a test case in lanczos_svds.cu that exercises
sparse_lanczos_svd with std::nullopt for U and/or Vt, since
sparse_lanczos_svd_config and the solver API support skipping these outputs.
Update the Run helper or add a new helper around sparse_lanczos_svd to call it
with nullopt, then verify the singular values still match expectations without
materializing the omitted outputs. Use the existing sparse_lanczos_svd,
sparse_lanczos_svd_config, and Run symbols to keep the new coverage aligned with
the current test structure.
- Around line 25-544: The test file repeats the same CSR setup,
sparse_lanczos_svd config construction, solve, host copy-back, and stream sync
flow in LanczosSvdsTest::Run, LanczosClusteredSpectrumTest::Run,
run_diagonal_hardening_case, and run_rotated_block_hardening_case. Extract that
sequence into a shared helper that accepts the CSR buffers and
sparse_lanczos_svd_config and returns the host S/U/Vt results, then update each
of those call sites to use it while keeping only the case-specific assertions
and matrix construction in place.
In `@python/pylibraft/pylibraft/tests/test_sparse.py`:
- Around line 228-405: The new Lanczos tests are implemented as methods on
TestSvds, but the test guideline prefers standalone test_foo_bar() functions.
Move these cases out of the class and make them module-level tests, preserving
the same assertions and parameters in functions like
test_lanczos_diagonal_spectrum, test_lanczos_return_singular_vectors_modes,
test_lanczos_rotated_clustered_spectrum,
test_lanczos_rank_deficient_and_edge_shapes,
test_lanczos_non_convergence_raises, and test_invalid_solver.
🪄 Autofix (Beta)
Fix all unresolved CodeRabbit comments on this PR:
- Push a commit to this branch (recommended)
- Create a new PR with the fixes
ℹ️ Review info
⚙️ Run configuration
Configuration used: Path: .coderabbit.yaml
Review profile: CHILL
Plan: Enterprise
Run ID: bf8109ed-b790-4078-9385-89d9febf71c1
📒 Files selected for processing (21)
cpp/CMakeLists.txtcpp/bench/prims/CMakeLists.txtcpp/bench/prims/sparse/svds.cucpp/include/raft/sparse/solver/detail/lanczos_svds.cuhcpp/include/raft/sparse/solver/detail/randomized_svds.cuhcpp/include/raft/sparse/solver/detail/svds_optional.hppcpp/include/raft/sparse/solver/lanczos_svds.cuhcpp/include/raft/sparse/solver/randomized_svds.cuhcpp/include/raft/sparse/solver/solver_types.hppcpp/include/raft/sparse/solver/svds_config.hppcpp/include/raft_runtime/solver/lanczos_svds.hppcpp/include/raft_runtime/solver/randomized_svds.hppcpp/src/raft_runtime/solver/lanczos_svds.cuhcpp/src/raft_runtime/solver/lanczos_svds_double.cucpp/src/raft_runtime/solver/lanczos_svds_float.cucpp/tests/CMakeLists.txtcpp/tests/sparse/solver/lanczos_svds.cudocs/source/cpp_api/sparse_solver.rstpython/pylibraft/benchmarks/sparse_svds.pypython/pylibraft/pylibraft/sparse/linalg/svds.pyxpython/pylibraft/pylibraft/tests/test_sparse.py
💤 Files with no reviewable changes (1)
- cpp/include/raft/sparse/solver/svds_config.hpp
| /* | ||
| * SPDX-FileCopyrightText: Copyright (c) 2026, NVIDIA CORPORATION. | ||
| * SPDX-License-Identifier: Apache-2.0 | ||
| */ | ||
|
|
||
| #pragma once | ||
|
|
||
| #include <raft/core/device_csr_matrix.hpp> | ||
| #include <raft/core/device_mdspan.hpp> | ||
| #include <raft/core/resources.hpp> | ||
| #include <raft/sparse/solver/detail/csr_linear_operator.cuh> | ||
| #include <raft/sparse/solver/detail/lanczos_svds.cuh> | ||
| #include <raft/sparse/solver/detail/svds_optional.hpp> | ||
|
|
||
| #include <optional> | ||
|
|
There was a problem hiding this comment.
🚀 Performance & Scalability | 🟠 Major | 🏗️ Heavy lift
🧩 Analysis chain
🏁 Script executed:
#!/bin/bash
# Check whether sibling solver headers (randomized_svds, eigensolver lanczos) follow the split-header pattern,
# to confirm whether this PR's lanczos_svds.cuh is a new deviation or matches existing convention.
fd -e cuh -e hpp . cpp/include/raft/sparse/solver cpp/include/raft/sparse/solver/detail
rg -n "RAFT_EXPLICIT_INSTANTIATE_ONLY|RAFT_COMPILED|extern template" cpp/include/raft/sparse/solver/randomized_svds.cuh cpp/include/raft/sparse/solver/detail/randomized_svds.cuh 2>/dev/nullRepository: NVIDIA/raft
Length of output: 1595
🏁 Script executed:
#!/bin/bash
set -euo pipefail
# Inspect the public and detail headers for Lanczos SVDS and nearby solver headers.
for f in \
cpp/include/raft/sparse/solver/lanczos_svds.cuh \
cpp/include/raft/sparse/solver/lanczos.cuh \
cpp/include/raft/sparse/solver/randomized_svds.cuh \
cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh \
cpp/include/raft/sparse/solver/detail/lanczos.cuh \
cpp/include/raft/sparse/solver/detail/randomized_svds.cuh
do
echo "===== $f ====="
wc -l "$f"
done
echo
echo "===== public header excerpts ====="
sed -n '1,220p' cpp/include/raft/sparse/solver/lanczos_svds.cuh
echo
sed -n '1,220p' cpp/include/raft/sparse/solver/lanczos.cuh
echo
sed -n '1,220p' cpp/include/raft/sparse/solver/randomized_svds.cuh
echo
echo "===== detail header excerpts ====="
sed -n '1,260p' cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh
echo
sed -n '1,260p' cpp/include/raft/sparse/solver/detail/lanczos.cuh
echo
sed -n '1,260p' cpp/include/raft/sparse/solver/detail/randomized_svds.cuh
echo
echo "===== explicit instantiation files ====="
fd -a -e cu -e cuh -e cpp 'lanczos_svds|randomized_svds|lanczos' cpp/src cpp/include | sort
rg -n "explicit instantiation|extern template|RAFT_EXPLICIT_INSTANTIATE_ONLY|RAFT_COMPILED|template class|template void|template struct" cpp/src cpp/include/raft/sparse/solver -g '*lanczos*' -g '*randomized_svds*'Repository: NVIDIA/raft
Length of output: 47320
🏁 Script executed:
#!/bin/bash
set -euo pipefail
# Inspect the files that were referenced by the initial results to identify how this area is organized.
fd -a -e cu -e cuh -e cpp 'lanczos_svds|randomized_svds|lanczos' cpp/src cpp/include | sort
rg -n "RAFT_EXPLICIT_INSTANTIATE_ONLY|RAFT_COMPILED|extern template|template class|template void|template struct" \
cpp/src cpp/include/raft/sparse/solver \
-g '*lanczos*' -g '*randomized_svds*'Repository: NVIDIA/raft
Length of output: 1355
🏁 Script executed:
#!/bin/bash
set -euo pipefail
# Map the relevant files and search for split-header patterns in this area.
printf "\n-- file list --\n"
git ls-files \
'cpp/include/raft/sparse/solver/*' \
'cpp/include/raft/sparse/solver/detail/*' \
'cpp/src/raft_runtime/solver/*' | \
rg 'lanczos|randomized_svds|svds'
printf "\n-- public headers --\n"
for f in \
cpp/include/raft/sparse/solver/lanczos_svds.cuh \
cpp/include/raft/sparse/solver/randomized_svds.cuh \
cpp/include/raft/sparse/solver/lanczos.cuh
do
echo "### $f"
sed -n '1,220p' "$f"
echo
done
printf "\n-- detail headers --\n"
for f in \
cpp/include/raft/sparse/solver/detail/lanczos_svds.cuh \
cpp/include/raft/sparse/solver/detail/randomized_svds.cuh \
cpp/include/raft/sparse/solver/detail/lanczos.cuh
do
echo "### $f"
sed -n '1,260p' "$f"
echo
done
printf "\n-- explicit instantiation files --\n"
for f in \
cpp/src/raft_runtime/solver/lanczos_svds_float.cu \
cpp/src/raft_runtime/solver/lanczos_svds_double.cu \
cpp/src/raft_runtime/solver/randomized_svds_float.cu \
cpp/src/raft_runtime/solver/randomized_svds_double.cu
do
if [ -f "$f" ]; then
echo "### $f"
sed -n '1,220p' "$f"
echo
fi
done
printf "\n-- search for RAFT_COMPILED / RAFT_EXPLICIT_INSTANTIATE_ONLY / extern template --\n"
rg -n "RAFT_COMPILED|RAFT_EXPLICIT_INSTANTIATE_ONLY|extern template" \
cpp/include/raft/sparse/solver cpp/src/raft_runtime/solverRepository: NVIDIA/raft
Length of output: 48062
🏁 Script executed:
#!/bin/bash
set -euo pipefail
# Inspect the runtime wrapper headers and any nearby instantiation plumbing.
for f in \
cpp/src/raft_runtime/solver/lanczos_svds.cuh \
cpp/src/raft_runtime/solver/randomized_svds.cuh \
cpp/src/raft_runtime/solver/lanczos_solver.cuh
do
echo "===== $f ====="
wc -l "$f"
sed -n '1,240p' "$f"
echo
done
echo "===== search for FUNC_DEF / explicit instantiation plumbing ====="
rg -n "FUNC_DEF|RAFT_EXPLICIT_INSTANTIATE_ONLY|RAFT_COMPILED|extern template|template class|template void" cpp/src/raft_runtime/solver cpp/include/raft/sparse/solverRepository: NVIDIA/raft
Length of output: 7292
Split Lanczos SVDS templates out of the public header. cpp/include/raft/sparse/solver/lanczos_svds.cuh still includes detail/lanczos_svds.cuh directly, so downstream TUs will instantiate the implementation even though the runtime .cu files only instantiate the wrapper overloads. Add the usual -inl.cuh/-ext.cuh split with RAFT_COMPILED / RAFT_EXPLICIT_INSTANTIATE_ONLY guards so the prebuilt float/double instantiations can be reused.
🤖 Prompt for AI Agents
Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.
In `@cpp/include/raft/sparse/solver/lanczos_svds.cuh` around lines 1 - 16, Move
the Lanczos SVDS implementation templates out of the public header by splitting
the existing wrapper in `lanczos_svds.cuh` into the usual
inline/explicit-instantiation pair. Keep only the public wrapper declarations in
the header, move the template definitions from `detail/lanczos_svds.cuh` into an
`-inl.cuh`/`-ext.cuh` layout, and add `RAFT_COMPILED` /
`RAFT_EXPLICIT_INSTANTIATE_ONLY` guards so downstream translation units do not
instantiate the implementation directly. Make sure the wrapper overloads remain
the only public entry points and that the float/double explicit instantiations
used by the `.cu` files are reused.
Source: Path instructions
There was a problem hiding this comment.
I dont see this done anywhere else.
There was a problem hiding this comment.
🧩 Analysis chain
🏁 Script executed:
#!/bin/bash
set -euo pipefail
echo "== search for -inl.cuh / -ext.cuh pattern anywhere in repo =="
fd -e cuh . cpp/include | rg -- '-inl\.cuh|-ext\.cuh' | head -50
echo
echo "== count how many such files exist =="
fd -e cuh . cpp/include | rg -c -- '-inl\.cuh|-ext\.cuh' || true
echo
echo "== check other sparse solver headers for RAFT_COMPILED / RAFT_EXPLICIT_INSTANTIATE_ONLY pattern =="
rg -n "RAFT_COMPILED|RAFT_EXPLICIT_INSTANTIATE_ONLY" cpp/include/raft/sparse -g '*.cuh' -g '*.hpp'
echo
echo "== check broader raft for this split convention, e.g. in linalg or distance =="
rg -n "RAFT_EXPLICIT_INSTANTIATE_ONLY" cpp/include/raft -l | head -30Length of output: 1326
@Intron7, you're right — I checked more broadly, and this split pattern is not the established convention here. Only 3 files in the entire cpp/include tree use the -inl.cuh/-ext.cuh split (linalg/detail/coalesced_reduction-inl.cuh, matrix/detail/select_k-inl.cuh, sparse/matrix/detail/select_k-inl.cuh), none of which are in cpp/include/raft/sparse/solver. There's also no use of RAFT_COMPILED/RAFT_EXPLICIT_INSTANTIATE_ONLY anywhere under cpp/include/raft/sparse, and sibling headers like randomized_svds.cuh follow the same direct-include-of-detail-header pattern used in this PR. Withdrawing this suggestion — lanczos_svds.cuh is consistent with the existing convention in this area of the codebase.
(๑•́ ‧̫ •̀๑) noted, thanks for pushing back with the check
Signed-off-by: Intron7 <severin.dicks@icloud.com>
divyegala
left a comment
divyegala
left a comment
There was a problem hiding this comment.
Thank you for the PR! For the logical bits, I used AI's help and tagged those as [AI-reviewed]. If you find the comments invalid, please feel free to discard them.
| } | ||
| } | ||
|
|
||
| compute_restart_vector( |
There was a problem hiding this comment.
[AI-reviewed] restart_coeffs are coordinates in the active Krylov V basis starting at the pre-lock n_locked, but compute_ritz_vectors overwrites that segment and n_locked is incremented before compute_restart_vector runs. On partial convergence this forms the restart vector from a shifted/overwritten basis. Could we compute and preserve the restart vector before locking, or retain the original active basis until after restart construction?
There was a problem hiding this comment.
Good catch — the shifted/overwritten basis read is real. I tested constructing the restart vector from the intact pre-lock basis, but it caused silent correctness regressions: even with the shipped k + 10 ncv floor, the single-direction version returned the wrong top-k spectrum on 64/128 seeds of a geometric spectrum, while the current restart was correct on 128/128. A mixed pre-lock restart also failed cases that the current form passed.
The investigation did uncover a missing detail from rapids-singlecell: it allocates U_full and V_full with cp.zeros. When d > 1 Ritz vectors are locked together, the shifted restart window extends d - 1 columns beyond the basis written by the current bidiagonalization. RSC therefore reads deterministic zeros there, whereas RAFT’s uninitialized allocation could read arbitrary pooled-memory contents. I’ve restored that zero-initialization while retaining the restart behavior that performed best and matches RSC.
The sweep also exposed a separate porting regression in ncv: rapids-singlecell enforces ncv = max(ncv, k + 10) even for user-supplied values and uses a 2.5 * k default for the large-matrix/small-k branch. RAFT had reduced these to max(ncv, k) and approximately 2.2 * k. I restored both reference settings and documented the effective floor. With that floor, the current restart completed every tested production-reachable family correctly: 128/128 geometric, 64/64 linear, 64/64 repeated-spectrum, and 32/32 larger-k runs. This does not itself fix the shifted-tail read—the zero-initialization does—but it restores the convergence margin assumed by the reference algorithm.
| void sparse_lanczos_svd( | ||
| raft::resources const& handle, | ||
| sparse_lanczos_svd_config<ValueTypeT> const& config, | ||
| raft::device_csr_matrix_view<const ValueTypeT, int, int, NNZTypeT> A, |
There was a problem hiding this comment.
The generic API hard-codes uint32_t for all output and temporary views, while the operator shape and CSR row/column indices are int. Is signed 32-bit indexing an intentional constraint? If so, could we document and validate it explicitly; otherwise, should the index type be templated as in the existing sparse Lanczos eigensolver?
There was a problem hiding this comment.
It's an intentional constraint, matching sparse_randomized_svd: the operator reports its shape as int (the CSR path is additionally bound to int indptr/indices by the cuSPARSE SpMM wrapper), and all device views use uint32_t extents. I've documented this explicitly on both overloads (m, n < 2^31). Templating the index type like the existing eigensolver could be a follow-up if a >2^31 use case shows up, but I'd rather keep the two SVD solvers' APIs consistent for now.
Signed-off-by: Intron7 <severin.dicks@icloud.com>
5 participants
As discussed previously with @cjnolet I'm also adding my Lanczos SVD solver for sparse CSR matrices.
This is the more precise sparse SVD path next to the existing randomized solver. The solver repeatedly applies
A @ vandA.T @ uto build Krylov bases, computes the SVD of the small bidiagonal problem, uses the resulting Ritz vectors to identify converged singular triplets, locks those vectors, and restarts on the remaining unconverged part. It also uses full reorthogonalization and a finalA @ Vrefinement step to improve singular values and left singular vectors.Compared with randomized SVD, this is aimed at quality: clustered spectra, slow singular-value decay, near-rank-deficient inputs, and PCA workloads where ARPACK-like accuracy matters.
Ran the #2999 -style row sweep on the same singlecell dataset, k=50, n_oversamples=10, n_power_iters=2, best-of-3 GPU timings. GPU: RTX PRO 6000 Blackwell.
Using the 2999 CPU baselines for context: