Use the sparse kernels in mul! of two sparse matrices into a caller-supplied dense matrix

[ViralBShah]

Rebased on #865, which merged the rest of what this PR first proposed for #35 (any AbstractVecOrMat dense factor on either side, and no more generic_matmatmul_wrapper! methods). What is left is the case #865 kept on the generic product: both factors sparse, destination dense. This is only about mul! into a dense matrix the caller supplies; A * B of two sparse matrices is untouched and still returns a sparse matrix.

using LinearAlgebra, SparseArrays
n = 500; A = sprandn(n, n, 5/n); B = sprandn(n, n, 5/n); C = zeros(n, n)
@time mul!(C, A', B, 2.0, 0.5)             # 636 ms; 865 ms with Symmetric(B)

mul!(C::StridedMatrix, tA, tB, A::sparse, B::sparse, α, β) called LinearAlgebra._generic_matmatmul!, which reads both factors elementwise; with a transposed left factor that is a sparse getindex per entry of an n × n sweep. It now calls spdensemul!, which walks the stored columns of the left factor and only indexes the right one. Two cases keep a cheaper route: A' * B' / transpose(A) * transpose(B) run as (B * A)' into the wrapped destination (real and complex eltypes only, since that needs commutativity), and symmetric/Hermitian × symmetric/Hermitian stays on the generic product, which measured faster. Results now follow the sparse kernels' convention of skipping structural zeros of the left factor, as sparse × dense already does. No new methods, so dispatch is unchanged; detect_ambiguities against LinearAlgebra and Base reports none.

mul!(C, fA, fB, 2.0, 0.5), n = 500, 5 nonzeros per column, nightly, ms main → PR, 0 allocations throughout:

	B	B'	Symmetric(B)	view(B, :, 1:n)
A	5.01 → 1.64	5.31 → 1.78	5.28 → 1.89	5.03 → 1.61
A'	636 → 7.07	1.72 → 1.73	865 → 10.8	637 → 7.13
Symmetric(A)	6.71 → 6.34	7.16 → 5.39	6.87 → 6.81	6.56 → 6.43

The benchmark tables of #865, #866 and #867 (n = 3000 sparse-dense products, dot on column views, dense × sparse triangular) rerun on main and on this branch agree within run-to-run noise, with identical allocations.

Tests extend the #865 testset: each sparse factor times S' into a dense destination, real and complex, and two stored-zero proofs that fail on main. Full Pkg.test() and test/ambiguous.jl pass on nightly.

Left out: a kernel that walks both patterns, which would avoid the remaining sparse getindex per stored entry of the left factor and output column; Symmetric(A) * SparseVector still reaches the generic fallback.

Written by Claude Code (Claude Fable 5.1).

🤖 Generated with Claude Code

[codecov]

Codecov Report

❌ Patch coverage is 85.71429% with 1 line in your changes missing coverage. Please review.
✅ Project coverage is 93.01%. Comparing base (72d9e9f) to head (f4740f8).

Files with missing lines	Patch %	Lines
src/linalg.jl	85.71%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #868      +/-   ##
==========================================
+ Coverage   92.97%   93.01%   +0.03%     
==========================================
  Files          12       12              
  Lines        9072     9076       +4     
==========================================
+ Hits         8435     8442       +7     
+ Misses        637      634       -3     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.