# Time dependent inputs to reaction network in Catalyst

**URL:** <https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078>\
**Category:** Chemistry\
**Tags:** catalyst\
**Created:** [October 9, 2024, 2:26am UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078 "2024-10-09T02:26:54Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![contejuz](https://avatars.discourse-cdn.com/v4/letter/c/e79b87/32.png) [@contejuz](https://discourse.julialang.org/u/contejuz)\
**Post date:** [October 9, 2024, 2:26am UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078/1 "2024-10-09T02:26:54Z")

</div>

I’m trying to create a simple reaction network that takes external input. Here is one example (adapted from the docs):

```julia
# Create model.
model = @reaction_network begin
    kB*I(t), S + E --> SE
    kD, SE --> S + E
    kP, SE --> P + E
end

```

Here I want `I(t)` to be an arbitrary external function that I can pass in along with my parameters. If I were to use `OrdinaryDiffEq` to solve this as an inhomogenous ODE, here’s what I would do (directly from [Getting Started with Differential Equations in Julia · DifferentialEquations.jl](https://docs.sciml.ai/DiffEqDocs/stable/getting_started)):

```julia
l = 1.0 # length [m]
m = 1.0 # mass [kg]
g = 9.81 # gravitational acceleration [m/s²]

function pendulum!(du, u, p, t)
    du[1] = u[2] # θ'(t) = ω(t)
    du[2] = -3g / (2l) * sin(u[1]) + 3 / (m * l^2) * p(t) # ω'(t) = -3g/(2l) sin θ(t) + 3/(ml^2)M(t)
end

θ₀ = 0.01 # initial angular deflection [rad]
ω₀ = 0.0 # initial angular velocity [rad/s]
u₀ = [θ₀, ω₀] # initial state vector
tspan = (0.0, 10.0) # time interval

M = t -> 0.1sin(t) # external torque [Nm]

prob = ODEProblem(pendulum!, u₀, tspan, M)
sol = solve(prob)

```

Here, the function `M` is passed as a parameter and called directly in `f(u,p,t)`. However, when I try to similary call an external function `I` from my reaction network, I get an undefined variable error. I was hoping that Catalyst would automatically recognize this as time dependent function. I could hardcode the function inside the reaction network, but I want the flexibility to change this input without having to recompile the model.

Is this possible in Catalyst?

---

<div class="post-metadata">

**Author:** ![johannesnauta](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/johannesnauta/32/47434_2.png) [@johannesnauta](https://discourse.julialang.org/u/johannesnauta)\
**Post date:** [October 9, 2024, 1:00pm UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078/2 "2024-10-09T13:00:29Z")

</div>

When you build the reactions programmatically, it seems to work:

```julia
@variables t
@species S(t) E(t) SE(t) P(t)
@parameters kB kD kP

I(t) = 0.1*sin(t)

rxB = Reaction(kB*I(t), [S,E], [SE])
rxD = Reaction(kD, [SE], [S, E])
rxP = Reaction(kP, [SE], [P,E])

@named rs = ReactionSystem([rxB, rxD, rxP], t)

```

and this seems also what you actually mean to do, i.e.

```julia
julia> rs
Model rs
Unknowns (4):
  S(t)
  E(t)
  SE(t)
  P(t)
Parameters (3):
  kB
  kD
  kP

julia> reactions(rs)
3-element Vector{Reaction}:
 0.1kB*sin(t), S + E --> SE
 kD, SE --> S + E
 kP, SE --> P + E

```

and you can create the `ODEProblem` simply by giving the initial conditions and the values of the parameters, as you did in your example. My hunch (but I might be wrong!) on why your reaction network did not work is that, at the time of creating the reaction, the function `I(t)` is unknown, so `Catalyst` does not know what to do. Doing it the programmatic way creates `I(t)` as a function in the scope, which can then be used as input to your equations. I have done this a few times in the past.

---

<div class="post-metadata">

**Author:** ![contejuz](https://avatars.discourse-cdn.com/v4/letter/c/e79b87/32.png) [@contejuz](https://discourse.julialang.org/u/contejuz)\
**Post date:** [October 9, 2024, 1:57pm UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078/3 "2024-10-09T13:57:27Z")

</div>

Thanks! That makes sense. I guess if I wanted to continue using the DSL, I could create a similar global function and call it within `@reaction_network` and recompile when I change the function.

I’m guessing that you would also have to recompile `rxB` and `rs` if you change `I(t)`? I’m trying it in Pluto but it automatically recompiles `rs` so I can’t tell. I guess performance-wise it doesn’t hurt to recompile `rs` if `I(t)` changes?

---

<div class="post-metadata">

**Author:** ![johannesnauta](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/johannesnauta/32/47434_2.png) [@johannesnauta](https://discourse.julialang.org/u/johannesnauta)\
**Post date:** [October 9, 2024, 2:12pm UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078/4 "2024-10-09T14:12:01Z")

</div>

Good question. While setting the parameters, and also changing them, is straightforward (check `ModelingToolkit.setp`, or `remake`, but I recommend the first), when the function changes in its entirety I believe the problem should be remade from scratch — but I am not sure of this. Others with more expertise on this manner should step in for that.

What you can do though (but I am not sure if this is what you want), is make your function be a function of a parameter. So something like (not tested),

```julia
@parameters f
I(f,t) = sin(f*t)

```

Then, you can relatively easily change the frequency and solve the ODE with a bunch of different frequencies when you change `f` in the `prob_func` of an `EnsembleProblem`, for example. But, again, if the function changes, i.e. from `sin(t)` to `cos(t)`, I am not sure you can avoid a full re-creation of the `ReactionSystem` and/or the `ODEProblem`.

---

<div class="post-metadata">

**Author:** ![contejuz](https://avatars.discourse-cdn.com/v4/letter/c/e79b87/32.png) [@contejuz](https://discourse.julialang.org/u/contejuz)\
**Post date:** [October 9, 2024, 8:52pm UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078/5 "2024-10-09T20:52:04Z")

</div>

Awesome thanks! That’s exactly what I have now. It doesn’t give me full control over the form of the function, but with enough parameters I can get enough flexibility.

---

<div class="post-metadata">

**Author:** ![Torkel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/torkel/32/5030_2.png) [@Torkel](https://discourse.julialang.org/u/Torkel)\
**Post date:** [February 11, 2025, 3:20pm UTC](https://discourse.julialang.org/t/time-dependent-inputs-to-reaction-network-in-catalyst/121078/6 "2025-02-11T15:20:25Z")

</div>

I am late to the party, but is there any reason that simply writing e.g.

```julia
using Catalyst
myfunc(t) = 1/t
model = @reaction_network begin
    kB*myfunc(t), S + E --> SE
    kD, SE --> S + E
    kP, SE --> P + E
end

```

doesn’t work? It works fine for me, and seems to be what you want to achieve?