Tests: growth_tol for LDLᵀ, seeds, SparseMatrixCSC constructors, blocked regime-B coverage, refinement bound

TASKS.md

@@ -2184,8 +2184,8 @@ repeat it once T21 has landed. The bar is cuDSS with `"matching_alg" =
 
 ### Report
 
-- Status: [!] (done on the CPU backend; the "100× on a badly scaled SPD matrix" assertion is `@test_broken`, see
-  deviations; CUDA from CI)
+- Status: [!] (done on the CPU backend; the "100× on a badly scaled SPD matrix" assertion is replaced by what holds,
+  see deviations; CUDA from CI)
 - What was built:
   - `src/solve/refinement.jl` (new): `RefinementWorkspace` (device CSR of the full matrix over the contributions
     of the user's `nzval`, signed source index = conjugated mirror entry for `"H"`/`"HPD"`, built from the existing
@@ -2251,7 +2251,9 @@ repeat it once T21 has landed. The bar is cuDSS with `"matching_alg" =
     is backward stable and invariant under symmetric diagonal scaling (the only scaling that keeps the matrix
     SPD). Measured on 3 generators × 3 scalings × 5 right-hand-side choices: the unrefined relres is either at
     rounding level (`b = A x`) or limited by `cond(A)` (random `b`), and one step gains ≤ 7× (once 42×). The SPD
-    assertion is `@test_broken` with that comment; the 100× reduction is asserted on a KKT matrix whose
+    assertion was `@test_broken` with that comment; since the post-T21 test hardening (#79 review: a backend with
+    other rounding could flip it to an unexpected pass) it asserts what holds instead, `r0 ≤ 100 eps` and
+    `r1 ≤ r0`; the 100× reduction is asserted on a KKT matrix whose
     factorization carries static pivot perturbations (the case refinement is meant for, issue #71), together with
     the `ir_tol`/early-exit/steps-performed checks (also run on the SPD matrix, where they pass).
   - `ir_tol` is the largest `‖Rₖ‖₂/‖Bₖ‖₂` over the right-hand sides; checking it costs one host synchronization

test/test_api.jl

@@ -63,6 +63,42 @@ end
     @test maximum(batch_relres([A, A], api_solve(backend, sb, bb), bb)) <= tol(T)
 end
 
+@testset "SparseMatrixCSC constructor and update! (CPU, $T, $INT)" for backend in filter(b -> b isa CPU, BACKENDS),
+                                                                       T in ELTYPES, INT in INTTYPES
+    # a SparseMatrixCSC is converted to host CSR arrays (index type kept) once, for the CPU backend
+    Random.seed!(666)
+    n = 60
+    A = random_spd(T, n, 0.05)
+    s = spd_structure(T)
+    for (view, index) in (('L', 'O'), ('U', 'Z'), ('F', 'O'))
+        solver = DirectSolver(SparseMatrixCSC{T, INT}(triangle_view(A, view)), s, view; index)
+        @test solver isa DirectSolver{T, INT}
+        @test solver.backend isa CPU
+        @test solver.A.index == SDS._index_base(index)
+        @test size(solver) == (n, n)
+        b = rand(T, n)
+        execute!("analysis", solver, nothing, nothing)
+        execute!("factorization", solver, nothing, nothing)
+        @test relres(A, api_solve(backend, solver, b), b) <= tol(T)
+        # new values of the same pattern, then refactorization
+        A2 = A + 2 * I
+        update!(solver, SparseMatrixCSC{T, INT}(triangle_view(A2, view)))
+        execute!("refactorization", solver, nothing, nothing)
+        @test getparam(solver, "info") == 0
+        @test relres(A2, api_solve(backend, solver, b), b) <= tol(T)
+        # a different size or number of stored entries
+        @test thrown(() -> update!(solver, SparseMatrixCSC{T, INT}(triangle_view(random_spd(T, n + 1, 0.05), view)))) isa
+              InvalidValueError
+        @test thrown(() -> update!(solver, SparseMatrixCSC{T, INT}(triangle_view(A2 + sprand(T, n, n, 0.2), view)))) isa
+              InvalidValueError
+        # another index type
+        @test thrown(() -> update!(solver, SparseMatrixCSC{T, INT == Int32 ? Int64 : Int32}(triangle_view(A2, view)))) isa
+              InvalidValueError
+    end
+    # a rectangular matrix
+    @test thrown(() -> DirectSolver(SparseMatrixCSC{T, INT}(sprand(T, n, n + 1, 0.05)), s, 'F')) isa InvalidValueError
+end
+
 @testset "phases ($(backend_name(backend)), $T, $INT)" for backend in BACKENDS, T in ELTYPES, INT in INTTYPES
     A = laplacian2d(T, 15, 12) + spdiagm(0 => rand(real(T), 180))
     n = size(A, 1)

test/test_fgmres.jl

@@ -37,6 +37,7 @@ end
 
 @testset "FGMRES-IR: element types, multiple right-hand sides, solve_mode ($(backend_name(backend)), $T)" for backend in BACKENDS,
                                                                                                                T in ELTYPES
+    Random.seed!(666)
     n = 200
     for (A, structure) in ((random_spd(T, n, 0.02), spd_structure(T)), (random_symindef(T, n, 0.02), sym_structure(T)))
         solver = ir_solver(backend, A, structure; params = (("ir_mode", "fgmres"), ("deterministic_mode", 1)))
@@ -74,6 +75,7 @@ end
 
 @testset "FGMRES-IR: uniform batch ($(backend_name(backend)), $T)" for backend in BACKENDS,
                                                                        T in filter(in((Float64, ComplexF32)), ELTYPES)
+    Random.seed!(666)
     n, nb, nrhs = 120, 3, 2
     members = batch_members(random_symindef(T, n, 0.03), nb)
     solver = DirectSolver(api_batch_matrix(backend, members, 'L'), sym_structure(T), 'L')

test/test_matching.jl

@@ -13,6 +13,7 @@ host_matching(A::SparseMatrixCSC, structure, alg; view = 'F') =
                                                           Options(matching_alg = alg); view))
 
 RUN_SHARED && @testset "MC64 jobs against brute force" begin
