# Is there a way to run differential equations to a specific point and then change input values, continue running, and graph it together?

**URL:** <https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232>\
**Category:** New to Julia\
**Tags:** question, diffeq\
**Created:** [January 3, 2019, 10:23am UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232 "2019-01-03T10:23:23Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![HelgavonLichtenstein](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/helgavonlichtenstein/32/5005_2.png) [@HelgavonLichtenstein](https://discourse.julialang.org/u/HelgavonLichtenstein)\
**Post date:** [January 3, 2019, 10:23am UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/1 "2019-01-03T10:23:23Z")

</div>

I have a simple 3-species food chain that I would like to run for 200 steps so the system stabilizes. I would then like to start harvesting and show how populations fluctuate. Is there a way to do this?

The natural/non-harvested system is:

> using DifferentialEquations  
> f = @ode\_def\_bare Pristine begin  
> dBB = r\*(1- BB/K)_BB -xI_yIB_BI_((BB^h/(B0^h + BB^h))/eIB)  
> dBI = -xI_BI + xI_yIB_BI_(BB^h/(B0^h + BB^h))-xT _yTI_BT\*((BI^h/(B0^h + BI^h))/eTI)  
> dBT = -xT_BT + xT_yTI_BT_(BI^h/(B0^h + BI^h))-q_BT_E  
> end K r yIB yTI eIB eTI B0 h xI xT q E  
> u0 = [155.0,107.0,93.0]  
> tspan = (0.0,1000.0)  
> p = (700, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0, 0)  
> prob = ODEProblem(f,u0,tspan,p)  
> sol=solve(prob,Rosenbrock23())

And it looks like this (_sorry at the moment I haven’t figured out how to add a main title to the plots and not for each subplot_):

 ![image](https://global.discourse-cdn.com/julialang/original/3X/8/e/8efe17a542f00668d10586c0851016e1dc5cd709.png)

Plot code is

> using Plots  
> plot(sol,xlabel = “Time” ,ylabel = “Density”, title = “Food Chain Pristine”, lw=0.5, layout = (3,1))

Then I would like to start harvesting after 200 steps so q = 0.01 and E = 50 (instead of both 0 in the pristine scenario).

Is this possible?

---

<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:** [January 3, 2019, 10:28am UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/2 "2019-01-03T10:28:28Z")

</div>

The easiest way is to just make it a single diffeq solve with a callback that implements the change.

---

<div class="post-metadata">

**Author:** ![HelgavonLichtenstein](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/helgavonlichtenstein/32/5005_2.png) [@HelgavonLichtenstein](https://discourse.julialang.org/u/HelgavonLichtenstein)\
**Post date:** [January 7, 2019, 10:46am UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/3 "2019-01-07T10:46:24Z")

</div>

Thank you for your response here and on my user\_cache question. I have tried to implement the callback using the [Callback example](http://docs.juliadiffeq.org/latest/features/diffeq_arrays.html#Example:-A-Control-Problem-1) but I don’t know how to have two parameters changed. Here is my attempt which generates `ERROR: LoadError: invalid redefinition of constant f`.

I probably need to put an f2 component here

```julia
using DifferentialEquations
mutable struct SimType{T} <: DEDataVector{T}
    x::Array{T,1}
    f1::T

```

and another error source is probably stating initial values of q and E as 0.0 while specifying 0.01 and 50 under parameters. q and E do not change over time in this model (they will in future ones) so perhaps they should not be put into u0?

```julia
u0 = [155.0,107.0,93.0,0.0,0.0]
tspan = (0.0,1000.0)
p = (450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50)

```

Whole code:

```julia
using DifferentialEquations
mutable struct SimType{T} <: DEDataVector{T}
    x::Array{T,1}
    f1::T
end
f = @ode_def_bare ManagedLinear begin
dBB = r*(1- BB/K)*BB -xI*yIB*BI*((BB^h/(B0^h + BB^h))/eIB)
dBI = -xI*BI + xI*yIB*BI*(BB^h/(B0^h + BB^h))-xT *yTI*BT*((BI^h/(B0^h + BI^h))/eTI)
dBT = -xT*BT + xT*yTI*BT*(BI^h/(B0^h + BI^h))-q*BT*E + u.f1 +u.f2
end K r yIB yTI eIB eTI B0 h xI xT q E
const tstop1 = [200.]

function condition(u,t,integrator)
  t in tstop1
end

function affect!(integrator)
  for c in full_cache(integrator)
    c.f1 = 0.1
  end
end

function affect2!(integrator)
  for c in full_cache(integrator)
    c.f2 = 50
  end
end
save_positions = (true,true)

cb = DiscreteCallback(condition, affect!, save_positions=save_positions)

save_positions = (false,true)

cb2 = DiscreteCallback(condition2, affect2!, save_positions=save_positions)

cbs = CallbackSet(cb,cb2)
u0 = [155.0,107.0,93.0,0.0,0.0]
tspan = (0.0,1000.0)
p = (450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50)
# K r yIB yTI eIB eTI B0 h xI xT q E
prob = ODEProblem(f,u0,tspan,p)
sol=solve(prob,Rosenbrock23())

using Plots

plot(sol, title = "Food Chain Managed Linear", xlabel = "Time" ,ylabel = "Density", lw=0.5, layout = (3,1))

```

---

<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:** [January 8, 2019, 5:42pm UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/4 "2019-01-08T17:42:43Z")

</div>

> [@HelgavonLichtenstein](#):
>
> Here is my attempt which generates `ERROR: LoadError: invalid redefinition of constant f` .

That happens when you define a function twice, the first time having it be a constant. Try your code in a fresh REPL session.

For your code, your issues are you haven’t defined `condition2` (which I assume needs a second tstop and f2 handling), and you need to add the callback to the solving: `sol=solve(prob,Rosenbrock23(),callback=cbs)`.

Note that you can do this without DEDataVectors though. This is probably an example of how you want to work:

```julia
using DifferentialEquations

f = @ode_def_bare ManagedLinear begin
  dBB = r*(1- BB/K)*BB -xI*yIB*BI*((BB^h/(B0^h + BB^h))/eIB)
  dBI = -xI*BI + xI*yIB*BI*(BB^h/(B0^h + BB^h))-xT *yTI*BT*((BI^h/(B0^h + BI^h))/eTI)
  dBT = -xT*BT + xT*yTI*BT*(BI^h/(B0^h + BI^h))-q*BT*E + f1
end K r yIB yTI eIB eTI B0 h xI xT q E f1

function condition(u,t,integrator)
  t == 200.
end
function affect!(integrator)
  for c in full_cache(integrator)
    integrator.p[13] = 0.1
  end
end
cb = DiscreteCallback(condition, affect!)

u0 = [155.0,107.0,93.0]
tspan = (0.0,1000.0)
p = [450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50, 0.0] # initial f1 = 0.0
# K r yIB yTI eIB eTI B0 h xI xT q E
prob = ODEProblem(f,u0,tspan,p)
sol=solve(prob,Rosenbrock23(),callback=cb)

using Plots

plot(sol, title = "Food Chain Managed Linear", xlabel = "Time" ,ylabel = "Density", lw=0.5)

```

---

<div class="post-metadata">

**Author:** ![HelgavonLichtenstein](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/helgavonlichtenstein/32/5005_2.png) [@HelgavonLichtenstein](https://discourse.julialang.org/u/HelgavonLichtenstein)\
**Post date:** [January 9, 2019, 12:13pm UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/5 "2019-01-09T12:13:41Z")

</div>

I have included a `condition2`, ‘affect2’ and grouped it in a ‘CallBackSet’, but when I run it, the plot generated does not differ from the original code without the callback. What am I missing?

```julia
using DifferentialEquations

f = @ode_def_bare ManagedLinear begin
  dBB = r*(1- BB/K)*BB -xI*yIB*BI*((BB^h/(B0^h + BB^h))/eIB)
  dBI = -xI*BI + xI*yIB*BI*(BB^h/(B0^h + BB^h))-xT *yTI*BT*((BI^h/(B0^h + BI^h))/eTI)
  dBT = -xT*BT + xT*yTI*BT*(BI^h/(B0^h + BI^h))-q*BT*E + f1 + f2
end K r yIB yTI eIB eTI B0 h xI xT q E f1 f2

function condition(u,t,integrator)
  t == 200.
end
function affect!(integrator)
  for c in full_cache(integrator)
    integrator.p[13] = 0.01
  end
end
cb = DiscreteCallback(condition, affect!)

function condition2(u,t,integrator)
  t == 200.
end
function affect2!(integrator)
  for c in full_cache(integrator)
    integrator.p[14] = 50
  end
end
cb2 = DiscreteCallback(condition2, affect2!)

cbs = CallbackSet(cb,cb2)

u0 = [500.0,200.0,100.0]
tspan = (0.0,1000.0)
p = [450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50, 0, 0] # initial f1 = 0.0 f2 = 0.0
# K r yIB yTI eIB eTI B0 h xI xT q E
prob = ODEProblem(f,u0,tspan,p)
sol=solve(prob,Rosenbrock23(),callbacks=cbs)

using Plots

plot(sol, title = "Food Chain Managed Linear Callback", xlabel = "Time" ,ylabel = "Density", lw=0.5)

```

---

<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:** [January 9, 2019, 3:52pm UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/6 "2019-01-09T15:52:16Z")

</div>

> [@HelgavonLichtenstein](#):
>
> sol=solve(prob,Rosenbrock23(),callbacks=cbs)

My bad here. I forgot to set the `tstops` so that it hits `t=200.0` exactly to trigger the callbacks.

```julia
sol=solve(prob,Rosenbrock23(),callbacks=cbs,tstops=[200.0])

```

Hope that makes sense.

---

<div class="post-metadata">

**Author:** ![HelgavonLichtenstein](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/helgavonlichtenstein/32/5005_2.png) [@HelgavonLichtenstein](https://discourse.julialang.org/u/HelgavonLichtenstein)\
**Post date:** [January 10, 2019, 1:41pm UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/7 "2019-01-10T13:41:58Z")

</div>

Hello again,

Even with the `tstops` it generates the same graph. Perhaps the problem is with…

```julia
p = [450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50, 0, 0] # initial f1 = 0.0 f2 = 0.0
# K r yIB yTI eIB eTI B0 h xI xT q E f1 f2

```

…since q and E should both be 0 until `t=200.0` and then become 0.01 and 50?

When I set them both to 0 it generates the Pristine graph, which is what I want, but then nothing happens once `t=200.0`.

```julia
using DifferentialEquations

f = @ode_def_bare ManagedLinear begin
  dBB = r*(1- BB/K)*BB -xI*yIB*BI*((BB^h/(B0^h + BB^h))/eIB)
  dBI = -xI*BI + xI*yIB*BI*(BB^h/(B0^h + BB^h))-xT *yTI*BT*((BI^h/(B0^h + BI^h))/eTI)
  dBT = -xT*BT + xT*yTI*BT*(BI^h/(B0^h + BI^h))-q*BT*E + f1 + f2
end K r yIB yTI eIB eTI B0 h xI xT q E f1 f2

function condition(u,t,integrator)
  t == 200.
end
function affect!(integrator)
  for c in full_cache(integrator)
    integrator.p[13] = 0.01
  end
end
cb = DiscreteCallback(condition, affect!)

function condition2(u,t,integrator)
  t == 200.
end
function affect2!(integrator)
  for c in full_cache(integrator)
    integrator.p[14] = 50
  end
end
cb2 = DiscreteCallback(condition2, affect2!)

cbs = CallbackSet(cb,cb2)

u0 = [500.0,200.0,100.0]
tspan = (0.0,1000.0)
p = [450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50, 0, 0] # initial f1 = 0.0 f2 = 0.0
# K r yIB yTI eIB eTI B0 h xI xT q E
prob = ODEProblem(f,u0,tspan,p)
sol=solve(prob,Rosenbrock23(),callbacks=cbs,tstops=[200.0])

using Plots

plot(sol, title = "Food Chain Managed Linear ", xlabel = "Time" ,ylabel = "Density", lw=0.5)

```

---

<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:** [January 10, 2019, 7:24pm UTC](https://discourse.julialang.org/t/is-there-a-way-to-run-differential-equations-to-a-specific-point-and-then-change-input-values-continue-running-and-graph-it-together/19232/8 "2019-01-10T19:24:46Z")

</div>

The issue is that you used the argument `callbacks` instead of `callback`. That was my bad (and we really need to error for that). A corrected code is:

```julia
using DifferentialEquations

f = @ode_def_bare ManagedLinear begin
  dBB = r*(1- BB/K)*BB -xI*yIB*BI*((BB^h/(B0^h + BB^h))/eIB)
  dBI = -xI*BI + xI*yIB*BI*(BB^h/(B0^h + BB^h))-xT *yTI*BT*((BI^h/(B0^h + BI^h))/eTI)
  dBT = -xT*BT + xT*yTI*BT*(BI^h/(B0^h + BI^h))-q*BT*E + f1 + f2
end K r yIB yTI eIB eTI B0 h xI xT q E f1 f2

function condition(u,t,integrator)
  t == 200.
end
function affect!(integrator)
  integrator.p[13] = 0.01
end
cb = DiscreteCallback(condition, affect!)

function condition2(u,t,integrator)
  t == 200.
end
function affect2!(integrator)
  integrator.p[14] = 50
end
cb2 = DiscreteCallback(condition2, affect2!)

cbs = CallbackSet(cb,cb2)

u0 = [500.0,200.0,100.0]
tspan = (0.0,1000.0)
p = [450, 1.1, 10, 10, 0.66, 0.85, 80, 1.2, 0.15, 0.06, 0.01, 50, 0, 0] # initial f1 = 0.0 f2 = 0.0
# K r yIB yTI eIB eTI B0 h xI xT q E
prob = ODEProblem(f,u0,tspan,p)
sol=solve(prob,Rosenbrock23(),callback=cbs,tstops=[200.0],reltol=1e-6)

using Plots

plot(sol, title = "Food Chain Managed Linear ", xlabel = "Time" ,ylabel = "Density", lw=0.5)

```

![plot](https://global.discourse-cdn.com/julialang/original/3X/b/6/b6c7402043fb674ccc2954b12fd877d41665cda9.png)