From 3547c6af6a3bd152a12340bede2d2e94309f263f Mon Sep 17 00:00:00 2001 From: says1117 Date: Mon, 20 Jul 2026 02:55:03 -0400 Subject: [PATCH 1/7] Replace hardcoded Lanczos eigenvalue goldens with property-based checks Golden-value comparisons via devArrMatch/CompareApprox were flaky across CUDA versions/architectures because FP reduction order isn't deterministic. Replace with residual, orthogonality, unit-norm, and ordering checks that verify the defining properties of an eigenpair instead of exact values. Fixes #2519 --- cpp/tests/sparse/solver/lanczos.cu | 270 ++++++++++++++++++++--------- 1 file changed, 187 insertions(+), 83 deletions(-) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index cd10cee71a..f803502b3e 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -8,9 +8,12 @@ #include #include #include +#include #include #include #include +#include +#include #include #include #include @@ -35,7 +38,9 @@ #include #include +#include #include +#include namespace raft::sparse::solver { @@ -52,7 +57,6 @@ struct lanczos_inputs { std::vector rows; // indptr std::vector cols; // indices std::vector vals; // data - std::vector expected_eigenvalues; }; template @@ -68,9 +72,143 @@ struct rmat_lanczos_inputs { int r_scale; int c_scale; float sparsity; - std::vector expected_eigenvalues; }; +// Frobenius norm of the sparse matrix, used to scale residual/orthogonality +// tolerances to the magnitude of the problem instead of an absolute value. +template +ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* vals, NnzType nnz) +{ + auto vals_view = + raft::make_device_vector_view(vals, static_cast(nnz)); + ValueType sum_sq{}; + raft::linalg::dot(handle, vals_view, vals_view, raft::make_host_scalar_view(&sum_sq)); + return std::sqrt(sum_sq); +} + +// Validates that (eigenvalues[i], eigenvectors[:, i]) are genuine eigenpairs +// of `A`, without relying on hardcoded golden eigenvalues. This checks the +// defining characteristics of a converged real-symmetric eigensolve: +// - ||A*v_i - lambda_i*v_i|| is within the noise floor implied by the +// solver's own convergence criterion (see lanczos_solver_config::tolerance, +// which is documented as governing exactly this residual). +// - lambda_i matches the Rayleigh quotient of v_i. This pins the eigenvalue +// down more sharply than the residual alone on matrices with a large +// ||A||, where a small eigenvalue's residual budget is dominated by the +// matrix scale. +// - each eigenvector is unit-normalized. +// - distinct eigenvectors are mutually orthogonal (catches "ghost" +// duplicate eigenpairs from loss of orthogonality, which a residual-only +// check would miss). +// - eigenvalues are returned in ascending algebraic order (verified +// against every previously-hardcoded golden array in this file, which +// were all ascending regardless of `which`). +template +void expect_valid_eigenpairs(raft::resources const& handle, + raft::spectral::matrix::sparse_matrix_t const& A, + ValueType const* eigenvalues, + ValueType* eigenvectors, + IndexType n, + int n_components, + ValueType frobenius_norm) +{ + auto stream = resource::get_cuda_stream(handle); + auto eps = std::numeric_limits::epsilon(); + + // Tolerances are expressed as multiples of the floating-point noise floor + // for each kind of quantity, since the exact low-order bits depend on + // cuSPARSE/cuBLAS reduction order, which varies by CUDA version and GPU + // architecture (the root cause of issue #2519). + // + // Quantities carrying the magnitude of A (the residual ||A*v - lambda*v|| + // and the Rayleigh quotient) have a noise floor proportional to + // ||A||_F * eps. Measured across all fixtures in this file on sm_75 / + // CUDA 13, these stay within 15 * ||A||_F * eps and show no dependence on + // n (n=100 and n=4096 both land in that range), so n is deliberately not a + // factor here -- including it would inflate the tolerance ~40x on the RMAT + // fixture and let a 1% eigenvalue error pass undetected. + // + // Dimensionless quantities (deviation from unit norm, mutual + // orthogonality) have a noise floor of a small multiple of eps; measured + // values stay within 100 * eps. + const ValueType matrix_tol = ValueType(500) * std::max(frobenius_norm, ValueType(1)) * eps; + const ValueType unit_tol = ValueType(1000) * eps; + const ValueType norm_tol = unit_tol; + const ValueType ortho_tol = unit_tol; + const ValueType residual_tol = matrix_tol; + + std::vector host_eigenvalues(n_components); + raft::update_host(host_eigenvalues.data(), eigenvalues, n_components, stream); + RAFT_CUDA_TRY(cudaStreamSynchronize(stream)); + + for (int i = 1; i < n_components; ++i) { + EXPECT_LE(host_eigenvalues[i - 1], host_eigenvalues[i] + matrix_tol) + << "eigenvalues not returned in ascending algebraic order at index " << i; + } + + auto residual = raft::make_device_vector(handle, n); + + for (int i = 0; i < n_components; ++i) { + ValueType* v_i = eigenvectors + static_cast(i) * n; + auto v_i_view = raft::make_device_vector_view(v_i, n); + + // residual = A * v_i + A.mv(ValueType{1}, v_i, ValueType{0}, residual.data_handle()); + + ValueType rayleigh{}; + raft::linalg::dot(handle, + v_i_view, + raft::make_const_mdspan(residual.view()), + raft::make_host_scalar_view(&rayleigh)); + + // residual = A*v_i - lambda_i*v_i + ValueType neg_lambda = -host_eigenvalues[i]; + raft::linalg::axpy( + handle, raft::make_host_scalar_view(&neg_lambda), v_i_view, residual.view()); + + ValueType residual_sq{}; + raft::linalg::dot(handle, + raft::make_const_mdspan(residual.view()), + raft::make_const_mdspan(residual.view()), + raft::make_host_scalar_view(&residual_sq)); + + ValueType v_norm_sq{}; + raft::linalg::dot( + handle, v_i_view, v_i_view, raft::make_host_scalar_view(&v_norm_sq)); + + ValueType residual_norm = std::sqrt(residual_sq); + ValueType v_norm = std::sqrt(v_norm_sq); + + EXPECT_NEAR(v_norm, ValueType(1), norm_tol) + << "eigenvector " << i << " is not unit-normalized: ||v||=" << v_norm; + EXPECT_LE(residual_norm, residual_tol) + << "eigenpair " << i + << " failed ||A*v - lambda*v|| residual check: residual=" << residual_norm + << " lambda=" << host_eigenvalues[i] << " tol=" << residual_tol; + + // The Rayleigh quotient (v^T A v) / (v^T v) is the eigenvalue implied by + // the returned eigenvector. Comparing it against the returned eigenvalue + // pins down lambda independently of the residual check, which on a + // large-||A|| matrix is too coarse to notice a small relative error in a + // small-magnitude eigenvalue. + EXPECT_NEAR(host_eigenvalues[i], rayleigh / v_norm_sq, matrix_tol) + << "eigenvalue " << i + << " disagrees with the Rayleigh quotient of its eigenvector: lambda=" << host_eigenvalues[i] + << " rayleigh=" << rayleigh / v_norm_sq; + + for (int j = 0; j < i; ++j) { + ValueType const* v_j = eigenvectors + static_cast(j) * n; + auto v_j_view = raft::make_device_vector_view(v_j, n); + ValueType dot_ij{}; + raft::linalg::dot( + handle, v_i_view, v_j_view, raft::make_host_scalar_view(&dot_ij)); + EXPECT_LE(std::abs(dot_ij), ortho_tol) + << "eigenvectors " << i << " and " << j + << " are not orthogonal: |v_i . v_j|=" << std::abs(dot_ij); + } + } +} + template class rmat_lanczos_tests : public ::testing::TestWithParam> { @@ -79,8 +217,6 @@ class rmat_lanczos_tests : params(::testing::TestWithParam>::GetParam()), stream(resource::get_cuda_stream(handle)), rng(params.seed), - expected_eigenvalues(raft::make_device_vector( - handle, params.n_components)), r_scale(params.r_scale), c_scale(params.c_scale), sparsity(params.sparsity) @@ -88,13 +224,7 @@ class rmat_lanczos_tests } protected: - void SetUp() override - { - raft::copy(expected_eigenvalues.data_handle(), - params.expected_eigenvalues.data(), - params.n_components, - stream); - } + void SetUp() override {} void TearDown() override {} @@ -200,11 +330,15 @@ class rmat_lanczos_tests eigenvalues.view(), eigenvectors.view()); - ASSERT_TRUE(raft::devArrMatch(eigenvalues.data_handle(), - expected_eigenvalues.data_handle(), - n_components, - raft::CompareApprox(1e-5), - stream)); + ValueType frobenius_norm = + compute_frobenius_norm(handle, symmetric_coo.vals(), symmetric_coo.nnz); + expect_valid_eigenpairs(handle, + csr_m, + eigenvalues.data_handle(), + eigenvectors.data_handle(), + (IndexType)symmetric_coo.n_rows, + n_components, + frobenius_norm); // Reproducibility test - run again with same seed and verify exact match raft::device_vector eigenvalues2 = @@ -254,11 +388,13 @@ class rmat_lanczos_tests eigenvalues_coo.view(), eigenvectors_coo.view()); - ASSERT_TRUE(raft::devArrMatch(eigenvalues_coo.data_handle(), - expected_eigenvalues.data_handle(), - n_components, - raft::CompareApprox(1e-4), - stream)); + expect_valid_eigenpairs(handle, + csr_m, + eigenvalues_coo.data_handle(), + eigenvectors_coo.data_handle(), + (IndexType)symmetric_coo.n_rows, + n_components, + frobenius_norm); } protected: @@ -269,7 +405,6 @@ class rmat_lanczos_tests int r_scale; int c_scale; float sparsity; - raft::device_vector expected_eigenvalues; }; template @@ -288,9 +423,7 @@ class lanczos_tests : public ::testing::TestWithParam( handle, params.n_components)), eigenvectors(raft::make_device_matrix( - handle, n, params.n_components)), - expected_eigenvalues( - raft::make_device_vector(handle, params.n_components)) + handle, n, params.n_components)) { } @@ -300,10 +433,6 @@ class lanczos_tests : public ::testing::TestWithParam( const_cast(vals.data_handle()), csr_structure); + raft::spectral::matrix::sparse_matrix_t const A_check{ + handle, rows.data_handle(), cols.data_handle(), vals.data_handle(), n, (uint64_t)nnz}; + ValueType frobenius_norm = compute_frobenius_norm(handle, vals.data_handle(), nnz); + std::get<0>(stats) = raft::sparse::solver::lanczos_compute_eigenpairs( handle, config, @@ -348,13 +481,13 @@ class lanczos_tests : public ::testing::TestWithParam( - eigenvalues.data_handle(), - expected_eigenvalues.data_handle(), - params.n_components, - raft::CompareApprox( - params.which == raft::sparse::solver::LANCZOS_WHICH::SM ? 5e-5 : 1e-5), - stream)); + expect_valid_eigenpairs(handle, + A_check, + eigenvalues.data_handle(), + eigenvectors.data_handle(), + (IndexType)n, + params.n_components, + frobenius_norm); // Reproducibility test - run again with same seed and verify exact match raft::device_vector eigenvalues2 = @@ -403,13 +536,13 @@ class lanczos_tests : public ::testing::TestWithParam( - eigenvalues_coo.data_handle(), - expected_eigenvalues.data_handle(), - params.n_components, - raft::CompareApprox( - params.which == raft::sparse::solver::LANCZOS_WHICH::SM ? 5e-5 : 1e-5), - stream)); + expect_valid_eigenpairs(handle, + A_check, + eigenvalues_coo.data_handle(), + eigenvectors_coo.data_handle(), + (IndexType)n, + params.n_components, + frobenius_norm); } protected: @@ -425,7 +558,6 @@ class lanczos_tests : public ::testing::TestWithParam v0; raft::device_vector eigenvalues; raft::device_matrix eigenvectors; - raft::device_vector expected_eigenvalues; }; template @@ -695,8 +827,7 @@ const std::vector> inputsf = { 0.4350971, 0.6997072, 0.4320931, 0.3315690, 0.0844443, 0.1445242, 0.3059566, 0.6594226, 0.8961608, 0.6498466, 0.9585592, 0.7827352, 0.6498466, 0.2812338, 0.1767728, 0.5810611, 0.7269946, 0.6997072, 0.1705930, 0.1792683, 0.1077409, 0.9368132, 0.4823034, 0.8311127, - 0.7194629, 0.6273088, 0.2909178, 0.5188584, 0.5876446, 0.2812338}, - {-2.0369630, -1.7673520}}}; + 0.7194629, 0.6273088, 0.2909178, 0.5188584, 0.5876446, 0.2812338}}}; const std::vector> inputsd = { {2, @@ -746,8 +877,7 @@ const std::vector> inputsd = { 0.4350971, 0.6997072, 0.4320931, 0.3315690, 0.0844443, 0.1445242, 0.3059566, 0.6594226, 0.8961608, 0.6498466, 0.9585592, 0.7827352, 0.6498466, 0.2812338, 0.1767728, 0.5810611, 0.7269946, 0.6997072, 0.1705930, 0.1792683, 0.1077409, 0.9368132, 0.4823034, 0.8311127, - 0.7194629, 0.6273088, 0.2909178, 0.5188584, 0.5876446, 0.2812338}, - {-2.0369630, -1.7673520}}}; + 0.7194629, 0.6273088, 0.2909178, 0.5188584, 0.5876446, 0.2812338}}}; const std::vector> inputsd_SM = { {2, @@ -760,8 +890,7 @@ const std::vector> inputsd_SM = { 42, rows(), cols(), - vals(), - {-0.03944135, 0.01367824}}}; + vals()}}; const std::vector> inputsd_LM = { {2, @@ -774,8 +903,7 @@ const std::vector> inputsd_LM = { 42, rows(), cols(), - vals(), - {-2.00968758, 3.45939575}}}; + vals()}}; const std::vector> inputsd_LA = { {2, @@ -788,8 +916,7 @@ const std::vector> inputsd_LA = { 42, rows(), cols(), - vals(), - {1.99482678, 3.45939575}}}; + vals()}}; const std::vector> inputsd_SA = { {2, @@ -802,8 +929,7 @@ const std::vector> inputsd_SA = { 42, rows(), cols(), - vals(), - {-2.00968758, -1.83487935}}}; + vals()}}; const std::vector> inputsf_SM = { {2, @@ -816,8 +942,7 @@ const std::vector> inputsf_SM = { 42, rows(), cols(), - vals(), - {-0.03944135, 0.01367824}}}; + vals()}}; const std::vector> inputsf_LM = { {2, @@ -830,8 +955,7 @@ const std::vector> inputsf_LM = { 42, rows(), cols(), - vals(), - {-2.00968758, 3.45939575}}}; + vals()}}; const std::vector> inputsf_LA = { {2, @@ -844,8 +968,7 @@ const std::vector> inputsf_LA = { 42, rows(), cols(), - vals(), - {1.99482678, 3.45939575}}}; + vals()}}; const std::vector> inputsf_SA = { {2, @@ -858,29 +981,10 @@ const std::vector> inputsf_SA = { 42, rows(), cols(), - vals(), - {-2.00968758, -1.83487935}}}; + vals()}}; const std::vector> rmat_inputsf = { - {50, - 100, - 10000, - 0, - 0, - 1e-9, - raft::sparse::solver::LANCZOS_WHICH::SA, - 42, - 12, - 12, - 1, - {-122.526794, -74.00686, -59.698284, -54.68617, -49.686813, -34.02644, -32.130703, - -31.26906, -30.32097, -22.946098, -20.497862, -20.23817, -19.269697, -18.42496, - -17.675667, -17.013401, -16.734581, -15.820215, -15.73925, -15.448187, -15.044634, - -14.692028, -14.127425, -13.967386, -13.6237755, -13.469393, -13.181225, -12.777589, - -12.623185, -12.55508, -12.2874565, -12.053391, -11.677346, -11.558279, -11.163732, - -10.922034, -10.7936945, -10.558049, -10.205776, -10.005316, -9.559181, -9.491834, - -9.242631, -8.883637, -8.765364, -8.688508, -8.458255, -8.385196, -8.217982, - -8.0442095}}}; + {50, 100, 10000, 0, 0, 1e-9, raft::sparse::solver::LANCZOS_WHICH::SA, 42, 12, 12, 1}}; using LanczosTestF = lanczos_tests; TEST_P(LanczosTestF, Result) { Run(); } From 38d9c6947347e0eefe9ce3ac867f4c5f71a73cdd Mon Sep 17 00:00:00 2001 From: says1117 Date: Mon, 20 Jul 2026 05:02:33 -0400 Subject: [PATCH 2/7] Address CodeRabbit review: sync before host-scalar dot reads, add Doxygen docstrings --- cpp/tests/sparse/solver/lanczos.cu | 73 ++++++++++++++++++++++-------- 1 file changed, 53 insertions(+), 20 deletions(-) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index f803502b3e..c66aa2f8e3 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -11,6 +11,7 @@ #include #include #include +#include #include #include #include @@ -74,8 +75,19 @@ struct rmat_lanczos_inputs { float sparsity; }; -// Frobenius norm of the sparse matrix, used to scale residual/orthogonality -// tolerances to the magnitude of the problem instead of an absolute value. +/** + * @brief Compute the Frobenius norm of a sparse matrix from its stored values. + * + * Used to scale residual/orthogonality tolerances to the magnitude of the + * problem instead of an absolute value. + * + * @tparam ValueType data type of the matrix values (float or double) + * @tparam NnzType integral type of the nonzero count + * @param[in] handle raft resources + * @param[in] vals device pointer to the nnz stored values of the matrix + * @param[in] nnz number of stored values + * @return the Frobenius norm sqrt(sum(vals[i]^2)) + */ template ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* vals, NnzType nnz) { @@ -83,26 +95,42 @@ ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* raft::make_device_vector_view(vals, static_cast(nnz)); ValueType sum_sq{}; raft::linalg::dot(handle, vals_view, vals_view, raft::make_host_scalar_view(&sum_sq)); + raft::resource::sync_stream(handle); return std::sqrt(sum_sq); } -// Validates that (eigenvalues[i], eigenvectors[:, i]) are genuine eigenpairs -// of `A`, without relying on hardcoded golden eigenvalues. This checks the -// defining characteristics of a converged real-symmetric eigensolve: -// - ||A*v_i - lambda_i*v_i|| is within the noise floor implied by the -// solver's own convergence criterion (see lanczos_solver_config::tolerance, -// which is documented as governing exactly this residual). -// - lambda_i matches the Rayleigh quotient of v_i. This pins the eigenvalue -// down more sharply than the residual alone on matrices with a large -// ||A||, where a small eigenvalue's residual budget is dominated by the -// matrix scale. -// - each eigenvector is unit-normalized. -// - distinct eigenvectors are mutually orthogonal (catches "ghost" -// duplicate eigenpairs from loss of orthogonality, which a residual-only -// check would miss). -// - eigenvalues are returned in ascending algebraic order (verified -// against every previously-hardcoded golden array in this file, which -// were all ascending regardless of `which`). +/** + * @brief Validate that (eigenvalues[i], eigenvectors[:, i]) are genuine + * eigenpairs of `A`, without relying on hardcoded golden eigenvalues. + * + * Checks the defining characteristics of a converged real-symmetric + * eigensolve: + * - ||A*v_i - lambda_i*v_i|| is within the noise floor implied by the + * solver's own convergence criterion (see lanczos_solver_config::tolerance, + * which is documented as governing exactly this residual). + * - lambda_i matches the Rayleigh quotient of v_i. This pins the eigenvalue + * down more sharply than the residual alone on matrices with a large + * ||A||, where a small eigenvalue's residual budget is dominated by the + * matrix scale. + * - each eigenvector is unit-normalized. + * - distinct eigenvectors are mutually orthogonal (catches "ghost" + * duplicate eigenpairs from loss of orthogonality, which a residual-only + * check would miss). + * - eigenvalues are returned in ascending algebraic order (verified + * against every previously-hardcoded golden array in this file, which + * were all ascending regardless of `which`). + * + * @tparam IndexType integral type of the matrix indices + * @tparam ValueType data type of matrix values and eigenpairs (float or double) + * @param[in] handle raft resources + * @param[in] A sparse matrix wrapper providing the A*v product + * @param[in] eigenvalues device pointer to n_components computed eigenvalues + * @param[in] eigenvectors device pointer to the column-major n x n_components + * eigenvector matrix + * @param[in] n dimension of the (square) matrix + * @param[in] n_components number of eigenpairs to validate + * @param[in] frobenius_norm Frobenius norm of A, used to scale tolerances + */ template void expect_valid_eigenpairs(raft::resources const& handle, raft::spectral::matrix::sparse_matrix_t const& A, @@ -139,7 +167,7 @@ void expect_valid_eigenpairs(raft::resources const& handle, std::vector host_eigenvalues(n_components); raft::update_host(host_eigenvalues.data(), eigenvalues, n_components, stream); - RAFT_CUDA_TRY(cudaStreamSynchronize(stream)); + raft::resource::sync_stream(handle); for (int i = 1; i < n_components; ++i) { EXPECT_LE(host_eigenvalues[i - 1], host_eigenvalues[i] + matrix_tol) @@ -176,6 +204,10 @@ void expect_valid_eigenpairs(raft::resources const& handle, raft::linalg::dot( handle, v_i_view, v_i_view, raft::make_host_scalar_view(&v_norm_sq)); + // dot with a host-scalar output is not guaranteed to have synchronized; + // make the host-side reads below safe. + raft::resource::sync_stream(handle); + ValueType residual_norm = std::sqrt(residual_sq); ValueType v_norm = std::sqrt(v_norm_sq); @@ -202,6 +234,7 @@ void expect_valid_eigenpairs(raft::resources const& handle, ValueType dot_ij{}; raft::linalg::dot( handle, v_i_view, v_j_view, raft::make_host_scalar_view(&dot_ij)); + raft::resource::sync_stream(handle); EXPECT_LE(std::abs(dot_ij), ortho_tol) << "eigenvectors " << i << " and " << j << " are not orthogonal: |v_i . v_j|=" << std::abs(dot_ij); From 908c8641a75d05572a22a11f2461e4367329d10d Mon Sep 17 00:00:00 2001 From: says1117 Date: Wed, 22 Jul 2026 21:48:30 -0400 Subject: [PATCH 3/7] Relax SM eigenvector tolerances; keep SM eigenvalues verified via Rayleigh quotient CI showed LanczosTestD_SM failing on H100 and V100: the smallest-magnitude eigenVECTOR residual is ill-conditioned and varies across CUDA versions and GPU architectures, while the eigenVALUE stays correct on every platform (Rayleigh quotient passes everywhere). Split the checks into eigenvalue- accuracy (ordering, Rayleigh -- held tight for all modes) and eigenvector- accuracy (residual, norm, orthogonality -- relaxed for SM only). Non-SM and all eigenvalue tolerances are unchanged, so no previously-passing test can regress. See #2705, #2758. --- cpp/tests/sparse/solver/lanczos.cu | 116 ++++++++++++++++++----------- 1 file changed, 73 insertions(+), 43 deletions(-) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index c66aa2f8e3..c0369038fb 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -104,21 +104,32 @@ ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* * eigenpairs of `A`, without relying on hardcoded golden eigenvalues. * * Checks the defining characteristics of a converged real-symmetric - * eigensolve: - * - ||A*v_i - lambda_i*v_i|| is within the noise floor implied by the - * solver's own convergence criterion (see lanczos_solver_config::tolerance, - * which is documented as governing exactly this residual). - * - lambda_i matches the Rayleigh quotient of v_i. This pins the eigenvalue - * down more sharply than the residual alone on matrices with a large - * ||A||, where a small eigenvalue's residual budget is dominated by the - * matrix scale. - * - each eigenvector is unit-normalized. - * - distinct eigenvectors are mutually orthogonal (catches "ghost" - * duplicate eigenpairs from loss of orthogonality, which a residual-only - * check would miss). + * eigensolve, split into eigenVALUE-accuracy checks (held tight for every + * mode) and eigenVECTOR-accuracy checks (relaxed for smallest-magnitude + * problems, see below): + * - lambda_i matches the Rayleigh quotient (v_i^T A v_i)/(v_i^T v_i). This + * is the primary correctness check: it confirms the returned eigenvalue + * is the true eigenvalue for the returned eigenvector direction, and it + * is stationary at eigenvectors, so its error is second-order in the + * eigenvector error while the residual's is first-order. * - eigenvalues are returned in ascending algebraic order (verified * against every previously-hardcoded golden array in this file, which * were all ascending regardless of `which`). + * - ||A*v_i - lambda_i*v_i|| is small (see lanczos_solver_config::tolerance, + * documented as governing exactly this residual). + * - each eigenvector is unit-normalized. + * - distinct eigenvectors are mutually orthogonal (catches "ghost" + * duplicate eigenpairs from loss of orthogonality). + * + * The three eigenVECTOR checks (residual, normalization, orthogonality) are + * loosened for LANCZOS_WHICH::SM. Smallest-magnitude eigenvectors are + * ill-conditioned -- a tiny eigenvalue buried deep in the spectrum is highly + * sensitive to rounding -- so the eigenvector residual grows by orders of + * magnitude and varies across CUDA versions and GPU architectures, even + * though the eigenVALUE stays accurate (see issues #2705 and #2758). For SM + * we therefore assert eigenvalue correctness tightly via the Rayleigh + * quotient and only sanity-bound the eigenvector; for all other modes the + * eigenvector is well conditioned and held to the tight noise floor. * * @tparam IndexType integral type of the matrix indices * @tparam ValueType data type of matrix values and eigenpairs (float or double) @@ -130,6 +141,8 @@ ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* * @param[in] n dimension of the (square) matrix * @param[in] n_components number of eigenpairs to validate * @param[in] frobenius_norm Frobenius norm of A, used to scale tolerances + * @param[in] which which eigenvalues were requested; selects the tight or + * relaxed eigenvector tolerance (SM is relaxed) */ template void expect_valid_eigenpairs(raft::resources const& handle, @@ -138,39 +151,51 @@ void expect_valid_eigenpairs(raft::resources const& handle, ValueType* eigenvectors, IndexType n, int n_components, - ValueType frobenius_norm) + ValueType frobenius_norm, + raft::sparse::solver::LANCZOS_WHICH which) { auto stream = resource::get_cuda_stream(handle); auto eps = std::numeric_limits::epsilon(); - // Tolerances are expressed as multiples of the floating-point noise floor - // for each kind of quantity, since the exact low-order bits depend on - // cuSPARSE/cuBLAS reduction order, which varies by CUDA version and GPU - // architecture (the root cause of issue #2519). - // - // Quantities carrying the magnitude of A (the residual ||A*v - lambda*v|| - // and the Rayleigh quotient) have a noise floor proportional to - // ||A||_F * eps. Measured across all fixtures in this file on sm_75 / - // CUDA 13, these stay within 15 * ||A||_F * eps and show no dependence on - // n (n=100 and n=4096 both land in that range), so n is deliberately not a - // factor here -- including it would inflate the tolerance ~40x on the RMAT - // fixture and let a 1% eigenvalue error pass undetected. + // Eigenvalue-accuracy tolerance -- applied to the ascending-order check and + // the Rayleigh-quotient check. Both carry the magnitude of A, so the noise + // floor is proportional to ||A||_F * eps. Measured across all fixtures on + // sm_75/CUDA 13 these stay within ~15 * ||A||_F * eps with no dependence on + // n, so n is deliberately not a factor (including it would inflate the RMAT + // tolerance ~40x and let a 1% eigenvalue error pass). This is held tight for + // every mode: Lanczos returns accurate eigenVALUES even when the + // eigenVECTOR is poorly conditioned. + const ValueType eigenvalue_tol = ValueType(500) * std::max(frobenius_norm, ValueType(1)) * eps; + + // Eigenvector-accuracy tolerances -- applied to the residual, unit-norm and + // orthogonality checks. For smallest-magnitude (SM) eigenpairs the returned + // eigenVECTOR is ill-conditioned: its residual is observed to grow to ~5e-6 + // and to vary across CUDA versions / GPU architectures (H100, V100), while + // the eigenVALUE itself stays correct (see #2705, #2758). We therefore relax + // the eigenvector bounds for SM only -- correctness for SM is asserted + // through the tight Rayleigh-quotient check above -- and keep the tight + // floating-point noise floor for every other mode. // - // Dimensionless quantities (deviation from unit norm, mutual - // orthogonality) have a noise floor of a small multiple of eps; measured - // values stay within 100 * eps. - const ValueType matrix_tol = ValueType(500) * std::max(frobenius_norm, ValueType(1)) * eps; - const ValueType unit_tol = ValueType(1000) * eps; - const ValueType norm_tol = unit_tol; - const ValueType ortho_tol = unit_tol; - const ValueType residual_tol = matrix_tol; + // The SM residual reflects conditioning rather than rounding, so it is + // bounded by an absolute empirical fraction of ||A||_F (well above the worst + // observed CI residual, ~1e-5 for both float and double, yet far below any + // gross error). The unit-norm and orthogonality deviations are ordinary + // dot-product noise and scale with eps; the SM values seen in CI reach + // ~3000 * eps (double) and ~100 * eps (float), so 1e5 * eps covers both + // precisions with margin while still catching gross eigenvector corruption. + const bool is_sm = (which == raft::sparse::solver::LANCZOS_WHICH::SM); + const ValueType residual_tol = is_sm + ? ValueType(1e-4) * std::max(frobenius_norm, ValueType(1)) + : ValueType(500) * std::max(frobenius_norm, ValueType(1)) * eps; + const ValueType norm_tol = is_sm ? ValueType(1e5) * eps : ValueType(1000) * eps; + const ValueType ortho_tol = is_sm ? ValueType(1e5) * eps : ValueType(1000) * eps; std::vector host_eigenvalues(n_components); raft::update_host(host_eigenvalues.data(), eigenvalues, n_components, stream); raft::resource::sync_stream(handle); for (int i = 1; i < n_components; ++i) { - EXPECT_LE(host_eigenvalues[i - 1], host_eigenvalues[i] + matrix_tol) + EXPECT_LE(host_eigenvalues[i - 1], host_eigenvalues[i] + eigenvalue_tol) << "eigenvalues not returned in ascending algebraic order at index " << i; } @@ -219,11 +244,12 @@ void expect_valid_eigenpairs(raft::resources const& handle, << " lambda=" << host_eigenvalues[i] << " tol=" << residual_tol; // The Rayleigh quotient (v^T A v) / (v^T v) is the eigenvalue implied by - // the returned eigenvector. Comparing it against the returned eigenvalue - // pins down lambda independently of the residual check, which on a - // large-||A|| matrix is too coarse to notice a small relative error in a - // small-magnitude eigenvalue. - EXPECT_NEAR(host_eigenvalues[i], rayleigh / v_norm_sq, matrix_tol) + // the returned eigenvector. This is the primary eigenvalue-correctness + // check and is held tight for every mode (including SM, where it is the + // sole tight guarantee since the eigenvector residual is relaxed): a + // wrong eigenvalue -- e.g. a 1% perturbation -- fails here by many orders + // of magnitude regardless of eigenvector conditioning. + EXPECT_NEAR(host_eigenvalues[i], rayleigh / v_norm_sq, eigenvalue_tol) << "eigenvalue " << i << " disagrees with the Rayleigh quotient of its eigenvector: lambda=" << host_eigenvalues[i] << " rayleigh=" << rayleigh / v_norm_sq; @@ -371,7 +397,8 @@ class rmat_lanczos_tests eigenvectors.data_handle(), (IndexType)symmetric_coo.n_rows, n_components, - frobenius_norm); + frobenius_norm, + params.which); // Reproducibility test - run again with same seed and verify exact match raft::device_vector eigenvalues2 = @@ -427,7 +454,8 @@ class rmat_lanczos_tests eigenvectors_coo.data_handle(), (IndexType)symmetric_coo.n_rows, n_components, - frobenius_norm); + frobenius_norm, + params.which); } protected: @@ -520,7 +548,8 @@ class lanczos_tests : public ::testing::TestWithParam eigenvalues2 = @@ -575,7 +604,8 @@ class lanczos_tests : public ::testing::TestWithParam Date: Wed, 12 Aug 2026 00:04:40 -0400 Subject: [PATCH 4/7] Validate which-selection via independent dense eigensolver; use device_csr_matrix_view Address review: expect_valid_eigenpairs only confirmed returned pairs were *some* genuine eigenpairs of A, not that they were the correct SA/SM/LA/LM selection. Add compute_full_spectrum() (dense CSR->eigDC, independent of Lanczos) and expect_correct_selection() to check against it. Also switch the sparse-matrix parameter to device_csr_matrix_view, matching the type already used at the solve call sites, instead of the legacy sparse_matrix_t wrapper. --- cpp/tests/sparse/solver/lanczos.cu | 221 ++++++++++++++++++++++++----- 1 file changed, 189 insertions(+), 32 deletions(-) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index c0369038fb..e2940488ec 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -6,21 +6,25 @@ #include "../../test_utils.cuh" #include +#include #include #include #include #include #include #include +#include #include #include #include +#include #include #include #include #include #include #include +#include #include #include #include @@ -99,6 +103,23 @@ ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* return std::sqrt(sum_sq); } +/** + * @brief Absolute tolerance for eigenvalue-accuracy comparisons (Rayleigh + * quotient, ascending order, and selection correctness). Shared between + * expect_valid_eigenpairs() and expect_correct_selection() so both checks + * agree on how tightly an eigenvalue must be pinned down. See + * expect_valid_eigenpairs() for the derivation. + * + * @tparam ValueType data type of the eigenvalues (float or double) + * @param[in] frobenius_norm Frobenius norm of A + */ +template +ValueType eigenvalue_tolerance(ValueType frobenius_norm) +{ + auto eps = std::numeric_limits::epsilon(); + return ValueType(500) * std::max(frobenius_norm, ValueType(1)) * eps; +} + /** * @brief Validate that (eigenvalues[i], eigenvectors[:, i]) are genuine * eigenpairs of `A`, without relying on hardcoded golden eigenvalues. @@ -134,28 +155,40 @@ ValueType compute_frobenius_norm(raft::resources const& handle, ValueType const* * @tparam IndexType integral type of the matrix indices * @tparam ValueType data type of matrix values and eigenpairs (float or double) * @param[in] handle raft resources - * @param[in] A sparse matrix wrapper providing the A*v product + * @param[in] A CSR view of the (symmetric) matrix the eigenpairs belong to * @param[in] eigenvalues device pointer to n_components computed eigenvalues * @param[in] eigenvectors device pointer to the column-major n x n_components * eigenvector matrix - * @param[in] n dimension of the (square) matrix * @param[in] n_components number of eigenpairs to validate * @param[in] frobenius_norm Frobenius norm of A, used to scale tolerances * @param[in] which which eigenvalues were requested; selects the tight or * relaxed eigenvector tolerance (SM is relaxed) */ template -void expect_valid_eigenpairs(raft::resources const& handle, - raft::spectral::matrix::sparse_matrix_t const& A, - ValueType const* eigenvalues, - ValueType* eigenvectors, - IndexType n, - int n_components, - ValueType frobenius_norm, - raft::sparse::solver::LANCZOS_WHICH which) +void expect_valid_eigenpairs( + raft::resources const& handle, + raft::device_csr_matrix_view A, + ValueType const* eigenvalues, + ValueType* eigenvectors, + int n_components, + ValueType frobenius_norm, + raft::sparse::solver::LANCZOS_WHICH which) { - auto stream = resource::get_cuda_stream(handle); - auto eps = std::numeric_limits::epsilon(); + auto stream = resource::get_cuda_stream(handle); + auto eps = std::numeric_limits::epsilon(); + auto structure = A.structure_view(); + IndexType n = structure.get_n_rows(); + + // sparse_matrix_t provides the A*v product (mv()) used below; built here + // from the CSR view's own storage rather than requiring the caller to + // construct one separately. + raft::spectral::matrix::sparse_matrix_t const A_mv{ + handle, + structure.get_indptr().data(), + structure.get_indices().data(), + A.get_elements().data(), + n, + static_cast(structure.get_nnz())}; // Eigenvalue-accuracy tolerance -- applied to the ascending-order check and // the Rayleigh-quotient check. Both carry the magnitude of A, so the noise @@ -165,7 +198,7 @@ void expect_valid_eigenpairs(raft::resources const& handle, // tolerance ~40x and let a 1% eigenvalue error pass). This is held tight for // every mode: Lanczos returns accurate eigenVALUES even when the // eigenVECTOR is poorly conditioned. - const ValueType eigenvalue_tol = ValueType(500) * std::max(frobenius_norm, ValueType(1)) * eps; + const ValueType eigenvalue_tol = eigenvalue_tolerance(frobenius_norm); // Eigenvector-accuracy tolerances -- applied to the residual, unit-norm and // orthogonality checks. For smallest-magnitude (SM) eigenpairs the returned @@ -206,7 +239,7 @@ void expect_valid_eigenpairs(raft::resources const& handle, auto v_i_view = raft::make_device_vector_view(v_i, n); // residual = A * v_i - A.mv(ValueType{1}, v_i, ValueType{0}, residual.data_handle()); + A_mv.mv(ValueType{1}, v_i, ValueType{0}, residual.data_handle()); ValueType rayleigh{}; raft::linalg::dot(handle, @@ -268,6 +301,118 @@ void expect_valid_eigenpairs(raft::resources const& handle, } } +/** + * @brief Compute the full eigenvalue spectrum of a sparse symmetric matrix + * using an independent dense eigensolver (cuSOLVER syevd via + * raft::linalg::eigDC), for validating that Lanczos selected the correct + * extremal subset of the spectrum (see expect_correct_selection()). + * + * expect_valid_eigenpairs() only confirms each returned (lambda_i, v_i) is + * *some* genuine eigenpair of A; it cannot distinguish a correctly-selected + * eigenvalue from a genuine-but-wrong one (e.g. an interior eigenvalue + * returned when the smallest algebraic eigenvalues were requested). This + * spectrum, computed by a completely different algorithm than Lanczos and + * not subject to the same restart/subspace dynamics, is the independent + * reference that answers that question. + * + * @tparam IndexType integral type of the matrix indices + * @tparam ValueType data type of matrix values (float or double) + * @param[in] handle raft resources + * @param[in] A CSR view of the (symmetric) matrix + * @return ascending host vector of all n eigenvalues of A + */ +template +std::vector compute_full_spectrum( + raft::resources const& handle, + raft::device_csr_matrix_view A) +{ + auto stream = resource::get_cuda_stream(handle); + auto structure = A.structure_view(); + IndexType n = structure.get_n_rows(); + + auto dense = raft::make_device_matrix(handle, n, n); + RAFT_CUDA_TRY(cudaMemsetAsync(dense.data_handle(), 0, dense.size() * sizeof(ValueType), stream)); + raft::sparse::convert::csr_to_dense( + resource::get_cusparse_handle(handle), + n, + n, + static_cast(structure.get_nnz()), + structure.get_indptr().data(), + structure.get_indices().data(), + A.get_elements().data(), + n, + dense.data_handle(), + stream, + false); // column-major output; A is symmetric so row/col-major coincide anyway + + auto ref_vectors = raft::make_device_matrix(handle, n, n); + auto ref_values = raft::make_device_vector(handle, n); + raft::linalg::eigDC( + handle, dense.data_handle(), n, n, ref_vectors.data_handle(), ref_values.data_handle(), stream); + + std::vector host_spectrum(n); + raft::update_host(host_spectrum.data(), ref_values.data_handle(), n, stream); + raft::resource::sync_stream(handle); + return host_spectrum; +} + +/** + * @brief Validate that the returned eigenvalues are the requested extremal + * subset (SA/SM/LA/LM) of A's spectrum, not merely *some* genuine eigenpairs + * of A. See compute_full_spectrum() for why this needs an independent + * reference rather than reusing Lanczos's own output. + * + * @tparam ValueType data type of the eigenvalues (float or double) + * @param[in] handle raft resources + * @param[in] eigenvalues device pointer to n_components computed eigenvalues + * @param[in] n_components number of eigenpairs + * @param[in] which which eigenvalues were requested + * @param[in] full_spectrum ascending reference spectrum from compute_full_spectrum() + * @param[in] tol absolute tolerance for the comparison + */ +template +void expect_correct_selection(raft::resources const& handle, + ValueType const* eigenvalues, + int n_components, + raft::sparse::solver::LANCZOS_WHICH which, + std::vector const& full_spectrum, + ValueType tol) +{ + auto stream = resource::get_cuda_stream(handle); + std::vector host_eigenvalues(n_components); + raft::update_host(host_eigenvalues.data(), eigenvalues, n_components, stream); + raft::resource::sync_stream(handle); + + std::vector expected(n_components); + using raft::sparse::solver::LANCZOS_WHICH; + if (which == LANCZOS_WHICH::SA) { + std::copy(full_spectrum.begin(), full_spectrum.begin() + n_components, expected.begin()); + } else if (which == LANCZOS_WHICH::LA) { + std::copy(full_spectrum.end() - n_components, full_spectrum.end(), expected.begin()); + } else { + std::vector by_magnitude(full_spectrum); + std::stable_sort(by_magnitude.begin(), by_magnitude.end(), [](ValueType a, ValueType b) { + return std::abs(a) < std::abs(b); + }); + if (which == LANCZOS_WHICH::SM) { + std::copy(by_magnitude.begin(), by_magnitude.begin() + n_components, expected.begin()); + } else { // LM + std::copy(by_magnitude.end() - n_components, by_magnitude.end(), expected.begin()); + } + // Lanczos always returns its selection sorted ascending, regardless of + // `which`; sort the expected subset the same way for a like-for-like + // comparison. + std::sort(expected.begin(), expected.end()); + } + + for (int i = 0; i < n_components; ++i) { + EXPECT_NEAR(host_eigenvalues[i], expected[i], tol) + << "eigenvalue " << i << " (which=" << static_cast(which) + << ") is not part of the requested selection: got=" << host_eigenvalues[i] + << " expected (from independent dense reference spectrum)=" << expected[i]; + } +} + template class rmat_lanczos_tests : public ::testing::TestWithParam> { @@ -360,13 +505,6 @@ class rmat_lanczos_tests raft::make_device_matrix( handle, symmetric_coo.n_rows, n_components); - raft::spectral::matrix::sparse_matrix_t const csr_m{ - handle, - row_indices.data_handle(), - symmetric_coo.cols(), - symmetric_coo.vals(), - symmetric_coo.n_rows, - (uint64_t)symmetric_coo.nnz}; raft::sparse::solver::lanczos_solver_config config{ n_components, params.maxiter, params.restartiter, params.tol, params.which, rng.seed}; @@ -391,14 +529,18 @@ class rmat_lanczos_tests ValueType frobenius_norm = compute_frobenius_norm(handle, symmetric_coo.vals(), symmetric_coo.nnz); + std::vector full_spectrum = compute_full_spectrum(handle, csr_matrix); + ValueType eigenvalue_tol = eigenvalue_tolerance(frobenius_norm); + expect_valid_eigenpairs(handle, - csr_m, + csr_matrix, eigenvalues.data_handle(), eigenvectors.data_handle(), - (IndexType)symmetric_coo.n_rows, n_components, frobenius_norm, params.which); + expect_correct_selection( + handle, eigenvalues.data_handle(), n_components, params.which, full_spectrum, eigenvalue_tol); // Reproducibility test - run again with same seed and verify exact match raft::device_vector eigenvalues2 = @@ -449,13 +591,18 @@ class rmat_lanczos_tests eigenvectors_coo.view()); expect_valid_eigenpairs(handle, - csr_m, + csr_matrix, eigenvalues_coo.data_handle(), eigenvectors_coo.data_handle(), - (IndexType)symmetric_coo.n_rows, n_components, frobenius_norm, params.which); + expect_correct_selection(handle, + eigenvalues_coo.data_handle(), + n_components, + params.which, + full_spectrum, + eigenvalue_tol); } protected: @@ -530,9 +677,9 @@ class lanczos_tests : public ::testing::TestWithParam( const_cast(vals.data_handle()), csr_structure); - raft::spectral::matrix::sparse_matrix_t const A_check{ - handle, rows.data_handle(), cols.data_handle(), vals.data_handle(), n, (uint64_t)nnz}; - ValueType frobenius_norm = compute_frobenius_norm(handle, vals.data_handle(), nnz); + ValueType frobenius_norm = compute_frobenius_norm(handle, vals.data_handle(), nnz); + std::vector full_spectrum = compute_full_spectrum(handle, csr_matrix); + ValueType eigenvalue_tol = eigenvalue_tolerance(frobenius_norm); std::get<0>(stats) = raft::sparse::solver::lanczos_compute_eigenpairs( handle, @@ -543,13 +690,18 @@ class lanczos_tests : public ::testing::TestWithParam eigenvalues2 = @@ -599,13 +751,18 @@ class lanczos_tests : public ::testing::TestWithParam Date: Wed, 12 Aug 2026 01:14:52 -0400 Subject: [PATCH 5/7] Address nitpicks: document RMAT dense-check cost, use eig_dc mdspan API --- cpp/tests/sparse/solver/lanczos.cu | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index e2940488ec..867a744cc0 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -304,7 +304,7 @@ void expect_valid_eigenpairs( /** * @brief Compute the full eigenvalue spectrum of a sparse symmetric matrix * using an independent dense eigensolver (cuSOLVER syevd via - * raft::linalg::eigDC), for validating that Lanczos selected the correct + * raft::linalg::eig_dc), for validating that Lanczos selected the correct * extremal subset of the spectrum (see expect_correct_selection()). * * expect_valid_eigenpairs() only confirms each returned (lambda_i, v_i) is @@ -315,6 +315,12 @@ void expect_valid_eigenpairs( * not subject to the same restart/subspace dynamics, is the independent * reference that answers that question. * + * Densifies A and runs an O(n^3) dense eigensolve, so this scales poorly for + * large n: at the RMAT fixture's current r_scale = c_scale = 12 (n = 4096) + * this allocates two ~64 MiB dense matrices and adds a few seconds per test. + * If RMAT test scales grow significantly, prefer bounding the selection + * check to a smaller fixture over raising this one's n. + * * @tparam IndexType integral type of the matrix indices * @tparam ValueType data type of matrix values (float or double) * @param[in] handle raft resources @@ -347,8 +353,8 @@ std::vector compute_full_spectrum( auto ref_vectors = raft::make_device_matrix(handle, n, n); auto ref_values = raft::make_device_vector(handle, n); - raft::linalg::eigDC( - handle, dense.data_handle(), n, n, ref_vectors.data_handle(), ref_values.data_handle(), stream); + raft::linalg::eig_dc( + handle, raft::make_const_mdspan(dense.view()), ref_vectors.view(), ref_values.view()); std::vector host_spectrum(n); raft::update_host(host_spectrum.data(), ref_values.data_handle(), n, stream); From bac2b7711d037b3823ef15576275b94a298429cc Mon Sep 17 00:00:00 2001 From: says1117 Date: Wed, 12 Aug 2026 02:26:22 -0400 Subject: [PATCH 6/7] Guard subset extraction against n_components >= full_spectrum.size() --- cpp/tests/sparse/solver/lanczos.cu | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index 867a744cc0..53ab7e4086 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -384,6 +384,16 @@ void expect_correct_selection(raft::resources const& handle, std::vector const& full_spectrum, ValueType tol) { + // The subset-extraction branches below advance iterators by n_components + // in both directions (begin()+n_components, end()-n_components); with + // n_components >= full_spectrum.size() that arithmetic goes out of bounds, + // which is undefined behavior rather than a clean test failure. No current + // fixture violates this (n_components is always well under n), but assert + // it so a future fixture that does fails clearly instead of silently. + ASSERT_LE(n_components, static_cast(full_spectrum.size())) + << "n_components (" << n_components << ") must not exceed the matrix dimension (" + << full_spectrum.size() << ")"; + auto stream = resource::get_cuda_stream(handle); std::vector host_eigenvalues(n_components); raft::update_host(host_eigenvalues.data(), eigenvalues, n_components, stream); From 45e4085baf198d68b31537a74cc0889923dd9f29 Mon Sep 17 00:00:00 2001 From: says1117 Date: Wed, 12 Aug 2026 02:40:47 -0400 Subject: [PATCH 7/7] Document RMAT scale limit next to rmat_inputsf --- cpp/tests/sparse/solver/lanczos.cu | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/cpp/tests/sparse/solver/lanczos.cu b/cpp/tests/sparse/solver/lanczos.cu index 53ab7e4086..0b005ab922 100644 --- a/cpp/tests/sparse/solver/lanczos.cu +++ b/cpp/tests/sparse/solver/lanczos.cu @@ -1219,6 +1219,11 @@ const std::vector> inputsf_SA = { cols(), vals()}}; +// NOTE: r_scale/c_scale = 12 gives n = 4096. expect_correct_selection() validates +// against a dense reference spectrum (see compute_full_spectrum()), which costs +// O(n^2) memory and O(n^3) runtime -- at n = 4096 that is two ~64 MiB matrices and +// a few seconds. Raising these scales grows that cost cubically; prefer adding a +// separate smaller fixture for the selection check over increasing them here. const std::vector> rmat_inputsf = { {50, 100, 10000, 0, 0, 1e-9, raft::sparse::solver::LANCZOS_WHICH::SA, 42, 12, 12, 1}};