Reference

PoincareInvariants

A Julia library for the computation of Poincaré integral invariants.

compute!(pinv::AbstractPoincareInvariant, points::AbstractMatrix, t::Real=NaN, p=nothing)

computes a Poincaré invariant using setup object pinv, where points represents the image of the curve or surface parameterisation evaluated on some set of points, e.g. a grid.

t is the time at which the invariant is evaluated p are any user supplied optional arguments. Both t and p are passed directly to the differential form.

The differential form is an in-place function following the convention form(out, t, z, p): it writes its value at the phase space point z (one point, a length-D vector) into the preallocated output out — a length-D vector for a one form θ, a D×D matrix for a two form ω. This argument order matches the one- and two-forms exported by GeometricEquations/GeometricProblems (ϑ(Θ, t, q, params), ω(Ω, t, q, params)), so a problem's own forms can be passed directly. A constant ω::AbstractMatrix (e.g. a CanonicalSymplecticMatrix) is used as-is and is not called.

Plan implementations should define a method compute!(pinv, t::Real, p), which acts on the internal points storage pinv.points.

compute!(pinv::AbstractPoincareInvariant, points::AbstractVector{<:AbstractVector},
    times::Union{AbstractVector{<:Real}, Real}=NaN, p=nothing)

computes a Poincaré invariant, using the setup object pinv, at each time for a set of trajectories given by points.

points represents an AbstractVector of trajectories. Each element of points is a trajectory, which is itself an AbstractVector of some kind of iterable or vector. times is either a constant, like NaN, or a vector of the times at which the trajectories have been evaluated, i.e. the i-th phase space position in each trajectory was evaluated at time times[i]. p is an arbitrary optional parameter which is passed to the differential form, like the the time.

compute!(pinv::AbstractPoincareInvariant, sol::EnsembleSolution, p=nothing)

computes the Poincaré invariant at every saved time step of an EnsembleSolution returned by GeometricIntegrators.integrate. Each ensemble member provides the trajectory of one point on the curve or surface. Returns a Vector holding one invariant value per time step.

The optional parameter p is passed to the differential form, just like the time.

FirstPoincareInvariant{T, D, θT, P} <: AbstractPoincareInvariant{T, D}
FirstPI{T, D, θT, P} <: AbstractPoincareInvariant{T, D}

setup object to compute the first invariant.

FirstPoincareInvariant{T, D}(θ::θT, N::Integer, plan::P)

constructs a FirstPoincareInvariant setup object to calculate integral invariants, given a numeric type T, a phase space dimension D, a differential form θ, a number of points N and a plan. Note that the number of points may not be exactly N. The true number used depends on the implementation and is guaranteed to be no smaller than N and not too much larger.

FirstPoincareInvariant{T, D}(θ::θT, N::Integer, P::Type=DEFAULT_FIRST_PLAN)

constructs a FirstPoincareInvariant setup object to calculate integral invariants, given a numeric type T, a phase space dimension D, a differential form θ, a number of points N and a plan type P, to be initialised and used for future computation.

The plan type defaults to DEFAULT_FIRST_PLAN, which is currently set to FirstFourierPlan.

FirstPoincareInvariant{T, D, typeof(canonical_one_form!)}(N::Integer, P=DEFAULT_FIRST_PLAN)
CanonicalFirstPI{T, D}(N::Integer, P=DEFAULT_FIRST_PLAN)

creates a setup object to compute the first integral invariant using the canonical one form in phase space of dimension D with numeric type T, N points and plan type P.

SecondPoincareInvariant{T, D, ωT, PS, P} <: AbstractPoincareInvariant{T, D}
SecondPI{T, D, ωT, PS, P} <: AbstractPoincareInvariant{T, D}

setup object used to compute the second invariant.

SecondPoincareInvariant{T, D}(ω, N, plan::P)

constructs a SecondPoincareInvariant setup object to calculate integral invariants, given a numeric type T, a phase space dimension D, a differential form ω, a point specification N, usually a number of points, and a plan.

SecondPoincareInvariant{T, D}(ω, N, P::Type=DEFAULT_SECOND_PLAN)

constructs a SecondPoincareInvariant setup object to calculate integral invariants, given a numeric type T, a phase space dimension D, a differential form ω, a point specification N and a plan type P, to be initialised and used for future computation.

