flatten and reduction on flatten

Author: shashi

src/systems/abstractsystem.jl

@@ -194,7 +194,7 @@ pins(sys::AbstractSystem) = isempty(sys.systems) ? sys.pins : [sys.pins;reduce(v
 function observed(sys::AbstractSystem)
     [sys.observed;
      reduce(vcat,
-            (namespace_equation.(s.observed, s.name, s.iv.name) for s in sys.systems),
+            (namespace_equation.(observed(s), s.name, s.iv.name) for s in sys.systems),
             init=Equation[])]
 end
 

src/systems/reduction.jl

@@ -1,11 +1,29 @@
 export alias_elimination
 
+function flatten(sys::ODESystem)
+    ODESystem(equations(sys),
+              independent_variable(sys),
+              states(sys),
+              parameters(sys),
+              observed=observed(sys))
+end
+
+
 function substitute_aliases(diffeqs, outputs)
     lhss(diffeqs) .~ substitute.(rhss(diffeqs), (Dict(lhss(outputs) .=> rhss(outputs)),))
 end
 
+function make_lhs_0(eq)
+    if eq.lhs isa Constant && iszero(eq.lhs)
+        return eq
+    else
+        0 ~ eq.lhs - eq.rhs
+    end
+end
+
 function alias_elimination(sys::ODESystem)
-    eqs = equations(sys)
+    eqs = vcat(equations(sys),
+               make_lhs_0.(observed(sys)))
 
     new_stateops = map(eqs) do eq
         if eq.lhs isa Operation && eq.lhs.op isa Differential

test/reduction.jl

@@ -1,7 +1,7 @@
 using ModelingToolkit, OrdinaryDiffEq, Test
 
 @parameters t σ ρ β F(t)
-@variables x(t) y(t) z(t) a(t)
+@variables x(t) y(t) z(t) a(t) u(t)
 @derivatives D'~t
 
 eqs = [D(x) ~ σ*(y-x),
@@ -21,57 +21,64 @@ eqs = [D(x) ~ σ*(y-x),
 
 @test lorenz1_aliased == ODESystem(eqs,t,[x,y,z],[σ,ρ,β],observed=[a ~ x],name=:lorenz1)
 
-#=
 # Multi-System Reduction
 
 eqs1 = [D(x) ~ σ*(y-x) + F,
        D(y) ~ x*(ρ-z)-u,
        D(z) ~ x*y - β*z]
+
 aliases = [u ~ x + y - z]
+
 lorenz1 = ODESystem(eqs1,pins=[F],observed=aliases,name=:lorenz1)
 
 eqs2 = [D(x) ~ F,
        D(y) ~ x*(ρ-z)-x,
        D(z) ~ x*y - β*z]
+
 aliases2 = [u ~ x - y - z]
+
 lorenz2 = ODESystem(eqs2,pins=[F],observed=aliases2,name=:lorenz2)
 
 connections = [lorenz1.F ~ lorenz2.u,
                lorenz2.F ~ lorenz1.u]
+
 connected = ODESystem([0 ~ a + lorenz1.x - lorenz2.y],t,[a],[],observed=connections,systems=[lorenz1,lorenz2])
 
 # Reduced Unflattened System
+#=
 
 connections2 = [lorenz1.F ~ lorenz2.u,
                 lorenz2.F ~ lorenz1.u,
                 a ~ -lorenz1.x + lorenz2.y]
 connected = ODESystem(Equation[],t,[],[],observed=connections2,systems=[lorenz1,lorenz2])
+=#
 
 # Reduced Flattened System
 
-flattened_system = flatten(connected)
-flatten(sys::AbstractSystem) = ODESystem(equations(sys),states(sys),parameters(sys),independent_variable(sys))
+flattened_system = ModelingToolkit.flatten(connected)
 
 aliased_flattened_system = alias_elimination(flattened_system)
 
-states(reduced_flattened_system) == [
+@test states(aliased_flattened_system) == convert.(Variable, [
         lorenz1.x
         lorenz1.y
         lorenz1.z
         lorenz2.x
         lorenz2.y
         lorenz2.z
-]
+       ])
 
-parameters(reduced_flattened_system) == [
+@test setdiff(parameters(aliased_flattened_system), convert.(Variable, [
         lorenz1.σ
         lorenz1.ρ
         lorenz1.β
-        lorenz2.σ
+        lorenz1.F
+        lorenz2.F
         lorenz2.ρ
         lorenz2.β
-]
+       ])) |> isempty
 
+#=
 equations(reduced_flattened_system) == [
        D(lorenz1.x) ~ lorenz1.σ*(lorenz1.y-lorenz1.x) + lorenz2.x - lorenz2.y - lorenz2.z,
        D(lorenz1.y) ~ lorenz1.x*(ρ-z)-lorenz1.x - lorenz1.y + lorenz1.z,