OptimalControl.jl

Appearance

Poisson bracket ​

For two functions    , the Poisson bracket is

julia
using OptimalControl
using OrdinaryDiffEqTsit5

Computing a bracket ​

With    and   ,  ,  ,   and  , so   , which is at  ,  :

julia
H(x, p) = p[1] * x[2] + p[2] * x[1]
G(x, p) = x[1]^2 + p[2]^2

x, p = [1.0, 2.0], [3.0, 4.0]
Poisson(H, G)(x, p)
-20.0

On Hamiltonians, Poisson returns a Hamiltonian, which can be bracketed again:

julia
HT = Hamiltonian(H)
GT = Hamiltonian(G)

B = Poisson(HT, GT)
B isa Hamiltonian, Poisson(B, GT)(x, p)   # {{H, G}, G}
(true, -32.0)

For functions of , pass is_autonomous=false, and for functions of , is_variable=true. The derivatives are taken with respect to and only, so adding terms in does not change the bracket:

julia
Ht(t, x, p) = t + H(x, p)
Gt(t, x, p) = t^2 + G(x, p)

Poisson(Ht, Gt; is_autonomous=false)(1.0, x, p)
-20.0

Properties ​

The Poisson bracket is

  • bilinear:     for scalars ;

  • antisymmetric:   ;

  • a derivation (Leibniz rule):   ;

  • and it satisfies the Jacobi identity:    .

Checked at the point above, with a third function :

julia
K(x, p) = x[1] * p[1] * p[2] + x[2]^3

antisymmetry = Poisson(H, G)(x, p) + Poisson(G, H)(x, p)
leibniz = Poisson((x, p) -> H(x, p) * G(x, p), K)(x, p) -
          (H(x, p) * Poisson(G, K)(x, p) + G(x, p) * Poisson(H, K)(x, p))
jacobi = Poisson(Poisson(H, G), K)(x, p) + Poisson(Poisson(K, H), G)(x, p) +
         Poisson(Poisson(G, K), H)(x, p)
antisymmetry, leibniz, jacobi
(0.0, 0.0, 0.0)

The bracket of two lifts is the lift of the Lie bracket,   (see The bridge identity).

Time derivative along an extremal ​

The Hamiltonian vector field of is   . The Poisson bracket is the Lie derivative along it:

So along an integral curve of , that is, an extremal  ,   , a function varies as

the partial time derivative, for the explicit dependence on , plus the Poisson bracket, for the motion along the flow. To check it, integrate the flow of a Hamiltonian that depends on time, and compare the derivative of along it, by finite differences, with  :

julia
Hn(t, x, p) = p[1] * x[2] + p[2]^2 / 2 + t * x[1]
Gn(t, x, p) = t * p[1] + x[1] * x[2]

flow = Flow(Hamiltonian(Hn; is_autonomous=false); reltol=1e-12, abstol=1e-12)
G_along(t) = Gn(t, flow(0.0, x, p, t)...)

t, h = 0.5, 1e-4
xt, pt = flow(0.0, x, p, t)
lhs = (G_along(t + h) - G_along(t - h)) / 2h
rhs = ∂ₜ(Gn)(t, xt, pt) + Poisson(Hn, Gn; is_autonomous=false)(t, xt, pt)
lhs, rhs
(21.929699666296543, 21.929699707031244)

In particular, an autonomous is constant along its extremals, since  .

Application: singular controls ​

Take a system   , whose pseudo-Hamiltonian    is linear in , with  , the lifts of . The control maximises , so it is given by the sign of the switching function , except on an arc where vanishes: a singular arc. There, the time derivatives of vanish too. By the formula above, with autonomous,

where  ,   and  . On the arc,   , and the second derivative gives the singular control

For example, take the vehicle of Singular control:  ,   ,  , so     and  . Then

On the singular arc,   and  , and the formula gives  :

julia
H0(q, p) = p[1] * cos(q[3]) + p[2] * (sin(q[3]) + q[1])
H1(q, p) = p[3]

H01 = Poisson(H0, H1)
H001 = Poisson(H0, H01)
H101 = Poisson(H1, H01)

q = [0.3, -0.1, 0.5]
p_s = [1.7, 1.7 * tan(q[3]), 0]   # a point where H₁ = H₀₁ = 0
H1(q, p_s), H01(q, p_s), -H001(q, p_s) / H101(q, p_s), sin(q[3])^2
(0.0, 0.0, 0.22984884706593012, 0.22984884706593015)

The example page integrates the extremal with this control and checks it against the direct solution.

Vector fields are not Hamiltonians ​

Poisson rejects a VectorField: lift it first, or use ad for the Lie bracket.

julia
julia> XV = VectorField(x -> [x[2], -x[1]]);

julia> Poisson(XV, x -> x[1])
IncorrectArgument → top-level scope, REPL[2]:2
│
│  Poisson is not defined for AbstractVectorField operands
│
│  Context  Poisson on AbstractVectorField
│  Hint     Use ad(X, Y) for the Lie bracket of VectorFields
└─

See also ​