AlgebraicJulia / GATlab.jl

GATlab: a computer algebra system based on generalized algebraic theories (GATs)
https://algebraicjulia.github.io/GATlab.jl/
MIT License
24 stars 2 forks source link

Lens model of periodic parameters for another model #48

Closed jpfairbanks closed 1 year ago

jpfairbanks commented 1 year ago

@olynch, I was trying to use the new lens stuff to build an SIR model where beta and gamma oscillate.

sir = @lens ThRing begin
  dom = [s, i, r] | [ds, di, dr, dβ, dγ]
  codom = [s, i, r] | [dβ, dγ, β, γ]
  expose = begin
    s = s
    i = i
    r = r
  end
  update = begin
    ds = -β * (s * i)
    di = β * (s * i) + (- γ) * i
    dr = γ * i
    dβ = dβ
    dγ = dγ
  end
end

periodic_params = @lens ThRing begin
  dom = [s, i, r] | [dβ, dγ, β, γ]
  codom = [] | [β₀, kᵦ, γ₀, kᵧ]
  expose = begin
  end
  update = begin
    dβ = -kᵦ*(β - β₀)
    dγ = -kᵧ*(γ - γ₀)
    β = β
    γ = γ
  end
end

I get the error ERROR: KeyError: key SymLit(:β) not found when executing the definition of periodic_params.

jpfairbanks commented 1 year ago

Did I do it right?

sir = @lens ThRing begin
  dom = [s, i, r, β, γ] | [ds, di, dr, dβ, dγ]
  codom = [s, i, r, β, γ] | [dβ, dγ]
  expose = begin
    s = s
    i = i
    r = r
    β = β
    γ = γ
  end
  update = begin
    ds = -β * (s * i)
    di = β * (s * i) + (- γ) * i
    dr = γ * i
    dβ = dβ
    dγ = dγ
  end
end
periodic_params = @lens ThRing begin
  dom = [s, i, r, β, γ] | [dβ, dγ]
  codom = [i, β, γ] | [β₀, kᵦ, γ₀, kᵧ]
  expose = begin
    i = i
    β = β
    γ = γ
  end
  update = begin
    dβ = -kᵦ*(β - β₀)
    dγ = -kᵧ*(γ - γ₀)
  end
end

composed = compose(sir, periodic_params)
jpfairbanks commented 1 year ago

This is the output. It looks good.

dom: Gatlab.Stdlib.StdModels.SimpleLenses.Impl.Arena{ThRing}([s: 1, i: 1, r: 1, β: 1, γ: 1], [ds: 1, di: 1, dr: 1, dβ: 1, dγ: 1])
codom: Gatlab.Stdlib.StdModels.SimpleLenses.Impl.Arena{ThRing}([i: 1, β: 1, γ: 1], [β₀: 1, kᵦ: 1, γ₀: 1, kᵧ: 1])
expose: 
i = i
β = β
γ = γ

update: 
ds = (-(β) * (s * i))
di = ((β * (s * i)) + (-(γ) * i))
dr = (γ * i)
dβ = (-(kᵦ) * (β - β₀))
dγ = (-(kᵧ) * (γ - γ₀))
olynch commented 1 year ago

The way that I was planning on doing this was with 3 lenses. One is normal sir, one is the oscillating system, and then the last one is like a wiring diagram. I compose sir and oscillating in parallel, and then compose in series with the wiring diagram.

I'm not exactly sure what went wrong in yours, but I think something is off in that the domain of sir has three state variables and five update variables; the dimension of the tangent space should be the same as the dimension of the space.

Edit: I see the specific problem. In your update function, \beta isn't in the list of variables given by the domain positions plus the codomain directions.

olynch commented 1 year ago

Did I do it right?

sir = @lens ThRing begin
  dom = [s, i, r, β, γ] | [ds, di, dr, dβ, dγ]
  codom = [s, i, r, β, γ] | [dβ, dγ]
  expose = begin
    s = s
    i = i
    r = r
    β = β
    γ = γ
  end
  update = begin
    ds = -β * (s * i)
    di = β * (s * i) + (- γ) * i
    dr = γ * i
    dβ = dβ
    dγ = dγ
  end
end
periodic_params = @lens ThRing begin
  dom = [s, i, r, β, γ] | [dβ, dγ]
  codom = [i, β, γ] | [β₀, kᵦ, γ₀, kᵧ]
  expose = begin
    i = i
    β = β
    γ = γ
  end
  update = begin
    dβ = -kᵦ*(β - β₀)
    dγ = -kᵧ*(γ - γ₀)
  end
end

composed = compose(sir, periodic_params)

Yeah, this looks right!

jpfairbanks commented 1 year ago

I was trying to do it with three lenses too! But then I gave up on that. Should I make a PR with this in it?

olynch commented 1 year ago

Yeah, sure! Maybe this can go in an examples folder?

jpfairbanks commented 1 year ago

By the way, what are we going to do about solving?

olynch commented 1 year ago

I'm writing code that converts a lens to an ODESystem, which then can be used with any diffeq solvers.

olynch commented 1 year ago

Incidentally, there's a great Lawvere paper which mentions that integration of ODEs is given by left and right adjoints to the functor which takes a flow (i.e. R x X -> X) to the vector field on X given by taking the derivative at 0.

Page 11-12 of https://raw.githubusercontent.com/mattearnshaw/lawvere/master/pdfs/1978-categorical-dynamics.pdf