# Define an ODE which depends on an external function with the @ode\_def macro

**URL:** <https://discourse.julialang.org/t/define-an-ode-which-depends-on-an-external-function-with-the-ode-def-macro/98067>\
**Category:** Numerics\
**Tags:** question, differentialequation\
**Created:** [April 28, 2023, 10:42pm UTC](https://discourse.julialang.org/t/define-an-ode-which-depends-on-an-external-function-with-the-ode-def-macro/98067 "2023-04-28T22:42:26Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![atmyo4](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/atmyo4/32/49587_2.png) [@atmyo4](https://discourse.julialang.org/u/atmyo4)\
**Post date:** [April 28, 2023, 10:42pm UTC](https://discourse.julialang.org/t/define-an-ode-which-depends-on-an-external-function-with-the-ode-def-macro/98067/1 "2023-04-28T22:42:26Z")

</div>

Hello! I want to solve a system of ODEs that includes a parameter that varies over time. To obtain the value of the parameter I use an external function that uses the solution of another system of ODEs that has already been solved. I obtained the desired result by defining the second system with a generic function, but I want to see if it is possible to achieve it by defining it with the macro @ode\_def of ParameterizedFunctions.jl.

I have thought of two ways to do this:

1. Just define the second system with @ode\_def and include the function in the respective equation, but I think this is not possible.
2. Use a callback. This does not throw me any error but the result is not correct and I don’t know why.

Thank you for your time. I hope I have been clear enough. I attach the code below.

> **Code that generates the desired result using a generic function to define the second system of ODEs**
>
> ```julia
> using DifferentialEquations
> using ParameterizedFunctions
> 
> # Solve the first system of ODEs 
> odeFB = @ode_def begin 
> dGp = ((Vmax * Ge * Gp)/(Km + Gp)) - (k4 * Gp) - (k2 * A * Gp) + (km2 * A) + (km4 * Gi)
> dGi = (k4 * Gp) - (k5 * Gi) - (km4 * Gi) + (km5 * A)
> dA = - (k3 * A) + (k5 * Gi) - (km5 * A)
> end Vmax Ge Km k2 k3 k4 k5 km2 km4 km5
> 
> x0 = ones(length(odeFB.syms))
> p = (1400, 24.15, 1100, 5.3, 4, 2, 2.3, 5.3, 2, 2.3)
> tstops = 80000
> tspan = (0,81074)
> solFB = solve(ODEProblem(odeFB,x0,tspan,p); reltol = 1e-12, tstops = tstops);
> 
> # External function to calculate a constant of the equation du[2] of the second system 
> function f(t)
> if t > 0
> i = findfirst(times .>= t)
> return ((((sol[i] - sol[i-1]) / (times[i] - times[i-1])) * (t .- times[i-1])) .+ sol[i-1])
> else 
> return x0[1]
> end
> end 
> 
> # Solve the second system of ODEs 
> function odeNF(du, u, p, t)
> Vmax, Ge, Km, k2, k3, k4, k5, km2, km4, km5 = p 
> du[1] = ((Vmax * Ge * u[1])/(Km + u[1])) - (k4 * u[1]) - (k2 * u[3] * u[1]) + (km2 * u[3]) + (km4 * u[2])
> du[2] = (k4 * f(t))- (k5 * u[2]) - (km4 * u[2]) + (km5 * u[3])
> du[3] = - (k3 * u[3]) + (k5 * u[2]) - (km5 * u[3])
> end 
> 
> k = findfirst(solFB.t .== tstops)
> x0 = solFB[k]
> sol = solFB[1, k:end]
> times = solFB.t[k:end] .- solFB.t[k]
> tspan = times[end]
> solNF = solve(ODEProblem(odeNF, x0, tspan, p); reltol = 1e-12, saveat = times);
> 
> ```

> **Just define the second system with @ode\_def**
>
> The code would be exactly the same as the previous one, only that the second system is defined in this way:
> 
> ```julia
> odeNF = @ode_def begin 
> dGp = ((Vmax * Ge * Gp)/(Km + Gp)) - (k4 * Gp) - (k2 * A * Gp) + (km2 * A) + (km4 * Gi)
> dGi = (k4 * f(t)) - (k5 * Gi) - (km4 * Gi) + (km5 * A)
> dA = - (k3 * A) + (k5 * Gi) - (km5 * A)
> end Vmax Ge Km k2 k3 k4 k5 km2 km4 km5 
> 
> ```
> 
> However, doing this throws an error → LoadError: TypeError: non-boolean (Term{Bool, Nothing}) used in boolean context in expression starting at Untitled-1.ipynb:28

> **Callback**
>
> ```julia
> using DifferentialEquations
> using ParameterizedFunctions
> 
> # Solve the first system of ODEs 
> odeFB = @ode_def begin 
> dGp = ((Vmax * Ge * Gp)/(Km + Gp)) - (k4 * Gp) - (k2 * A * Gp) + (km2 * A) + (km4 * Gi)
> dGi = (k4 * Gp) - (k5 * Gi) - (km4 * Gi) + (km5 * A)
> dA = - (k3 * A) + (k5 * Gi) - (km5 * A)
> end Vmax Ge Km k2 k3 k4 k5 k6 km2 km4 km5
> 
> x0 = ones(length(odeFB.syms))
> p = [1400, 24.15, 1100, 5.3, 4, 2, 2.3, NaN, 5.3, 2, 2.3]
> tstops = 80000
> tspan = (0,81074)
> solFB = solve(ODEProblem(odeFB,x0,tspan,p); reltol = 1e-12, tstops = tstops);
> 
> # Set callback 
> condition(u, t, integrator) = t > 0
> 
> function affect!(integrator)
> i = findfirst(times .>= integrator.t)
> integrator.p[8] = integrator.p[6] * ((((sol[i] - sol[i-1]) / (times[i] - times[i-1])) * (integrator.t .- times[i-1])) .+ sol[i-1])
> end 
> 
> cb = DiscreteCallback(condition, affect!, save_positions=(false,false))
> 
> # Solve the second system of ODEs 
> odeNF = @ode_def begin 
> dGp = ((Vmax * Ge * Gp)/(Km + Gp)) - (k4 * Gp) - (k2 * A * Gp) + (km2 * A) + (km4 * Gi)
> dGi = (k4 * k6) - (k5 * Gi) - (km4 * Gi) + (km5 * A)
> dA = - (k3 * A) + (k5 * Gi) - (km5 * A)
> end Vmax Ge Km k2 k3 k4 k5 k6 km2 km4 km5 
> 
> k = findfirst(solFB.t .== tstops)
> x0 = solFB[k]
> p[8] = x0[1]
> sol = solFB[1, k:end]
> times = solFB.t[k:end] .- solFB.t[k]
> tspan = times[end]
> solNF = solve(ODEProblem(odeNF, x0, tspan, p); reltol = 1e-12, saveat = times, callback = cb);
> 
> ```
> 
> This does not throw an error but the result is not correct.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [May 6, 2023, 12:21am UTC](https://discourse.julialang.org/t/define-an-ode-which-depends-on-an-external-function-with-the-ode-def-macro/98067/2 "2023-05-06T00:21:28Z")

</div>

This is a duplicate of another thread which has the answer:

> [@How to define a system of ODEs with @ode\_def and update the value of a variable using a function that depends on time?](https://discourse.julialang.org/t/how-to-define-a-system-of-odes-with-ode-def-and-update-the-value-of-a-variable-using-a-function-that-depends-on-time/98404/2):
>
> You need to @register\_symbolic it so that way it’s kept as a node in the symbolic graph. using DifferentialEquations using ParameterizedFunctions # Solve odeFB odeFB = @ode\_def begin dGp = ((Vmax \* Ge \* Gp)/(Km + Gp)) - (k4 \* Gp) - (k2 \* A \* Gp) + (km2 \* A) + (km4 \* Gi) dGi = (k4 \* Gp) - (k5 \* Gi) - (km4 \* Gi) + (km5 \* A) dA = - (k3 \* A) + (k5 \* Gi) - (km5 \* A) end Vmax Ge Km k2 k3 k4 k5 km2 km4 km5 x0 = ones(length(odeFB.syms)) p = (1400, 24.15, 1100, 5.3, 4, 2, 2.3, 5.3, 2, 2.…