@@ -32,10 +32,9 @@ test_equal(a, b) = @test isequal(simplify(a, polynorm=true), simplify(b, polynor
32
32
eqs = [
33
33
D (x) ~ σ* (y- x)
34
34
D (y) ~ x* (ρ- z)- y + β
35
- 0 ~ sin (z) - x + y
36
- sin (u) ~ x + y
37
- 2 x ~ 3 a
38
- 2 u ~ 3 x
35
+ 0 ~ z - x + y
36
+ 0 ~ a + z
37
+ u ~ z + a
39
38
]
40
39
41
40
lorenz1 = ODESystem (eqs,t,name= :lorenz1 )
@@ -46,14 +45,13 @@ io = IOBuffer(); show(io, lorenz1_aliased); str = String(take!(io))
46
45
reduced_eqs = [
47
46
D (x) ~ σ* (y - x)
48
47
D (y) ~ β + x* (ρ - z) - y
49
- 0 ~ y + sin (z) - x
50
- 0 ~ x + y - sin (1.5 x)
48
+ 0 ~ y + z - x
51
49
]
52
50
test_equal .(equations (lorenz1_aliased), reduced_eqs)
53
- @test isempty (setdiff (states (lorenz1_aliased), [u, x, y, z]))
51
+ @test isempty (setdiff (states (lorenz1_aliased), [x, y, z]))
54
52
test_equal .(observed (lorenz1_aliased), [
55
- a ~ 2 / 3 * x ,
56
- u ~ 3 / 2 * x ,
53
+ u ~ 0 ,
54
+ a ~ - z ,
57
55
])
58
56
59
57
# Multi-System Reduction
@@ -70,8 +68,6 @@ lorenz = name -> ODESystem(eqs1,t,name=name)
70
68
lorenz1 = lorenz (:lorenz1 )
71
69
ss = ModelingToolkit. get_structure (initialize_system_structure (lorenz1))
72
70
@test isequal (ss. fullvars, [x, y, z, D (x), D (y), D (z), F, u])
73
- @test ss. solvable_graph. fadjlist == [[7 ], [8 ], [], [8 ]]
74
- @test ss. solvable_graph. metadata == [[1 ], [1 ], [], [1 ]]
75
71
lorenz2 = lorenz (:lorenz2 )
76
72
77
73
connected = ODESystem ([s ~ a + lorenz1. x
@@ -81,16 +77,19 @@ connected = ODESystem([s ~ a + lorenz1.x
81
77
82
78
# Reduced Flattened System
83
79
84
- reduced_system = alias_elimination (connected; conservative = false )
80
+ reduced_system = alias_elimination (connected)
85
81
86
82
@test setdiff (states (reduced_system), [
83
+ s
87
84
a
88
85
lorenz1. x
89
86
lorenz1. y
90
87
lorenz1. z
88
+ lorenz1. u
91
89
lorenz2. x
92
90
lorenz2. y
93
91
lorenz2. z
92
+ lorenz2. u
94
93
]) |> isempty
95
94
96
95
@test setdiff (parameters (reduced_system), [
@@ -103,21 +102,21 @@ reduced_system = alias_elimination(connected; conservative=false)
103
102
]) |> isempty
104
103
105
104
reduced_eqs = [
106
- 0 ~ a + lorenz1. x - (lorenz2. y)
107
- D (lorenz1. x) ~ lorenz2. x + lorenz2. y + lorenz1. σ* ((lorenz1. y) - (lorenz1. x)) - (lorenz2. z)
108
- D (lorenz1. y) ~ lorenz1. x* (lorenz1. ρ - (lorenz1. z)) - ((lorenz1. x) + (lorenz1. y) - (lorenz1. z))
105
+ 0 ~ a + lorenz1. x - (s)
106
+ 0 ~ s - (lorenz2. y)
107
+ D (lorenz1. x) ~ lorenz2. u + lorenz1. σ* ((lorenz1. y) - (lorenz1. x))
108
+ D (lorenz1. y) ~ lorenz1. x* (lorenz1. ρ - (lorenz1. z)) - (lorenz1. u)
109
109
D (lorenz1. z) ~ lorenz1. x* lorenz1. y - (lorenz1. β* (lorenz1. z))
110
- D (lorenz2. x) ~ lorenz1. x + lorenz1. y + lorenz2. σ* ((lorenz2. y) - (lorenz2. x)) - (lorenz1. z)
111
- D (lorenz2. y) ~ lorenz2. x* (lorenz2. ρ - (lorenz2. z)) - ((lorenz2. x) + (lorenz2. y) - (lorenz2. z))
110
+ 0 ~ lorenz1. x + lorenz1. y - (lorenz1. u) - (lorenz1. z)
111
+ D (lorenz2. x) ~ lorenz1. u + lorenz2. σ* ((lorenz2. y) - (lorenz2. x))
112
+ D (lorenz2. y) ~ lorenz2. x* (lorenz2. ρ - (lorenz2. z)) - (lorenz2. u)
112
113
D (lorenz2. z) ~ lorenz2. x* lorenz2. y - (lorenz2. β* (lorenz2. z))
114
+ 0 ~ lorenz2. x + lorenz2. y - (lorenz2. u) - (lorenz2. z)
113
115
]
114
116
115
117
test_equal .(equations (reduced_system), reduced_eqs)
116
118
117
119
observed_eqs = [
118
- s ~ a + lorenz1. x
119
- lorenz1. u ~ lorenz1. x + lorenz1. y - lorenz1. z
120
- lorenz2. u ~ lorenz2. x + lorenz2. y - lorenz2. z
121
120
lorenz2. F ~ lorenz1. u
122
121
lorenz1. F ~ lorenz2. u
123
122
]
@@ -133,32 +132,19 @@ pp = [
133
132
]
134
133
u0 = [
135
134
a => 1.0
135
+ s => 1.0
136
136
lorenz1. x => 1.0
137
137
lorenz1. y => 0.0
138
138
lorenz1. z => 0.0
139
+ lorenz1. u => 0.0
139
140
lorenz2. x => 1.0
140
141
lorenz2. y => 0.0
141
142
lorenz2. z => 0.0
143
+ lorenz2. u => 0.0
142
144
]
143
145
prob1 = ODEProblem (reduced_system, u0, (0.0 , 100.0 ), pp)
144
146
solve (prob1, Rodas5 ())
145
147
146
- # issue #578
147
-
148
- let
149
- @variables t x (t) y (t) z (t)
150
- D = Differential (t)
151
- eqs = [
152
- D (x) ~ x + y
153
- x ~ y
154
- ]
155
- sys = ODESystem (eqs, t)
156
- asys = alias_elimination (sys)
157
-
158
- test_equal .(equations (asys), [D (x) ~ 2 x])
159
- test_equal .(observed (asys), [y ~ x])
160
- end
161
-
162
148
# issue #724 and #716
163
149
let
164
150
@parameters t
@@ -186,29 +172,25 @@ end
186
172
187
173
# NonlinearSystem
188
174
@parameters t
189
- @variables u1 (t) u2 (t) u3 (t) u4 (t) u5 (t)
175
+ @variables u1 (t) u2 (t) u3 (t)
190
176
@parameters p
191
177
eqs = [
192
- 2 u1 ~ 3 u5
193
- u2 ~ u1
194
- u3 ~ 2 u1 - u2
195
- u4 ~ u2 + u3^ p
196
- u5 ~ u4^ 2 - u1
178
+ u1 ~ u2
179
+ u3 ~ u1 + u2 + p
180
+ u3 ~ hypot (u1, u2) * p
197
181
]
198
- sys = NonlinearSystem (eqs, [u1, u2, u3, u4, u5 ], [p])
182
+ sys = NonlinearSystem (eqs, [u1, u2, u3], [p])
199
183
reducedsys = alias_elimination (sys)
200
- @test observed (reducedsys) == [u1 ~ 3 / 2 * u5 ]
184
+ @test observed (reducedsys) == [u1 ~ u2 ]
201
185
202
186
u0 = [
203
187
u1 => 1
204
188
u2 => 1
205
189
u3 => 0.3
206
- u4 => 0.6
207
- u5 => 2 / 3
208
190
]
209
191
pp = [2 ]
210
192
nlprob = NonlinearProblem (reducedsys, u0, pp)
211
193
reducedsol = solve (nlprob, NewtonRaphson ())
212
- residual = fill (100.0 , 4 )
194
+ residual = fill (100.0 , length ( states (reducedsys)) )
213
195
nlprob. f (residual, reducedsol. u, pp)
214
196
@test all (x-> abs (x) < 1e-5 , residual)
0 commit comments