Fig 7.21 - MMSB

repressilator model

using OrdinaryDiffEq
using ComponentArrays: ComponentArray
using SimpleUnPack
using CairoMakie

function model721!(D, u, p, t)
    @unpack α, α0, n, β = p
    @unpack mA, mB, mC, A, B, C = u
    D.mA = α / (1 + C^n) + α0 - mA
    D.mB = α / (1 + A^n) + α0 - mB
    D.mC = α / (1 + B^n) + α0 - mC
    D.A = β * (mA - A)
    D.B = β * (mB - B)
    D.C = β * (mC - C)
    nothing
end

model721! (generic function with 1 method)

ps721 = ComponentArray(
    α = 298.2,
    α0 = 0.03,
    n = 2.0,
    β = 0.2
)
u0721 = ComponentArray(
    mA = 0.2,
    mB = 0.3,
    mC = 0.4,
    A = 0.1,
    B = 0.1,
    C = 0.5
)

tspan = (0.0, 400.0)

prob721 = ODEProblem(model721!, u0721, tspan, ps721)

ODEProblem with uType ComponentArrays.ComponentVector{Float64, Vector{Float64}, Tuple{ComponentArrays.Axis{(mA = 1, mB = 2, mC = 3, A = 4, B = 5, C = 6)}}} and tType Float64. In-place: true
Non-trivial mass matrix: false
timespan: (0.0, 400.0)
u0: ComponentVector{Float64}(mA = 0.2, mB = 0.3, mC = 0.4, A = 0.1, B = 0.1, C = 0.5)

@time sol721 = solve(prob721, KenCarp47())

  2.718222 seconds (9.60 M allocations: 647.257 MiB, 16.08% gc time, 102.21% compilation time)

retcode: Success
Interpolation: 3rd order Hermite
t: 90-element Vector{Float64}:
   0.0
   0.0014286712480576914
   0.015715383728634603
   0.11569483191195365
   0.27128377991731645
   0.4855260770621106
   0.8019632936409467
   1.3268305727420127
   2.175436806177097
   3.152367394116522
   ⋮
 359.3461786716006
 364.3621423574526
 369.50657030688245
 373.4523114874227
 379.2473603223246
 385.0424091572265
 390.83745799212835
 397.74010422928177
 400.0
u: 90-element Vector{ComponentArrays.ComponentVector{Float64, Vector{Float64}, Tuple{ComponentArrays.Axis{(mA = 1, mB = 2, mC = 3, A = 4, B = 5, C = 6)}}}}:
 ComponentVector{Float64}(mA = 0.2, mB = 0.3, mC = 0.4, A = 0.1, B = 0.1, C = 0.5)
 ComponentVector{Float64}(mA = 0.5403361928101721, mB = 0.7211224652974476, mC = 0.8209781853695932, A = 0.1000771993468159, B = 0.1001173119895755, C = 0.5000315834241665)
 ComponentVector{Float64}(mA = 3.910351841334096, mB = 4.897484476883872, mC = 4.995334228504263, A = 0.1061593599856239, B = 0.1078655992669062, C = 0.5069212046571355)
 ComponentVector{Float64}(mA = 23.581351528532206, mB = 31.368439235693298, mC = 31.003082686509487, A = 0.39137876640500624, B = 0.47564522130856657, C = 0.8664356684728488)
 ComponentVector{Float64}(mA = 34.750272262056036, mB = 52.535701548167204, mC = 47.11815083412874, A = 1.3246358053503131, B = 1.8193687372272824, C = 2.1168450640636887)
 ComponentVector{Float64}(mA = 33.904530160876085, mB = 54.182870954272396, mC = 44.62085070376357, A = 2.7334531048276616, B = 4.033840650382787, C = 3.984656060392383)
 ComponentVector{Float64}(mA = 27.72505709091729, mB = 45.32876740377409, mC = 35.21078971393252, A = 4.4586645888267595, B = 6.851879185989259, C = 6.186154409260148)
 ComponentVector{Float64}(mA = 18.58156493323217, mB = 30.766358317677685, mC = 22.496415363730996, A = 6.287681024849206, B = 9.912382812174412, C = 8.39150281192476)
 ComponentVector{Float64}(mA = 9.958522250551779, mB = 16.55046694321977, mC = 10.991041772867716, A = 7.430020771011377, B = 11.892150724001857, C = 9.553518970362447)
 ComponentVector{Float64}(mA = 5.831838074676068, mB = 9.49887806438225, mC = 5.428967612857784, A = 7.438417643050988, B = 11.97224737632704, C = 9.220160987249487)
 ⋮
 ComponentVector{Float64}(mA = 25.671258390712744, mB = 8.329502479094478, mC = 0.07681505537493012, A = 9.148639164537085, B = 69.55081599384637, C = 2.7501302597177313)
 ComponentVector{Float64}(mA = 110.13920048909831, mB = 0.45834205127732935, mC = 0.3269175950302204, A = 47.0877042960301, B = 26.812022980924297, C = 1.1245435300035869)
 ComponentVector{Float64}(mA = 146.92722376624346, mB = 0.06733386191283046, mC = 2.287864128854908, A = 112.63170977504758, B = 9.670364960717244, C = 1.150198672495098)
 ComponentVector{Float64}(mA = 39.01340189300651, mB = 0.05650766533658095, mC = 10.52077524526735, A = 98.81633696700739, B = 4.423406747205617, C = 3.7750874595222865)
 ComponentVector{Float64}(mA = 1.173546751359415, mB = 0.1866077958695485, mC = 74.70953376104684, A = 37.050492672689124, B = 1.4636695304563394, C = 29.16254536526912)
 ComponentVector{Float64}(mA = 0.09284365926452581, mB = 1.5500055798408008, mC = 159.49176346091508, A = 11.798470615812858, B = 0.9798564400428811, C = 103.64544117845003)
 ComponentVector{Float64}(mA = 0.061107759644824594, mB = 14.434285136833898, mC = 18.84789076331902, A = 3.744049043792873, B = 5.152084679147567, C = 88.86373150826076)
 ComponentVector{Float64}(mA = 0.36851732213213606, mB = 116.22311198062266, mC = 0.3237712737717123, A = 1.081185213474111, B = 50.98267819674458, C = 25.118809895610323)
 ComponentVector{Float64}(mA = 0.855931648104186, mB = 154.2511765678168, mC = 0.11437575112284777, A = 0.9042973000745865, B = 82.92122176892073, C = 16.052561049712704)

fig = Figure()
ax = Axis(fig[1, 1],
    xlabel = "Time",
    ylabel = "Concentration",
    title = "Fig 7.21"
)
lines!(ax, 0..400, t-> sol721(t).A, label = "A")
lines!(ax, 0..400, t-> sol721(t).B, label = "B")
lines!(ax, 0..400, t-> sol721(t).C, label = "C")
axislegend(ax, position = :lt)
fig

This notebook was generated using Literate.jl.