Skip to content

Commit

Permalink
Convert LabelledArrays to ComponentArrays (#144)
Browse files Browse the repository at this point in the history
* Convert LabelledArrays to ComponentArrays
  • Loading branch information
jboubarcelo authored Jul 30, 2024
1 parent fb52920 commit 62c1a66
Show file tree
Hide file tree
Showing 16 changed files with 616 additions and 55 deletions.
3 changes: 1 addition & 2 deletions Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -5,8 +5,8 @@ version = "0.2.2"

[deps]
Catlab = "134e5e36-593f-5add-ad60-77f754baafbe"
ComponentArrays = "b0b7db55-cfe3-40fc-9ded-d10e2dbeff66"
Compose = "a81c6b42-2e10-5240-aca2-a61377ecd94b"
LabelledArrays = "2ee39098-c373-598a-b85f-a56591580800"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
LinearMaps = "7a12625a-238d-50fd-b39a-03d52299707e"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
Expand All @@ -32,7 +32,6 @@ AlgebraicPetri = "0.9"
Catlab = "0.15, 0.16"
Compose = "^0.9.1"
DelayDiffEq = "5"
LabelledArrays = "1"
LinearAlgebra = "1.9"
LinearMaps = "^3.1.0"
OrdinaryDiffEq = "^5.49.0, 6"
Expand Down
2 changes: 1 addition & 1 deletion docs/Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -2,11 +2,11 @@
AlgebraicDynamics = "5fd6ff03-a254-427e-8840-ba658f502e32"
AlgebraicPetri = "4f99eebe-17bf-4e98-b6a1-2c4f205a959b"
Catlab = "134e5e36-593f-5add-ad60-77f754baafbe"
ComponentArrays = "b0b7db55-cfe3-40fc-9ded-d10e2dbeff66"
DelayDiffEq = "bcd4f6db-9728-5f36-b5f7-82caef46ccdb"
DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
DynamicalSystems = "61744808-ddfa-5f27-97ff-6e42cc95d634"
LabelledArrays = "2ee39098-c373-598a-b85f-a56591580800"
Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306"
NLSolvers = "337daf1e-9722-11e9-073e-8b9effe078ba"
NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
Expand Down
10 changes: 4 additions & 6 deletions docs/src/AlgebraicPetri.md
Original file line number Diff line number Diff line change
Expand Up @@ -58,15 +58,13 @@ open_bruss = Open([:A, :D], Brusselator, [:A, :B])
rs = ContinuousResourceSharer{Float64}(open_bruss)
````

!!! note
@example AlgPetri
SciML integration for LabelledArrays is broken due to this [issue](https://github.com/SciML/LabelledArrays.jl/issues/162) on LabelledArrays.jl. Once that issue is closed, you should be able to solve LabelledPetriNets with SciML Solvers.

````julia
using LabelledArrays
using ComponentArrays
tspan = (0.0,100.0)
params = LVector(t1=1.0, t2=1.2, t3=3.14, t4=0.1)
u0 = LVector(A=1.0, B=3.17, D=0.0, E=0.0, X=1.0, Y=1.9)
params = ComponentArray(t1=1.0, t2=1.2, t3=3.14, t4=0.1)
u0 = ComponentArray(A=1.0, B=3.17, D=0.0, E=0.0, X=1.0, Y=1.9)
eval_dynamics(rs, u0, params, 0.0)
prob = ODEProblem((u,p,t) -> eval_dynamics(rs, u, p, t), u0, tspan, params)
sol = solve(prob, Tsit5())
plot(sol, idxs=[:X, :Y])
Expand Down
13 changes: 7 additions & 6 deletions examples/Cyber-Physical.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,7 @@

using AlgebraicDynamics
using Catlab

using LabelledArrays
using ComponentArrays
using DifferentialEquations
using Plots

Expand Down Expand Up @@ -58,10 +57,11 @@ to_graphviz(UAV)


# Then we assign behaviors to inhabit the boxes.

function 𝗟(𝐖)
𝐿(u, x, p, t) = [ -p.𝓐l * (u[1] - x[1] - x[2]) ] # sc
𝐶(u, x, p, t) = [ -p.𝓐c * (u[1] + p.𝓑c*x[1] - x[2]) ] # sl
𝐷(u, x, p, t) = LVector= -0.313*u[1] + 56.7*u[2] + 0.232*x[1],
𝐷(u, x, p, t) = ComponentArray= -0.313*u[1] + 56.7*u[2] + 0.232*x[1],
q = -0.013*u[1] - 0.426*u[2] + 0.0203*x[1],
θ = 56.7*u[2] )

Expand All @@ -81,17 +81,18 @@ end
# Lastly, we compute and plot the solution.

## initial values
xₒ = LVector( e = 0.01, # [e, d] -> [θ offset, 𝛿 control input]

x₀ = ComponentArray( e = 0.01, # [e, d] -> [θ offset, 𝛿 control input]
d = 0.05);

uₒ = [0.0, 0, 0, 0, 0]
u₀ = [0.0, 0, 0, 0, 0]
tspan = (0, 20.0)

params = (𝓐l = 100, # decay constant of sensor
𝓐c = 100, # decay constant of controller
𝓑c = 0) # ratio of velocity to reference velocity

prob = ODEProblem(𝑢ᵤₐᵥ, uₒ, xₒ, tspan, params)
prob = ODEProblem(𝑢ᵤₐᵥ, u₀, x₀, tspan, params)
sol = solve(prob, alg_hints=[:stiff]);

#-
Expand Down
5 changes: 3 additions & 2 deletions examples/Ecosystem.jl
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@
using AlgebraicDynamics
using Catlab

using LabelledArrays
using ComponentArrays
using OrdinaryDiffEq
using Plots, Plots.PlotMeasures

Expand All @@ -16,7 +16,8 @@ using Plots, Plots.PlotMeasures
# - parameters $\gamma$ represent the rate at which a population of predators grows in a predation interaction
# - parameters $\delta$ represent the rate at with a species population declines

params = LVector(αr=0.3, βrf=0.015, γrf=0.015, δf=0.7,

params = ComponentArray(αr=0.3, βrf=0.015, γrf=0.015, δf=0.7,
βrh=0.01, γrh=0.01, δh=0.5,
γfishh=0.001, βfishh=0.003,
αfish=0.35, βfishF=0.015, γfishF=0.015,
Expand Down
10 changes: 5 additions & 5 deletions examples/Lotka-Volterra.jl
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,7 @@
using AlgebraicDynamics
using Catlab.WiringDiagrams, Catlab.Programs

using LabelledArrays
using ComponentArrays
using OrdinaryDiffEq, Plots, Plots.PlotMeasures

const UWD = UndirectedWiringDiagram
Expand All @@ -47,10 +47,11 @@ end

## Compose
rabbitfox_system = oapply(rf, [rabbit_growth, rabbitfox_predation, fox_decline])


## Solve and plot
u0 = [10.0, 100.0]
params = LVector=.3, β=0.015, γ=0.015, δ=0.7)
params = ComponentArray=.3, β=0.015, γ=0.015, δ=0.7)
tspan = (0.0, 100.0)

prob = ODEProblem(rabbitfox_system, u0, tspan, params)
Expand Down Expand Up @@ -100,7 +101,7 @@ rabbitfox_system = oapply(rabbitfox_pattern, [rabbit, fox])

## Solve and plot
u0 = [10.0, 100.0]
params = LVector=.3, β=0.015, γ=0.015, δ=0.7)
params = ComponentArray=.3, β=0.015, γ=0.015, δ=0.7)
tspan = (0.0, 100.0)

