|
| 1 | +# Modeling Discrete Systems |
| 2 | + |
| 3 | +In this example, we will use the new [`DiscreteSystem`](@ref) API |
| 4 | +to create an SIR model. |
| 5 | + |
| 6 | +```@example discrete |
| 7 | +using ModelingToolkit |
| 8 | +using ModelingToolkit: t_nounits as t |
| 9 | +using OrdinaryDiffEq: FunctionMap |
| 10 | +
|
| 11 | +@inline function rate_to_proportion(r, t) |
| 12 | + 1 - exp(-r * t) |
| 13 | +end |
| 14 | +@parameters c δt β γ |
| 15 | +@constants h = 1 |
| 16 | +@variables S(t) I(t) R(t) |
| 17 | +k = ShiftIndex(t) |
| 18 | +infection = rate_to_proportion(β * c * I / (S * h + I + R), δt * h) * S |
| 19 | +recovery = rate_to_proportion(γ * h, δt) * I |
| 20 | +
|
| 21 | +# Equations |
| 22 | +eqs = [S(k + 1) ~ S - infection * h, |
| 23 | + I(k + 1) ~ I + infection - recovery, |
| 24 | + R(k + 1) ~ R + recovery] |
| 25 | +@mtkbuild sys = DiscreteSystem(eqs, t) |
| 26 | +
|
| 27 | +u0 = [S => 990.0, I => 10.0, R => 0.0] |
| 28 | +p = [β => 0.05, c => 10.0, γ => 0.25, δt => 0.1] |
| 29 | +tspan = (0.0, 100.0) # value function (from Symbolics) is used to convert a Num to Float64 |
| 30 | +prob = DiscreteProblem(sys, u0, tspan, p) |
| 31 | +sol = solve(prob, FunctionMap()) |
| 32 | +``` |
| 33 | + |
| 34 | +Currently, all updates happen at integer intervals (`t = 1, 2, 3...`). All shifts |
| 35 | +in the equations must be positive. For example, instead of specifying the update equation |
| 36 | +for the Fibonacci series as `x ~ x(k-1) + x(k-2)` use `x(k+2) ~ x(k+1) + x`. |
0 commit comments