DifferentialEquations: Confusion regarding time points and values in the solution array

Specific Domains Biology, Health, and Medicine
1

I am a bit confused about the timepoints that the solution array represents.

Since the algorithms used to solve are adaptive, the size of solutions might be different. Now, if I say that all my growth rates of species have the units hour^(-1) and all other constants have appropriate units and specify the time span to be 0.0 to 100.0 and I wish to track the population size in each hour. Given the model equations, if the length of solution comes out be 1000, then:

1) What time points do these 1000 values represent?    
2) How to get the values of population size in each hour i.e., get the 101 values which would represent population size at each hour? Can this array be obtained by using `solve(problem,saveat = 1.0)`?

Thank you.

2 
2) 
ys = [sol(t) for t in 0.0:100]

You get some interpolated values which should be with the specified tolerances.

3

See here: https://docs.sciml.ai/stable/basics/solution/#Interpolations-1

4

Yes, if you want outputs saved by 1, then saveat=1.0

1 Like
5

So from I understand, if I set up the time span to be 0.0 to 10.0 and my units are still hour^(-1) and the algorithm somhow calculates stepsize as 0.5 (assume, I dont know if its possible or not), I must get the solution vector to be of length 10.0/0.5 i.e., 20, right?

And,would each of the 20 values would represent values for time points 0 hr, 0.5 hr, 1.0 hr, 1.5 hr…,19.5 hr,20 hr?

6

If you don’t use saveat, the solver will typically pick how many points it needs to represent the solution well enough that the interpolation will be accurate in between. I’ve almost exclusively used Tsit5() and this seems to work really well. As a result, I don’t worry about what points the solver actually finds and use the interpolation interface to find values at times of interest. The one time I got away from that was when doing a large and long simulation in which I was running low on RAM, but that should be the exception.

7

The time points the solver will use will probably not be evenly spaced. They will depend on how fast things are changing in your simulation. I would encourage you to try using the interpolation - it is a good feature.

8

I am quiet new to solving ODEs, rather this is my first time.
Can you shed some light upon the interpolating function as I did not find any good material on it online?
Is it anything like the integrated function of the ODE/s?

9

https://docs.sciml.ai/latest/basics/solution/

Just add a line like Intu' = u

10

Here’s a very simple example:

using DifferentialEquations

function ode!(du,u,p,t)
    du[1] =1- u[1]
end

problem=ODEProblem(ode!,[0.0],[0.0,5.0])
sol=solve(problem,Tsit5())

If you examine the sol object after running that code, you see this:

julia> sol
retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 15-element Array{Float64,1}:
 0.0
 9.999999999999999e-5
 0.0010999999999999998
 0.011099999999999997
 0.07674209860034185
 0.2256220209761885
 0.4455760494689467
 0.7272505875899613
 1.08996450301086
 1.5331654123549805
 2.0697241857623103
 2.705750702626476
 3.4562447902280993
 4.337840756338852
 5.0
u: 15-element Array{Array{Float64,1},1}:
 [0.0]
 [9.999500016666247e-5]
 [0.001099395221772342]
 [0.011038622307372232]
 [0.07387132730531631]
 [0.20198031990033405]
 [0.35954475580746376]
 [0.5167641828239731]
 [0.6637714219615871]
 [0.7841482550875905]
 [0.8737784281094817]
 [0.9331779992336123]
 [0.9684489722571047]
 [0.9869310694285471]
 [0.9932596937930706]

It only stored 15 points over the 5 seconds of simulation. However, you can get the value at any time by executing something like sol(1.0) (note I’m using parentheses, not brackets) and it will interpolate and give an accurate answer.
Or you can do this:

julia> sol(0.0:0.5:5.0)
t: 0.0:0.5:5.0
u: 11-element Array{Array{Float64,1},1}:
 [0.0]
 [0.39346948341153787]
 [0.6321206591004334]
 [0.7768694671005957]
 [0.8646641097009187]
 [0.9179166359464174]
 [0.9502181305661904]
 [0.9698001352429271]
 [0.9816881279747315]
 [0.9888878525731617]
 [0.9932596937930706]

to get interpolation of results on a regular time basis.

Docs are here: https://docs.sciml.ai/stable/basics/solution/#Interpolations-1

1 Like
11

Thanks a ton!
These makes sense to me now.

1 Like