Tearing tests pass

Author: YingboMa

src/structural_transformation/codegen.jl

@@ -197,7 +197,7 @@ function build_torn_function(
     )
 
     s = structure(sys)
-    states = s.fullvars[diffvars_range(s)]
+    states = map(i->s.fullvars[i], diffvars_range(s))
     syms = map(Symbol, states)
 
     expr = SymbolicUtils.Code.toexpr(
@@ -266,8 +266,8 @@ function build_observed_function(
     syms_set = Set(syms)
     s = structure(sys)
     @unpack partitions, fullvars, graph = s
-    diffvars = fullvars[diffvars_range(s)]
-    algvars = fullvars[algvars_range(s)]
+    diffvars = map(i->fullvars[i], diffvars_range(s))
+    algvars = map(i->fullvars[i], algvars_range(s))
 
     required_algvars = Set(intersect(algvars, syms_set))
     obs = observed(sys)
@@ -337,7 +337,7 @@ function ODAEProblem{iip}(
                          ) where {iip}
     s = structure(sys)
     @unpack fullvars = s
-    dvs = fullvars[diffvars_range(s)]
+    dvs = map(i->fullvars[i], diffvars_range(s))
     ps = parameters(sys)
     defs = defaults(sys)
 

src/structural_transformation/tearing.jl

@@ -177,11 +177,14 @@ function tearing_reassemble(sys; simplify=false)
 
     obseqs = solvars .~ rhss
 
-    @set! sys.structure.graph = newgraph
-    @set! sys.structure.scc = newscc
-    @set! sys.structure.fullvars = fullvars[active_vars]
-    @set! sys.structure.partitions = newpartitions
-    @set! sys.structure.algeqs = newalgeqs
+    @set! s.graph = newgraph
+    @set! s.scc = newscc
+    @set! s.fullvars = fullvars[active_vars]
+    @set! s.vartype = s.vartype[active_vars]
+    @set! s.partitions = newpartitions
+    @set! s.algeqs = newalgeqs
+
+    @set! sys.structure = s
     @set! sys.eqs = neweqs
     @set! sys.states = newstates
     @set! sys.reduced_states = [get_reduced_states(sys); solvars]

src/systems/systemstructure.jl

@@ -61,9 +61,9 @@ isdervar(s::SystemStructure, var::Integer) = s.vartype[var] === DERIVATIVE_VARIA
 isdiffvar(s::SystemStructure, var::Integer) = s.vartype[var] === DIFFERENTIAL_VARIABLE
 isalgvar(s::SystemStructure, var::Integer) = s.vartype[var] === ALGEBRAIC_VARIABLE
 
-dervars_range(s::SystemStructure) = Iterators.filter(Base.Fix1(s, isdervar), eachindex(s.vartype))
-diffvars_range(s::SystemStructure) = Iterators.filter(Base.Fix1(s, isdiffvar), eachindex(s.vartype))
-algvars_range(s::SystemStructure) = Iterators.filter(Base.Fix1(s, isalgeq), eachindex(s.vartype))
+dervars_range(s::SystemStructure) = Iterators.filter(Base.Fix1(isdervar, s), eachindex(s.vartype))
+diffvars_range(s::SystemStructure) = Iterators.filter(Base.Fix1(isdiffvar, s), eachindex(s.vartype))
+algvars_range(s::SystemStructure) = Iterators.filter(Base.Fix1(isalgvar, s), eachindex(s.vartype))
 
 isalgeq(s::SystemStructure, eq::Integer) = s.algeqs[eq]
 isdiffeq(s::SystemStructure, eq::Integer) = !isalgeq(s, eq)

test/structural_transformation/tearing.jl

@@ -29,15 +29,12 @@ show(io, sss)
 prt = String(take!(io))
 
 if VERSION >= v"1.6"
-@test prt == """xvars: Any[]
-dxvars: Any[]
-algvars: Any[u1(t), u2(t), u3(t), u4(t), u5(t)]
-Incidence matrix:
- ×  ⋅  ⋅  ⋅  ×
+@test prt == "Incidence matrix:
  ×  ×  ⋅  ⋅  ⋅
- ×  ×  ×  ⋅  ⋅
- ⋅  ×  ×  ×  ⋅
- ×  ⋅  ⋅  ×  ×"""
+ ×  ⋅  ×  ⋅  ⋅
+ ×  ⋅  ×  ×  ⋅
+ ⋅  ⋅  ×  ×  ×
+ ×  ×  ⋅  ⋅  ×"
 end
 
 # u1 = f1(u5)
@@ -49,8 +46,8 @@ sys = initialize_system_structure(sys)
 find_solvables!(sys)
 sss = structure(sys)
 @unpack graph, solvable_graph, assign, partitions = sss
-@test graph.fadjlist == [[1, 5], [1, 2], [1, 2, 3], [2, 3, 4], [1, 4, 5]]
-@test solvable_graph.fadjlist == map(x->[x], 1:5)
+@test graph.fadjlist == [[1, 2], [1, 3], [1, 3, 4], [3, 4, 5], [1, 2, 5]]
+@test solvable_graph.fadjlist == map(x->[x], [1, 3, 4, 5, 2])
 
 tornsys = tearing(sys)
 sss = structure(tornsys)
@@ -148,7 +145,7 @@ eqs = [
 daesys = ODESystem(eqs, t)
 newdaesys = tearing(daesys)
 @test equations(newdaesys) == [D(x) ~ z; 0 ~ x + sin(z) - p*t]
-@test isequal(states(newdaesys), [x, z])
+@test isequal(states(newdaesys), [z, x])
 @test_throws ArgumentError ODAEProblem(newdaesys, [x=>1.0], (0, 1.0))
 prob = ODAEProblem(newdaesys, [x=>1.0], (0, 1.0), [p=>0.2])
 du = [0.0]; u = [1.0]; pr = 0.2; tt = 0.1