# Solving ODE based on experimental data (DifferentialEquations.jl)

**URL:** <https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002>\
**Category:** Numerics\
**Tags:** package\
**Created:** [November 30, 2020, 5:11pm UTC](https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002 "2020-11-30T17:11:36Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![maajdl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maajdl/32/22838_2.png) [@maajdl](https://discourse.julialang.org/u/maajdl)\
**Post date:** [November 30, 2020, 5:11pm UTC](https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002/1 "2020-11-30T17:11:36Z")

</div>

Hello,

I am solving a simple ODE system with 5 state variables.  
It is essentially made of 5 coupled RC sytems with a few (current) excitations.  
These excitations are 6000 experimental data measured every 10 minutes.  
I am using DifferentialEquations.jl and wanted to test the various solvers.  
The experimental data are provided to the model by LinearInterpolation.

To my big surprise, almost all solvers failed.  
For example:

- for Euler I got this warning: “Warning: Instability detected. Aborting”

- for Tsit5: “Warning: Interrupted. Larger maxiters is needed.”  
I checked that the model was accessed 6000001 times !  
Even stranger, the solver only explored the 80 first time steps of the data.

- the same for RK4, …

However, TRBDF2 solves the problem very easily and accessed the model only 134187 times, and of course accessed all the data.

Any suggestion?

Thanks,

Michel

---

<div class="post-metadata">

**Author:** ![jlchan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jlchan/32/10958_2.png) [@jlchan](https://discourse.julialang.org/u/jlchan)\
**Post date:** [November 30, 2020, 5:31pm UTC](https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002/2 "2020-11-30T17:31:41Z")

</div>

Seems like all the solvers that are failing are based on explicit time-stepping methods. Is the system of ODEs fairly stiff?

---

<div class="post-metadata">

**Author:** ![maajdl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maajdl/32/22838_2.png) [@maajdl](https://discourse.julialang.org/u/maajdl)\
**Post date:** [November 30, 2020, 8:13pm UTC](https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002/3 "2020-11-30T20:13:34Z")

</div>

These equations should not be stiff.  
It is a very simple model of heat conduction in a building.  
Below are the equations.  
The heat sources (i1, …) are measurements averaged over the time step (10 min).  
These sources are a bit noisy, but that should not matter a lot.  
The time constants are all rather large, except maybe for one of them.  
Thanks!

```julia
function building(dT,T,p,t)
    T1, T2, T3, T4, T5 = T
    τ1, τ2, τ3, τ4, τ5, U1, U2, U3, U4, U5, k1, AgS = abs.(p)
    T0 = te(t)
    i2 = i21(t) + AgS*i22(t)
    U1k = 1/(1/U1 + k1*sh(t))
    dT[1] = (i1(t) + U1k*(T0 - T1) + U2*(T2 - T1) ) /(τ1*U1)
    dT[2] = (i2 + U2*(T1 - T2) + U3*(T3 - T2) + U4*(T4 - T2) + U5*(T5 - T2)) /(τ2*U2)
    dT[3] = (i3(t) + U3*(T2 - T3) ) /(τ3*U3)
    dT[4] = (i4(t) + U4*(T2 - T4) ) /(τ4*U4)
    dT[5] = (i5(t) + U5*(T2 - T5) ) /(τ5*U5)
end

```

---

<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:** [November 30, 2020, 9:05pm UTC](https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002/4 "2020-11-30T21:05:40Z")

</div>

> [@maajdl](#):
>
> These sources are a bit noisy, but that should not matter a lot.

That would likely lead to a form of stiffness in the definition. You might want to use DataInterpolations to generate a regression spline instead of doing a direct interpolation if that’s the case.

> [@maajdl](#):
>
> It is a very simple model of heat conduction in a building.

Those models are usually stiff.

> [@maajdl](#):
>
> - for Euler I got this warning: “Warning: Instability detected. Aborting”
> - for Tsit5: “Warning: Interrupted. Larger maxiters is needed.”  
> I checked that the model was accessed 6000001 times !  
> Even stranger, the solver only explored the 80 first time steps of the data.

> [@maajdl](#):
>
> However, TRBDF2 solves the problem very easily and accessed the model only 134187 times

All of this is pointing to the equation either being stiff, or, the other issue is that

> [@maajdl](#):
>
> The experimental data are provided to the model by LinearInterpolation.

So if you have a lot of data points, then your `f` is discontinuous at every data point, in which case the error estimators will have a very difficult time handling your problem. You probably want to use a higher order interpolation to smooth it out.

---

<div class="post-metadata">

**Author:** ![maajdl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/maajdl/32/22838_2.png) [@maajdl](https://discourse.julialang.org/u/maajdl)\
**Post date:** [November 30, 2020, 9:09pm UTC](https://discourse.julialang.org/t/solving-ode-based-on-experimental-data-differentialequations-jl/51002/5 "2020-11-30T21:09:51Z")

</div>

Just observed that `ImplicitEuler` worked very well.  
I overlooked the fact that there is a small time constant in this problem.

Thanks ChrisRackauckas  
Thanks jlchan

Further suggestions welcome …
