strides method not implemented for NotIPIV

[autocorr]

Hi, first of all let me say thank you for developing this excellent package. It sped up the LU factorization in my project Jadex.jl considerably. I've started to update the aforementioned project to newer dependencies and have run into an issue on Julia 1.8 and the versions of RecursiveFactorization.jl of 0.2.15 and later.

The issue comes up when using lu! without pivoting, and strides function complains that the NotIPIV type has no method. It seems that part of the issue is that the LU factorization type stores the ipiv vector not as a Vector{Int} but as NotIPIV (which is an instance of AbstractArray so I guess it works). When calling ldiv! using the LU instance, it runs into several method errors

Version information below and running latest v2.18 of RecursiveFactorization:

julia> versioninfo()
Julia Version 1.8.5
Commit 17cfb8e65ea (2023-01-08 06:45 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin21.4.0)
  CPU: 8 × Intel(R) Core(TM) i5-1038NG7 CPU @ 2.00GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-13.0.1 (ORCJIT, icelake-client)
  Threads: 1 on 8 virtual cores

A minimum working example:

using LinearAlgebra
using RecursiveFactorization
A = rand(5,5)
a = rand(5)
b = rand(5)
Q = RecursiveFactorization.lu!(A, NoPivot())
ldiv!(a, Q, b)

gives the results for Q, note the type signature:

LU{Float64, Matrix{Float64}, RecursiveFactorization.NotIPIV}
L factor:
5×5 Matrix{Float64}:
 1.0       0.0          0.0       0.0       0.0
 3.30285   1.0          0.0       0.0       0.0
 1.12176   0.0460152    1.0       0.0       0.0
 1.31507   1.20737      3.37898   1.0       0.0
 4.80136  -0.56637    -10.4022   -0.925234  1.0
U factor:
5×5 Matrix{Float64}:
 0.194241  0.0821204   0.589417   0.433037   0.234671
 0.0       0.60333    -1.32537   -1.32246   -0.0335227
 0.0       0.0         0.330024   0.210092   0.42771
 0.0       0.0         0.0        0.44405   -1.19933
 0.0       0.0         0.0        0.0        2.57653

which results in an error when passing to ldiv!(a, Q, b):

ERROR: MethodError: no method matching strides(::RecursiveFactorization.NotIPIV)
Closest candidates are:
  strides(::LoopVectorization.LowDimArray) at ~/.julia/packages/LoopVectorization/IkdFM/src/broadcast.jl:47
  strides(::StaticArraysCore.MArray) at ~/.julia/packages/StaticArrays/VLqRb/src/abstractarray.jl:22
  strides(::VectorizationBase.OffsetPrecalc) at ~/.julia/packages/VectorizationBase/0dXyA/src/strided_pointers/cse_stridemultiples.jl:27
  ...
Stacktrace:
 [1] stride(A::RecursiveFactorization.NotIPIV, k::Int64)
   @ Base ./abstractarray.jl:547
 [2] stride1(x::RecursiveFactorization.NotIPIV)
   @ LinearAlgebra ~/.julia/juliaup/julia-1.8.5+0.x64.apple.darwin14/share/julia/stdlib/v1.8/LinearAlgebra/src/LinearAlgebra.jl:208
 [3] _chkstride1 (repeats 3 times)
   @ ~/.julia/juliaup/julia-1.8.5+0.x64.apple.darwin14/share/julia/stdlib/v1.8/LinearAlgebra/src/LinearAlgebra.jl:214 [inlined]
 [4] chkstride1
   @ ~/.julia/juliaup/julia-1.8.5+0.x64.apple.darwin14/share/julia/stdlib/v1.8/LinearAlgebra/src/LinearAlgebra.jl:212 [inlined]
 [5] getrs!(trans::Char, A::Matrix{Float64}, ipiv::RecursiveFactorization.NotIPIV, B::Vector{Float64})
   @ LinearAlgebra.LAPACK ~/.julia/juliaup/julia-1.8.5+0.x64.apple.darwin14/share/julia/stdlib/v1.8/LinearAlgebra/src/lapack.jl:1008
 [6] ldiv!(A::LU{Float64, Matrix{Float64}, RecursiveFactorization.NotIPIV}, B::Vector{Float64})
   @ LinearAlgebra ~/.julia/juliaup/julia-1.8.5+0.x64.apple.darwin14/share/julia/stdlib/v1.8/LinearAlgebra/src/lu.jl:404
 [7] ldiv!(Y::Vector{Float64}, A::LU{Float64, Matrix{Float64}, RecursiveFactorization.NotIPIV}, B::Vector{Float64})
   @ LinearAlgebra ~/.julia/juliaup/julia-1.8.5+0.x64.apple.darwin14/share/julia/stdlib/v1.8/LinearAlgebra/src/factorization.jl:126
 [8] top-level scope
   @ REPL[8]:1

I tried hacking the stride method to return (1,), but it seems that there are other functions in LinearAlgebra that complain. I'm not a hundred percent sure if the solution is appropriate, but I think that the fix could be as simple as converting the ipiv vector into a Vector{BlasInt} before constructing the LU instance at the end of lu! here. So changing:

LU(A, ipiv, info)
to
LU(A, convert(Vector{BlasInt}, ipiv), info)

or probably ideally something that didn't allocate.

I've found that Julia 1.8.5 and RecursiveFactorization v0.2.14 work, although there is a deprecation warning issued that LU{T, S, typeof(ipiv)}(fatcors, ipiv, info) is the preferred constructor.

[YingboMa]

Ah, this looks like a missing overload for ldiv! with no pivoting.

[autocorr]

Great! Thank you very much for taking a look into this. :) Do you have a guess when the next release will be? My package is currently pinned to v0.2.14.

[chriselrod]

Both 0.2.19 and 0.2.20 have been released since #83.