|
321 | 321 |
|
322 | 322 | To achieve even better results for this kind of problematic systems than with previuosly introduced `Irregular Grids` we provide a functionality to construct `Subdivision Based Grids` in which |
323 | 323 | one can obtain more coarse or dense structure not only along some axis but for a specific regions where the state space flow has |
324 | | -significantly different speed. [`subdivided_based_grid`](@ref) enables automatic evaluation of velocity vectors for regions of originally user specified |
| 324 | +significantly different speed. [`subdivision_based_grid`](@ref) enables automatic evaluation of velocity vectors for regions of originally user specified |
325 | 325 | grid to further treat those areas as having more dense or coarse structure than others. |
326 | 326 |
|
327 | 327 | ```@example MAIN |
@@ -536,7 +536,7 @@ plot_basins_curves(aggregated_fractions, prange; |
536 | 536 |
|
537 | 537 | ## Trivial featurizing and grouping for basins fractions |
538 | 538 |
|
539 | | -This is a rather trivial example showcasing the usage of [`AttractorsViaFeaturizing`](@ref). Let us use once again the magnetic pendulum example. For it, we have a really good idea of what features will uniquely describe each attractor: the last points of a trajectory (which should be very close to the magnetic the trajectory converged to). To provide this information to the [`AttractorsviaFeaturizing`](@ref) we just create a julia function that returns this last point |
| 539 | +This is a rather trivial example showcasing the usage of [`AttractorsViaFeaturizing`](@ref). Let us use once again the magnetic pendulum example. For it, we have a really good idea of what features will uniquely describe each attractor: the last points of a trajectory (which should be very close to the magnetic the trajectory converged to). To provide this information to the [`AttractorsViaFeaturizing`](@ref) we just create a julia function that returns this last point |
540 | 540 |
|
541 | 541 | ```@example MAIN |
542 | 542 | using Attractors |
@@ -871,3 +871,92 @@ animate_attractors_via_recurrences(mapper, u0s) |
871 | 871 | <source src="https://raw.githubusercontent.com/JuliaDynamics/JuliaDynamics/master/videos/attractors/recurrence_algorithm.mp4?raw=true" type="video/mp4"> |
872 | 872 | </video> |
873 | 873 | ``` |
| 874 | + |
| 875 | +## Edge tracking |
| 876 | + |
| 877 | +To showcase how to run the [`edgetracking`](@ref) algorithm, let us use it to find the |
| 878 | +saddle point of the bistable FitzHugh-Nagumo (FHN) model, a two-dimensional ODE system |
| 879 | +originally conceived to represent a spiking neuron. |
| 880 | +We define the system in the following form: |
| 881 | + |
| 882 | +```@example MAIN |
| 883 | +function fitzhugh_nagumo(u,p,t) |
| 884 | + x, y = u |
| 885 | + eps, beta = p |
| 886 | + dx = (x - x^3 - y)/eps |
| 887 | + dy = -beta*y + x |
| 888 | + return SVector{2}([dx, dy]) |
| 889 | +end |
| 890 | +
|
| 891 | +params = [0.1, 3.0] |
| 892 | +ds = CoupledODEs(fitzhugh_nagumo, ones(2), params, diffeq=(;alg = Vern9(), reltol=1e-11)) |
| 893 | +``` |
| 894 | + |
| 895 | +Now, we can use Attractors.jl to compute the fixed points and basins of attraction of the |
| 896 | +FHN model. |
| 897 | + |
| 898 | +```@example MAIN |
| 899 | +xg = yg = range(-1.5, 1.5; length = 201) |
| 900 | +grid = (xg, yg) |
| 901 | +mapper = AttractorsViaRecurrences(ds, grid; sparse=false) |
| 902 | +basins, attractors = basins_of_attraction(mapper) |
| 903 | +
|
| 904 | +for i in 1:length(attractors) |
| 905 | + println(attractors[i][1]) |
| 906 | +end |
| 907 | +``` |
| 908 | + |
| 909 | +The `basins_of_attraction` function found three fixed points: the two stable nodes of the |
| 910 | +system (labeled A and B) and the saddle point at the origin. The saddle is an unstable |
| 911 | +equilibrium and can therefore not be found by simulation, but we can find it using the |
| 912 | +[`edgetracking`](@ref) algorithm. For illustration, let us initialize the algorithm from |
| 913 | +two initial conditions `init1` and `init2` (which must belong to different basins |
| 914 | +of attraction, see figure below). |
| 915 | + |
| 916 | +```julia |
| 917 | +attractors_AB = Dict(1 => attractors[1], 2 => attractors[2]) |
| 918 | +init1, init2 = [-1.0, -1.0], [-1.0, 0.2] |
| 919 | +``` |
| 920 | + |
| 921 | +Now, we run the edge tracking algorithm: |
| 922 | + |
| 923 | +```@example MAIN |
| 924 | +bisect_thresh, diverge_thresh, Δt, abstol = 1e-3, 2e-3, 1e-5, 1e-3 |
| 925 | +et = edgetracking(ds, attractors_AB; u1=init1, u2=init2, |
| 926 | + bisect_thresh, diverge_thresh, Δt, abstol) |
| 927 | +
|
| 928 | +et.edge[end] |
| 929 | +``` |
| 930 | + |
| 931 | +The algorithm has converged to the origin (up to the specified accuracy) where the saddle |
| 932 | +is located. The figure below shows how the algorithm has iteratively tracked along the basin |
| 933 | +boundary from the two initial conditions (red points) to the saddle (green square). Points |
| 934 | +of the edge track (orange) at which a re-bisection occured are marked with a white border. |
| 935 | +The figure also depicts two trajectories (blue) intialized on either side of the basin |
| 936 | +boundary at the first bisection point. We see that these trajectories follow the basin |
| 937 | +boundary for a while but then relax to either attractor before reaching the saddle. By |
| 938 | +counteracting the instability of the saddle, the edge tracking algorithm instead allows to |
| 939 | +track the basin boundary all the way to the saddle, or edge state. |
| 940 | + |
| 941 | +```@example MAIN |
| 942 | +traj1 = trajectory(ds, 2, et.track1[et.bisect_idx[1]], Δt=1e-5) |
| 943 | +traj2 = trajectory(ds, 2, et.track2[et.bisect_idx[1]], Δt=1e-5) |
| 944 | +
|
| 945 | +fig = Figure() |
| 946 | +ax = Axis(fig[1,1], xlabel="x", ylabel="y") |
| 947 | +heatmap_basins_attractors!(ax, grid, basins, attractors, add_legend=false, labels=Dict(1=>"Attractor A", 2=>"Attractor B", 3=>"Saddle")) |
| 948 | +lines!(ax, traj1[1][:,1], traj1[1][:,2], color=:dodgerblue, linewidth=2, label="Trajectories") |
| 949 | +lines!(ax, traj2[1][:,1], traj2[1][:,2], color=:dodgerblue, linewidth=2) |
| 950 | +lines!(ax, et.edge[:,1], et.edge[:,2], color=:orange, linestyle=:dash) |
| 951 | +scatter!(ax, et.edge[et.bisect_idx,1], et.edge[et.bisect_idx,2], color=:white, markersize=15, marker=:circle, zorder=10) |
| 952 | +scatter!(ax, et.edge[:,1], et.edge[:,2], color=:orange, markersize=11, marker=:circle, zorder=10, label="Edge track") |
| 953 | +scatter!(ax, [-1.0,-1.0], [-1.0, 0.2], color=:red, markersize=15, label="Initial conditions") |
| 954 | +xlims!(ax, -1.2, 1.1); ylims!(ax, -1.3, 0.8) |
| 955 | +axislegend(ax, position=:rb) |
| 956 | +fig |
| 957 | +``` |
| 958 | + |
| 959 | +In this simple two-dimensional model, we could of course have found the saddle directly by |
| 960 | +computing the zeroes of the ODE system. However, the edge tracking algorithm allows finding |
| 961 | +edge states also in high-dimensional and chaotic systems where a simple computation of |
| 962 | +unstable equilibria becomes infeasible. |
0 commit comments