Merge pull request #982 from SciML/myb/promote

Add system promotion (toodesystem)

Author: YingboMa

src/ModelingToolkit.jl

@@ -136,7 +136,7 @@ include("systems/alias_elimination.jl")
 include("structural_transformation/StructuralTransformations.jl")
 @reexport using .StructuralTransformations
 
-export ODESystem, ODEFunction, ODEFunctionExpr, ODEProblemExpr
+export ODESystem, ODEFunction, ODEFunctionExpr, ODEProblemExpr, convert_system
 export DAEFunctionExpr, DAEProblemExpr
 export SDESystem, SDEFunction, SDEFunctionExpr, SDESystemExpr
 export SystemStructure

src/systems/diffeqs/odesystem.jl

@@ -300,3 +300,42 @@ function collect_differential_variables(sys::ODESystem)
     end
     return diffvars
 end
+
+# We have a stand-alone function to convert a `NonlinearSystem` or `ODESystem`
+# to an `ODESystem` to connect systems, and we later can reply on
+# `structural_simplify` to convert `ODESystem`s to `NonlinearSystem`s.
+"""
+$(TYPEDSIGNATURES)
+
+Convert a `NonlinearSystem` to an `ODESystem` or converts an `ODESystem` to a
+new `ODESystem` with a different independent variable.
+"""
+function convert_system(::Type{<:ODESystem}, sys, t; name=nameof(sys))
+    isempty(observed(sys)) || throw(ArgumentError("`convert_system` cannot handle reduced model (i.e. observed(sys) is non-empty)."))
+    t = value(t)
+    varmap = Dict()
+    sts = states(sys)
+    newsts = similar(sts, Any)
+    for (i, s) in enumerate(sts)
+        if istree(s)
+            args = arguments(s)
+            length(args) == 1 || throw(InvalidSystemException("Illegal state: $s. The state can have at most one argument like `x(t)`."))
+            arg = args[1]
+            if isequal(arg, t)
+                newsts[i] = s
+                continue
+            end
+            ns = operation(s)(t)
+            newsts[i] = ns
+            varmap[s] = ns
+        else
+            ns = indepvar2depvar(s, t)
+            newsts[i] = ns
+            varmap[s] = ns
+        end
+    end
+    sub = Base.Fix2(substitute, varmap)
+    neweqs = map(sub, equations(sys))
+    defs = Dict(sub(k) => sub(v) for (k, v) in defaults(sys))
+    return ODESystem(neweqs, t, newsts, parameters(sys); defaults=defs, name=name)
+end

src/utils.jl

@@ -78,3 +78,9 @@ function todict(d)
 end
 
 _merge(d1, d2) = merge(todict(d1), todict(d2))
+
+function indepvar2depvar(s::Sym, args...)
+    T = FnType{NTuple{length(args)}, symtype(s)}
+    ns = Sym{T}(nameof(s))(args...)
+    @set! ns.metadata = s.metadata
+end

test/nonlinearsystem.jl

@@ -95,3 +95,15 @@ connected = NonlinearSystem([s ~ a + lorenz1.x
                              lorenz1.F ~ lorenz2.u
                              lorenz2.F ~ lorenz1.u], [s, a], [], systems=[lorenz1,lorenz2])
 @test_nowarn alias_elimination(connected)
+
+# system promotion
+using OrdinaryDiffEq
+@variables t
+D = Differential(t)
+@named subsys = convert_system(ODESystem, lorenz1, t)
+@named sys = ODESystem([D(subsys.x) ~ subsys.x + subsys.x], t, systems=[subsys])
+sys = structural_simplify(sys)
+prob = ODEProblem(sys, [subsys.x => 1, subsys.z => 2.0], (0, 1.0), [subsys.σ=>1,subsys.ρ=>2,subsys.β=>3])
+sol = solve(prob, Rodas5())
+@test sol[subsys.x] + sol[subsys.y] - sol[subsys.z] ≈ sol[subsys.u]
+@test_throws ArgumentError convert_system(ODESystem, sys, t)