@@ -5,19 +5,17 @@ One can directly use symbolic variables to index into SciML solution objects.
5
5
Moreover, observables can also be evaluated in this way. For example,
6
6
consider the system
7
7
``` @example faq1
8
- using Catalyst, DifferentialEquations , Plots
8
+ using Catalyst, OrdinaryDiffEq , Plots
9
9
rn = @reaction_network ABtoC begin
10
10
(k₊,k₋), A + B <--> C
11
11
end
12
-
13
- # initial condition and parameter values
14
- setdefaults!(rn, [:A => 1.0, :B => 2.0, :C => 0.0, :k₊ => 1.0, :k₋ => 1.0])
15
12
nothing # hide
16
13
```
17
14
Let's convert it to a system of ODEs, using the conservation laws of the system
18
15
to eliminate two of the species:
19
16
``` @example faq1
20
17
osys = convert(ODESystem, rn; remove_conserved = true)
18
+ osys = complete(osys)
21
19
```
22
20
Notice the resulting ODE system has just one ODE, while algebraic observables
23
21
have been added for the two removed species (in terms of the conservation law
@@ -28,7 +26,9 @@ observed(osys)
28
26
Let's solve the system and see how to index the solution using our symbolic
29
27
variables
30
28
``` @example faq1
31
- oprob = ODEProblem(osys, [], (0.0, 10.0), [])
29
+ u0 = [osys.A => 1.0, osys.B => 2.0, osys.C => 0.0]
30
+ ps = [osys.k₊ => 1.0, osys.k₋ => 1.0]
31
+ oprob = ODEProblem(osys, u0, (0.0, 10.0), ps)
32
32
sol = solve(oprob, Tsit5())
33
33
```
34
34
Suppose we want to plot just species ` C ` , without having to know its integer
44
44
```
45
45
To evaluate ` C ` at specific times and plot it we can just do
46
46
``` @example faq1
47
- t = range(0.0, 10.0, length= 101)
48
- plot(t, sol(t, idxs = C), label = "C(t)", xlabel = "t")
47
+ t = range(0.0, 10.0, length = 101)
48
+ plot(sol(t, idxs = C), label = "C(t)", xlabel = "t")
49
49
```
50
50
If we want to get multiple variables we can just do
51
51
``` @example faq1
@@ -65,7 +65,7 @@ constant, giving `k*X^2/2` instead of `k*X^2` for ODEs and `k*X*(X-1)/2` instead
65
65
of ` k*X*(X-1) ` for jumps. This can be disabled when directly ` convert ` ing a
66
66
[ ` ReactionSystem ` ] ( @ref ) . If ` rn ` is a generated [ ` ReactionSystem ` ] ( @ref ) , we can
67
67
do
68
- ``` julia
68
+ ``` @example faq1
69
69
osys = convert(ODESystem, rn; combinatoric_ratelaws=false)
70
70
```
71
71
Disabling these rescalings should work for all conversions of ` ReactionSystem ` s
@@ -87,7 +87,9 @@ rx1 = Reaction(k,[B,C],[B,D], [2.5,1],[3.5, 2.5])
87
87
rx2 = Reaction(2*k, [B], [D], [1], [2.5])
88
88
rx3 = Reaction(2*k, [B], [D], [2.5], [2])
89
89
@named mixedsys = ReactionSystem([rx1, rx2, rx3], t, [A, B, C, D], [k, b])
90
+ mixedsys = complete(mixedsys)
90
91
osys = convert(ODESystem, mixedsys; combinatoric_ratelaws = false)
92
+ osys = complete(osys)
91
93
```
92
94
Note, when using ` convert(ODESystem, mixedsys; combinatoric_ratelaws=false) ` the
93
95
` combinatoric_ratelaws=false ` parameter must be passed. This is also true when
@@ -127,6 +129,7 @@ t = default_t()
127
129
rx1 = Reaction(β, [S, I], [I], [1,1], [2])
128
130
rx2 = Reaction(ν, [I], [R])
129
131
@named sir = ReactionSystem([rx1, rx2], t)
132
+ sir = complete(sir)
130
133
oprob = ODEProblem(sir, [], (0.0, 250.0))
131
134
sol = solve(oprob, Tsit5())
132
135
plot(sol)
@@ -162,7 +165,7 @@ Julia `Symbol`s corresponding to each variable/parameter to their values, or
162
165
from ModelingToolkit symbolic variables/parameters to their values. Using
163
166
` Symbol ` s we have
164
167
``` @example faq4
165
- using Catalyst, DifferentialEquations
168
+ using Catalyst, OrdinaryDiffEq
166
169
rn = @reaction_network begin
167
170
α, S + I --> 2I
168
171
β, I --> R
@@ -199,6 +202,7 @@ the second example, or one can use the `symmap_to_varmap` function to convert th
199
202
` Symbol ` mapping to a symbolic mapping. I.e. this works
200
203
``` @example faq4
201
204
osys = convert(ODESystem, rn)
205
+ osys = complete(osys)
202
206
203
207
# this works
204
208
u0 = symmap_to_varmap(rn, [:S => 999.0, :I => 1.0, :R => 0.0])
@@ -221,6 +225,7 @@ rx1 = @reaction k, A --> 0
221
225
rx2 = @reaction $f, 0 --> A
222
226
eq = f ~ (1 + sin(t))
223
227
@named rs = ReactionSystem([rx1, rx2, eq], t)
228
+ rs = complete(rs)
224
229
osys = convert(ODESystem, rs)
225
230
```
226
231
In the final ODE model, ` f ` can be eliminated by using
0 commit comments