Skip to content

Commit d96b8c1

Browse files
committed
Add up to SD6
1 parent fe402c7 commit d96b8c1

File tree

1 file changed

+128
-6
lines changed

1 file changed

+128
-6
lines changed

src/ode/enright_pryce.jl

Lines changed: 128 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -2,14 +2,15 @@ using ModelingToolkit
22
using DiffEqBase, StaticArrays
33

44
@variables t y(t)[1:10]
5+
y = collect(y)
56
D = Differential(t)
67

78
sa1sys = let
89
sa1eqs = [D(y[1]) ~ -0.5 * y[1], D(y[2]) ~ -y[2], D(y[3]) ~ -100*y[3], D(y[4]) ~ -90*y[4]]
910
ODESystem(sa1eqs, t, name = :sa1)
1011
end
1112

12-
sa1prob = ODEProblem{false}(sa1sys, collect(y[1:4]) .=> 1.0, (0, 20.0), dt = 1e-2)
13+
sa1prob = ODEProblem{false}(sa1sys, y[1:4] .=> 1.0, (0, 20.0), dt = 1e-2)
1314

1415
sa2sys = let
1516
sa2eqs = [D(y[1]) ~ -1800 * y[1] + 900 * y[2]]
@@ -20,7 +21,7 @@ sa2sys = let
2021
ODESystem(sa2eqs, t, name = :sa2)
2122
end
2223

23-
sa2prob = ODEProblem{false}(sa2sys, collect(y[1:9]) .=> 0.0, (0, 120.0), dt = 5e-4)
24+
sa2prob = ODEProblem{false}(sa2sys, y[1:9] .=> 0.0, (0, 120.0), dt = 5e-4)
2425

