17 August, 2026

Notes

2026-10-30

To do!

On period orbits

Other…

Questions

  • what dos it mean to be a strobocopic Map \(T\)?

Finding periodic Orbits

The algorithm I use for finding periodic orbits initializes the integrator. Then the periodic orbit gets sought for, after which the actual periodicity is found and the trajectory in configuration space is generated.

p       =   (1.0, 1.0, 1.0)
E0      =   0.1127
nmax_search = 40
nmax        =   1
tmax    =   20_000.0
y0, py0, c0 =    0.0, 0.2, C_TEAL
y1, py1, c1 =   -0.35, 0.0, C_ORANGE
y2, py2, c2 =    0.27, 0.3, C_GREEN
y3, py3, c3 =    0.3, 0.0, C_PURPLE
v_init      =   [[y0, py0], [y1, py1], [y2, py2], [y3, py3]] #, 
color_scheme= [c0, c1, c2, c3]

bff = create_integrators(p; nmax=Ref(nmax_search))

prm = SectionParams(bff.integ_fast, bff.yf, bff.pyf, bff.tsf,
                         bff.integ_dense, bff.yd, bff.pyd, bff.tsd,
                         Ref(nmax), Ref(nmax_search), E0, p, tmax)

orbits = orbit_table()
for (i,v0) in enumerate(v_init)
    orbit = analyse_seed(v0, nmax, prm; id=i)
    orbit === nothing || push!(orbits, orbit)
end

With the function creating the integrator for the time evolution \(\dot{\x}=\mathbf{f}(\x)\) which in turn defines the poincaré map \(T\) by the ‘condition(u, t, integ)=u[1]’ i.e. \(x=0, p_x\geq0\) creates the ‘affect!(integ)’ of collecting \((y, p_y, t)\), the integration gets terminated if there are more than ‘nmax’ itersections with poincaré section. This is embedded problem to be solved (‘ODEProblem’) by the ‘ContinuousCallback’ function. The integrator gets initialized with the integration method (Runge-Kutter method, \(9^{\text{th}}\) order) and can easily be reinitilized for new initial conditions \(\x_0\).


"""
    Two integrators, one for finding periodic orbit (integ_fast) and another (integ_dense) of for finding minimal periodicity and finding the trajectories of the orbit.

return (; integ_fast, yf, pyf, tsf, integ_dense, yd, pyd, tsd)
"""
function create_integrators(p; tmax = 20_000.0, nmax = Ref(1))
    # fast: for residual evaluations, no trajectory saved
    yf, pyf, tsf = Float64[], Float64[], Float64[]
    condition(u, t, integ) = u[1]
    affect!(integ) = begin
        push!(yf, integ.u[2]); push!(pyf, integ.u[4]); push!(tsf, integ.t)
        length(yf) ≥ nmax[] && terminate!(integ)
    end
    cbf = ContinuousCallback(condition, affect!, nothing; abstol = 1e-13)
    integ_fast = init(ODEProblem(HenonHeiles.equations!, zeros(4), (0.0, tmax), p),
                      Vern9(); abstol = 1e-14, reltol = 1e-14,
                      save_everystep = false, save_start = false, callback = cbf)

    # dense: for plotting
    yd, pyd, tsd = Float64[], Float64[], Float64[]
    affect2!(integ) = begin
        push!(yd, integ.u[2]); push!(pyd, integ.u[4]); push!(tsd, integ.t)
        length(yd) ≥ nmax[] && terminate!(integ)
    end
    cbd = ContinuousCallback(condition, affect2!, nothing; abstol = 1e-13)
    integ_dense = init(ODEProblem(HenonHeiles.equations!, zeros(4), (0.0, tmax), p),
                       Vern9(); abstol = 1e-14, reltol = 1e-14,
                       saveat = dt, callback = cbd)

    return (; integ_fast, yf, pyf, tsf, integ_dense, yd, pyd, tsd)
end

The key analysis occurs in ‘analyse_seed’. This contains the function that findes the orbit by locating the root (Newtons Method or similar). If this is converged, the root is being checked to find the minimum number of reoccurences of the orbit, (i.e. \(n=6\) may find orbits with periods \(p=1,2,3\) which will be reffered to as primes of the orbit). Then the orbit gets characterized by analysing the linear approximation of the Hamiltonian flow of the poincaré Map \(DT\), which is a \(2\times2\)–Matrix and characterized by its eigenvalues, hence by the trace of the Matrix \(DT\). Due to the symplectic behaviour \(det(DT)\overset{!}{=}1\).

