Skip to content

Commit

Permalink
update for [email protected]
Browse files Browse the repository at this point in the history
  • Loading branch information
rveltz committed Nov 27, 2024
1 parent e58fa80 commit 4da7e33
Show file tree
Hide file tree
Showing 9 changed files with 95 additions and 96 deletions.
4 changes: 2 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "HclinicBifurcationKit"
uuid = "0e751702-a6d5-4e4e-bf0f-194619fbec56"
authors = ["romain veltz <[email protected]>"]
version = "0.1.1"
version = "0.1.2"

[deps]
BifurcationKit = "0f109fa4-8a5d-4b75-95aa-f515264e7665"
Expand All @@ -14,7 +14,7 @@ SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"

[compat]
BifurcationKit = "^0.3.6"
BifurcationKit = "^0.3.6, ^0.4"
DocStringExtensions = "0.9"
ForwardDiff = "^0.10"
Parameters = "^0.12"
Expand Down
37 changes: 18 additions & 19 deletions docs/src/tutorials/ode/OPL.md
Original file line number Diff line number Diff line change
Expand Up @@ -20,12 +20,11 @@ It is easy to encode the ODE as follows

```@example TUTOPL
using Revise, Plots
using Setfield
using BifurcationKit
using HclinicBifurcationKit
const BK = BifurcationKit
record_from_solution(x, p) = (D₂₃ = x[6], β = x[1],)
record_from_solution(x, p; k...) = (D₂₃ = x[6], β = x[1],)
####################################################################################################
function OPL!(dz, u, p, t = 0)
(;b, σ, g, a, D₂₁⁰, D₂₃⁰) = p
Expand All @@ -41,7 +40,7 @@ end
par_OPL = (b = 1.2, σ = 2.0, g=50., a = 1., D₂₁⁰ = -1., D₂₃⁰ = 0.)
z0 = zeros(6)
prob = BK.BifurcationProblem(OPL!, z0, par_OPL, (@lens _.a); record_from_solution)
prob = BK.BifurcationProblem(OPL!, z0, par_OPL, (@optic _.a); record_from_solution)
nothing #hide
```
Expand All @@ -62,7 +61,7 @@ scene = plot(br, br2)


```@example TUTOPL
sn_br = continuation(br, 1, (@lens _.b), ContinuationPar(opts_br, detect_bifurcation = 1, max_steps = 80) ;
sn_br = continuation(br, 1, (@optic _.b), ContinuationPar(opts_br, detect_bifurcation = 1, max_steps = 80) ;
detect_codim2_bifurcation = 2,
start_with_eigen = true,
bothside = true,
Expand All @@ -71,12 +70,12 @@ show(sn_br)
```

