Add a component system

Author: ChrisRackauckas

Project.toml

@@ -16,6 +16,7 @@ NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
+TreeViews = "a2a6695c-b41b-5b7d-aed9-dbfdeacea5d7"
 
 [compat]
 DiffEqBase = "6"
@@ -31,8 +32,8 @@ julia = "1.2"
 
 [extras]
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
-Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
 test = ["OrdinaryDiffEq", "Test", "StochasticDiffEq"]

src/ModelingToolkit.jl

@@ -11,6 +11,8 @@ import GeneralizedGenerated
 using DocStringExtensions
 using Base: RefValue
 
+import TreeViews
+
 """
 $(TYPEDEF)
 
@@ -109,10 +111,9 @@ export simplify_constants
 export Operation, Expression
 export calculate_jacobian, generate_jacobian, generate_function
 export calculate_massmatrix, generate_diffusion_function
-export independent_variables, dependent_variables, parameters
+export independent_variable, states, parameters, equations
 export simplified_expr, eval_function
 export @register, @I
 export modelingtoolkitize
 export Variable, @variables, @parameters
-
 end # module

src/systems/diffeqs/abstractodesystem.jl

@@ -159,13 +159,54 @@ function DiffEqBase.ODEFunction{iip}(sys::AbstractODESystem, dvs = sys.dvs, ps =
                       syms = Symbol.(sys.dvs))
 end
 
+renamespace(namespace,name) = Symbol(string(namespace)*"′"*string(name))
+
 function DiffEqBase.ODEFunction(sys::AbstractODESystem, args...; kwargs...)
     ODEFunction{true}(sys, args...; kwargs...)
 end
 
-independent_variables(sys::AbstractODESystem) = Set{Variable}([sys.iv])
-dependent_variables(sys::AbstractODESystem) = Set{Variable}(sys.dvs)
-parameters(sys::AbstractODESystem) = Set{Variable}(sys.ps)
+function namespace_variables(sys::AbstractODESystem)
+    [rename(x,renamespace(sys.name,x.name)) for x in sys.dvs]
+end
+
+function namespace_parameters(sys::AbstractODESystem)
+    [rename(x,renamespace(sys.name,x.name)) for x in sys.ps]
+end
+
+namespace_equations(sys::AbstractODESystem) = namespace_equation.(sys.eqs,sys.name,sys.iv.name)
+
+function namespace_equation(eq::Equation,name,ivname)
+    _lhs = namespace_operation(eq.lhs,name,ivname)
+    _rhs = namespace_operation(eq.rhs,name,ivname)
+    _lhs ~ _rhs
+end
+
+function namespace_operation(O::Operation,name,ivname)
+    if O.op isa Variable && O.op.name != ivname
+        Operation(rename(O.op,renamespace(name,O.op.name)),namespace_operation.(O.args,name,ivname))
+    else
+        Operation(O.op,namespace_operation.(O.args,name,ivname))
+    end
+end
+namespace_operation(O::Constant,name,ivname) = O
+
+independent_variable(sys::AbstractODESystem) = sys.iv
+states(sys::AbstractODESystem) = [sys.dvs;reduce(vcat,namespace_variables.(sys.systems))]
+parameters(sys::AbstractODESystem) = [sys.ps;reduce(vcat,namespace_parameters.(sys.systems))]
+
+function equations(sys::AbstractODESystem)
+    [sys.eqs;reduce(vcat,namespace_equations.(sys.systems))]
+end
+
+function states(sys::AbstractODESystem,name::Symbol)
+    x = sys.dvs[findfirst(x->x.name==name,sys.dvs)]
+    Variable(Symbol(string(sys.name)*"′"*string(x.name)),known=x.known)(sys.iv())
+end
+
+function parameters(sys::AbstractODESystem,name::Symbol)
+    x = sys.ps[findfirst(x->x.name==name,sys.ps)]
+    Variable(Symbol(string(sys.name)*"′"*string(x.name)),known=x.known)(sys.iv())
+end
 
 function _eq_unordered(a, b)
     length(a) === length(b) || return false

src/systems/diffeqs/odesystem.jl

@@ -55,9 +55,14 @@ struct ODESystem <: AbstractODESystem
     Name: the name of the system
     """
     name::Symbol
+    """
+    systems: The internal systems
+    """
+    systems::Vector{ODESystem}
 end
 
 function ODESystem(deqs::AbstractVector{<:Equation}, iv, dvs, ps;
+                   systems = ODESystem[],
                    name=gensym(:ODESystem))
     iv′ = clean(iv)
     dvs′ = [clean(dv) for dv ∈ dvs]
