# Comparison of integrators DifferentialEquations.jl
**URL:** https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271
**Category:** Modelling & Simulations
**Created:** [August 24, 2022, 4:27pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271 "2022-08-24T16:27:13Z")
**Posts on this page:** 16
**Page:** 1
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 24, 2022, 4:27pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/1 "2022-08-24T16:27:13Z")
I get two different trajectoriesf or identical parametres, time, and initial condition.
for comparison i using integrator from Julia and python. I have used several integrators from Julia and python and get completely different results. Even when I imported an integrator with scipy to Julia, the result did not change.
Code from python
```julia
Gamma = lambda x_i: 1 / (1 + np.exp(-lambd * (x_i - Theta)))
def two_HRs(t, X, k1, k2):
x1, y1, z1, x2, y2, z2 = X
dx1 = y1 - a*x1 **3 + b*x1** 2 -z1 + I - k1*(x1 - v_s) * Gamma(x2)
dy1 = c - d*x1**2 - y1
dz1 = r*(s*(x1 - x_r) - z1)
dx2 = y2 - a*x2 **3 + b*x2** 2 -z2 + I - k2*(x2 - v_s) * Gamma(x1)
dy2 = c - d*x2**2 - y2
dz2 = r*(s*(x2 - x_r) - z2)
return [dx1, dy1, dz1, dx2, dy2, dz2]
a = 1
b = 3
c = 1
d = 5
x_r = -1.6
r = 0.01
s = 5
I = 4
v_s = 2
lambd = 10
Theta = -0.25
k1 = -0.17
k2 = -0.17
initials = [0.1, 0.1, 0.1, 0.1, 0.1, 0.1]
sol = solve_ivp(two_HRs, [0, 100000], initials, rtol = 1e-11, atol = 1e-11, args = (-0.17, -0.17), method="LSODA")
ts = sol.t
x1s, y1s, z1s, x2s, y2s, z2s = sol.y
```
i try use RK45, LSODA.
Code from Julia:
```julia
function sigma(x)
return 1.0 / ( 1.0 + exp( -10.0 * ( x - ( - 0.25 ) ) ) )
end
function HR(u, p, t)
a, b, c, d, s, xr, r, I, vs, k1, k2, el_link = p
x1, y1, z1, x2, y2, z2 = u
du1 = y1 + b * x1 ^ 2 - a * x1 ^3 - z1 + I - k1 * ( x1 - vs ) * sigma(x2) + el_link * ( x2 - x1 )
du2 = c - d * x1 ^2 - y1
du3 = r * ( s * ( x1 - xr ) - z1 )
du4 = y2 + b * x2 ^ 2 - a * x2 ^3 - z2 + I - k2 * ( x2 - vs ) * sigma(x1) + el_link * ( x1 - x2 )
du5 = c - d * x2 ^2 - y2
du6 = r * ( s * ( x2 - xr ) - z2 )
return SVector(du1, du2, du3,
du4, du5, du6)
end
a = 1.0; b = 3.0; c = 1.0; d = 5.0
xr = -1.6; r = 0.01; s = 5.0; I = 4.0; xv = 2.0
k12 = -0.17
k1= k12; k2 = k12
k = 0.0
condition = SA [0.1, 0.1, 0.1, 0.1, 0.1, 0.1]
tspan = (0.0, 100000.0)
p = SA[a, b, c, d, s, xr, r, I, xv, k1, k2, k];
probsa = ODEProblem(HR, condition, tspan, p)
sol = solve(probsa, AutoVern9(Rodas5()), abstol = 1e-11, reltol = 1e-11, dense=false, maxiters = 50000000);
```
Result from python
and Julia
in time series i using sum x1 + x2.
Also i try use `dense = false` and `dense = true`
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [August 24, 2022, 4:37pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/2 "2022-08-24T16:37:37Z")
> [@Sergey\_Novak](#):
>
> for comparison i using integrator from Julia and python. I have used several integrators from Julia and python and get completely different results.
Is the system chaotic?
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 24, 2022, 4:38pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/3 "2022-08-24T16:38:15Z")
Yes, system with such parameters has multistability, that is, there is a regular chaotic attractor.
Using completely identical initial conditions Julia show regular attractor and python chaotic.
I also calculated the Lyapunov spectrum using two different algorithms jitcode and DynamicalSystems.jl, and the result indicates hyperhaos.
In addition, I identify the attractors and their basins. In result i have two attractor this regular and chaotic.
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 24, 2022, 4:39pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/4 "2022-08-24T16:39:27Z")
python integrators have too low accuracy and therefore fly into chaotic attractor?
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [August 24, 2022, 4:42pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/5 "2022-08-24T16:42:34Z")
No, this is just literally a property of chaos. Changes of 1e-16 will give O(1) differences in the trajectory. The difference of using `muladd` vs `x*y + z` is then enough to, over time after a Lyopunov time, give O(1) differences. You can only get exactly the same result out of chaotic systems if you have exactly, floating point exact, the same implementation.
Take a look at [https://diffeqcallbacks.sciml.ai/stable/uncertainty\_quantification/](https://diffeqcallbacks.sciml.ai/stable/uncertainty_quantification/) to explore the uncertainty of such systems.
> [@Sergey\_Novak](#):
>
> python integrators have too low accuracy and therefore fly into chaotic attractor?
Everything that approximates a chaotic attractor is incorrect. You only have the shadowing lemma: the fact that there’s a trajectory with an epsilon different initial condition on the attractor that you happen to have pulled, but it can be O(1) different from the real trajectory. Numerical approximations to chaotic systems are only trustworthy in a probabilistic sense.
Here’s a nice post that shows some of the details:
> **[Shadowing Methods for Forward and Adjoint Sensitivity Analysis of Chaotic...](https://frankschae.github.io/post/shadowing/)**
>
> In this post, we dig into sensitivity analysis of chaotic systems. Chaotic systems are dynamical, deterministic systems that are extremely sensitive to small changes in the initial state or the system parameters.
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 24, 2022, 4:49pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/6 "2022-08-24T16:49:41Z")
Thank you. Is it possible to somehow accurately determine and say that this trajectory with an initial condition is a regular mode? it exists only in one single case if all the initial conditions are the same
---
### Author:  [@josuagrw](https://discourse.julialang.org/u/josuagrw)
#### Post date: [August 24, 2022, 5:05pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/7 "2022-08-24T17:05:19Z")
Can you define the regularity (i.e. the rules) of the regular mode?
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 24, 2022, 5:14pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/8 "2022-08-24T17:14:09Z")
The regular regime can be determined by a period or lyapunov spectrum or largest lyapunov exponent, also power spectrum
---
### Author:  [@josuagrw](https://discourse.julialang.org/u/josuagrw)
#### Post date: [August 24, 2022, 5:22pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/9 "2022-08-24T17:22:07Z")
So you are saying a regular mode in your case means a periodic orbit
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 24, 2022, 5:25pm UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/10 "2022-08-24T17:25:40Z")
Yes, in articles people often write a regular mode, either a periodic orbit or a periodic trajectory
---
### Author:  [@josuagrw](https://discourse.julialang.org/u/josuagrw)
#### Post date: [August 25, 2022, 8:34am UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/11 "2022-08-25T08:34:09Z")
I would personally use `ApproxFun.jl` combined with an optimization package to find periodic solutions. (Of course it’s quite hard to show non-existence of a periodic solution this way.)
EDIT: Especially since strange attractors can have unstable periodic orbits which would show up in this approach.
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 25, 2022, 9:29am UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/12 "2022-08-25T09:29:15Z")
Thank you, I’ll try to use this packages. I started writing code on Julia relatively recently and don’t know many packages
---
### Author:  [@josuagrw](https://discourse.julialang.org/u/josuagrw)
#### Post date: [August 25, 2022, 9:41am UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/13 "2022-08-25T09:41:59Z")
The gist of this approach is explained in [section 2.3.4 of this article](https://iopscience.iop.org/article/10.1088/1367-2630/aae998/meta#njpaae998s2-3-4) (it’s definitely not original; just a methods section)
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 25, 2022, 10:02am UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/14 "2022-08-25T10:02:36Z")
Thank you. Do you happen to have a code example for using this with a continuous system?
---
### Author:  [@josuagrw](https://discourse.julialang.org/u/josuagrw)
#### Post date: [August 25, 2022, 10:50am UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/15 "2022-08-25T10:50:01Z")
Yes, have a look at this—it’s doesn’t have much documentation, unfortunately. And since it’s a few years old the code might not work completely with current versions of the packages.
> **[Jupyter Notebook Viewer](https://nbviewer.org/gist/gwater/a070be4869bc0e3d650c487d3133f4b5)**
>
> Check out this Jupyter notebook!
---
### Author:  [@Sergey\_Novak](https://discourse.julialang.org/u/Sergey_Novak)
#### Post date: [August 25, 2022, 11:11am UTC](https://discourse.julialang.org/t/comparison-of-integrators-differentialequations-jl/86271/16 "2022-08-25T11:11:56Z")
Thank you