BlockArray errors with generic linear algebra functions

[imciner2]

When using the generic fallback functions provided by LinearAlgebra, BlockArray type matrices don't work properly. I'm not sure what the best way of handling this would be, but I guess we could theoretically generate versions of the linear algebra functions that fail and have them take the block array and do a conversion to a normal array when calling the standard linear algebra function.

For example, the cond function fails due to a missing bidiagonalize_tall! function. Of course, explicitly converting the block array to a normal array allows the computation to succeed (since it no longer needs the specialization of bidiagonalize_tall!).

julia> D = BlockArray( randn(4,4), [2, 2], [2, 2] )                                                                                                                                                                                                       
2×2-blocked 4×4 BlockArray{Float64,2}:
  0.923874   -0.0829803  │  0.286211   0.46231 
  0.195939    0.14953    │  1.8477    -0.321445
 ────────────────────────┼─────────────────────
 -0.0989605  -1.26688    │  1.42871   -0.659864
  0.403813   -0.86728    │  0.72237   -0.621602

julia> cond( D )                                                                                                                                                                                                                                          
ERROR: MethodError: no method matching bidiagonalize_tall!(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}})
Closest candidates are:
  bidiagonalize_tall!(::Array{T,2}) where T at /home/imcinerney/.julia/packages/GenericSVD/cT5Cu/src/bidiagonalize.jl:51
  bidiagonalize_tall!(::Array{T,2}, ::Bidiagonal) where T at /home/imcinerney/.julia/packages/GenericSVD/cT5Cu/src/bidiagonalize.jl:19
  bidiagonalize_tall!(::Adjoint{T2,Array{T,2}}) where {T, T2} at /home/imcinerney/.julia/packages/GenericSVD/cT5Cu/src/bidiagonalize.jl:47
Stacktrace:
 [1] generic_svdvals!(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}}; sorted::Bool) at /home/imcinerney/.julia/packages/GenericSVD/cT5Cu/src/GenericSVD.jl:57
 [2] generic_svdvals!(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}}) at /home/imcinerney/.julia/packages/GenericSVD/cT5Cu/src/GenericSVD.jl:53
 [3] svdvals!(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}}) at /home/imcinerney/.julia/packages/GenericSVD/cT5Cu/src/GenericSVD.jl:21
 [4] svdvals(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}}) at /home/imcinerney/dev/julia/main/1.4/usr/share/julia/stdlib/v1.4/LinearAlgebra/src/svd.jl:194
 [5] cond(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}}, ::Int64) at /home/imcinerney/dev/julia/main/1.4/usr/share/julia/stdlib/v1.4/LinearAlgebra/src/dense.jl:1419
 [6] cond(::BlockArray{Float64,2,Array{Array{Float64,2},2},Tuple{BlockedUnitRange{Array{Int64,1}},BlockedUnitRange{Array{Int64,1}}}}) at /home/imcinerney/dev/julia/main/1.4/usr/share/julia/stdlib/v1.4/LinearAlgebra/src/dense.jl:1418
 [7] top-level scope at REPL[15]:1
 [8] eval(::Module, ::Any) at ./boot.jl:331
 [9] eval_user_input(::Any, ::REPL.REPLBackend) at /home/imcinerney/dev/julia/main/1.4/usr/share/julia/stdlib/v1.4/REPL/src/REPL.jl:86
 [10] run_backend(::REPL.REPLBackend) at /home/imcinerney/.julia/packages/Revise/XFtoQ/src/Revise.jl:1162
 [11] top-level scope at REPL[1]:0

julia> cond( Array( D ) )                                                                                                                                                                                                                                 
8.508044778284596

So far the functions that I have found fail are (but there are probably more):

• cond
• svd
• svdvals
• eigen

[dlfivefifty]

Note there is no guarantee such functions work for all arrays, including those in Base:

julia> cond(Tridiagonal(randn(5),randn(6),randn(5)))
ERROR: MethodError: no method matching svdvals!(::Tridiagonal{Float64,Array{Float64,1}})
Closest candidates are:
  svdvals!(::SymTridiagonal) at /Users/solver/Projects/julia-1.5/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/tridiag.jl:350
  svdvals!(::StridedArray{T, 2}) where T<:Union{Complex{Float32}, Complex{Float64}, Float32, Float64} at /Users/solver/Projects/julia-1.5/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/svd.jl:193
  svdvals!(::StridedArray{T, 2}, ::StridedArray{T, 2}) where T<:Union{Complex{Float32}, Complex{Float64}, Float32, Float64} at /Users/solver/Projects/julia-1.5/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/svd.jl:501
  ...
Stacktrace:
 [1] svdvals(::Tridiagonal{Float64,Array{Float64,1}}) at /Users/solver/Projects/julia-1.5/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/svd.jl:194
 [2] cond(::Tridiagonal{Float64,Array{Float64,1}}, ::Int64) at /Users/solver/Projects/julia-1.5/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/dense.jl:1430
 [3] cond(::Tridiagonal{Float64,Array{Float64,1}}) at /Users/solver/Projects/julia-1.5/usr/share/julia/stdlib/v1.5/LinearAlgebra/src/dense.jl:1429
 [4] top-level scope at REPL[3]:1

So while your suggestion is a good one, it isn't required.

That said, I'll be happy to merge a PR, though svd will require some thought on what the natural block structure is.