@@ -66,7 +71,7 @@ function ODESystem(deqs::AbstractVector{<:Equation}, iv, dvs, ps;
     jac = RefValue(Matrix{Expression}(undef, 0, 0))
     Wfact   = RefValue(Matrix{Expression}(undef, 0, 0))
     Wfact_t = RefValue(Matrix{Expression}(undef, 0, 0))
-    ODESystem(deqs, iv′, dvs′, ps′, tgrad, jac, Wfact, Wfact_t, name)
+    ODESystem(deqs, iv′, dvs′, ps′, tgrad, jac, Wfact, Wfact_t, name, systems)
 end
 
 var_from_nested_derivative(x) = var_from_nested_derivative(x,0)

src/systems/diffeqs/sdesystem.jl

@@ -33,17 +33,23 @@ struct SDESystem <: AbstractODESystem
     Name: the name of the system
     """
     name::Symbol
+    """
+    Systems: the internal systems
+    """
+    systems::Vector{SDESystem}
 end
 
-function SDESystem(deqs::AbstractVector{<:Equation}, neqs, iv, dvs, ps; name = gensym(:SDESystem))
+function SDESystem(deqs::AbstractVector{<:Equation}, neqs, iv, dvs, ps;
+                   systems = SDESystem[],
+                   name = gensym(:SDESystem))
     dvs′ = [clean(dv) for dv ∈ dvs]
     ps′ = [clean(p) for p ∈ ps]
     iv′ = clean(iv)
     tgrad = RefValue(Vector{Expression}(undef, 0))
     jac = RefValue(Matrix{Expression}(undef, 0, 0))
     Wfact   = RefValue(Matrix{Expression}(undef, 0, 0))
     Wfact_t = RefValue(Matrix{Expression}(undef, 0, 0))
-    SDESystem(deqs, neqs, iv′, dvs′, ps′, tgrad, jac, Wfact, Wfact_t, name)
+    SDESystem(deqs, neqs, iv′, dvs′, ps′, tgrad, jac, Wfact, Wfact_t, name, systems)
 end
 
 function generate_diffusion_function(sys::SDESystem, dvs = sys.dvs, ps = sys.ps, expression = Val{true}; kwargs...)

src/variables.jl

@@ -168,3 +168,13 @@ Define one or more known variables.
 macro parameters(xs...)
     esc(_parse_vars(:parameters, true, xs))
 end
+
+function rename(x::Variable,name::Symbol)
+    Variable(name,known=x.known)
+end
+
+TreeViews.hastreeview(x::Variable) = true
+function TreeViews.treelabel(io::IO,x::Variable,
+                             mime::MIME"text/plain" = MIME"text/plain"())
+  show(io,mime,Text(x.name))
+end

test/components.jl

@@ -9,29 +9,28 @@ using Test
 
 eqs = [D(x) ~ σ*(y-x),
        D(y) ~ x*(ρ-z)-y,
-       0 ~ x + y + z]
+       0 ~ x + y + β*z]
 
-de1 = ODESystem(eqs,t,[x,y,z],[σ,ρ,β],name=:lorenz1)
-de2 = ODESystem(eqs,t,[x,y,z],[σ,ρ,β],name=:lorenz2)
+lorenz1 = ODESystem(eqs,t,[x,y,z],[σ,ρ,β],name=:lorenz1)
+lorenz2 = ODESystem(eqs,t,[x,y,z],[σ,ρ,β],name=:lorenz2)
 
 @parameters α
 @variables a(t)
-for
+connnectedeqs = [D(a) ~ a*states(lorenz1,:x)]
 
-convert(DAESystem,sys)
-convert(ODESystem,sys)
+connected = ODESystem(connnectedeqs,t,[a],[α],systems=[lorenz1,lorenz2],name=:connectedlorenz)
 
-connected = ODESystems(connnectedeqs,[a],[α],systems=[de1,de2],name=:connectedlorenz)
+@variables lorenz1′x(t) lorenz1′y(t) lorenz1′z(t) lorenz2′x(t) lorenz2′y(t) lorenz2′z(t)
+@parameters lorenz1′σ lorenz1′ρ lorenz1′β lorenz2′σ lorenz2′ρ lorenz2′β
 
-deleteat!(de2.eqs,3)
-de2
-push!(de2.eqs,convert(ModelingToolkit.ODEExpr,D(z) ~ z^2)[2])
+eqs_flat = [D(a) ~ a*lorenz1′x,
+            D(lorenz1′x) ~ lorenz1′σ*(lorenz1′y-lorenz1′x),
+            D(lorenz1′y) ~ lorenz1′x*(lorenz1′ρ-lorenz1′z)-lorenz1′y,
+            0 ~ lorenz1′x + lorenz1′y + lorenz1′β*lorenz1′z,
+            D(lorenz2′x) ~ lorenz2′σ*(lorenz2′y-lorenz2′x),
+            D(lorenz2′y) ~ lorenz2′x*(lorenz2′ρ-lorenz2′z)-lorenz2′y,
+            0 ~ lorenz2′x + lorenz2′y + lorenz2′β*lorenz2′z]
 
-check_consistency(de2)
-
-b = a*z
-           [states(de1,:x) ~ α*states(de2,:z),
-            D(a) ~ b*states(de1,:x)]
-
-states(states) == [lorenz1′x,lorenz1′y,lorenz1′z,lorenz2′x,lorenz2′y,lorenz2′z,a]
-parameters(connected) == [lorenz1′σ,lorenz1′ρ,lorenz1′β,lorenz2′σ,lorenz2′ρ,lorenz2′β,α]
+@test [x.name for x in states(connected)] == [:a,:lorenz1′x,:lorenz1′y,:lorenz1′z,:lorenz2′x,:lorenz2′y,:lorenz2′z]
+@test [x.name for x in parameters(connected)] == [:α,:lorenz1′σ,:lorenz1′ρ,:lorenz1′β,:lorenz2′σ,:lorenz2′ρ,:lorenz2′β]
+@test eqs_flat == equations(connected)