```@example TUTOPL
hopf_br = continuation(br, 2, (@lens _.b), ContinuationPar(opts_br, detect_bifurcation = 1, max_steps = 140),
hopf_br = continuation(br, 2, (@optic _.b), ContinuationPar(opts_br, detect_bifurcation = 1, max_steps = 140),
detect_codim2_bifurcation = 2,
bothside = true,
)
hopf_br2 = continuation(br2, 1, (@lens _.b), ContinuationPar(opts_br, detect_bifurcation = 1, max_steps = 140),
hopf_br2 = continuation(br2, 1, (@optic _.b), ContinuationPar(opts_br, detect_bifurcation = 1, max_steps = 140),
detect_codim2_bifurcation = 2,
bothside = true,
)
Expand All @@ -103,7 +102,7 @@ function plotPO(x, p; k...)
end
# record function
function recordPO(x, p)
function recordPO(x, p; k...)
xtt = BK.get_periodic_orbit(p.prob, x, set(getparams(p.prob), BK.getlens(p.prob), p.p))
period = BK.getperiod(p.prob, x, p.p)
return (max = maximum(xtt[6,:]), min = minimum(xtt[6,:]), period = period, )
Expand Down Expand Up @@ -151,10 +150,10 @@ probhom, solh = generate_hom_problem(
# ϵ0 = 1e-7, ϵ1 = 1e-7, # WORK BEST
# ϵ0 = 1e-8, ϵ1 = 1e-5, # maxT = 70,
t0 = 0., t1 = 120.,
# freeparams = ((@lens _.T), (@lens _.ϵ1),)
# freeparams = ((@lens _.T), (@lens _.ϵ0)),
freeparams = ((@lens _.ϵ0), (@lens _.ϵ1)), # WORK BEST
# freeparams = ((@lens _.T),),
# freeparams = ((@optic _.T), (@optic _.ϵ1),)
# freeparams = ((@optic _.T), (@optic _.ϵ0)),
freeparams = ((@optic _.ϵ0), (@optic _.ϵ1)), # WORK BEST
# freeparams = ((@optic _.T),),
testOrbitFlip = false,
testInclinationFlip = false
)
Expand All @@ -169,15 +168,15 @@ optn_hom = NewtonPar(verbose = true, tol = 1e-10, max_iterations = 5)
optc_hom = ContinuationPar(newton_options = optn_hom, ds = -1e-4, dsmin = 1e-5, detect_bifurcation = 0, detect_event = 2, p_min = -1.01, dsmax = 4e-2, p_max = 1.5, max_steps = 100)
br_hom_c = continuation(
deepcopy(probhom), solh, (@lens _.b),
deepcopy(probhom), solh, (@optic _.b),
PALC(tangent = Bordered()),
optc_hom;
verbosity = 1, plot = true,
normC = norminf,
plot_solution = (x,p;k...) -> begin
𝐇𝐨𝐦 = p.prob
par0 = set(BK.getparams(𝐇𝐨𝐦), BK.getlens(𝐇𝐨𝐦), x.x[end][1])
par0 = set(par0, (@lens _.b), p.p)
par0 = set(par0, (@optic _.b), p.p)
sol = get_homoclinic_orbit(𝐇𝐨𝐦, x, par0)
m = (𝐇𝐨𝐦.bvp isa PeriodicOrbitOCollProblem && 𝐇𝐨𝐦.bvp.meshadapt) ? :d : :none
plot!(sol.t, sol[:,:]',subplot=3, markersize = 1, marker=m)
Expand Down Expand Up @@ -238,10 +237,10 @@ probhom, solh = generate_hom_problem(
update_every_step = 4,
# ϵ0 = 1e-6, ϵ1 = 1e-5,
t0 = 75, t1 = 25,
# freeparams = ((@lens _.T), (@lens _.ϵ1),)
# freeparams = ((@lens _.T), (@lens _.ϵ0)),
# freeparams = ((@lens _.ϵ0), (@lens _.ϵ1)), # WORK BEST
freeparams = ((@lens _.T),),
# freeparams = ((@optic _.T), (@optic _.ϵ1),)
# freeparams = ((@optic _.T), (@optic _.ϵ0)),
# freeparams = ((@optic _.ϵ0), (@optic _.ϵ1)), # WORK BEST
freeparams = ((@optic _.T),),
)
_sol = get_homoclinic_orbit(probhom, solh, BK.getparams(probhom); saveat=.1)
Expand All @@ -251,7 +250,7 @@ optn_hom = NewtonPar(verbose = true, tol = 1e-9, max_iterations = 7)
optc_hom = ContinuationPar(newton_options = optn_hom, ds = -1e-4, dsmin = 1e-6, max_steps = 300, detect_bifurcation = 0, dsmax = 12e-2, plot_every_step = 3, p_max = 7., detect_event = 2, a = 0.9)
br_hom_sh = continuation(
deepcopy(probhom), solh, (@lens _.b),
deepcopy(probhom), solh, (@optic _.b),
PALC(),
optc_hom;
verbosity = 3, plot = true,
Expand All @@ -260,7 +259,7 @@ br_hom_sh = continuation(
plot_solution = (x,p;k...) -> begin
𝐇𝐨𝐦 = p.prob
par0 = set(BK.getparams(𝐇𝐨𝐦), BK.getlens(𝐇𝐨𝐦), x.x[end][1])
par0 = set(par0, (@lens _.b), p.p)
par0 = set(par0, (@optic _.b), p.p)
sol = get_homoclinic_orbit(𝐇𝐨𝐦, x, par0)
m = (𝐇𝐨𝐦.bvp isa PeriodicOrbitOCollProblem && 𝐇𝐨𝐦.bvp.meshadapt) ? :d : :none
plot!(sol.t, sol[1:6,:]',subplot=3, markersize = 1, marker=m)
Expand Down
24 changes: 12 additions & 12 deletions docs/src/tutorials/ode/tutorialsFreire.md
Original file line number Diff line number Diff line change
Expand Up @@ -18,12 +18,12 @@ It is easy to encode the ODE as follows

```@example TUTFREIRE
using Revise, Plots
using Setfield, LinearAlgebra, Test, ForwardDiff
using LinearAlgebra, Test, ForwardDiff
using BifurcationKit, Test
using HclinicBifurcationKit
const BK = BifurcationKit
recordFromSolution(x, p) = (x = x[1], y = x[2])
recordFromSolution(x, p; k...) = (x = x[1], y = x[2])
function freire!(dz, u, p, t)
(;ν, β, A₃, B₃, r, ϵ) = p
Expand All @@ -38,7 +38,7 @@ freire(z, p) = freire!(similar(z), z, p, 0)
par_freire = (ν = -0.75, β = -0.1, A₃ = 0.328578, B₃ = 0.933578, r = 0.6, ϵ = 0.01)
z0 = [0.7,0.3,0.1]
z0 = zeros(3)
prob = BK.BifurcationProblem(freire, z0, par_freire, (@lens _.β); record_from_solution = recordFromSolution)
prob = BK.BifurcationProblem(freire, z0, par_freire, (@optic _.β); record_from_solution = recordFromSolution)
nothing #hide
```
Expand All @@ -63,15 +63,15 @@ br


```@example TUTFREIRE
sn_br = continuation(br, 2, (@lens _.ν), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, dsmax = 0.01, max_steps = 80) ;
sn_br = continuation(br, 2, (@optic _.ν), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, dsmax = 0.01, max_steps = 80) ;
alg = PALC(),
detect_codim2_bifurcation = 2,
start_with_eigen = true,
update_minaug_every_step = 1,
bothside = true,
)
hopf_br = continuation(br, 4, (@lens _.ν), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 140, dsmax = 0.02, n_inversion = 6),
hopf_br = continuation(br, 4, (@optic _.ν), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 140, dsmax = 0.02, n_inversion = 6),
detect_codim2_bifurcation = 2,
start_with_eigen = true,
update_minaug_every_step = 1,
Expand Down Expand Up @@ -112,9 +112,9 @@ br_hom_c = continuation(
setproperties(opts_br, max_steps = 30, save_sol_every_step = 1, dsmax = 1e-2, plot_every_step = 1, p_min = -1.01, ds = 0.001, detect_event = 2, detect_bifurcation = 0);
verbosity = 0, plot = false,
ϵ0 = 1e-5, amplitude = 2e-3,
# freeparams = ((@lens _.T), (@lens _.ϵ1),)
# freeparams = ((@lens _.T), (@lens _.ϵ0)),
freeparams = ((@lens _.ϵ0), (@lens _.ϵ1)),
# freeparams = ((@optic _.T), (@optic _.ϵ1),)
# freeparams = ((@optic _.T), (@optic _.ϵ0)),
freeparams = ((@optic _.ϵ0), (@optic _.ϵ1)),
normC = norminf,
update_every_step = 4,
)
Expand All @@ -140,10 +140,10 @@ br_hom_sh = continuation(
verbosity = 1, plot = true,
ϵ0 = 1e-6, amplitude = 2e-2,
update_every_step = 2,
# freeparams = ((@lens _.T), (@lens _.ϵ1),),
# freeparams = ((@lens _.T), (@lens _.ϵ0)),
# freeparams = ((@lens _.ϵ0), (@lens _.ϵ1)),
freeparams = ((@lens _.T),),
# freeparams = ((@optic _.T), (@optic _.ϵ1),),
# freeparams = ((@optic _.T), (@optic _.ϵ0)),
# freeparams = ((@optic _.ϵ0), (@optic _.ϵ1)),
freeparams = ((@optic _.T),),
normC = norminf,
plot_solution = plotHom,
maxT = 45.,
Expand Down
36 changes: 18 additions & 18 deletions examples/OPL.jl
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ using BifurcationKit, Test
using HclinicBifurcationKit
const BK = BifurcationKit

recordFromSolution(x, p) = (D₂₃ = x[6], β = x[1],)
recordFromSolution(x, p;k...) = (D₂₃ = x[6], β = x[1],)
####################################################################################################
function OPL!(dz, u, p, t)
(;b, σ, g, a, D₂₁⁰, D₂₃⁰) = p
Expand All @@ -21,7 +21,7 @@ end
OPL(z, p) = OPL!(similar(z), z, p, 0)
par_OPL = (b = 1.2, σ = 2.0, g=50., a = 1., D₂₁⁰ = -1., D₂₃⁰ = 0.)
z0 = zeros(6)
prob = BK.BifurcationProblem(OPL, z0, par_OPL, (@lens _.a); record_from_solution = recordFromSolution)
prob = BK.BifurcationProblem(OPL, z0, par_OPL, (@optic _.a); record_from_solution = recordFromSolution)

opts_br = ContinuationPar(p_min = -1., p_max = 8., ds = 0.001, dsmax = 0.06, n_inversion = 6, detect_bifurcation = 3, max_bisection_steps = 25, nev = 6, plot_every_step = 20, max_steps = 100, save_sol_every_step = 1, detect_fold = true)
br = continuation(prob, PALC(tangent = Secant()), opts_br;
Expand All @@ -33,7 +33,7 @@ br2 = continuation(br, 1)

plot(br, br2)
####################################################################################################
sn_br = continuation(br, 1, (@lens _.b), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 80) ;
sn_br = continuation(br, 1, (@optic _.b), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 80) ;
alg = PALC(),
verbosity = 0,
detect_codim2_bifurcation = 2,
Expand All @@ -45,14 +45,14 @@ sn_br = continuation(br, 1, (@lens _.b), ContinuationPar(opts_br, detect_bifurca
plot(sn_br)
bt = get_normal_form(sn_br, 2, verbose = true, detailed = true, autodiff = false)

hopf_br = continuation(br, 2, (@lens _.b), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 140),
hopf_br = continuation(br, 2, (@optic _.b), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 140),
detect_codim2_bifurcation = 2,
start_with_eigen = true,
update_minaug_every_step = 1,
bothside = true,
)

hopf_br2 = continuation(br2, 1, (@lens _.b), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 140),
hopf_br2 = continuation(br2, 1, (@optic _.b), ContinuationPar(opts_br, detect_bifurcation = 1, save_sol_every_step = 1, max_steps = 140),
detect_codim2_bifurcation = 2,
start_with_eigen = true,
update_minaug_every_step = 1,
Expand All @@ -74,7 +74,7 @@ function plotPO(x, p; k...)
end

# record function
function recordPO(x, p)
function recordPO(x, p; k...)
xtt = BK.get_periodic_orbit(p.prob, x, set(getparams(p.prob), BK.getlens(p.prob), p.p))
period = BK.getperiod(p.prob, x, p.p)
return (max = maximum(xtt[6,:]), min = minimum(xtt[6,:]), period = period, )
Expand Down Expand Up @@ -122,10 +122,10 @@ probhom, solh = generate_hom_problem(
# ϵ0 = 1e-7, ϵ1 = 1e-7, # WORK BEST
# ϵ0 = 1e-8, ϵ1 = 1e-5, # maxT = 70,
t0 = 0., t1 = 120.,
# freeparams = ((@lens _.T), (@lens _.ϵ1),)
# freeparams = ((@lens _.T), (@lens _.ϵ0)),
freeparams = ((@lens _.ϵ0), (@lens _.ϵ1)), # WORK BEST
# freeparams = ((@lens _.T),),
# freeparams = ((@optic _.T), (@optic _.ϵ1),)
# freeparams = ((@optic _.T), (@optic _.ϵ0)),
freeparams = ((@optic _.ϵ0), (@optic _.ϵ1)), # WORK BEST
# freeparams = ((@optic _.T),),
testOrbitFlip = false,
testInclinationFlip = false
)
Expand All @@ -139,7 +139,7 @@ optn_hom = NewtonPar(verbose = true, tol = 1e-10, max_iterations = 5)
optc_hom = ContinuationPar(newton_options = optn_hom, ds = -0.0001, dsmin = 1e-5, plot_every_step = 10, max_steps = 100, detect_bifurcation = 0, detect_event = 2, save_sol_every_step = 1, p_min = -1.01)

br_hom_c = continuation(
deepcopy(probhom), solh, (@lens _.b),
deepcopy(probhom), solh, (@optic _.b),
PALC(tangent = Bordered()),
# PALC(),
# MoorePenrose(),
Expand All @@ -151,7 +151,7 @@ br_hom_c = continuation(
plot_solution = (x,p;k...) -> begin
𝐇𝐨𝐦 = p.prob
par0 = set(BK.getparams(𝐇𝐨𝐦), BK.getlens(𝐇𝐨𝐦), x.x[end][1])
par0 = set(par0, (@lens _.b), p.p)
par0 = set(par0, (@optic _.b), p.p)
sol = get_homoclinic_orbit(𝐇𝐨𝐦, x, par0)
m = (𝐇𝐨𝐦.bvp isa PeriodicOrbitOCollProblem && 𝐇𝐨𝐦.bvp.meshadapt) ? :d : :none
plot!(sol.t, sol[:,:]',subplot=3, markersize = 1, marker=m)
Expand Down Expand Up @@ -212,10 +212,10 @@ probhom, solh = generate_hom_problem(
update_every_step = 4,
# ϵ0 = 1e-6, ϵ1 = 1e-5,
t0 = 75, t1 = 25,
# freeparams = ((@lens _.T), (@lens _.ϵ1),)
# freeparams = ((@lens _.T), (@lens _.ϵ0)),
# freeparams = ((@lens _.ϵ0), (@lens _.ϵ1)), # WORK BEST
freeparams = ((@lens _.T),),
# freeparams = ((@optic _.T), (@optic _.ϵ1),)
# freeparams = ((@optic _.T), (@optic _.ϵ0)),
# freeparams = ((@optic _.ϵ0), (@optic _.ϵ1)), # WORK BEST
freeparams = ((@optic _.T),),
)

_sol = get_homoclinic_orbit(probhom, solh, BK.getparams(probhom); saveat=.1)
Expand All @@ -225,7 +225,7 @@ optn_hom = NewtonPar(verbose = true, tol = 1e-9, max_iterations = 7)
optc_hom = ContinuationPar(newton_options = optn_hom, ds = -1e-4, dsmin = 1e-6, dsmax = 1e-3, plot_every_step = 1, max_steps = 10, detect_bifurcation = 0, save_sol_every_step = 1)

br_hom_sh = continuation(
deepcopy(probhom), solh, (@lens _.b),
deepcopy(probhom), solh, (@optic _.b),
# PALC(tangent = Bordered()),
PALC(),
# ANM(6, 1e-8)
Expand All @@ -237,7 +237,7 @@ br_hom_sh = continuation(
plot_solution = (x,p;k...) -> begin
𝐇𝐨𝐦 = p.prob
par0 = set(BK.getparams(𝐇𝐨𝐦), BK.getlens(𝐇𝐨𝐦), x.x[end][1])
par0 = set(par0, (@lens _.b), p.p)
par0 = set(par0, (@optic _.b), p.p)
sol = get_homoclinic_orbit(𝐇𝐨𝐦, x, par0)
m = (𝐇𝐨𝐦.bvp isa PeriodicOrbitOCollProblem && 𝐇𝐨𝐦.bvp.meshadapt) ? :d : :none
plot!(sol.t, sol[1:6,:]',subplot=3, markersize = 1, marker=m)
Expand Down
Loading

0 comments on commit 4da7e33

Please sign in to comment.