Forward pins and add pins to the states

Author: YingboMa

src/systems/abstractsystem.jl

@@ -215,7 +215,7 @@ end
 independent_variable(sys::AbstractSystem) = sys.iv
 function states(sys::AbstractSystem)
     unique(isempty(sys.systems) ?
-           setdiff(sys.states, value.(sys.pins)) :
+           sys.states :
            [sys.states;reduce(vcat,namespace_variables.(sys.systems))])
 end
 parameters(sys::AbstractSystem) = isempty(sys.systems) ? sys.ps : [sys.ps;reduce(vcat,namespace_parameters.(sys.systems))]

src/systems/diffeqs/odesystem.jl

@@ -71,11 +71,12 @@ function ODESystem(deqs::AbstractVector{<:Equation}, iv, dvs, ps;
     iv′ = value(iv)
     dvs′ = value.(dvs)
     ps′ = value.(ps)
+    pins′ = value.(pins)
     tgrad = RefValue(Vector{Num}(undef, 0))
     jac = RefValue{Any}(Matrix{Num}(undef, 0, 0))
     Wfact   = RefValue(Matrix{Num}(undef, 0, 0))
     Wfact_t = RefValue(Matrix{Num}(undef, 0, 0))
-    ODESystem(deqs, iv′, dvs′, ps′, pins, observed, tgrad, jac, Wfact, Wfact_t, name, systems)
+    ODESystem(deqs, iv′, dvs′, ps′, pins′, observed, tgrad, jac, Wfact, Wfact_t, name, systems)
 end
 
 var_from_nested_derivative(x, i=0) = (missing, missing)

src/systems/nonlinear/nonlinearsystem.jl

@@ -42,7 +42,7 @@ function NonlinearSystem(eqs, states, ps;
                          observed = [],
                          name = gensym(:NonlinearSystem),
                          systems = NonlinearSystem[])
-    NonlinearSystem(eqs, value.(states), value.(ps), value.(pins), observed, name, systems)
+    NonlinearSystem(value.(eqs), value.(states), value.(ps), value.(pins), observed, name, systems)
 end
 
 function calculate_jacobian(sys::NonlinearSystem;sparse=false,simplify=false)

src/systems/reduction.jl

@@ -8,6 +8,7 @@ function flatten(sys::ODESystem)
                          independent_variable(sys),
                          states(sys),
                          parameters(sys),
+                         pins=pins(sys),
                          observed=observed(sys))
     end
 end
@@ -56,8 +57,6 @@ end
 
 function alias_elimination(sys::ODESystem)
     eqs = vcat(equations(sys), observed(sys))
-    neweqs = Equation[]; sizehint!(neweqs, length(eqs))
-    subs = Pair[]
     diff_vars = filter(!isnothing, map(eqs) do eq
             if isdiffeq(eq)
                 arguments(eq.lhs)[1]
@@ -67,6 +66,9 @@ function alias_elimination(sys::ODESystem)
         end) |> Set
 
     deps = Set()
+    subs = Pair[]
+    neweqs = Equation[]; sizehint!(neweqs, length(eqs))
+
     for (i, eq) in enumerate(eqs)
         # only substitute when the variable is algebraic
         if isdiffeq(eq)
@@ -109,7 +111,7 @@ function alias_elimination(sys::ODESystem)
 
     alias_vars = first.(subs)
     sys_states = states(sys)
-    alias_eqs = alias_vars .~ last.(subs)
+    alias_eqs = topsort_equations(alias_vars .~ last.(subs), sys_states)
 
     newstates = setdiff(sys_states, alias_vars)
     ODESystem(neweqs, sys.iv, newstates, parameters(sys), observed=alias_eqs)