/ (rdiv) doesn't work with BunchKaufman types

[ojwoodford]

MWE:

using LinearAlgebra, Test, StaticArrays
@testset "bunchkaufman-rdiv" begin
    A = rand(3, 3)
    A = A' * A
    @test A / cholesky(A) ≈ I # No error
    @test A / bunchkaufman(A) ≈ I # Error
    A = SMatrix{3, 3}(A)
    @test A / cholesky(A) ≈ I # No error
    @test A / bunchkaufman(A) ≈ I # No error
end

Output:

MethodError: no method matching rdiv!(::Matrix{Float64}, ::BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}})

Versions

julia> versioninfo()
Julia Version 1.10.0-beta2
Commit a468aa198d0 (2023-08-17 06:27 UTC)
Build Info:
  Official https://julialang.org/ release
Platform Info:
  OS: macOS (arm64-apple-darwin22.4.0)
  CPU: 10 × Apple M1 Pro
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-15.0.7 (ORCJIT, apple-m1)
  Threads: 1 on 8 virtual cores
Environment:
  JULIA_NUM_THREADS = 0
  JULIA_EDITOR = code

[aravindh-krishnamoorthy]

I tried to reproduce this with the latest commit https://github.com/JuliaLang/julia/commit/eaef647957ca5e085eea299cfa9f699c6afe6d8f with updated packages, specifically, with StaticArrays v1.6.5. In my case, both tests with bunchkaufman fail, i.e., the test also fails for SMatrix.

Details

```julia julia> versioninfo() Julia Version 1.11.0-DEV.911 Commit eaef647957 (2023-11-15 02:42 UTC) Platform Info: OS: Linux (x86_64-linux-gnu) CPU: 8 × 11th Gen Intel(R) Core(TM) i7-1165G7 @ 2.80GHz WORD_SIZE: 64 LLVM: libLLVM-15.0.7 (ORCJIT, tigerlake) Threads: 1 on 8 virtual cores

julia> using LinearAlgebra, Test, StaticArrays

julia> @testset "bunchkaufman-rdiv" begin A = rand(3, 3) A = A' * A @test A / cholesky(A) ≈ I # No error @test A / bunchkaufman(A) ≈ I # Error A = SMatrix{3, 3}(A) @test A / cholesky(A) ≈ I # No error @test A / bunchkaufman(A) ≈ I # No error end bunchkaufman-rdiv: Error During Test at REPL[5]:5 Test threw exception Expression: A / bunchkaufman(A) ≈ I MethodError: no method matching rdiv!(::Matrix{Float64}, ::BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}})

Closest candidates are: rdiv!(::AbstractVecOrMat, ::Diagonal) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/diagonal.jl:444 rdiv!(::AbstractVecOrMat, ::LinearAlgebra.AbstractQ) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/abstractq.jl:233 rdiv!(::AbstractMatrix, ::UpperHessenberg; shift) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/hessenberg.jl:229 ...

Stacktrace: [1] /(B::Matrix{Float64}, F::BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}}) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/factorization.jl:187 [2] macro expansion @ ~/julia/julia/usr/share/julia/stdlib/v1.11/Test/src/Test.jl:676 [inlined] [3] macro expansion @ REPL[5]:5 [inlined] [4] macro expansion @ ~/julia/julia/usr/share/julia/stdlib/v1.11/Test/src/Test.jl:1598 [inlined] [5] top-level scope @ REPL[5]:2 bunchkaufman-rdiv: Error During Test at REPL[5]:8 Test threw exception Expression: A / bunchkaufman(A) ≈ I MethodError: no method matching rdiv!(::Matrix{Float64}, ::BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}})

Closest candidates are: rdiv!(::AbstractVecOrMat, ::Diagonal) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/diagonal.jl:444 rdiv!(::AbstractVecOrMat, ::LinearAlgebra.AbstractQ) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/abstractq.jl:233 rdiv!(::AbstractMatrix, ::UpperHessenberg; shift) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/hessenberg.jl:229 ...

Stacktrace: [1] /(B::SMatrix{3, 3, Float64, 9}, F::BunchKaufman{Float64, Matrix{Float64}, Vector{Int64}}) @ LinearAlgebra ~/julia/julia/usr/share/julia/stdlib/v1.11/LinearAlgebra/src/factorization.jl:187 [2] macro expansion @ ~/julia/julia/usr/share/julia/stdlib/v1.11/Test/src/Test.jl:676 [inlined] [3] macro expansion @ REPL[5]:8 [inlined] [4] macro expansion @ ~/julia/julia/usr/share/julia/stdlib/v1.11/Test/src/Test.jl:1598 [inlined] [5] top-level scope @ REPL[5]:2 Test Summary: | Pass Error Total Time bunchkaufman-rdiv | 2 2 4 3.1s ERROR: Some tests did not pass: 2 passed, 0 failed, 2 errored, 0 broken.

julia>

</details>

[dkarrasch]

ldiv! for BunchKaufman is implemented in LAPACK. One would need to check if LAPACK also provides right-solves. Alternatively, this may fix it:

 /(A::StridedMatrix, B::BunchKaufman) = copy((B' \ A')')

at the expense of quite a number of (seemingly unavoidable) allocations. First, A' is materialized when \ calls ldiv!, and finally the adjoint of the result needs to be materialized again to avoid an adjoint return type.

[aravindh-krishnamoorthy]

Unfortunately, so far as I can see, there does not seem to be a LAPACK routine for rdiv.

Quick note: The documentation for LAPACK has improved a LOT in the recent 3.12.0 release. Please click on "Explore LAPACK code" here: LAPACK online document browser and navigate to LAPACK->Routines->LAPACK->Linear solve, AX=B->LDL for the relevant routines.