fix v0.6 perf

Author: mohdibntarek

src/aggregation.jl

@@ -1,4 +1,5 @@
-function smoothed_aggregation(A::SparseMatrixCSC{T,V}, ::Type{Val{bs}}=Val{1},
+function smoothed_aggregation(A::TA, 
+                        ::Type{Val{bs}}=Val{1},
                         symmetry = HermitianSymmetry(),
                         strength = SymmetricStrength(),
                         aggregate = StandardAggregation(),
@@ -10,7 +11,7 @@ function smoothed_aggregation(A::SparseMatrixCSC{T,V}, ::Type{Val{bs}}=Val{1},
                         max_coarse = 10,
                         diagonal_dominance = false,
                         keep = false,
-                        coarse_solver = Pinv) where {T,V,bs}
+                        coarse_solver = Pinv) where {T,V,bs,TA<:SparseMatrixCSC{T,V}}
 
     n = size(A, 1)
     # B = kron(ones(n, 1), eye(1))
@@ -27,7 +28,11 @@ function smoothed_aggregation(A::SparseMatrixCSC{T,V}, ::Type{Val{bs}}=Val{1},
     # agg = [aggregate for _ in 1:max_levels - 1]
     # sm = [smooth for _ in 1:max_levels]
 
-    levels = Vector{Level{T,V}}()
+    @static if VERSION < v"0.7-"
+        levels = Vector{Level{TA, TA, TA}}()
+    else
+        levels = Vector{Level{TA, TA, Adjoint{T, TA}}}()
+    end
     bsr_flag = false
     w = MultiLevelWorkspace(Val{bs}, eltype(A))
 
@@ -57,7 +62,11 @@ function extend_hierarchy!(levels, strength, aggregate, smooth,
                             symmetry, bsr_flag)
 
     # Calculate strength of connection matrix
-    S = strength(A, bsr_flag)
+    if symmetry isa HermitianSymmetry
+        S, _T = strength(A, bsr_flag)
+    else
+        S, _T = strength(adjoint(A), bsr_flag)
+    end
 
     # Aggregation operator
     AggOp = aggregate(S)
@@ -80,11 +89,11 @@ function extend_hierarchy!(levels, strength, aggregate, smooth,
 
     A, B, bsr_flag
 end
-construct_R(::HermitianSymmetry, P) = copy(P')
+construct_R(::HermitianSymmetry, P) = P'
 
 function fit_candidates(AggOp, B, tol = 1e-10)
 
-    A = copy(AggOp')
+    A = adjoint(AggOp)
     n_fine, n_coarse = size(A)
     n_col = n_coarse
 

src/classical.jl

@@ -7,33 +7,56 @@ struct Solver{S,T,P,PS}
     max_coarse::Int64
 end
 
-function ruge_stuben(A::SparseMatrixCSC{Ti,Tv}, ::Type{Val{bs}}=Val{1};
+function ruge_stuben(_A::Union{TA, Symmetric{Ti, TA}, Hermitian{Ti, TA}}, 
+                ::Type{Val{bs}}=Val{1};
                 strength = Classical(0.25),
                 CF = RS(),
                 presmoother = GaussSeidel(),
                 postsmoother = GaussSeidel(),
                 max_levels = 10,
                 max_coarse = 10,
-                coarse_solver = Pinv) where {Ti,Tv,bs}
+                coarse_solver = Pinv) where {Ti,Tv,bs,TA<:SparseMatrixCSC{Ti,Tv}}
 
     s = Solver(strength, CF, presmoother,
                 postsmoother, max_levels, max_levels)
 
-    levels = Vector{Level{Ti,Tv}}()
+    if _A isa Symmetric && Ti <: Real || _A isa Hermitian
+        A = _A.data
+        At = A
+        symmetric = true
+        @static if VERSION < v"0.7-"
+            levels = Vector{Level{TA, TA}}()
+        else
+            levels = Vector{Level{TA, Adjoint{Ti, TA}, TA}}()
+        end
+    else
+        symmetric = false
+        A = _A
+        At = adjoint(A)
+        @static if VERSION < v"0.7-"
+            levels = Vector{Level{TA, TA, TA}}()
+        else
+            levels = Vector{Level{TA, Adjoint{Ti, TA}, TA}}()
+        end
+    end
     w = MultiLevelWorkspace(Val{bs}, eltype(A))
 
     while length(levels) + 1 < max_levels && size(A, 1) > max_coarse
         residual!(w, size(A, 1))
-        A = extend_heirarchy!(levels, strength, CF, A)
+        A = extend_heirarchy!(levels, strength, CF, A, symmetric)
         coarse_x!(w, size(A, 1))
         coarse_b!(w, size(A, 1))
     end
 
     MultiLevel(levels, A, coarse_solver(A), presmoother, postsmoother, w)
 end
 
-function extend_heirarchy!(levels::Vector{Level{Ti,Tv}}, strength, CF, A::SparseMatrixCSC{Ti,Tv}) where {Ti,Tv}
-    At = copy(A')
+function extend_heirarchy!(levels, strength, CF, A::SparseMatrixCSC{Ti,Tv}, symmetric) where {Ti,Tv}
+    if symmetric
+        At = A
+    else
+        At = adjoint(A)
+    end
     S, T = strength(At)
     splitting = CF(S)
     P, R = direct_interpolation(At, T, splitting)
@@ -48,7 +71,7 @@ function direct_interpolation(At, T, splitting)
     Pp = rs_direct_interpolation_pass1(T, splitting)
     Px, Pj, Pp = rs_direct_interpolation_pass2(At, T, splitting, Pp)
     R = SparseMatrixCSC(maximum(Pj), size(At, 1), Pp, Pj, Px)
-    P = copy(R')
+    P = R'
 
     P, R
 end

src/multilevel.jl

@@ -1,12 +1,12 @@
-struct Level{T,V}
-    A::SparseMatrixCSC{T,V}
-    P::SparseMatrixCSC{T,V}
-    R::SparseMatrixCSC{T,V}
+struct Level{TA, TP, TR}
+    A::TA
+    P::TP
+    R::TR
 end
 
-struct MultiLevel{S, Pre, Post, Ti, Tv, TW}
-    levels::Vector{Level{Ti,Tv}}
-    final_A::SparseMatrixCSC{Ti,Tv}
+struct MultiLevel{S, Pre, Post, TA, TP, TR, TW}
+    levels::Vector{Level{TA, TP, TR}}
+    final_A::TA
     coarse_solver::S
     presmoother::Pre
     postsmoother::Post

src/splitting.jl

@@ -19,7 +19,7 @@ end
 
 function (::RS)(S)
 	remove_diag!(S)
-	RS_CF_splitting(S, copy(S'))
+	RS_CF_splitting(S, adjoint(S))
 end
 
 function RS_CF_splitting(S::SparseMatrixCSC, T::SparseMatrixCSC)

src/strength.jl

@@ -33,7 +33,7 @@ function (c::Classical)(At::SparseMatrixCSC{Tv,Ti}) where {Ti,Tv}
 
     scale_cols_by_largest_entry!(T)
 
-    copy(T'), T
+    adjoint(T), T
 end
 
 function find_max_off_diag(A, i)
@@ -81,7 +81,7 @@ function (s::SymmetricStrength{T})(A, bsr_flag = false) where {T}
     if bsr_flag && θ == 0
         S = SparseMatrixCSC(size(A)...,
                     A.colptr, A.rowval, ones(eltype(A), size(A.rowval)))
-        return S
+        return S, S
     else
         S = deepcopy(A)
     end
@@ -118,5 +118,5 @@ function (s::SymmetricStrength{T})(A, bsr_flag = false) where {T}
     S.nzval .= abs.(S.nzval)
     scale_cols_by_largest_entry!(S)    
 
-    S
+    S, S
 end

src/utils.jl

@@ -1,3 +1,11 @@
+function adjoint(A)
+    @static if VERSION < v"0.7-"
+        A'
+    else
+        copy(A')
+    end
+end
+
 function approximate_spectral_radius(A, tol = 0.01,
                                         maxiter = 15, restart = 5)
 

test/gmg.jl

@@ -1,9 +1,9 @@
 import AlgebraicMultigrid: Level, MultiLevel, GaussSeidel
 
-function multigrid(A::SparseMatrixCSC{T,V}; max_levels = 10, max_coarse = 10,
-                    presmoother = GaussSeidel(), postsmoother = GaussSeidel()) where {T,V}
+function multigrid(A::TA; max_levels = 10, max_coarse = 10,
+                    presmoother = GaussSeidel(), postsmoother = GaussSeidel()) where {T,V,TA<:SparseMatrixCSC{T,V}}
 
-    levels = Vector{Level{T,V}}()
+    levels = Vector{Level{TA,TA,TA}}()
     w = AlgebraicMultigrid.MultiLevelWorkspace(Val{1}, eltype(A))
 
     while length(levels) + 1 < max_levels && size(A, 1) > max_coarse
@@ -43,7 +43,7 @@ function extend!(levels, A::SparseMatrixCSC{Ti,Tv}) where {Ti,Tv}
 
     P = sparse(I, J, V, size_F, size_C)
 
-    R = copy(P')
+    R = AlgebraicMultigrid.adjoint(P)
 
     push!(levels, Level(A, P, R))
 

test/sa_tests.jl

@@ -30,7 +30,7 @@ function test_symmetric_soc()
     for matrix in cases
         for θ in (0.0, 0.1, 0.5, 1., 10.)
             ref_matrix = symmetric_soc(matrix, θ)
-            calc_matrix = SymmetricStrength(θ)(matrix)
+            calc_matrix, _ = SymmetricStrength(θ)(matrix)
 
             @test sum(abs2, ref_matrix - calc_matrix) < 1e-6
         end