Re-write source code of lyapunov to use set_state!

[Datseris]

The source code of lyapunov was one of the very first functions of DynamicalSystems.jl. It does not use the set_state! function for parallel integrators. That is a mistake (and in fact made us miss a bug).

One should re-write the rescale! function of the code to simply generate the new reference state and the just call set_state!.

[Datseris]

(the solution of this issue is literally changing 10 lines of code, it is very straightforward)

[SudoMishra]

How to give two indices in the set_state function for this rescale function :
function rescale!(integ::AbstractODEIntegrator{Alg, IIP, M}, a) where {Alg, IIP, M<:Matrix}
I used CartesianIndex(i,2) however there's no set_state function that takes k as a Cartesian Index, it takes an Int.
The rest of the rescale functions now use set_state. I'll send a PR once I get this one.
Thank you.

[Datseris]

Don't give two indices, but rather call set_state! twice with the first and then the second index.

[SudoMishra]

I'm not sure I get what you mean. For the matrix case, I need to assign integ.u[i,2] = u. However, last parameter of set_state! takes an int. How do you call set_state! twice so that it assigns u at the matrix index(i,2)?

[Datseris]

I'm confused, can you please link the source code line you are referring to in the rescale function? :)

[SudoMishra]

Yeah sure. I'm talking about the rescale! function on line 339 in lyapunov.jl. Specifically, line 341. It currently uses
integ.u[i, 2] = integ.u[i,1] + (integ.u[i,2] - integ.u[i,1])/a. To rewrite this using set_state!, I'm trying to do it like this

du = integ.u[i,1] + (integ.u[i,2] - integ.u[i,1])/a
set_state!(integ,du,CartesianIndex(i,2))

which would be equivalent to integ.u[i,2] = du.
However, with this I get an error that no set_state! function matches with Cartesian_Index.
Thank you for responding.

[Datseris]

@SudoMishra yeah, you're confused because the point is to re-write the function. Have a read at the documentation string of set_state!. The point is, you need to create the entire rescaled state and the set it. Eg. do

function rescale!(integ::AbstractODEIntegrator{Alg, IIP, M}, a) where {Alg, IIP, M<:Matrix}
    r = get_state(integ, 1) .+ (get_state(integ, 2) .- get_state(integ, 1)) ./ a
    set_state!(integ, r, 2)
end

Of course, we can do some optimization so that we have a dummy preinitialized r that we carry around and modify its value in place so that we don't create a new r all the time, but that's the code. It's pretty simple. And as you can see it kind of unifies the different rescale functions. there is no longer a need for many rescale functions,

[Datseris]

if you open a pr I will give more advice