Merge pull request #273 from SciML/validation

add unit validation of equations

Author: ChrisRackauckas

Project.toml

@@ -17,6 +17,7 @@ SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 TreeViews = "a2a6695c-b41b-5b7d-aed9-dbfdeacea5d7"
+Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [compat]
 DiffEqBase = "6"

src/ModelingToolkit.jl

@@ -2,7 +2,7 @@ module ModelingToolkit
 
 using DiffEqBase, Distributed
 using StaticArrays, LinearAlgebra, SparseArrays
-using Latexify
+using Latexify, Unitful
 using MacroTools
 
 using MacroTools
@@ -92,6 +92,7 @@ include("systems/diffeqs/sdesystem.jl")
 include("systems/diffeqs/abstractodesystem.jl")
 include("systems/diffeqs/first_order_transform.jl")
 include("systems/diffeqs/modelingtoolkitize.jl")
+include("systems/diffeqs/validation.jl")
 include("systems/nonlinear/nonlinear_system.jl")
 include("systems/pde/pdesystem.jl")
 include("systems/reaction/reactionsystem.jl")

src/systems/diffeqs/validation.jl

@@ -1,17 +1,26 @@
-isintermediate(eq::Equation) = !(isa(eq.lhs, Operation) && isa(eq.lhs.op, Differential))
+Base.:*(x::Expression,y::Unitful.AbstractQuantity) = x * ModelingToolkit.Constant(y)
 
-struct ODEExpr  # dⁿx/dtⁿ = rhs
-    x::Variable
-    n::Int
-    rhs::Expression
+instantiate(x::ModelingToolkit.Constant) = x.value
+instantiate(x::ModelingToolkit.Variable{Number}) = 1.0
+instantiate(x::ModelingToolkit.Variable) = oneunit(1*ModelingToolkit.vartype(x))
+function instantiate(x::ModelingToolkit.Operation)
+    if x.op isa Variable
+        return instantiate(x.op)
+    elseif x.op isa Differential
+        instantiate(x.args[1])/instantiate(x.args[1].args[1])
+    else
+        x.op(instantiate.(x.args)...)
+    end
 end
-function Base.convert(::Type{ODEExpr},eq::Equation)
-    isintermediate(eq) && throw(ArgumentError("intermediate equation received"))
-    (x, t, n) = flatten_differential(eq.lhs)
-    (isa(t, Operation) && isa(t.op, Variable) && isempty(t.args)) ||
-        throw(ArgumentError("invalid independent variable $t"))
-    (isa(x, Operation) && isa(x.op, Variable) && length(x.args) == 1 && isequal(first(x.args), t)) ||
-        throw(ArgumentError("invalid dependent variable $x"))
-    return t.op, ODEExpr(x.op, n, eq.rhs)
+
+function validate(eq::ModelingToolkit.Equation)
+    try
+        return typeof(instantiate(eq.lhs)) == typeof(instantiate(eq.rhs))
+    catch
+        return false
+    end
+end
+
+function validate(sys::AbstractODESystem)
+    all(validate.(equations(sys)))
 end
-Base.:(==)(a::ODEExpr, b::ODEExpr) = isequal((a.x, a.n, a.rhs), (b.x, b.n, b.rhs))

test/units.jl

@@ -0,0 +1,25 @@
+using ModelingToolkit, Unitful
+using Test
+
+t = Variable{u"s"}(:t)()
+x = Variable{u"kg"}(:x)(t)
+y = Variable{u"kg"}(:y)(t)
+D = Differential(t)
+
+eq1 = x ~ y*t
+eq2 = x*10u"s" ~ y*t
+
+@test ModelingToolkit.instantiate(t) == 1u"s"
+@test ModelingToolkit.instantiate(x) == 1u"kg"
+@test ModelingToolkit.instantiate(y) == 1u"kg"
+
+@test !ModelingToolkit.validate(eq1)
+@test ModelingToolkit.validate(eq2)
+
+eqs = [
+        D(x) ~ y/t
+        D(y) ~ (x*y)/(t*10u"kg")
+]
+
+sys = ODESystem(eqs,t,[x,y],[])
+@test ModelingToolkit.validate(sys)