ComponentArrays support

[gdalle]

Similar to #766

julia> using Krylov, ComponentArrays

julia> A = rand(2, 2)
2×2 Matrix{Float64}:
 0.106809  0.817828
 0.506456  0.728449

julia> b = ComponentVector(b1=rand(), b2=rand())
ComponentVector{Float64}(b1 = 0.11865752943169328, b2 = 0.5098226106263868)

julia> A \ b
2-element Vector{Float64}:
 0.9825279457773433
 0.016769695792115385

julia> gmres(A, b)
ERROR: MethodError: Cannot `convert` an object of type UndefInitializer to an object of type Vector{Float64}

Closest candidates are:
  convert(::Type{T}, ::LinearAlgebra.Factorization) where T<:AbstractArray
   @ LinearAlgebra ~/.julia/juliaup/julia-1.9.2+0.x64.linux.gnu/share/julia/stdlib/v1.9/LinearAlgebra/src/factorization.jl:59
  convert(::Type{<:Array}, ::ComponentArray)
   @ ComponentArrays ~/.julia/packages/ComponentArrays/ehUwN/src/similar_convert_copy.jl:58
  convert(::Type{T}, ::AbstractArray) where T<:Array
   @ Base array.jl:613
  ...

Stacktrace:
 [1] ComponentVector{Float64, Vector{Float64}, Tuple{Axis{(b1 = 1, b2 = 2)}}}(data::UndefInitializer, axes::Int64)
   @ ComponentArrays ~/.julia/packages/ComponentArrays/ehUwN/src/componentarray.jl:35
 [2] GmresSolver(m::Int64, n::Int64, memory::Int64, S::Type{ComponentVector{Float64, Vector{Float64}, Tuple{Axis{(b1 = 1, b2 = 2)}}}})
   @ Krylov ~/.julia/packages/Krylov/Iuavd/src/krylov_solvers.jl:1549
 [3] GmresSolver(A::Matrix{Float64}, b::ComponentVector{Float64, Vector{Float64}, Tuple{Axis{(b1 = 1, b2 = 2)}}}, memory::Int64)
   @ Krylov ~/.julia/packages/Krylov/Iuavd/src/krylov_solvers.jl:1567
 [4] gmres(A::Matrix{Float64}, b::ComponentVector{Float64, Vector{Float64}, Tuple{Axis{(b1 = 1, b2 = 2)}}}; memory::Int64, M::LinearAlgebra.UniformScaling{Bool}, N::LinearAlgebra.UniformScaling{Bool}, ldiv::Bool, restart::Bool, reorthogonalization::Bool, atol::Float64, rtol::Float64, itmax::Int64, timemax::Float64, verbose::Int64, history::Bool, callback::Krylov.var"#164#171", iostream::Core.CoreSTDOUT)
   @ Krylov ~/.julia/packages/Krylov/Iuavd/src/gmres.jl:121
 [5] gmres(A::Matrix{Float64}, b::ComponentVector{Float64, Vector{Float64}, Tuple{Axis{(b1 = 1, b2 = 2)}}})
   @ Krylov ~/.julia/packages/Krylov/Iuavd/src/gmres.jl:119
 [6] top-level scope
   @ REPL[8]:1

[amontoison]

The issue of component array is that we don't have a constructor for similar(b, m) where m is different than the length of b.

[gdalle]

I just tried the same fix as for StaticArrays and it seems to work

julia> Krylov.ktypeof(b::ComponentVector{T,V}) where {T,V} = V

julia> gmres(A, b)
([0.9825279457773434, 0.016769695792115402], SimpleStats
 niter: 2
 solved: true
 inconsistent: false
 residuals: []
 Aresiduals: []
 κ₂(A): []
 timer: 9.33μs
 status: solution good enough given atol and rtol
)

Gonna add it to #767

[gdalle]

Done! Hope this is convincing enough to merge it ;)

[gdalle]

Actually this wasn't closed by #768 because, to my great surprise, ComponentVector is actually a subtype of DenseVector. I didn't even think to check, sorry. That means the dispatch misses the check S.name.name == :ComponentArray.

A problem we wouldn't have had with the extensions 😉 But at least it's not blocking my release

[gdalle]

See also #701

[amontoison]

@gdalle With the release 0.9.9, I added a new constructor for all Krylov workspaces.

using Krylov
using ComponentArrays

# Specify that we also want some ComponentVector
# in the workspaces.
Krylov.ktypeof(v::ComponentArray) = typeof(v)

A = rand(2, 2)
b = ComponentVector(b1=rand(), b2=rand())

kc = KrylovConstructor(b)
solver = GmresSolver(kc)
gmres!(solver, A, b, verbose=1)
x = solution(solver)

I don't know what the most efficient vector type in the Krylov workspaces for your package but you have more flexibility now.