# How to add extra variables when defining ODEProblem in DifferentialEquations?

**URL:** <https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126>\
**Category:** New to Julia\
**Tags:** ode, differentialequation\
**Created:** [August 1, 2022, 10:26pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126 "2022-08-01T22:26:43Z")\
**Posts on this page:** 9\
**Page:** 1

<div class="post-metadata">

**Author:** ![laurar1891](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/laurar1891/32/38443_2.png) [@laurar1891](https://discourse.julialang.org/u/laurar1891)\
**Post date:** [August 1, 2022, 10:26pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/1 "2022-08-01T22:26:43Z")

</div>

I’m trying to define a system of differential equations using matrix notation. The issue is that this matrix varies depending on the problem and in all the examples that I’ve seen the function _f_ that goes in the ODEProblem() function, only has the three variables _f(u,p,t)_.

In the following example, I want to define a function for which I can change the value that defines the size of the identity matrix (_A_):

```julia
using DifferentialEquations  
using LinearAlgebrape or paste code here

function f(u,p,t)
  A = I(4)
  p[1].*(A*u).*(1/t^(1+p[2]))
end

u0 = [0.1,0.3,0.5,0.8]
p = [5.5, 0.5];

```

As you can deduct, I’m using the matrix A to define 4 different state variables. Each state variable will be defined by the same parameters, but the initial state of each variable is different.

```julia
tspan = (3.0,30.0)
prob = ODEProblem(f,u0,tspan,p)

Age1 = collect(5:5:30)
sol = solve(prob, saveat=Age1); #you can plot this for clarity

```

Ideally, I want something like this:

```julia
n = 4
function f(u,p,t,n)
  A = I(n)
  p[1].*(A*u).*(1/t^(1+p[2]))
end

```

I know this doesn’t work, but any ideas of how I can modify the code so it admits different sizes for the matrix A?

Thanks!

---

<div class="post-metadata">

**Author:** ![ANasc](https://avatars.discourse-cdn.com/v4/letter/a/41988e/32.png) [@ANasc](https://discourse.julialang.org/u/ANasc)\
**Post date:** [August 1, 2022, 11:33pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/2 "2022-08-01T23:33:21Z")

</div>

Could you just make `p[3] = n` ?

---

<div class="post-metadata">

**Author:** ![mkitti](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkitti/32/12459_2.png) [@mkitti](https://discourse.julialang.org/u/mkitti)\
**Post date:** [August 2, 2022, 12:06am UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/3 "2022-08-02T00:06:06Z")

</div>

You could create a function that returns a function:

```julia
julia> function generate_f(n)
           A = I(n)
           f(u, p, t) = p[1].*(A*u).*(1/t^(1+p[2]))
           return f
       end
generate_f (generic function with 1 method)

julia> f4 = generate_f(4)
(::var"#f#8"{Diagonal{Bool, Vector{Bool}}}) (generic function with 1 method)

julia> f4(inits[1:4], inits[(nplot+1):(nplot+2)], 5)
4-element Vector{Float64}:
 2.726706751006412
 2.7139936889300844
 2.742306107133945
 2.6816291968428603

```

Another way is to create a functor.

```julia
julia> struct MyFunctor
           A::Diagonal{Bool, Vector{Bool}}
           MyFunctor(n::Int) = new(I(n))
       end

julia> (mf::MyFunctor)(u, p, t) = p[1].*(mf.A*u).*(1/t^(1+p[2]))

julia> f4(inits[1:4], inits[(nplot+1):(nplot+2)], 5)
4-element Vector{Float64}:
 2.7531945173074637
 2.74035795804277
 2.768945408652123
 2.7076790709064475

```

A third way would be to use a constant global `Ref`:

```julia
julia> const A_ref = Ref(I(4))
Base.RefValue{Diagonal{Bool, Vector{Bool}}}(Bool[1 0 0 0; 0 1 0 0; 0 0 1 0; 0 0 0 1])

julia> function f(u, p, t)
           p[1].*(A_ref[]*u).*(1/t^(1+p[2]))
       end
f (generic function with 1 method)

julia> f(inits[1:4], inits[(nplot+1):(nplot+2)], 5)
4-element Vector{Float64}:
 2.7531945173074637
 2.74035795804277
 2.768945408652123
 2.7076790709064475

```

---

<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:** [August 2, 2022, 10:53am UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/4 "2022-08-02T10:53:46Z")

</div>

It’s the identity matrix, so why not replace `A*u` with `u` and save the allocation?

---

<div class="post-metadata">

**Author:** ![laurar1891](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/laurar1891/32/38443_2.png) [@laurar1891](https://discourse.julialang.org/u/laurar1891)\
**Post date:** [August 2, 2022, 12:57pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/5 "2022-08-02T12:57:12Z")

</div>

You’re right. I thought I needed the matrix so it will recognize that it was a system of equations. But just by defining u as a vector, it recognizes it as a system.  
So, this is just the same, without needing to define the matrix:

```julia
function f(u,p,t)
  p[1].*(u).*(1/t^(1+p[2]))
end

```

Tank you 🙂

---

<div class="post-metadata">

**Author:** ![laurar1891](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/laurar1891/32/38443_2.png) [@laurar1891](https://discourse.julialang.org/u/laurar1891)\
**Post date:** [November 1, 2022, 9:08pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/6 "2022-11-01T21:08:46Z")

</div>

I found this question to be relevant again. Now for an external variable that changes with time for each one of the variables in the system. This is an example:

Let’s say I have this, which works fine:

```julia
using DifferentialEquations  
using LinearAlgebra

function f(u,p,t)
    (p[1]+p[2]).*(u).*(1/t^(1+p[3]))
  end

u0 = [0.1,0.3,0.5,0.8]
p = [1,3, 0.5];

tspan = (3.0,30.0)
Age = collect(5:5:30)
ndata = length(Age)*length(u0)

prob = ODEProblem(f,u0,tspan,p)
sol = solve(prob, saveat=Age)
plot(sol)

```

And I want to add the influence of another variable which varies for each one of the equations of the system and for each time for which we are solving the problem, something like this:

```julia
ext_var = rand(Float64, (length(u0),length(Age))) #matrix 

function g(u,p,t)
    (p[1]+p[2].*ext_var).*(u).*(1/t^(1+p[3])) #added extra var
  end

probg = ODEProblem(g,u0,tspan,p)
solg = solve(probg, saveat=Age)
plot(solg)

```

How should I define this problem. What I thought could work didn’t work.

---

<div class="post-metadata">

**Author:** ![stillyslalom](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stillyslalom/32/45687_2.png) [@stillyslalom](https://discourse.julialang.org/u/stillyslalom)\
**Post date:** [November 1, 2022, 9:43pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/7 "2022-11-01T21:43:32Z")

</div>

All external parameters and data can be squeezed into `p`. The parameter argument `p` doesn’t have to be a vector - try using a tuple of parameters (including `ext_var`) instead. Since your `ext_var` parameter varies in time, you’ll need to use a callback to modify it. Following [this example](https://diffeq.sciml.ai/stable/features/callback_functions/#Example-1:-Interventions-at-Preset-Times), try something like this:

```julia
p = (1, 3, 0.5,
    collect(eachcol(ext_var)), # time-varying parameter
    Ref(1)) # Ref container for ext_var column index

function g(u,p,t)
    ext_vec = p[4][p[5][]]
    @. (p[1]+p[2]*ext_vec)*(u)*(1/t^(1+p[3]))
end

# at the start of each age, increment the column index
affect!(integrator) = integrator.p[5][] += 1
cb = PresetTimeCallback(Age,affect!, save_positions=(false, false))

probg = ODEProblem(g,u0,tspan,p)
solg = solve(probg, callback=cb, saveat=Age)

```

I’m using the messy `Ref[]` indexing because tuples are immutable, but they can contain mutable elements, and `Ref` is a good single-element mutable container.

---

<div class="post-metadata">

**Author:** ![bgroenks](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bgroenks/32/21784_2.png) [@bgroenks](https://discourse.julialang.org/u/bgroenks)\
**Post date:** [November 3, 2022, 10:10am UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/8 "2022-11-03T10:10:39Z")

</div>

As already suggested, for numerical parameters you can/should use `p`. While non-array types are technically possible, this will cause problems with basically any and all parameter estimation and sensitivity analysis tools, so I would suggest [ComponentArrays](https://jonniedie.github.io/ComponentArrays.jl/stable/).

If you find yourself wanting to include other kinds of configuration options or static information, use a callable struct:

```julia
struct MyODE{T}
    opts::T
end

function (f::MyODE)(u,p,t)
    # do something with f.opts (or use it for dispatch)
    # ...
end

```

Edit: and for more complex use-cases with different model configurations, you might be interested in [ModelParameters](https://github.com/rafaqz/ModelParameters.jl)

---

<div class="post-metadata">

**Author:** ![bjarthur](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bjarthur/32/9638_2.png) [@bjarthur](https://discourse.julialang.org/u/bjarthur)\
**Post date:** [April 19, 2023, 5:10pm UTC](https://discourse.julialang.org/t/how-to-add-extra-variables-when-defining-odeproblem-in-differentialequations/85126/9 "2023-04-19T17:10:46Z")

</div>

an “external variable that changes with time” sounds like a forcing function, for which there is a [standard approach](https://docs.sciml.ai/ModelingToolkit/stable/tutorials/ode_modeling/#Specifying-a-time-variable-forcing-function).