Fig 7.17

Goodwin oscillator model: https://en.wikipedia.org/wiki/Goodwin_model_(biology)

using Catalyst
using ModelingToolkit
using OrdinaryDiffEq
using Plots
Plots.default(linewidth=2)

rn = @reaction_network begin
    (a / (k^n + Z^n), b), 0 <--> X
    (α * X, β), 0 <--> Y
    (γ * Y, δ), 0 <--> Z
end

\[\begin{split} \begin{align*} \varnothing &\xrightleftharpoons[b]{\frac{a}{k^{n} + Z^{n}}} \mathrm{X} \\ \varnothing &\xrightleftharpoons[\beta]{X \alpha} \mathrm{Y} \\ \varnothing &\xrightleftharpoons[\delta]{Y \gamma} \mathrm{Z} \end{align*} \end{split}\]

ps = [
    :a => 360,
    :k => 1.368,
    :b => 1,
    :α => 1,
    :β => 0.6,
    :γ => 1,
    :δ => 0.8,
    :n => 12
]

u0 = zeros(3)
tend = 35.0

35.0

prob = ODEProblem(rn, u0, tend, ps);
sol = solve(prob)

retcode: Success
Interpolation: 3rd order Hermite
t: 74-element Vector{Float64}:
  0.0
  9.999999999999999e-5
  0.0010999999999999998
  0.011099999999999997
  0.041561068228358915
  0.08505007717877648
  0.1400425786140992
  0.2146417543721738
  0.31201584322447784
  0.4390615514552524
  ⋮
 30.399216275051014
 31.019732084200392
 31.782434145750724
 32.3743621530088
 32.9240804229318
 33.47916779303648
 33.94226811022297
 34.51649669667564
 35.0
u: 74-element Vector{Vector{Float64}}:
 [0.0, 0.0, 0.0]
 [0.0008380069445083983, 4.1900207549215655e-8, 1.3966642737483395e-12]
 [0.009213468964138544, 5.067222002693897e-6, 1.8578451380267915e-9]
 [0.09250904681199522, 0.0005132339779722273, 1.8975593394034844e-6]
 [0.34116338821048725, 0.0070794913530376055, 9.780450237545616e-5]
 [0.683292323356918, 0.02897046980274338, 0.0008166279287941536]
 [1.095152006170422, 0.07629612549303964, 0.0035279811532653834]
 [1.6188553354727238, 0.1723356709595111, 0.012150996243756597]
 [2.246222895432165, 0.3461133927935709, 0.03523230985115964]
 [2.978081158465377, 0.6417560311991968, 0.09109022137495452]
 ⋮
 [0.5089918673039778, 1.766781968948389, 2.2228229156189716]
 [0.2880973416389677, 1.4141811569212317, 2.1240480319262973]
 [0.19686438454763613, 1.0306516147476084, 1.8357226708588616]
 [0.41078846058561214, 0.857256616858513, 1.5778546803847777]
 [1.2775009216741202, 0.9943227874496436, 1.414180979126931]
 [2.091840335870593, 1.5588431984569113, 1.4789312290822632]
 [1.7685358003391847, 1.98200785662198, 1.7201935243220485]
 [1.0772532733699713, 2.0788289937151023, 2.0395698694667495]
 [0.6803805235209356, 1.9141052038293493, 2.188474444402658]

plot(sol, title="Fig 7.17 (A)", xlabel="Time", ylabel="Concentration")

plot(sol, idxs=(rn.X, rn.Y, rn.Z), title="Fig 7.17 (B)", legend=false, size=(600, 600))

---

This notebook was generated using Literate.jl.