Edgetracking algorithm

[reykboerner]

Implementation of the edgetracking algorithm, as planned in issue #7

The implementation has been tested for the 2D bistable FitzHugh-Nagumo model (test added), a 5D ocean box model, and a 5D chaotic system of a 2D ocean coupled to a 3D atmosphere.

The code features some output control but this could be expanded. Also, currently the code has only been used for systems where one is interested in the edge between two adjacent basins of attraction. If there are more than two attractors with possibly fractal basin boundaries, the code should still work (raising a warning if a third basin is involved) but it becomes less clear which edge one actually wishes to track. I have some ideas how this could be addressed, happy to discuss.

[codecov-commenter]

Codecov Report

Attention: 8 lines in your changes are missing coverage. Please review.

Comparison is base (b049d3a) 84.15% compared to head (a32084a) 84.85%.

Files	Patch %	Lines
src/boundaries/edgetracking.jl	91.57%	8 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main      #61      +/-   ##
==========================================
+ Coverage   84.15%   84.85%   +0.69%     
==========================================
  Files          23       24       +1     
  Lines        1199     1294      +95     
==========================================
+ Hits         1009     1098      +89     
- Misses        190      196       +6     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[Datseris]

Great, thanks a lot! I am reviewing this! I hope to submit a review by the end of the weekend!

[awage]

Hi! This is quite extraordinary. I have an implementation of my own of this technique! Originally it was called saddle-straddle by Yorke in his book on the software Dynamics. We made a paper with this stuff, look at figure 1 in this paper:
https://arxiv.org/pdf/1901.10728.pdf

I had no idea that it was used in other context as well.

[Datseris]

Hi! This is quite extraordinary. I have an implementation of my own of this technique!

Hah, what a small world!!!

[reykboerner]

Thanks @awage, I wasn't aware that the method actually dates that far back! I should correct the references then and add the original Yorke papers. It seems that several people have implementations of this algorithm but usually for their system of interest.

I also read in your paper that you find the boundary between one basin and all other basins. To make the present code more general, one could add keyword arguments to specify two sets of attractors (by giving their respective Dict keys), so that the algorithm tracks the boundary separating these two sets.

[awage]

Hi! This is quite extraordinary. I have an implementation of my own of this technique!

Hah, what a small world!!!

This paper I wrote was my first serious work in Julia and the seed of Attractors.jl. So-... Cheers!

[awage]

Thanks @awage, I wasn't aware that the method actually dates that far back! I should correct the references then and add the original Yorke papers. It seems that several people have implementations of this algorithm but usually for their system of interest.

I also read in your paper that you find the boundary between one basin and all other basins. To make the present code more general, one could add keyword arguments to specify two sets of attractors (by giving their respective Dict keys), so that the algorithm tracks the boundary separating these two sets.

I will give you more detail about the original work. Also, they baked in the concept of "general basins", which is simply the combination of different basins. This is more general as you suggest, but I don't know if it is really useful unless in margin cases such as the one we used in the paper. We can discuss this.

I have a very simple trivial test for the algorithm, this is the equation of the cubic map:

$x_{n+1} = μ x_n − x_n^3$

they have two stable fixed point and and a saddle at x=0 for \mu > 1 . The algorithm should find the saddle straight away. We like simple predictable tests in DynamicalSystems.jl

[Datseris]

Great, thanks a lot for the contribution @reykboerner . To be honest, this is surprisingly good, and this is a complement. All basic checkpoints, like documentation, tests, legitimate code that runs are checked. So the PR is really good quality, something that I don't typically expect from a first time contributor.

I'm leaving s e v e r a l comments as part of my review in this PR. All comments are rather direct and to the point. This is a harsh style I hope you can get used to it quickly because it's by far the most efficient way to review PRs.

Feel free to address them all at once or step by step and ask questions whenever things are not clear. I think after you address the comments, we can have a second round of review.

---

There is one thing that bugs me: Let's try to be proactive here. Ideally, we would like to have code that finds the basin boundary. Is this possible to build upon this code base? Sure, one way to do that is via brute force based on the bijection. You pick all possible directions around a fixed point and then you biject until you find this boundary. THis has several problems: computationally unfeasible in high dim + it doesn't work in basins that are non-convex, such as those described in #60 .

So, given that your code finds one point on the boundary, what algorithms are there that "continue" this point to find the full boundary?

[reykboerner]

