Merge pull request #404 from caryan/iszero-for-operations

provide iszero for Operation to help sparse matrix addition and multiplication

Author: YingboMa

src/operations.jl

@@ -54,6 +54,10 @@ Base.isequal(::Variable , ::Operation) = false
 Base.isequal(::Operation, ::Constant ) = false
 Base.isequal(::Constant , ::Operation) = false
 
+# provide iszero for Operations to help sparse addition and multiplication
+# e.g. we want to tell the sparse library that iszero(zero(Operation) + zero(Operation)) == true
+Base.iszero(x::Operation) = (_x = simplify(x); _x isa Constant && iszero(_x.value))
+
 Base.show(io::IO, O::Operation) = print(io, convert(Expr, O))
 
 # For inv
@@ -73,12 +77,4 @@ Base.convert(::Type{Expr},x::Operation) = Expr(x)
 Base.promote_rule(::Type{<:Constant}, ::Type{<:Operation}) = Operation
 Base.promote_rule(::Type{<:Operation}, ::Type{<:Constant}) = Operation
 
-# Fix Sparse MatMul
-Base.:*(A::SparseMatrixCSC{Operation,S}, x::StridedVector{Operation}) where {S} =
-    (T = Operation; mul!(similar(x, T, A.m), A, x, true, false))
-Base.:*(A::SparseMatrixCSC{Tx,S}, x::StridedVector{Operation}) where {Tx,S} =
-    (T = LinearAlgebra.promote_op(LinearAlgebra.matprod, Operation, Tx); mul!(similar(x, T, A.m), A, x, true, false))
-Base.:*(A::SparseMatrixCSC{Operation,S}, x::StridedVector{Tx}) where {Tx,S} =
-    (T = LinearAlgebra.promote_op(LinearAlgebra.matprod, Operation, Tx); mul!(similar(x, T, A.m), A, x, true, false))
-
 LinearAlgebra.lu(O::AbstractMatrix{<:Operation};kwargs...) = lu(O,Val(false);kwargs...)

src/variables.jl

@@ -89,7 +89,8 @@ Get the value of a [`ModelingToolkit.Constant`](@ref).
 """
 Base.get(c::Constant) = c.value
 
-Base.iszero(ex::Expression) = isa(ex, Constant) && iszero(ex.value)
+Base.iszero(c::Constant) = iszero(c.value)
+
 Base.isone(ex::Expression)  = isa(ex, Constant) && isone(ex.value)
 
 # Variables use isequal for equality since == is an Operation

test/operation_overloads.jl

@@ -1,7 +1,34 @@
-using ModelingToolkit
-using LinearAlgebra
-@variables a,b,c,d
-X = [a b;c d]
-det(X)
-lu(X)
-inv(X)
+using ModelingToolkit
+using LinearAlgebra
+using SparseArrays: sparse
+using Test
+
+@variables a,b,c,d
+
+# test some matrix operations don't throw errors
+X = [a b;c d]
+det(X)
+lu(X)
+inv(X)
+
+# test operations with sparse arrays and Operations
+# note `isequal` instead of `==` because `==` would give another Operation
+
+# test that we can create a sparse array of Operation
+Oarray = zeros(Operation, 2,2)
+Oarray[2,2] = a
+@test isequal(sparse(Oarray), sparse([2], [2], [a]))
+
+# test Operation * sparse
+@test isequal(a * sparse([2], [2], [1]), sparse([2], [2], [a * 1]))
+
+# test sparse{Operation} + sparse{Operation}
+A = sparse([2], [2], [a])
+B = sparse([2], [2], [b])
+@test isequal(A + B, sparse([2], [2], [a+b]))
+
+# test sparse{Operation} * sparse{Operation}
+C = sparse([1, 2], [2, 1], [c, c])
+D = sparse([1, 2], [2, 1], [d, d])
+
+@test isequal(C * D, sparse([1,2], [1,2], [c * d, c * d]))