lu: pivot vector is always Vector{BlasInt} in generic and Tridiagonal paths

[devmotion]

lu returns an LU factorization whose ipiv field is hardcoded to Vector{BlasInt} in two generic code paths in src/lu.jl:

1. src/lu.jl#L211 — the default for the ipiv positional argument of generic_lufact!. Triggered by any AbstractMatrix with a non-BlasFloat element type (Rational, BigFloat, ForwardDiff.Dual, …).
2. src/lu.jl#L636 — inside lu!(::Tridiagonal). Triggered for all element types, since Tridiagonal lu has no LAPACK fast-path and always runs _lu_tridiag!.

This forces the resulting LU factorization into an asymmetric state: factors reflects the input array type (via _lucopy/copy_similar for dense, or stays as-is for Tridiagonal), while ipiv is always a Vector. By contrast, the dense BLAS path at src/lapack.jl#L787 already uses similar(A, BlasInt, min(m, n)) and produces a pivot of matching array type.

MWE 1 — Tridiagonal (any element type)

using LinearAlgebra, FixedSizeArrays  # FixedSizeArrays v1.3.0

dl = FixedSizeVector{Float64}(undef, 2); copyto!(dl, [1.0, 1.0])
d  = FixedSizeVector{Float64}(undef, 3); copyto!(d,  [4.0, 4.0, 4.0])
du = FixedSizeVector{Float64}(undef, 2); copyto!(du, [1.0, 1.0])
T  = Tridiagonal(dl, d, du)
# Tridiagonal{Float64, FixedSizeArray{Float64, 1, Memory{Float64}}}

F = lu(T)
typeof(F.factors.dl)  # FixedSizeArray{Float64, 1, Memory{Float64}}  ✓
typeof(F.ipiv)        # Vector{Int64}                                ✗

MWE 2 — dense, non-BlasFloat element

A = FixedSizeMatrix{Rational{Int}}(undef, 3, 3)
copyto!(A, Rational{Int}[2 1 0; 1 2 1; 0 1 2])
F = lu(A)
typeof(F.factors)   # FixedSizeArray{Rational{Int64}, 2, Memory{Rational{Int64}}}  ✓
typeof(F.ipiv)      # Vector{Int64}                                                ✗

For dense FixedSizeMatrix{Float64} the BLAS path runs and returns ipiv::FixedSizeVector{Int64} as expected, confirming the inconsistency is purely in the two generic code paths above.

Proposed fix

--- a/src/lu.jl
+++ b/src/lu.jl
@@ generic_lufact! signature @@
-                         ipiv::AbstractVector{BlasInt} = Vector{BlasInt}(undef, min(size(A)...));
+                         ipiv::AbstractVector{BlasInt} = similar(A, BlasInt, min(size(A)...));
@@ lu!(::Tridiagonal) @@
-    _lu_tridiag!(A.dl, A.d, A.du, du2, Vector{BlasInt}(undef, n), pivot, check, allowsingular)
+    _lu_tridiag!(A.dl, A.d, A.du, du2, similar(A.d, BlasInt, n), pivot, check, allowsingular)

Choosing similar(A.d, BlasInt, n) (rather than similar(A, …)) for the Tridiagonal case keeps the pivot in the same container family as the diagonal vectors, which is the natural choice given that Tridiagonal itself is a struct, not an array of A's container type.

Non-breaking: for Matrix and Tridiagonal{T, Vector{T}} (the defaults), similar returns Vector{BlasInt} as before. Custom AbstractMatrix types (or Tridiagonal with a custom diagonal vector type) automatically get a consistent pivot.

Happy to send a PR with the fix plus a regression test along the following lines (zero external dependencies):

@testset "ipiv storage respects array type" begin
    struct TaggedVec{T} <: AbstractVector{T}; v::Vector{T}; end
    Base.size(v::TaggedVec) = size(v.v)
    Base.getindex(v::TaggedVec, i::Int) = v.v[i]
    Base.setindex!(v::TaggedVec, x, i::Int) = (v.v[i] = x; v)
    Base.similar(::TaggedVec, ::Type{T}, dims::Dims{1}) where {T} =
        TaggedVec(Vector{T}(undef, dims))

    struct TaggedMat{T} <: AbstractMatrix{T}; m::Matrix{T}; end
    Base.size(A::TaggedMat) = size(A.m)
    Base.getindex(A::TaggedMat, i::Int, j::Int) = A.m[i,j]
    Base.setindex!(A::TaggedMat, v, i::Int, j::Int) = (A.m[i,j] = v; A)
    Base.similar(::TaggedMat, ::Type{T}, dims::Dims{1}) where {T} =
        TaggedVec(Vector{T}(undef, dims))
    Base.similar(::TaggedMat, ::Type{T}, dims::Dims{2}) where {T} =
        TaggedMat(Matrix{T}(undef, dims))

    A = TaggedMat(Rational{Int}[2 1 0; 1 2 1; 0 1 2])
    @test lu(A).ipiv isa TaggedVec{LinearAlgebra.BlasInt}

    T = Tridiagonal(TaggedVec([1.0, 1.0]),
                    TaggedVec([4.0, 4.0, 4.0]),
                    TaggedVec([1.0, 1.0]))
    @test lu(T).ipiv isa TaggedVec{LinearAlgebra.BlasInt}
end

---

_{Issue drafted with the help of Claude Code.}