Thanks a lot @Datseris for the feedback and detailed reviews! I'm completely on board with this style of reviewing. I'll respond to each point in the next days.

[Datseris]

bump!

[reykboerner]

bump!

Sorry, I will get back to this asap!

[Datseris]

bump @reykboerner

[reykboerner]

Thanks @awage, I wasn't aware that the method actually dates that far back! I should correct the references then and add the original Yorke papers. It seems that several people have implementations of this algorithm but usually for their system of interest.
I also read in your paper that you find the boundary between one basin and all other basins. To make the present code more general, one could add keyword arguments to specify two sets of attractors (by giving their respective Dict keys), so that the algorithm tracks the boundary separating these two sets.

I will give you more detail about the original work. Also, they baked in the concept of "general basins", which is simply the combination of different basins. This is more general as you suggest, but I don't know if it is really useful unless in margin cases such as the one we used in the paper. We can discuss this.

I have a very simple trivial test for the algorithm, this is the equation of the cubic map:

xn+1=μxn−xn3

they have two stable fixed point and and a saddle at x=0 for \mu > 1 . The algorithm should find the saddle straight away. We like simple predictable tests in DynamicalSystems.jl

[reykboerner]

I have a very simple trivial test for the algorithm, this is the equation of the cubic map:
...
they have two stable fixed point and and a saddle at x=0 for \mu > 1 . The algorithm should find the saddle straight away. We like simple predictable tests in DynamicalSystems.jl

I wanted to add this test but had an issue with AttractorsViaProximity.

The code below produces the error BoundsError: attempt to access 0-element Vector{Float64} at index [1], which has something to do with mapper.dist being an empty vector (line 118 of src/mapping/attractor_mapping_proximity.jl:

using Attractors
cubicmap(u, p, n) = SVector{1}(p[1]*u[1] - u[1]^3)
ds = DeterministicIteratedMap(cubicmap, [1.0], [2.0])
fp = [sqrt(2)]
attrs = Dict(1 => StateSpaceSet([fp]), 2 => StateSpaceSet([-fp]))
mapper = AttractorsViaProximity(ds, attrs)
label = mapper([2.0])

This error also comes up when trying to run edgetracking for this case, since mapper is called repeatedly inside edgetracking. Could you help with this?

[reykboerner]

Hi @Datseris, I made all changes we talked about and checked that the three testsets pass. However, I cannot request merging the PR since the docs fail to build. It seems like this is also the case on the origin main branch. Could you help please?

[Datseris]

ok, thanks, i'll review the PR and have a look at the documentation build.

[Datseris]

Yes !!! SO CLOSE! ITS HAPPENING

I have some comments left.

And, I have another comment: we also discussed that you would write a small example in the Examples in the documentation. Can you please copy one of the tests, hopefully the most visually appealing one, and make it into an example as well?

Meanwhike, I will try to fix why the test builds fail.

[reykboerner]

Thanks for your review, I have implemented all changes.

And, I have another comment: we also discussed that you would write a small example in the Examples in the documentation. Can you please copy one of the tests, hopefully the most visually appealing one, and make it into an example as well?

Sure, I would suggest taking the FitzHugh-Nagumo example as it is simple and easy to visualize. The Thomas' rule is maybe a bit more spectacular since it is chaotic but with 3D and multiple attractors and edge states, it becomes difficult to visualize. I'll prepare the example.

[Datseris]

I've fixed the failing docuentation build.

[reykboerner]

Hi @Datseris, thanks for fixing the docs! I've added the example - any further changes needed before merging?

[Datseris]

Great, thanks a lot @reykboerner . I am merging this in. I will tag a new release along with #120 !

[reykboerner]

Great, thanks a lot @reykboerner . I am merging this in. I will tag a new release along with #120 !

Thanks @Datseris, I'm glad I could contribute and learned a lot! There still seems to be an issue with the docs - I think running the code in my example requires the missing using OrdinaryDiffEq because I am specifying alg=Vern9() in the CoupledODEs definition of fitzhugh_nagumo.

How does using additional packages work in the docs examples? Should OrdinaryDiffEq be added to the deps in docs/Project.toml and then the line using OrdinaryDiffEq added in my example in docs/src/examples.md?

If this is correct I can make these changes and make another PR to include them.

[Datseris]

yes thats correct. Ordinary DIffeq should already be in docs/Project.toml so just the using should be missing

[Datseris]

(i was going to check docs after merging all other PRs and before a release, but thanks for catching this early)