More generally handle matrix operator

Author: ChrisRackauckas

Project.toml

@@ -1,14 +1,15 @@
 name = "LinearSolve"
 uuid = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae"
 authors = ["SciML"]
-version = "2.1.6"
+version = "2.1.7"
 
 [deps]
 ArrayInterface = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"
 DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 EnumX = "4e289a0a-7415-4d19-859d-a7e5c4648b56"
 FastLapackInterface = "29a986be-02c6-4525-aec4-84b980013641"
 GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527"
+InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
 KLU = "ef3ab10e-7fda-4108-b977-705223b18434"
 Krylov = "ba0b0d4f-ebba-5204-a429-3ac8c609bfb7"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

src/factorization.jl

@@ -103,16 +103,10 @@ function init_cacheval(alg::Union{LUFactorization, GenericLUFactorization}, A, b
     ArrayInterface.lu_instance(convert(AbstractMatrix, A))
 end
 
-function init_cacheval(alg::Union{LUFactorization, GenericLUFactorization}, A::MatrixOperator, b, u, Pl, Pr,
-    maxiters::Int, abstol, reltol, verbose::Bool,
-    assumptions::OperatorAssumptions)
-    ArrayInterface.lu_instance(convert(AbstractMatrix, A))
-end
-
 const PREALLOCATED_LU = ArrayInterface.lu_instance(rand(1, 1))
 
 function init_cacheval(alg::Union{LUFactorization, GenericLUFactorization},
-                       A::Union{Matrix{Float64},MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl, Pr,
+                       A::Matrix{Float64}, b, u, Pl, Pr,
                        maxiters::Int, abstol, reltol, verbose::Bool,
                        assumptions::OperatorAssumptions)
     PREALLOCATED_LU
@@ -181,7 +175,7 @@ end
 
 const PREALLOCATED_QR = ArrayInterface.qr_instance(rand(1, 1))
 
-function init_cacheval(alg::QRFactorization, A::Union{Matrix{Float64}, MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl, Pr,
+function init_cacheval(alg::QRFactorization, A::Matrix{Float64}, b, u, Pl, Pr,
                        maxiters::Int, abstol, reltol, verbose::Bool,
                        assumptions::OperatorAssumptions)
     PREALLOCATED_QR
@@ -260,7 +254,7 @@ end
 
 const PREALLOCATED_CHOLESKY = ArrayInterface.cholesky_instance(rand(1, 1), cholpivot)
 
-function init_cacheval(alg::CholeskyFactorization, A::Union{Matrix{Float64}, MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl, Pr,
+function init_cacheval(alg::CholeskyFactorization, A::Matrix{Float64}, b, u, Pl, Pr,
                        maxiters::Int, abstol, reltol, verbose::Bool,
                        assumptions::OperatorAssumptions)
     PREALLOCATED_CHOLESKY
@@ -352,7 +346,7 @@ end
 
 const PREALLOCATED_SVD = ArrayInterface.svd_instance(rand(1, 1))
 
-function init_cacheval(alg::SVDFactorization, A::Union{Matrix{Float64}, MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl, Pr,
+function init_cacheval(alg::SVDFactorization, A::Matrix{Float64}, b, u, Pl, Pr,
                        maxiters::Int, abstol, reltol, verbose::Bool,
                        assumptions::OperatorAssumptions)
     PREALLOCATED_SVD
@@ -677,7 +671,7 @@ function init_cacheval(alg::UMFPACKFactorization,
     nothing
 end
 
-function init_cacheval(alg::UMFPACKFactorization, A::Union{SparseMatrixCSC{Float64, Int}, MatrixOperator{Float64, Matrix{Float64}}}, b, u,
+function init_cacheval(alg::UMFPACKFactorization, A::SparseMatrixCSC{Float64, Int}, b, u,
                        Pl, Pr,
                        maxiters::Int, abstol, reltol,
                        verbose::Bool, assumptions::OperatorAssumptions)
@@ -774,7 +768,7 @@ function init_cacheval(alg::KLUFactorization,
     nothing
 end
 
-function init_cacheval(alg::KLUFactorization, A::Union{SparseMatrixCSC{Float64, Int}, MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl,
+function init_cacheval(alg::KLUFactorization, A::SparseMatrixCSC{Float64, Int}, b, u, Pl,
                        Pr,
                        maxiters::Int, abstol, reltol,
                        verbose::Bool, assumptions::OperatorAssumptions)
@@ -810,7 +804,7 @@ function SciMLBase.solve!(cache::LinearCache, alg::KLUFactorization; kwargs...)
                 # This won't recompute if it does.
                 KLU.klu_analyze!(cacheval)
                 copyto!(cacheval.nzval, nonzeros(A))
-                if cache.cacheval._numeric === C_NULL # We MUST have a numeric factorization for reuse, unlike UMFPACK.
+                if cacheval._numeric === C_NULL # We MUST have a numeric factorization for reuse, unlike UMFPACK.
                     KLU.klu_factor!(cacheval)
                 end
                 fact = KLU.klu!(cacheval,
@@ -861,7 +855,7 @@ function init_cacheval(alg::CHOLMODFactorization,
     nothing
 end
 
-function init_cacheval(alg::CHOLMODFactorization, A::Union{SparseMatrixCSC{Float64, Int}, MatrixOperator{Float64, Matrix{Float64}}}, b, u,
+function init_cacheval(alg::CHOLMODFactorization, A::SparseMatrixCSC{Float64, Int}, b, u,
                        Pl, Pr,
                        maxiters::Int, abstol, reltol,
                        verbose::Bool, assumptions::OperatorAssumptions)
@@ -917,13 +911,7 @@ function init_cacheval(alg::RFLUFactorization, A, b, u, Pl, Pr, maxiters::Int,
     ArrayInterface.lu_instance(convert(AbstractMatrix, A)), ipiv
 end
 
-function init_cacheval(alg::RFLUFactorization, A::MatrixOperator, b, u, Pl, Pr, maxiters::Int,
-                       abstol, reltol, verbose::Bool, assumptions::OperatorAssumptions)
-    ipiv = Vector{LinearAlgebra.BlasInt}(undef, min(size(A)...))
-    ArrayInterface.lu_instance(convert(AbstractMatrix, A)), ipiv
-end
-
-function init_cacheval(alg::RFLUFactorization, A::Union{Matrix{Float64}, MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl, Pr,
+function init_cacheval(alg::RFLUFactorization, A::Matrix{Float64}, b, u, Pl, Pr,
                        maxiters::Int,
                        abstol, reltol, verbose::Bool, assumptions::OperatorAssumptions)
     ipiv = Vector{LinearAlgebra.BlasInt}(undef, 0)
@@ -1278,7 +1266,7 @@ function init_cacheval(alg::SparspakFactorization,
     nothing
 end
 
-function init_cacheval(::SparspakFactorization, A::Union{SparseMatrixCSC{Float64, Int}, MatrixOperator{Float64, Matrix{Float64}}}, b, u, Pl,
+function init_cacheval(::SparspakFactorization, A::SparseMatrixCSC{Float64, Int}, b, u, Pl,
                        Pr, maxiters::Int, abstol,
                        reltol,
                        verbose::Bool, assumptions::OperatorAssumptions)
@@ -1316,3 +1304,13 @@ function SciMLBase.solve!(cache::LinearCache, alg::SparspakFactorization; kwargs
     y = ldiv!(cache.u, @get_cacheval(cache, :SparspakFactorization), cache.b)
     SciMLBase.build_linear_solution(alg, y, nothing, cache)
 end
+
+for alg in InteractiveUtils.subtypes(AbstractFactorization)
+    @eval function init_cacheval(alg::$alg, A::MatrixOperator, b, u, Pl, Pr,
+        maxiters::Int, abstol, reltol, verbose::Bool,
+        assumptions::OperatorAssumptions)
+        init_cacheval(alg, A.A, b, u, Pl, Pr,
+        maxiters::Int, abstol, reltol, verbose::Bool,
+        assumptions::OperatorAssumptions)
+    end
+end