+    Random.seed!(666)
     for trial in 1:12
         n = 6
         A = sprand(n, n, 0.4) + sparse(randperm(n), 1:n, rand(n) .+ 0.5)   # structurally nonsingular
@@ -63,6 +64,7 @@ RUN_SHARED && @testset "MC64 jobs against brute force" begin
 end
 
 RUN_SHARED && @testset "symmetric scaling and matching pairs" begin
+    Random.seed!(666)
     nh, nj = 40, 15
     K = kkt_matrix(Float64, nh, nj, 1.0e-10)
     for view in ('L', 'U', 'F')
@@ -101,6 +103,7 @@ RUN_SHARED && @testset "symmetric scaling and matching pairs" begin
 end
 
 @testset "LU with matching: badly scaled ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     n = 80
     INT = Int32
     A = badly_scaled_general(T, n, 0.05)
@@ -209,6 +212,7 @@ end
 end
 
 @testset "inertia with matching ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     INT = Int32
     structure = sym_structure(T)
     for (name, A) in (("random_symindef(100)", random_symindef(T, 100, 0.05)),
@@ -260,6 +264,7 @@ end
 end
 
 @testset "uniform batch with matching ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     n = 50
     A = badly_scaled_general(T, n, 0.08)
     members = [A, 2 * A, 3 * A]

test/test_numeric_cholesky_a.jl

@@ -31,6 +31,7 @@ function subtree_move_rounds(S, v)
 end
 
 RUN_SHARED && @testset "plan: regime-A groups" begin
+    Random.seed!(666)
     for (name, A, opts) in numeric_a_matrices(Float64)
         S, _, _, _, Nd, _ = numeric_setup(CPU(), A; opts)
         sc, plan = S.schedule, Nd.plan
@@ -69,6 +70,7 @@ RUN_SHARED && @testset "plan: regime-A groups" begin
 end
 
 @testset "panels, solves, determinism ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     for (name, A, opts) in numeric_a_matrices(T), INT in (name == "laplacian2d(40,40) A+B+C" ? INTTYPES : (Int32,))
         @testset "$name $INT" begin
             S, Nr, info_ref, Sd, Nd, nz = numeric_setup(backend, A, INT; opts)
@@ -95,6 +97,7 @@ end
 end
 
 @testset "budget classes and kernel sizes ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     # the same matrix under different budgets: other subtrees, other local sizes, same factor
     A = laplacian3d(T, 8, 8, 8)
     for opts in (Options(subtree_budgets = [8192], subtree_parallelism = 0),
@@ -114,6 +117,7 @@ end
 
 @testset "nlaunches on laplacian2d(100, 100) with AMD ($(backend_name(backend)), $T)" for backend in BACKENDS,
                                                                                             T in ELTYPES
+    Random.seed!(666)
     A = laplacian2d(T, 100, 100)
     b = rand(T, size(A, 1), 2)
     x = map((Options(reordering_alg = "algo3", subtree_parallelism = 0),
@@ -137,6 +141,7 @@ end
 end
 
 @testset "views, index bases, refactorization ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     A = random_spd(T, 300, 0.02)
     for opts in (NUMERIC_A_OPTS, NUMERIC_ABC_OPTS)
         S, Nr, _, Sd, Nd, nz = numeric_setup(backend, A; opts)
@@ -172,6 +177,7 @@ end
 end
 
 @testset "info: first non-positive pivot ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     n, j = 60, 23
     for opts in (NUMERIC_A_OPTS, Options(reordering_alg = "algo5", subtree_parallelism = 0),
                  Options(use_superpanels = 0, subtree_parallelism = 0), NUMERIC_ABC_OPTS)

test/test_numeric_cholesky_b.jl

@@ -49,6 +49,7 @@ end
 
 # synthetic assembled fronts through the dense part of the kernel, against LAPACK on the host
 @testset "front_cholesky! on synthetic fronts ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     for W in SDS.REGIME_B_WIDTHS
         shapes = [(w, m, cb) for w in unique((1, W ÷ 2 + 1, W)) for (m, cb) in ((0, false), (1, true), (37, true),
                                                                                (5, false))]
@@ -115,6 +116,7 @@ end
 end
 
 @testset "panels, solves, determinism ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     for (name, A, opts) in numeric_b_matrices(T), INT in (name == "laplacian2d(40,40)" ? INTTYPES : (Int32,))
         @testset "$name $INT" begin
             S, Nr, info_ref, Sd, Nd, nz = numeric_setup(backend, A, INT; opts)
@@ -139,6 +141,7 @@ end
 end
 
 @testset "factorization_alg algo1 vs algo2 ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     for A in (laplacian2d(T, 40, 40), random_spd(T, 500, 0.01))
         b = rand(T, size(A, 1), 2)
         x = map(("algo1", "algo2")) do alg
@@ -156,6 +159,7 @@ end
 end
 
 @testset "views, index bases, refactorization ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     A = random_spd(T, 300, 0.02)
     for opts in (NUMERIC_B_OPTS, NUMERIC_BC_OPTS)
         S, Nr, _, Sd, Nd, nz = numeric_setup(backend, A; opts)
@@ -191,6 +195,7 @@ end
 end
 
 @testset "info: first non-positive pivot ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     n, j = 60, 23
     for opts in (NUMERIC_B_OPTS, Options(subtree_budgets = Int[], reordering_alg = "algo5"),
                  Options(subtree_budgets = Int[], use_superpanels = 0), NUMERIC_BC_OPTS)

test/test_numeric_cholesky_c.jl

@@ -46,6 +46,7 @@ RUN_SHARED && @testset "storage and plan" begin
 end
 
 @testset "panels, solves, determinism ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     for (name, A) in numeric_c_matrices(T), INT in (name == "laplacian2d(40,40)" ? INTTYPES : (Int32,))
         @testset "$name $INT" begin
             S, Nr, info_ref, Sd, Nd, nz = numeric_c_setup(backend, A, INT)
@@ -70,6 +71,7 @@ end
 end
 
 @testset "dense implementations ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     A = laplacian2d(T, 25, 25)
     S, Nr, _, Sd, Nd, nz = numeric_c_setup(backend, A)
     b = rand(T, size(A, 1), 2)
@@ -85,6 +87,7 @@ end
 end
 
 @testset "views, index bases, refactorization ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     A = random_spd(T, 300, 0.02)
     S, Nr, _, Sd, Nd, nz = numeric_c_setup(backend, A)
     @test SDS.factorize!(Nd, Sd, nz) == 0
@@ -118,6 +121,7 @@ end
 end
 
 @testset "info: first non-positive pivot ($(backend_name(backend)), $T)" for backend in BACKENDS, T in ELTYPES
+    Random.seed!(666)
     n, j = 60, 23
     for opts in (NUMERIC_C_OPTS,
                  Options(subtree_budgets = Int[], factorization_alg = "algo2", reordering_alg = "algo5"),