2526
sa3sys = let
2627
sa3eqs = [
@@ -32,11 +33,132 @@ sa3sys = let
3233
ODESystem(sa3eqs, t, name = :sa3)
3334
end
3435

35-
sa3prob = ODEProblem{false}(sa3sys, collect(y[1:4]) .=> 1.0, (0, 20.0), dt = 1e-5)
36+
sa3prob = ODEProblem{false}(sa3sys, y[1:4] .=> 1.0, (0, 20.0), dt = 1e-5)
3637

3738
sa4sys = let
38-
sa3eqs = [D(y[i]) ~ -i^5 * y[i] for i in 1:10]
39-
ODESystem(sa3eqs, t, name = :sa3)
39+
sa4eqs = [D(y[i]) ~ -i^5 * y[i] for i in 1:10]
40+
ODESystem(sa4eqs, t, name = :sa4)
41+
end
42+
43+
sa4prob = ODEProblem{false}(sa4sys, y[1:10] .=> 1.0, (0, 1.0), dt = 1e-5)
44+
45+
#const SA_PROBLEMS = [sa1prob, sa2prob, sa3prob, sa4prob]
46+
47+
sb1sys = let
48+
sb1eqs = [D(y[1]) ~ -y[1] + y[2],
49+
D(y[2]) ~ -100y[1] - y[2],
50+
D(y[3]) ~ -100y[3] + y[4],
51+
D(y[4]) ~ -10_000y[3] - 100y[4]]
52+
53+
ODESystem(sb1eqs, t, name = :sb1)
54+
end
55+
56+
sb1prob = ODEProblem{false}(sb1sys, [y[[1, 3]] .=> 1.0; y[[2, 4]] .=> 0.0;], (0, 20.0), dt = 7e-3)
57+
58+
59+
@parameters α
60+
sb2sys = let
61+
sb2eqs = [D(y[1]) ~ -10y[1] + α * y[2],
62+
D(y[2]) ~ -α * y[1] - 10 * y[2],
63+
D(y[3]) ~ -4y[3],
64+
D(y[4]) ~ -y[4],
65+
D(y[5]) ~ -0.5y[5],
66+
D(y[6]) ~ -0.1y[6]]
67+
68+
ODESystem(sb2eqs, t, name = :sb2)
69+
end
70+
71+
sb2prob = ODEProblem{false}(sb2sys, y .=> 1.0, (0, 20.0), [α => 3], dt = 1e-2)
72+
sb3prob = ODEProblem{false}(sb2sys, y .=> 1.0, (0, 20.0), [α => 8], dt = 1e-2)
73+
sb4prob = ODEProblem{false}(sb2sys, y .=> 1.0, (0, 20.0), [α => 25], dt = 1e-2)
74+
sb5prob = ODEProblem{false}(sb2sys, y .=> 1.0, (0, 20.0), [α => 100], dt = 1e-2)
75+
76+
77+
sc1sys = let
78+
sc1eqs = [
79+
D(y[1]) ~ -y[1] + y[2]^2 + y[3]^2 + y[4]^2,
80+
D(y[2]) ~ -10y[2] + 10*(y[3]^2 + y[4]^2),
81+
D(y[3]) ~ -40y[3] + 40 * y[4]^2,
82+
D(y[4]) ~ -100y[4] + 2]
83+
84+
ODESystem(sc1eqs, t, name = :sc1)
85+
end
86+
87+
sc1prob = ODEProblem{false}(sc1sys, y .=> 1.0, (0, 20.0), dt = 1e-2)
88+
89+
90+
@parameters β
91+
sc2sys = let
92+
sc2eqs = [D(y[1]) ~ -y[1] + 2,
93+
D(y[2]) ~ -10y[2] + β * y[1]^2,
94+
D(y[3]) ~ -40y[3] + 4β * (y[1]^2 + y[2]^2),
95+
D(y[4]) ~ -100y[4] + 10β * (y[1]^2 + y[2]^2 + y[3]^2)]
96+
97+
ODESystem(sc2eqs, t, name = :sc2)
98+
end
99+
100+
sc2prob = ODEProblem{false}(sc2sys, y .=> 1.0, (0, 20.0), [β => 0.1], dt = 1e-2, cse = true)
101+
sc3prob = ODEProblem{false}(sc2sys, y .=> 1.0, (0, 20.0), [β => 1.0], dt = 1e-2, cse = true)
102+
sc3prob = ODEProblem{false}(sc2sys, y .=> 1.0, (0, 20.0), [β => 10.0], dt = 1e-2, cse = true)
103+
sc3prob = ODEProblem{false}(sc2sys, y .=> 1.0, (0, 20.0), [β => 20.0], dt = 1e-2, cse = true)
104+
105+
sd1sys = let
106+
sd1eqs = [D(y[1]) ~ 0.2 * (y[2] - y[1]),
107+
D(y[2]) ~ 10y[1] - (60 - 0.125y[3])*y[2] + 0.125y[3],
108+
D(y[3]) ~ 1]
109+
110+
ODESystem(sd1eqs, t, name = :sd1)
111+
end
112+
113+
sd1prob = ODEProblem{false}(sd1sys, y .=> 0.0, (0, 400.0), [β => 0.1], dt = 1.7e-2)
114+
115+
sd2sys = let
116+
sd2eqs = [D(y[1]) ~ -0.04y[1] + 0.01 * (y[2] * y[3]),
117+
D(y[2]) ~ 400y[1] - 100 * (y[2] * y[3]) - 3000 * y[2]^2,
118+
D(y[3]) ~ 30y[2]^2]
119+
120+
ODESystem(sd2eqs, t, name = :sd2)
121+
end
122+
123+
sd2prob = ODEProblem{false}(sd2sys, [y[1] => 1.0; y[2:3] .=> 0.0], (0, 40.0), dt = 1e-5, cse = true)
124+
125+
sd3sys = let
126+
sd3eqs = [D(y[1]) ~ y[3] - 100 * (y[1] * y[2]),
127+
D(y[2]) ~ y[3] + 2y[4] - 100 * (y[1] * y[2]) - 2e4 * y[2]^2,
128+
D(y[3]) ~ -y[3] + 100 * (y[1] * y[2]),
129+
D(y[4]) ~ -y[4] + 1e4*y[2]^2]
130+
131+
ODESystem(sd3eqs, t, name = :sd3)
132+
end
133+
134+
sd3prob = ODEProblem{false}(sd3sys, [y[1:2] .=> 1; y[3:4] .=> 0.0], (0, 20.0), dt = 2.5e-5, cse = true)
135+
136+
sd4sys = let
137+
sd4eqs = [D(y[1]) ~ -0.013y[1] - 1000 * (y[1] * y[3]),
138+
D(y[2]) ~ -2500*(y[2] * y[3]),
139+
D(y[3]) ~ -0.013y[1] - 1000 * (y[1] * y[3]) - 2500 * (y[2] * y[3])]
140+
141+
ODESystem(sd4eqs, t, name = :sd4)
142+
end
143+
144+
sd4prob = ODEProblem{false}(sd4sys, [y[1:2] .=> 1; y[3] => 0.0], (0, 50.0), dt = 2.9e-4, cse = true)
145+
146+
sd5sys = let
147+
sd5eqs = [D(y[1]) ~ 0.01 - (1 + (y[1] + 1000) * (y[1] + 1)) * (0.01 + y[1] + y[2]),
148+
D(y[2]) ~ 0.01 - (1 + y[2]^2) * (0.01 + y[1] + y[2])]
149+
150+
ODESystem(sd5eqs, t, name = :sd5)
151+
end
152+
153+
sd5prob = ODEProblem{false}(sd5sys, y[1:2] .=> 0.0, (0, 100.0), dt = 1e-4, cse = true)
154+
155+
sd6sys = let
156+
sd6eqs = [D(y[1]) ~ -y[1] + 10^8 * y[3] * (1 - y[1]),
157+
D(y[2]) ~ -10y[2] + 3e7 * y[3] * (1 - y[2]),
158+
D(y[3]) ~ -(-y[1] + 10^8 * y[3] * (1 - y[1])) - (-10y[2] + 3e7 * y[3] * (1 - y[2])),
159+
]
160+
161+
ODESystem(sd6eqs, t, name = :sd6)
40162
end
41163

42-
sa4prob = ODEProblem{false}(sa4sys, collect(y[1:10]) .=> 1.0, (0, 1.0), dt = 1e-5)
164+
sd6prob = ODEProblem{false}(sd6sys, [y[1] => 1.0; y[2:3] .=> 0.0], (0, 1.0), dt = 3.3e-8, cse = true)

0 commit comments

Comments
 (0)