prob = ODEProblem(rabbitfox_system, u0, tspan, params)
Expand All @@ -123,9 +124,8 @@ rabbitfox_pattern = barbell(1)
## Compose
rabbitfox_system = oapply(rabbitfox_pattern, [rabbit, fox])

## Solve and plot
u0 = [10.0, 100.0]
params = LVector=.3, β=0.015, γ=0.015, δ=0.7)
params = ComponentArray=.3, β=0.015, γ=0.015, δ=0.7)
tspan = (0.0, 100.0)

prob = ODEProblem(rabbitfox_system, u0, tspan, params)
Expand Down
8 changes: 4 additions & 4 deletions examples/Ross-Macdonald.jl
Original file line number Diff line number Diff line change
Expand Up @@ -7,8 +7,7 @@
using AlgebraicDynamics
using Catlab.WiringDiagrams, Catlab.Graphics

using LabelledArrays

using ComponentArrays
using DelayDiffEq, DifferentialEquations
using Plots

Expand Down Expand Up @@ -106,9 +105,10 @@ malaria_model = oapply(rm,
)

#-
params = LVector(a = 0.3, b = 0.55, c = 0.15,
params = ComponentArray(a = 0.3, b = 0.55, c = 0.15,
g = 0.1, n = 10, r = 1.0/200, m = 0.5)


u0 = [0.1, 0.3]
tspan = (0.0, 365.0*2)

Expand Down Expand Up @@ -254,7 +254,7 @@ malaria_delay_model = oapply(rm,
Dict(:humans => human_delay_model, :mosquitos => mosquito_delay_model)
)
#-
params = LVector(a = 0.3, b = 0.55, c = 0.15,
params = ComponentArray(a = 0.3, b = 0.55, c = 0.15,
g = 0.1, n = 10, r = 1.0/200, m = 0.5)

u0_delay = [0.09, .01, 0.3]
Expand Down
7 changes: 4 additions & 3 deletions ext/AlgebraicDynamicsAlgebraicPetriExt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@ module AlgebraicDynamicsAlgebraicPetriExt
using AlgebraicPetri
using AlgebraicDynamics
using Catlab
using LabelledArrays
using ComponentArrays
import AlgebraicDynamics.UWDDynam: ContinuousResourceSharer

export OpenNet
Expand All @@ -21,16 +21,17 @@ function dynamics(pn::OpenPetriNet,T::Type,ns::Int64)
vf(u, p, t) = begin f!(storage, u, p, t); return storage end
end


function dynamics(pn::OpenLabelledPetriNet,T::Type,ns::Int64)
f! = vectorfield(apex(pn))
storage = LVector(NamedTuple{tuple(snames(apex(pn))...)}(zeros(T,ns)))
storage = ComponentArray(NamedTuple{tuple(snames(apex(pn))...)}(zeros(T,ns)))
vf(u, p, t) = begin f!(storage, u, p, t); return storage end
return vf
end

function dynamics(pn::OpenLabelledReactionNet,T::Type,ns::Int64)
f! = vectorfield(apex(pn))
storage = LVector(NamedTuple{tuple(snames(apex(pn))...)}(zeros(T,ns)))
storage = ComponentArray(NamedTuple{tuple(snames(apex(pn))...)}(zeros(T,ns)))
rt = rates(apex(pn))
vf(u, p, t) = begin f!(storage, u, rt, t); return storage end
return vf
Expand Down
1 change: 0 additions & 1 deletion notebooks/Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,6 @@ Animations = "27a7e980-b3e6-11e9-2bcd-0b925532e340"
Catlab = "134e5e36-593f-5add-ad60-77f754baafbe"
Javis = "78b212ba-a7f9-42d4-b726-60726080707e"
JavisNB = "92afb270-2599-44f6-96a1-44c6efb1daf1"
LabelledArrays = "2ee39098-c373-598a-b85f-a56591580800"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
Expand Down
2 changes: 1 addition & 1 deletion notebooks/multi-city SIR.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@
"\n",
"using AlgebraicDynamics, AlgebraicDynamics.DWDDynam\n",
"\n",
"using LabelledArrays\n",
"using ComponentArrays\n",
"\n",
"using OrdinaryDiffEq\n",
"using Plots"
Expand Down
573 changes: 567 additions & 6 deletions notebooks/springs/springs.ipynb

Large diffs are not rendered by default.

4 changes: 2 additions & 2 deletions test/Project.toml
Original file line number Diff line number Diff line change
@@ -1,9 +1,9 @@
[deps]
AlgebraicDynamics = "5fd6ff03-a254-427e-8840-ba658f502e32"
AlgebraicPetri = "4f99eebe-17bf-4e98-b6a1-2c4f205a959b"
Catlab = "134e5e36-593f-5add-ad60-77f754baafbe"
ComponentArrays = "b0b7db55-cfe3-40fc-9ded-d10e2dbeff66"
DelayDiffEq = "bcd4f6db-9728-5f36-b5f7-82caef46ccdb"
LabelledArrays = "2ee39098-c373-598a-b85f-a56591580800"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
PrettyTables = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d"
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
AlgebraicPetri = "4f99eebe-17bf-4e98-b6a1-2c4f205a959b"
2 changes: 1 addition & 1 deletion test/cpg_dynam.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@ using AlgebraicDynamics.DWDDynam
using AlgebraicDynamics.CPortGraphDynam: draw, barbell, gridpath, grid, meshpath
using Catlab

using LabelledArrays
using ComponentArrays

using Test
using PrettyTables
Expand Down
2 changes: 1 addition & 1 deletion test/dwd_dynam.jl
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
using AlgebraicDynamics.DWDDynam
using Catlab

using LabelledArrays
using ComponentArrays
using DelayDiffEq
using Test

Expand Down
15 changes: 8 additions & 7 deletions test/ext/AlgebraicPetri.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@ using AlgebraicPetri
using Catlab
using Test
using AlgebraicDynamics.UWDDynam
using LabelledArrays
using ComponentArrays


@testset "AlgebraicPetri" begin
Expand All @@ -21,18 +21,19 @@ using LabelledArrays
@test nstates(rs) == 1
@test nports(rs) == 1
@test portmap(rs) == [1]
u = LVector(foo=1.0)
p = LVector(bar=2.0)
@test eval_dynamics(rs, u, p, 0) == [2.0]
u = ComponentArray(foo=1.0)
p = ComponentArray(bar=2.0)
@test eval_dynamics(rs, u, p, 0) == ComponentArray(foo=2.0)


labelled_birth_react = Open(LabelledReactionNet{Float64}{Float64}([:foo => 1.0], ((:bar => 2.0), (:foo => (:foo, :foo)))))
rs = ContinuousResourceSharer{Float64}(labelled_birth_react)
@test nstates(rs) == 1
@test nports(rs) == 1
@test portmap(rs) == [1]
u = LVector(NamedTuple(zip(snames(apex(labelled_birth_react)), Vector(subpart(apex(labelled_birth_react), :concentration)))))
p = LVector(bar=2.0)
@test eval_dynamics(rs, u, p, 0) == [2.0]
u = ComponentArray(NamedTuple(zip(snames(apex(labelled_birth_react)), Vector(subpart(apex(labelled_birth_react), :concentration)))))
p = ComponentArray(bar=2.0)
@test eval_dynamics(rs, u, p, 0) == ComponentArray(foo=2.0)

Brusselator = LabelledPetriNet([:A, :B, :D, :E, :X, :Y],
:t1 => (:A => (:X, :A)),
Expand Down
14 changes: 7 additions & 7 deletions test/ext/DelayDiffEq.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ using AlgebraicDynamics
using Catlab
using DelayDiffEq
using Test

using ComponentArrays
# DWDDynam Integration
######################

Expand Down Expand Up @@ -56,8 +56,8 @@ add_wires!(d_big, Pair[
u0 = [2.7]
x0 = [10.0]
τ = 10.0
p = LVector= τ)

p = ComponentArray= τ)
@test eval_dynamics(delay, u0, x0, hist, p) == [0.0]
@test eval_dynamics(delay_copy, u0, x0, hist, p) == [0.0]

Expand Down Expand Up @@ -109,14 +109,14 @@ add_wires!(d_big, Pair[
mult_delay = DelayMachine{Float64}(1,1,1, df, delay ? delay_rf : rf)
add_delay = DelayMachine{Float64}(1,1,1, ef, delay ? delay_rf : rf)

prob = DDEProblem(oapply(d, [mult_delay, add_delay]), u0, x0, hist, (0, 4.0), LVector= 2.0))
prob = DDEProblem(oapply(d, [mult_delay, add_delay]), u0, x0, hist, (0, 4.0), ComponentArray= 2.0))
sol1 = solve(prob, alg; dtmax = 0.1)
if delay
f = (u,h,p,t) -> [ (h(p,t - p.τ)[2] + x0[1]) * h(p,t - p.τ)[1], h(p, t - p.τ)[1] + h(p,t - p.τ)[2] ]
else
f = (u,h,p,t) -> [ (u[2] + x0[1]) * h(p,t - p.τ)[1], u[1] + h(p,t - p.τ)[2] ]
end
prob = DDEProblem(f, u0, hist, (0, 4.0), LVector= 2.0))
prob = DDEProblem(f, u0, hist, (0, 4.0), ComponentArray= 2.0))
sol2 = solve(prob, alg; dtmax = 0.1)
@test last(sol1) == last(sol2)
end
Expand Down Expand Up @@ -324,8 +324,8 @@ const UWD = UndirectedWiringDiagram
1, 1, 1, dxdt_delay, (u,h,p,t) -> [[u[1], h(p, t - p.n)[1]]])
rm_model = oapply(c, [mosquito_delay_model, human_delay_model])

params = LVector(a = 0.3, b = 0.55, c = 0.15,
g = 0.1, n = 10, r = 1.0/200, m = 0.5)
params = ComponentArray(a = 0.3, b = 0.55, c = 0.15,
g = 0.1, n = 10, r = 1.0/200, m = 0.5)

u0_delay = [0.09, .01, 0.3]
tspan = (0.0, 365.0*5)
Expand Down

0 comments on commit 62c1a66

Please sign in to comment.