The plan type defaults to DEFAULT_SECOND_PLAN, which is currently set to SecondChebyshevPlan. Note that the number of points may not be exactly N. The true number used depends on the implementation and is guaranteed to be no smaller than N and not too much larger. For the SecondFinDiffPlan, the grid of points may also be specified as a tuple (Nx, Ny), if non-square grids are sought.

SecondPoincareInvariant{T, D, CanonicalSymplecticMatrix{T}}(N, P=DEFAULT_FIRST_PLAN)
CanonicalSecondPI{T, D}(N, P=DEFAULT_FIRST_PLAN)

creates a setup object to compute the second integral invariant using the canonical two form in phase space of dimension D with numeric type T, point specification N and plan type P.

AbstractPoincareInvariant{T, D}

represents a Poincare integral invariant in a phase space of dimension Dusing numeric typeT` for calculations.

getpoints(pinv::AbstractPoincareInvariant)

returns points on which to evaluate the phase space line or surface parameterisation so as to compute! pinv.

getpointnum(pinv::AbstractPoincareInvariant)

returns number of points to sample in phase space to compute! pinv.

getpointspec(pinv::AbstractPoincareInvariant)

get point specification, which may, for example, be a tuple specifying a grid or a number giving the number of points used to sample in phase space.

getdim(pinv::AbstractPoincareInvariant)

returns dimension of phase space to compute! pinv in.

getform(pinv::AbstractPoincareInvariant)

get invariant one- or two-form.

getplan(pinv::AbstractPoincareInvariant)

returns the plan that will be used to compute! pinv.

PIEnsembleProblem(prob, pinv, init)

builds a GeometricEquations.EnsembleProblem from a base problem prob, a Poincaré invariant setup object pinv, and a curve or surface parameterisation init.

prob is a geometric EquationProblem (e.g. an ODEProblem, HODEProblem, IODEProblem or LODEProblem) specifying the dynamics, time span and time step. The parameterisation init is sampled with getpoints, and each sampled phase space point becomes the initial condition of one member of the ensemble. For HODEProblem/PODEProblem the point is split into the position and momentum halves (q, p); for IODEProblem/LODEProblem the momentum is initialised from the equation's one-form as p₀ = ϑ(t₀, q₀).

The returned ensemble can be integrated with GeometricIntegrators.integrate and the invariant computed from the resulting solution with compute!.

CanonicalSymplecticMatrix{T}(n::Integer)

constructs a canonical symplectic matrix of size (n, n) with eltype T. n must be even and positive. See the examples to see the form of the canonical symplectic matrix as defined here.

Examples

julia> CanonicalSymplecticMatrix(4)
4×4 CanonicalSymplecticMatrix{Int64}:
 0  0  -1   0
 0  0   0  -1
 1  0   0   0
 0  1   0   0

julia> CanonicalSymplecticMatrix{Int32}(6)
6×6 CanonicalSymplecticMatrix{Int32}:
 0  0  0  -1   0   0
 0  0  0   0  -1   0
 0  0  0   0   0  -1
 1  0  0   0   0   0
 0  1  0   0   0   0
 0  0  1   0   0   0

canonical_one_form!(out, t, z, p)

writes the canonical one form $\vartheta$ at the phase space point z into out, following the in-place form(out, t, z, p) convention. For a phase space point z = (q, p) of even length n = 2m, the canonical one form is $\vartheta = (p, 0)$, i.e. out[1:m] .= z[m+1:n] and out[m+1:n] .= 0. The time t and parameters p are ignored.

plot_invariant(pinv, sol; p=nothing, title="")
plot_invariant(pinv, label => sol...; p=nothing, title="", position=:lb)
plot_loop(sol; xlabel, ylabel, title="Advected Loop", nsteps=10)
plot_surface(grid, sol; xlabel, ylabel, title="Advected Surface", nsteps=10)

visualise Poincaré invariants and the advected curve/surface of an EnsembleSolution.

plot_invariant plots the relative error (I(t) - I(0)) / I(0) of the first or second invariant over time (a legend-labelled comparison when several label => sol pairs are given). plot_loop and plot_surface show the advected curve and surface in phase space, stacked along the time axis. Mutating plot_*! variants draw into an existing axis.

These functions are provided by a package extension; load a Makie backend (e.g. using CairoMakie) to enable them.