# Using DifferentialEquations.jl for a simple food-chain

**URL:** <https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576>\
**Category:** New to Julia\
**Tags:** question, diffeq\
**Created:** [August 16, 2018, 3:42pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576 "2018-08-16T15:42:29Z")\
**Posts on this page:** 13\
**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:** [August 16, 2018, 3:42pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/1 "2018-08-16T15:42:29Z")

</div>

Hello, I am a complete Julia beginner and am working through [Example 2: Solving Systems of Equations](http://docs.juliadiffeq.org/latest/tutorials/ode_example.html#Example-2:-Solving-Systems-of-Equations-1) and [Defining Parameterized Functions](http://docs.juliadiffeq.org/latest/tutorials/ode_example.html#Defining-Parameterized-Functions-1), trying to adapt my food-chain formulas to Julia. I have also watched [Intro to solving differential equations in Julia](https://www.youtube.com/watch?v=KPEqYtEd-zY).

I am using Julia 1.0 on Windows 10, and have updated all packages and have DifferentialEquations (also updated) installed.

The start equations and parameters are:

```
   du[1] = u1*(1-u1)-(1/e)*y2*x2*u1^(1-q)/(u012^(1+q)+u1^(q+1))*u2
   du[2] = y2*x2*u1^(1+q)/(u012^(1+q)+u1^(1+q))*u2-(1/e)*y3*x3*u2^(1+q)/(u023^(1+q)+u2^(1+q))*u3-x2*u2
   du[3] = y3*x3*u2^(1+q)/(u023^(1+q)+u2^(1+q))*u3-x3*u3

```

Initial values = u1= 0.01, u2 = 0.01, u3 = 0.01  
e = 1  
y2 = 2  
x2 = 3  
q = 4  
u012 = 5  
y3 = 6  
x3 = 7  
u023 = 8  
t-max = 10000  
t-interval = 0.05

* * *

Based on the tutorial my code currently looks like this:

using DifferentialEquations  
function Bio(du,u,p,t)  
du[1] = u[1]_(1-u[1])-(1/[1])_[2]\*[3]\*u[1]^(1-[4])/([5]^(1+[4])u[1]^([4]+1))_u[2]  
du[2] = [2]_[3]_u[1]^(1+[4])/([5]^(1+[4])+u[1]^(1+[4]))u[2]-(1/[1])[6]_[7]\*u[2]^(1+[4])/([8]^(1+[4])+u[2]^(1+[4]))\*u[3]-[3]_u[2]  
du[3] = [6]_[7]\*u[2]^(1+[4])/([8]^(1+[4])+u[2]^(1+[4]))\*u[3]-[7]\*u[3]  
end  
u0 = [0.01,0.01,0.01]  
p=(1,2.009,0.4,5,0.08,0,0.16129,0.5)  
tspan = (0.01, 10000)  
prob=ODEProblem(Bio,u0,tspan,p)  
sol=solve(prob)

* * *

I generate an error once _sol=solve(prob)_ runs.

ERROR: MethodError: no method matching -(::Int32, ::Array{Int32,1})  
Closest candidates are:  
-(::Union{Int128, Int16, Int32, Int64, Int8, UInt128, UInt16, UInt32, UInt64, UInt8}) at int.jl:51  
-(::T\<:Union{Int128, Int16, Int32, Int64, Int8, UInt128, UInt16, UInt32, UInt64, UInt8}, ::T\<:Union{Int128, Int16, Int32, Int64, Int8, UInt128, UInt16, UInt32, UInt64, UInt8}) where T\<:Union{Int128, Int16, Int32, Int64, Int8, UInt128, UInt16, UInt32, UInt64, UInt8} at int.jl:52  
-(::Union{Int16, Int32, Int8}, ::BigInt) at gmp.jl:458  
…  
Stacktrace:  
[1] Bio(::Array{Float64,1}, ::Array{Float64,1}, ::Tuple{Int32,Float64,Float64,Int32,Float64,Int32,Float64,Float64}, ::Float64) at .\REPL[12]:2  
[2] (::ODEFunction{true,typeof(Bio),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing})(::Array{Float64,1}, ::Array{Float64,1}, ::Vararg{Any,N} where N) at C:\Users\Helga.julia\packages\DiffEqBase\mcoTV\src\diffeqfunction.jl:106  
[3] initialize!..

… followed by 5 pages of text.

What am I doing wrong? And how may I incorporate time intervals of 0.05? Any assistant/advice would be greatly appreciated.

---

<div class="post-metadata">

**Author:** ![ElOceanografo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/eloceanografo/32/624_2.png) [@ElOceanografo](https://discourse.julialang.org/u/ElOceanografo)\
**Post date:** [August 16, 2018, 4:09pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/2 "2018-08-16T16:09:29Z")

</div>

In Julia 1.0, you need to use `.-` (i.e. broadcasting `-` by adding a dot) to subtract a scalar and an array. Julia 0.6 would let you do this with `-`, 0.7 added a deprecation warning, and 1.0 throws an error. Looks like the DifferentialEquations example problems haven’t been updated yet. For a deeper dive on Julia’s dot-broadcasting syntax, check out this [blog post](https://julialang.org/blog/2017/01/moredots).

Also, I just noticed your code has a bunch of numbers in brackets–not sure what your intention is, but each of these is going to evaluate to a single-element vector, which is probably not what you want?

```julia
 du[2] = [2][3]u[1]^(1+[4])/([5]^(1+[4])+u[1]^(1+[4]))u[2]-(1/[1])[6][7]*u[2]^(1+[4])/([8]^(1+[4])+u[2]^(1+[4]))*u[3]-[3]u[2]

```

(Also, make sure you use triple backticks (```) to quote your code blocks so they format correctly.)

---

<div class="post-metadata">

**Author:** ![platawiec](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/platawiec/32/31914_2.png) [@platawiec](https://discourse.julialang.org/u/platawiec)\
**Post date:** [August 16, 2018, 4:10pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/3 "2018-08-16T16:10:38Z")

</div>

This error:

```julia
ERROR: MethodError: no method matching -(::Int32, ::Array{Int32,1})

```

means that you should be “broadcasting” one of operations, or there is a typo where you expect a scalar but are currently getting a vector. It seems to me most likely that you have a typo in your code, if you can re-quote it in backticks it may be more clear. Based on how your code currently reads, are the `[1]` and so on supposed to be `p[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:** [August 17, 2018, 1:49pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/4 "2018-08-17T13:49:36Z")

</div>

@ElOceanografo @platawiec

I thought adjustable variables need to be replaced with `[1]` and to be specified in p, so `p=(1,2.009...)` would correspond with `[1][2]...` ect.

```julia
e =1 e=[1]
y2 = 2.009 y2=[2]
x2 = 0.4 x2=[3]
q = 0 q=[4]
u012 = 0.16129 u012=[5]
y3 =5 y3=[6]
x3 = 0.08 x3=[7]
u023 = 0.5 u023=[8]

```

My code block is (no change from previous)

```julia
using DifferentialEquations
function Bio(du,u,p,t)
du[1] = u[1]*(1-u[1])-(1/[1])*[2]*[3]*u[1]^(1-[4])/([5]^(1+[4])u[1]^([4]+1))*u[2]
du[2] = [2]*[3]*u[1]^(1+[4])/([5]^(1+[4])+u[1]^(1+[4]))*u[2]-(1/[1])*[6]*[7]*u[2]^(1+[4])/([8]^(1+[4])+u[2]^(1+[4]))*u[3]-[3]*u[2]
du[3] = [6]*[7]*u[2]^(1+[4])/([8]^(1+[4])+u[2]^(1+[4]))*u[3]-[7]*u[3]
end
u0 = [0.01,0.01,0.01]
p=(1,2.009,0.4,5,0.08,0,0.16129,0.5)
tspan = (0.01, 10000)
prob=ODEProblem(Bio,u0,tspan,p)
sol=solve(prob)

```

@ElOceanografo So I should replace each `-` with `.-` whenever I am subtracting a substituted number `[2]` such as

`du[1] = u[1]*(1-u[1])-(1/[1])*[2]*[3]*u[1]^(1.-[4])/([5]^(1+[4])u[1]^([4]+1))*u[2]` ?

@platawiec So I should place `p[#]` in front of each number that is an adjustable variable like

```julia
du[1] = u[1]*(1-u[1])-(1/p[1])*p[2]*p[3]*u[1]^(1-p[4])/(p[5]^(1+p[4])u[1]^(p[4]+1))*u[2]

```

I am sorry if these are very basic questions, but I am still a very basic person. Thank you for your patience.

---

<div class="post-metadata">

**Author:** ![platawiec](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/platawiec/32/31914_2.png) [@platawiec](https://discourse.julialang.org/u/platawiec)\
**Post date:** [August 17, 2018, 2:35pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/5 "2018-08-17T14:35:16Z")

</div>

> [@HelgavonLichtenstein](#):
>
> @platawiec So I should place `p[#]` in front of each number that is an adjustable variable like
> 
> ```julia
> du[1] = u[1]*(1-u[1])-(1/p[1])*p[2]*p[3]*u[1]^(1-p[4])/(p[5]^(1+p[4])u[1]^(p[4]+1))*u[2]
> 
> ```

Yes, this is correct. Have you tried running the code with this change?

Just like you write `u[1]` to get the first entry of `u`, you write `p[1]` to get the first entry of `p` (in this case `p[1] = 1`, etc.)

---

<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:** [August 21, 2018, 12:41pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/6 "2018-08-21T12:41:08Z")

</div>

I have put the p’s in and it now generates:

```julia
julia> function Bio(du,u,p,t)
       du[1] = u[1]*(1-u[1])-(1/p[1])*p[2]*p[3]*u[1]^(1-p[4])/p([5]^(1+p[4])u[1]^(p[4]+1))*u[2]
       du[2] = p[2]*p[3]*u[1]^(1+p[4])/(p[5]^(1+p[4])+u[1]^(1+p[4]))*u[2]-(1/p[1])*p[6]*p[7]*u[2]^(1+p[4])/(p[8]^(1+p[4])+u[2]^(1+p[4]))*u[3]-p[3]*u[2]
       du[3] = p[6]*p[7]*u[2]^(1+p[4])/(p[8]^(1+p[4])+u[2]^(1+p[4]))*u[3]-p[7]*u[3]
       end
Bio (generic function with 1 method)

julia> u0 = [0.01,0.01,0.01]
3-element Array{Float64,1}:
 0.01
 0.01
 0.01

julia> p=(1,2.009,0.4,5,0.08,0,0.16129,0.5)
(1, 2.009, 0.4, 5, 0.08, 0, 0.16129, 0.5)

julia> tspan = (0.01, 10000)
(0.01, 10000)

julia> prob=ODEProblem(Bio,u0,tspan,p)
ODEProblem with uType Array{Float64,1} and tType Float64. In-place: true
timespan: (0.01, 10000.0)
u0: [0.01, 0.01, 0.01]

julia> sol=solve(prob)
ERROR: MethodError: no method matching ^(::Array{Int32,1}, ::Float64)
Closest candidates are:
  ^(::Missing, ::Number) at missing.jl:93
  ^(::Float64, ::Float64) at math.jl:767
  ^(::Irrational{:ℯ}, ::Number) at mathconstants.jl:91
  ...

```

What is the next step?

---

<div class="post-metadata">

**Author:** ![rsrock](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rsrock/32/115_2.png) [@rsrock](https://discourse.julialang.org/u/rsrock)\
**Post date:** [August 21, 2018, 1:15pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/7 "2018-08-21T13:15:42Z")

</div>

You have one typo in there. The error message is giving you a hint where to look for it, so read it carefully.

---

<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 21, 2018, 1:22pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/8 "2018-08-21T13:22:51Z")

</div>

I am joining in this late so maybe it does more harm than good, but this probably would’ve been a good use-case for the `@ode_def` macro.

```julia
f = @ode_def FoodChain begin
   du1 = u1*(1-u1)-(1/e)*y2*x2*u1^(1-q)/(u012^(1+q)+u1^(q+1))*u2
   du2 = y2*x2*u1^(1+q)/(u012^(1+q)+u1^(1+q))*u2-(1/e)*y3*x3*u2^(1+q)/(u023^(1+q)+u2^(1+q))*u3-x2*u2
   du3 = y3*x3*u2^(1+q)/(u023^(1+q)+u2^(1+q))*u3-x3*u3
end e y2 x2 q u012 y3 x3 u023
u0 = [0.01,0.01,0.01]
tspan = (0.0,10000.0)
p = (1,2,3,4,5,6,7,8)
prob = ODEProblem(f,u0,tspan,p)
solve(prob,saveat=0.05)

```

---

<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:** [August 21, 2018, 2:01pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/9 "2018-08-21T14:01:08Z")

</div>

This might be a silly question, but does the error message mean that there is an error somewhere in the formula where there is a `^(`?

---

<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 21, 2018, 2:02pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/10 "2018-08-21T14:02:19Z")

</div>

> [@HelgavonLichtenstein](#):
>
> This might be a silly question, but does the error message mean that there is an error somewhere in the formula where there is a `^(` ?

Yes, somewhere you took the exponential of a vector of ints, likely `p^`. If you shared the error message it would be easier to find out though.

---

<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:** [August 21, 2018, 2:28pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/11 "2018-08-21T14:28:13Z")

</div>

Oh my.

I found the error of

```julia
p([5]^(1

```

instead of `(p[5]^(1`

@ChrisRackauckas the error message is

```julia
ERROR: MethodError: no method matching ^(::Array{Int32,1}, ::Float64)
Closest candidates are:
  ^(::Missing, ::Number) at missing.jl:93
  ^(::Float64, ::Float64) at math.jl:767
  ^(::Irrational{:ℯ}, ::Number) at mathconstants.jl:91
  ...

```

Unless you mean the accompanying 5 pages of text?

`retcode: Success` (Which I must say I really enjoy this about Julia, that and `Everything is okay` when installing)

I will try using `@ode_def` from now on to make things simpler.

Thank you everyone for your help!

---

<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:** [August 21, 2018, 2:39pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/12 "2018-08-21T14:39:12Z")

</div>

@ChrisRackauckas Sorry to bother you again, but when I run the `@ode_def` macro I get:

```julia
julia> f = @ode_def FoodChain begin
          du1 = u1*(1-u1)-(1/e)*y2*x2*u1^(1+q)/(u012^(1+q)+u1^(q+1))*u2
             du2 = y2*x2*u1^(1+q)/(u012^(1+q)+u1^(1+q))*u2-(1/e)*y3*x3*u2^(1+q)/(u023^(1+q)+u2^(1+q))*u3-x2*u2
                du3 = y3*x3*u2^(1+q)/(u023^(1+q)+u2^(1+q))*u3-x3*u3
                end e y2 x2 q u012 y3 x3 u023
┌ Warning: Failed to build the parameter derivatives.
└ @ ParameterizedFunctions C:\Users\Helga\.julia\packages\ParameterizedFunctions\47E7u\src\ode_def_opts.jl:226
(::FoodChain{getfield(Main, Symbol("##3#8")),getfield(Main, Symbol("##4#9")),getfield(Main, Symbol("##5#10")),getfield(Main, Symbol("##6#11")),getfield(Main, Symbol("##7#12")),Nothing,Expr,Expr}) (generic function with 2 methods)

```

What is the reason?

---

<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 21, 2018, 2:59pm UTC](https://discourse.julialang.org/t/using-differentialequations-jl-for-a-simple-food-chain/13576/13 "2018-08-21T14:59:27Z")

</div>

It just means that the Jacobian of the parameters doesn’t seem to be buildable, not sure why it isn’t, but that shouldn’t hurt your use case? You can use `@ode_def_bare` to turn off all extra computations.