Closed gdalle closed 2 months ago
I think we should only overload mul!, transpose and adjoint for now. We support the other ones later.
Okay, how about the 5-argument mul!
? From what I can tell, it uses pre-allocated storage that is specific to operator-vector multiplication?
We can just support the 3-args mul! for now with operator-matrix products. It's only what we need for Krylov.jl.
I added the minimal amount of changes and the minimal tests, can you approve the workflows?
@amontoison at your disposal to add more code, tests or docs depending on what you see fit in the review
Is it working with block_gmres
in Krylov.jl?
Yes!
julia> using Krylov, LinearOperators
julia> A, B = rand(2, 2), rand(2, 2);
julia> block_gmres(A, B)
([0.711998925881128 -0.31790146670357355; 1.0257506949284791 3.3662342660807827], SimpleStats
niter: 1
solved: true
inconsistent: false
residuals: []
Aresiduals: []
κ₂(A): []
timer: 189.24μs
status: solution good enough given atol and rtol
)
julia> block_gmres(LinearOperator(A), B)
([0.711998925881128 -0.31790146670357355; 1.0257506949284791 3.3662342660807827], SimpleStats
niter: 1
solved: true
inconsistent: false
residuals: []
Aresiduals: []
κ₂(A): []
timer: 521.58ms
status: solution good enough given atol and rtol
)
Goal: to fix #322 by adding
mul!(res::AbstractMatrix, op::LinearOperator, m::AbstractMatrix)
In practice, there are lots of things to overload and I don't know if you want to do all of them. I also don't know how you want to test this new functionality. But anyway, here's a first draft (without tests) to get the discussion started