"returns a row to the orbit"
function analyse_seed(v0, n, prm; id=0)::Union{Row,Nothing}
    res =   find_orbit(v0, n, prm)
    res.converged || return nothing
    orb =   minPeriodicity(res.v, prm; tol = 1e-8)
    orb.Nperiod === nothing && return nothing
    DT=     get_DT(res.v, orb.Nperiod, prm)
    trace=  tr(DT); detDT=det(DT)
    index = kind_index(trace; ε = 1e-6)
    return (; E = prm.E, n = n, seed_y = v0[1], seed_py = v0[2],
                     y = res.v[1], py = res.v[2],
                     prime = orb.Nperiod, T = orb.Tperiod,
                     trace = trace, detDT = detDT,
                     resnorm = res.resnorm, iters = size(res.history, 2),
                     sec_y = orb.pMap[1,:], sec_py=orb.pMap[2,:],
                     history_y = res.history[1,:], history_py = res.history[2,:],
                     traj_x = [u[1] for u in orb.traj], traj_y = [u[2] for u in orb.traj],
                     class = KIND_LABEL[index], index=index, id=id)
end

So I am currently wandering about the rotation, does the rotation work the way I think it does? I seem to find examples where it does not work the way I expect, there seems to be a extra reflection around the \(x\)-axis.

Configuration space.

Poincaré Map

I ran a few other parameters, here

Configuration space.

Poincaré Map

Configuration space.

Poincaré Map

Configuration space.

Poincaré Map

I also found an error, because I only reflected the periodic point and tried to ‘minPeriodicity’, the trajectory which was hypebollic (unstable) turned out not to find the periodic poit, note here, the orbit had periodicity of \(n=9\), so an unstable periodic orbit \(T^9 \x=\x\), for \(\x\in\Sigma\), need a high accuracy to the starting position, to not to diverge. The consequence I take, is that I reflect \(S(y,p_y)=(y,-p_y)\) and then find the periodic orbit ‘find_orbit’/‘solve_orbit’.

Configuration space.

Poincaré Map, we see the trajectory can not lie on a periodic orbit.

Here I think we can get a full fledged understanding of the effects of rotations of configuration space. We use orbits with periodicity of \(n=1\), which corresponds to orbits that return with one poincaré map.

Here we see the periodic orbits that are initially taken to see the effects of the symmetry transformations. The periodic orbits are found at \[o_i\in \{(0.0,0.237), (-0.354,0.0), (0.276, 0.296), (0.290, 0.0)\}\]

Configuration space, of the base selected periodic orbits.

Poincaré Map, we see the trajectory of the base selected periodic orbits.

Configuration space, rotations.

Poincaré Map, rotated states, periodic orbits.

Configuration space, rotations.

Poincaré Map, rotated states, periodic orbits.

Configuration space, rotations.

Poincaré Map, all orbits gained by rotations of states.

Configuration space, all orbits through one run with reflections.

Poincaré Map, all orbits obtained through one run, with reflection.

On Poincaré sections

A special kind of poincaré section are the surface of sections. These take a \(N\)-dimensional phase space \(\x\) and cut it by a \((N-1)\)-dimensional plane or manifold \(\Sigma\). The first return of to \(\Sigma\) defines a diffeomorphism, \(T:\Sigma \rightarrow \Sigma\), \(T(\x_0)=\x(\x_0,\tau)\), where \(\tau\) is the return time.

Periodic orbits are then defined by the condition \[T^n(\x_0)=\x_0\] they are fixed points to the diffeomorphism \(T^n\). We can linearize \(T^n \rightarrow DT^n\) at a given \(\x_0\) approximating the diffeomorphism in the small, open area \(U\subset\Sigma\) around \(\x_0\in U\). This gives the possibility to identify a stability matrix \(M\) to \(DT^n\) which is associated to the flow \(M=\exp(At)\) or \(A\ln(M)\,|t=1\) at \(t=1\); The eigenvalues \(\lambda\) of M correspond through \(\gamma=\ln(\lambda)\) to the eigenvalues of \(A\) as well as shared eigenvectors.

With the Hartman - Grobman theorem, we know that the linear approximation \(DT^n\) around the stable points \(\x_e\) are topologically equivalent to the \(T\) around \(\x_e\) itself

Lost Wonderers Thoughts: