Skip to content

Make lu/qr of ComponentMatrix preserve the wrapper on factors - #391

Merged
ChrisRackauckas merged 1 commit into
SciML:mainfrom
ChrisRackauckas-Claude:fix-lu-componentmatrix
Jul 11, 2026
Merged

Make lu/qr of ComponentMatrix preserve the wrapper on factors#391
ChrisRackauckas merged 1 commit into
SciML:mainfrom
ChrisRackauckas-Claude:fix-lu-componentmatrix

Conversation

@ChrisRackauckas-Claude

Copy link
Copy Markdown
Member

Note

Please ignore this PR until it has been reviewed by @ChrisRackauckas.

Fixes the Downstream CI failure seen in https://github.com/SciML/ComponentArrays.jl/actions/runs/29140401332/job/86512903689 (Issue 36 test in test/Downstream/diffeq_tests.jl).

Root cause

LinearSolve >= v3.87 (and v5) types its LUFactorization cache slot from ArrayInterface.lu_instance(A) but stores lu(A, check = false) when not aliasing A. For a ComponentMatrix:

  • ArrayInterface.lu_instance(::ComponentMatrix) (defined here) returns an LU wrapping a ComponentMatrix,
  • lu!(A) also keeps the wrapper (ComponentMatrix is strided, so it hits LAPACK directly),
  • but lu(A) copies through similar(A, T, dims), which drops the axes and demotes to a plain Matrix.

Storing the LU{..., Matrix, ...} into the LU{..., ComponentMatrix, ...}-typed slot throws TypeError: in setfield!, which is exactly the downstream error:

TypeError: in setfield!, expected LinearAlgebra.LU{Float64, ComponentArrays.ComponentMatrix{...}, Vector{Int64}}, got a value of type LinearAlgebra.LU{Float64, Matrix{Float64}, Vector{Int64}}

This reproduces without OrdinaryDiffEq via a plain solve(LinearProblem(A, b)) on a 2000x2000 ComponentMatrix.

Fix

Make lu(::ComponentMatrix) factorize getdata(A) (keeping the LAPACK fast path) and rewrap the factors, so lu, lu!, and lu_instance all return the same type. The same treatment is applied to qr/ArrayInterface.qr_instance, because LinearSolve's default solver falls back to QRFactorization on singular matrices (safetyfallback) and do_factorization uses the in-place qr!(A, pivot) there, which hits the identical mismatch (slot typed QRCompactWY{..., Matrix, ...} from qr_instance, value QRCompactWY{..., ComponentMatrix, ...}) — i.e. any ODE solve whose W matrix goes singular would crash instead of falling back.

Note GenericLUFactorization (chosen by the default solver for sizes <= 10) stores generic_lufact!(A), which wraps A itself, so its slot genuinely must be LU{..., ComponentMatrix, ...} — that is why the fix preserves the wrapper everywhere rather than demoting lu_instance to Matrix (which would just move the mismatch to the small-size path that currently passes).

Verification (all run locally)

  • MWE solve(LinearProblem(A, b)) on a 2000x2000 ComponentMatrix: fails with the CI TypeError before, passes after (LinearSolve v5.0.0, Julia 1.12.6).
  • Default solver sweep over sizes 4-2000 with alias_A in (false, true): all pass.
  • Explicit LUFactorization, GenericLUFactorization, QRFactorization, KrylovJL_GMRES: all pass.
  • Singular-matrix default solve (LU -> QR fallback): no longer throws, returns Success with the least-squares fallback.
  • Type consistency typeof(f(A)) === typeof(f!(copy(A))) === typeof(instance(A)) for lu and qr (NoPivot and ColumnNorm): verified on Julia 1.12.6 and 1.10.11.
  • GROUP=Core Pkg.test() (Julia 1.12): pass, including new regression tests in test/math_tests.jl.
  • GROUP=Downstream Pkg.test() (Julia 1.12): 6 passed, 1 broken (the broken one is the pre-existing @test_broken Sundials Issue 53), versus 4 passed, 1 errored, 1 broken on CI before this change.

Runic has been run on all changed files. Version bumped to 0.15.42.

🤖 Generated with Claude Code

https://claude.ai/code/session_01PZG7J62sLGGMUmvNsKHVMA

LinearSolve >= 3.87 types its LUFactorization cache slot from
ArrayInterface.lu_instance(A) but stores lu(A, check=false) when not
aliasing A. For a ComponentMatrix, lu_instance returned an LU wrapping a
ComponentMatrix while lu(A) demotes to Matrix (its copy goes through
similar(A, T, dims), which drops the axes), so storing the factorization
threw "TypeError: in setfield!" and broke the Downstream diffeq tests.

Define lu(::ComponentMatrix) to factorize getdata(A) (LAPACK fast path)
and rewrap the factors, so lu, lu! (already wrapper-preserving since
ComponentMatrix is strided), and lu_instance all return the same type.
Apply the same treatment to qr/qr_instance since LinearSolve's default
solver falls back to QR on singular matrices and hits the identical
mismatch there.

Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01PZG7J62sLGGMUmvNsKHVMA
@ChrisRackauckas
ChrisRackauckas marked this pull request as ready for review July 11, 2026 22:40
@ChrisRackauckas
ChrisRackauckas merged commit 652fd7c into SciML:main Jul 11, 2026
19 of 20 checks passed
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.

2 participants