diff --git a/src/minimal_period.jl b/src/minimal_period.jl index c267ae1..176e86c 100644 --- a/src/minimal_period.jl +++ b/src/minimal_period.jl @@ -15,6 +15,9 @@ period is also called prime, principal or fundamental period. * `maxiter = 40` : Maximum number of Poincare map iterations. Continuous-time systems only. If the number of Poincare map iterations exceeds `maxiter`, but the point `u0` has not returned to `atol` neighborhood of itself, the original period `po.T` is returned. +* `Δt = missing` : The time step between points in the trajectory `minT_po.points`. If `Δt` + is `missing`, then `Δt=minT_po.T/$(default_Δt_partition)` is used. Continuous-time + systems only. ## Description @@ -31,25 +34,35 @@ this hyperplane. Using the Poincare map, the hyperplane crossings are checked. T first crossing that is within `atol` distance of the initial point `u0` is the minimal period. At most `maxiter` crossings are checked. """ -function minimal_period(ds::DynamicalSystem, po::PeriodicOrbit; kwargs...) +function minimal_period( + ds::DynamicalSystem, + po::PeriodicOrbit; + Δt=missing, + kwargs... + ) type1 = isdiscretetime(ds) type2 = isdiscretetime(po) if type1 == type2 newT = _minimal_period(ds, po.points[1], po.T; kwargs...) - return _set_period(ds, po, newT) + return _set_period(ds, po, newT, Δt) else throw(ArgumentError("Both the periodic orbit and the dynamical system have to be either discrete or continuous.")) end end -function _set_period(ds::DynamicalSystem, po, newT) - if newT == po.T - return po - else - # ensure that continuous po has the same amount of points - isdiscretetime(ds) ? Δt = 1 : Δt = newT/default_Δt_partition - return PeriodicOrbit(complete_orbit(ds, po.points[1], newT; Δt=Δt), newT, po.stable) +function _set_period(ds::DynamicalSystem, po, newT, Δt) + newT == po.T && return po + + Dt = 1 + if !isdiscretetime(ds) + if ismissing(Δt) + # ensure that continuous po has the same amount of points + Dt = newT/default_Δt_partition + else + Dt = Δt + end end + return PeriodicOrbit(complete_orbit(ds, po.points[1], newT; Δt=Dt), newT, po.stable) end function _minimal_period(ds::DiscreteTimeDynamicalSystem, u0, T; atol=1e-4) diff --git a/test/minimal_period.jl b/test/minimal_period.jl index 006eab0..4a73199 100644 --- a/test/minimal_period.jl +++ b/test/minimal_period.jl @@ -31,9 +31,13 @@ end n = 10 po = PeriodicOrbit(ds, u0, n*T, 0.01) minT_po = minimal_period(ds, po) - @test length(minT_po.points) == minT_po.T / (minT_po.T / PeriodicOrbits.default_Δt_partition) + @test length(minT_po.points) == PeriodicOrbits.default_Δt_partition @test (ismissing(po.stable) && ismissing(minT_po.stable)) || (po.stable == minT_po.stable) @test isapprox(T, minT_po.T; atol=1e-4) + + Dt = 0.01 + minT_po = minimal_period(ds, po; Δt=Dt) + @test length(minT_po.points) == floor(minT_po.T / Dt) end function normalhopf(u, p, t)