Skip to content

Fix sparse CSC matrix multiplication - #602

Merged
michel2323 merged 2 commits into
mainfrom
fix-sparse-csc-mul
Aug 2, 2026
Merged

Fix sparse CSC matrix multiplication#602
michel2323 merged 2 commits into
mainfrom
fix-sparse-csc-mul

Conversation

@michel2323

@michel2323 michel2323 commented Aug 2, 2026

Copy link
Copy Markdown
Member

Two fixes to sparse CSC multiplication through the LinearAlgebra interface.

Transposed result. dA * x returned Aᵀx for a oneSparseMatrixCSC (same for dA * B and mul!). A CSC handle stores CSR(Aᵀ), so sparse_gemv!/sparse_gemm! already apply flip_trans internally. generic_matvecmul!/generic_matmatmul! applied it a second time, so the two cancelled and oneMKL computed S*x = Aᵀx. CSR was unaffected.

Complex element types. The CSC methods were constrained to BlasReal, so complex CSC never reached them and failed with CanonicalIndexError: getindex not defined. The underlying wrappers already have complex CSC specializations handling the 'C' conjugation identities, so widening to BlasFloat to match CSR is sufficient.

`dA * x` and `dA * B` returned Aᵀx and AᵀB for a oneSparseMatrixCSC.

A CSC handle stores CSR(Aᵀ), so sparse_gemv! and sparse_gemm! apply
flip_trans when calling into oneMKL and already present op(A) semantics
to their callers. generic_matvecmul! and generic_matmatmul! applied
flip_trans a second time, so 'N' reached oneMKL unchanged and it
computed S*x = Aᵀx. Pass the transpose flag through instead.

The existing sparse tests call sparse_gemv!/sparse_gemm! directly and
never exercised `*` or `mul!` on a sparse matrix, which is why this went
unnoticed; add a testset covering that path for both CSR and CSC.
The CSC generic_matvecmul!/generic_matmatmul! methods were constrained to
BlasReal, so `dA * x` on a complex oneSparseMatrixCSC never reached them
and fell through to the AbstractArray fallback, failing with

    CanonicalIndexError: getindex not defined for oneSparseMatrixCSC{ComplexF32, Int64}

The underlying sparse_gemv!/sparse_gemm! wrappers already carry complex
CSC specializations that handle the conjugation identities for the 'C'
case, so widening the constraint to BlasFloat, matching the CSR methods,
is all that is needed.

Extend the multiplication testset over identity/transpose/adjoint for
every supported element type, which covers the complex-only distinction
between transpose and adjoint.
@michel2323 michel2323 changed the title Fix transposed result from sparse CSC matrix multiplication Fix sparse CSC matrix multiplication Aug 2, 2026
@michel2323
michel2323 enabled auto-merge (squash) August 2, 2026 02:43
@michel2323
michel2323 merged commit 0a71d96 into main Aug 2, 2026
5 checks passed
@michel2323
michel2323 deleted the fix-sparse-csc-mul branch August 2, 2026 06:33
@codecov

codecov Bot commented Aug 2, 2026

Copy link
Copy Markdown

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 79.76%. Comparing base (54637d7) to head (1f20d19).
⚠️ Report is 2 commits behind head on main.

Additional details and impacted files
@@            Coverage Diff             @@
##             main     #602      +/-   ##
==========================================
+ Coverage   79.15%   79.76%   +0.61%     
==========================================
  Files          50       50              
  Lines        3391     3391              
==========================================
+ Hits         2684     2705      +21     
+ Misses        707      686      -21     

☔ 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.

michel2323 referenced this pull request Aug 3, 2026
Aurora and other large deployments run Intel's long-term-servicing branch of
the Compute Runtime rather than the rolling release oneAPI.jl targets. That
branch predates a number of driver and IGC fixes, several of which corrupt
results silently rather than raise an error. Add the workarounds behind a
single opt-in switch, oneL0.LTS[], resolved from ONEAPI_LTS at load time and
defaulting to off, so the rolling stack behaves exactly as before.

Gated on that switch:

* SPIR-V codegen goes through the Khronos translator instead of LLVM's
  SPIR-V back-end, whose output the LTS NEO/IGC runtime does not accept
  (#491). Both JLLs stay dependencies; GPUCompiler resolves the tool from
  the target's `backend` field, so the choice is made at compile time.

* BFloat16 is forced off. The LTS SPIR-V stack cannot codegen native bfloat
  in generic kernels -- a kernel that merely keeps a bfloat value fails with
  InvalidIRError, and declaring SPV_KHR_bfloat16 crashes the runtime.
  _device_supports_bfloat16() reports hardware capability and does not
  capture this, so the test suite and examples/bfloat16.jl gate on the
  switch instead.

* Reductions avoid the strided-read miscompile that silently breaks
  `sum(transpose(x))`, `a == transpose(b)` and `ishermitian`. Non-dense
  inputs are materialized before reducing, and reductions that keep the
  contiguous leading dimension use a coalesced one-work-item-per-slice
  kernel. Reductions that also reduce dim 1 keep a contiguous innermost
  axis and are left alone; an Int32 regression test covers that boundary.

* Buffers are freed only after draining the queues that may still reference
  them. LTS NEO advertises ZE_extension_memory_free_policies but ignores
  BLOCKING_FREE, so a GC-driven free of in-flight work faults and bans the
  context, surfacing later as a ZE_RESULT_ERROR_UNKNOWN at an unrelated
  call. A per-(context, device) registry tracks every queue -- including the
  replacement KA.priority! installs -- and the queue finalizer drains before
  destroying, with a bounded wait so a task that dies mid-submission cannot
  hang finalization.

Independently of that switch, ONEAPI_SYNC_EACH_SUBMISSION=1 synchronizes
after every command-list submission, working around a dropped-tail
corruption seen when several processes oversubscribe a single tile. It costs
roughly 3x throughput and is off by default, with getter, setter and scoped
forms available at runtime.

Also, not specific to the LTS stack: dlopen the NEO driver by full path
during __init__ so libsycl's bundled Level Zero loader finds it in-process.
Setting LD_LIBRARY_PATH there never could -- glibc captures it once at
process startup -- so it only ever served spawned worker processes.

The self-hosted Aurora runner is configured to exercise the LTS path, with
buildkite continuing to cover the rolling stack, and docs/src/lts.md
documents both switches, what they change, and what they cost.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant