fix newly introduced type piracies

Author: j-fu

ext/ExtendableSparseLinearSolveExt.jl

@@ -4,20 +4,23 @@ import ExtendableSparse: LinearSolvePreconBuilder
 import LinearAlgebra
 using SparseArrays: AbstractSparseMatrixCSC
 
-# Harrr!
-# Avoid type piracy by adding a wrapper struct
-function (method::LinearSolve.AbstractFactorization)(A)
+
+struct LinearSolvePrecon{T}
+    cache::T
+end
+
+function LinearSolvePrecon(A,method::LinearSolve.AbstractFactorization)
     pr = LinearProblem(A, zeros(eltype(A), size(A, 1)))
-    init(pr, method)
+    LinearSolvePrecon(init(pr, method))
 end
 
-function LinearAlgebra.ldiv!(u, cache::LinearSolve.LinearCache, b)
-    cache.b = b
-    sol = solve!(cache)
+function LinearAlgebra.ldiv!(u, P::LinearSolvePrecon, b)
+    P.cache.b = b
+    sol = solve!(P.cache)
     copyto!(u, sol.u)
 end
 
-(b::LinearSolvePreconBuilder)(A::AbstractSparseMatrixCSC,p) = (b.method(A), LinearAlgebra.I)
+(b::LinearSolvePreconBuilder)(A::AbstractSparseMatrixCSC,p) = (LinearSolvePrecon(A,b.method), LinearAlgebra.I)
 
 end
 

src/preconbuilders.jl

@@ -30,26 +30,29 @@ Base.@kwdef struct ILUZeroPreconBuilder
     blocksize::Int = 1
 end
 
-function (b::ILUZeroPreconBuilder)(A0,p)
-    A=SparseMatrixCSC(size(A0)..., getcolptr(A0), rowvals(A0),nonzeros(A0))
-    if b.blocksize==1
-        (ILUZero.ilu0(A),LinearAlgebra.I)
-    else
-        (ILUZero.ilu0(pointblock(A,b.blocksize),SVector{b.blocksize,eltype(A)}),LinearAlgebra.I)
-    end
+struct ILUBlockPrecon{N,NN,Tv,Ti}
+    ilu0::ILUZero.ILU0Precon{SMatrix{N, N, Tv, NN}, Ti, SVector{N, Tv}}
 end
 
-# Harrr!!! ☠
-# We could resolve this piracy by introducing a wrapper type for the block case.
 function LinearAlgebra.ldiv!(Y::Vector{Tv},
-                             A::ILUZero.ILU0Precon{SMatrix{N, N, Tv, NN}, Ti, SVector{N, Tv}},
+                             A::ILUBlockPrecon{N,NN,Tv,Ti},
                              B::Vector{Tv}) where {N,NN,Tv,Ti}
     BY=reinterpret(SVector{N,Tv},Y)
     BB=reinterpret(SVector{N,Tv},B)
-    ldiv!(BY,A,BB)
+    ldiv!(BY,A.ilu0,BB)
     Y
 end
 
+function (b::ILUZeroPreconBuilder)(A0,p)
+    A=SparseMatrixCSC(size(A0)..., getcolptr(A0), rowvals(A0),nonzeros(A0))
+    if b.blocksize==1
+        (ILUZero.ilu0(A),LinearAlgebra.I)
+    else
+        (ILUBlockPrecon(ILUZero.ilu0(pointblock(A,b.blocksize),SVector{b.blocksize,eltype(A)})),LinearAlgebra.I)
+    end
+end
+
+
 """
     ILUTPreconBuilder(; droptol=0.1)