search

JuliaDynamics / DynamicalSystems.jl

Award winning software library for nonlinear dynamics and nonlinear timeseries analysis
https://juliadynamics.github.io/DynamicalSystemsDocs.jl/dynamicalsystems/dev/
Other
847 stars 95 forks source link

Use `integrator` and `step` interfaces in `lyapunov` and `lyapunovs` functions for Continuous systems #9

Closed Datseris closed 7 years ago

Datseris commented 7 years ago

The methods that calculate lyapunovs of continuous systems are extremely slow. So far I am not sure whether this is due to DifferentialEquations, the qr-decomposition or that I Have simply written bad code. Profiling them will reveal issues in my code (which probably sucks).

In any case, the have to be made faster!

One way that this can be improved is by taking advantage of the integrator and step! interfaces of DiffEq.

ChrisRackauckas commented 7 years ago

Do you have a test code to profile?

Datseris commented 7 years ago

I'll cook one up in 30 minutes; I'll let you know then. The methods that are slow are lyapunovs(ds::ContinuousDS) and lyapunov(ds::ContinuousDS). The first one is really slow. For the second one, I think the reason is that it is slow for DifferentialEquations to initialize solve(), which for me stops and restarts all the time.

ChrisRackauckas commented 7 years ago

For the second one, I think the reason is that it is slow for DifferentialEquations to initialize solve(), which for me stops and restarts all the time.

What do you mean by this?

ChrisRackauckas commented 7 years ago

I see what you mean. I just took a look at

https://github.com/Datseris/DynamicalSystems.jl/blob/62bcf1a8a7007135ae5d5c8da41f3a8082d7072a/src/lyapunovs.jl#L171 https://github.com/Datseris/DynamicalSystems.jl/blob/76221d4ed9853d32827d7d9745f3c3c11ef9c744/src/continuous.jl#L120

Yes, solve will re-allocate the worker arrays each time. This probably isn't what you want to do. Instead, you likely want to use the integrator interface which is a dynamical iterator over a DiffEq. Only exists in the *DiffEq packages, but it works like this:

http://docs.juliadiffeq.org/latest/basics/integrator.html https://github.com/JuliaDiffEq/OrdinaryDiffEq.jl/blob/master/test/iterator_tests.jl

Essentially once you have the integrator you can step! whenever you want, and change integrator.u. If you change integrator.u, then you should u_modified!(integrator,true) so the solver knows FSAL is violated. But this will give you full control for the step-by-step behavior without having to reinitialize all of the time.

http://docs.juliadiffeq.org/latest/basics/integrator.html#Stepping-Controls-1

Datseris commented 7 years ago

Okay, I was writing a long comment but no reason anymore :D

EDIT: Renaming issue to "use integrator and step! interface of OrdinaryDiffEq.

Datseris commented 7 years ago

@ChrisRackauckas from your docs I see that u_modified!(integrator, true) should be used if I want to change u discontinuously.

This is not the case here, I was simply setting the last state as the "current" because I created a (pretty bad) DIY stepper. If I use the integrator that won't be necessary at all, since it happens automatically internally?

ChrisRackauckas commented 7 years ago

This is not the case here, I was simply setting the last state as the "current" because I created a (pretty bad) DIY stepper. If I use the integrator that won't be necessary at all, since it happens automatically internally?

Yes.

Datseris commented 7 years ago

This was solved many moons ago...