Welcome to Cadabra Q&A, where you can ask questions and receive answers from other members of the community.
+1 vote

A follow on question to my previous one:

[Not sure why the Latex typesetting is not displayed throughout. Looks fine in the preview.]

Context. I would like to now extend the coordinate computation from the configuration space ${\mathbf q }= (x, y)$ to phase space $({\mathbf q, \mathbf p}) = (x, y, p{x}, p\{y} )$ to compute the Hamiltonian vector field a.k.a. the Hamiltonian geodesic flow:

$Xk = g^{i j} p_j \frac{\partial } {\partial x^i} + g^{ kj} \Gamma^l\{j i} p_k p_l \frac{\partial } {\partial p_i}$

from the geometric Hamiltonian

$K = g^{i j}(q) p_i p_j$

with

$ \frac {dx^i} {ds} = \frac{\partial K} {\partial p_i} = g^{i j} p_j$

$ \frac {dp_i} {ds} = - \frac{\partial K} {\partial x^i} = g^{k j} \Gamma^l_{j i} p_k p_l$

where the metric is, as before,

$g_{i j} = 2[H_{o} - V(q)] \delta_v{i j}$

with a specific $V(q) = V(x, y)$ to be substituted at the end of the computation.

The code cell below is my current attempt implement $X_k$ using a separate set of indices for $q$ and $p$ so as to distinguish between the partial derivatives $\frac{\partial } {\partial x^i}$ and $\frac{\partial } {\partial p_{\nu}}$. However, will mixing the two sets of indices in the metric and its inverse work?

The code fragment

substitute(D1, h_flow);
substitute(D1, Christoffel);
evaluate(D1, metric, rhsonly=True);

returns $\nabla{q^{a}}$ for each line. To sum up, how would I go about implementing this?

# Version 2.0
# Geometric Hamiltonian dynamics
# Ref. "Geometric approach to Lyapunov analysis in Hamiltonian dynamics"
# Most recent version [V3]: https://arxiv.org/abs/nlin/0104005
# Published version: https://journals.aps.org/pre/abstract/10.1103/PhysRevE.64.066206

{x, y}::Coordinate;
{i, j, k, l, m, n, h#}::Integer(0..1);
{i, j, k, l, m, n, h#}::Indices(values={x, y}, position=fixed);

(px, py)::Coordinate;
{\iota, \kappa, \lambda, \mu, \nu, \xi, \eta#}::Integer(0..1);
{\iota, \kappa, \lambda, \mu, \nu, \xi, \eta#}::Indices(values={px, py}, position=fixed);

\nabla::Derivative;
\partial_{#}::PartialDerivative;

g_{i j}::Metric;
g_{i j}::Depends(x, y);

g^{i j}::InverseMetric;
g^{i j}::Depends(x, y);

# Hénon-Heiles Hamiltonian potential: V
V::Depends(x, y);
V := (x**2 + y**2)/2 + (x**2)*y - (y**3)/3;

# Ho is a constant of the motion: conservation of energy
metric := { g_{x x} = 2*(Ho - V), 
            g_{y y} = 2*(Ho - V) };
complete(metric, $g^{i j}$);

Christoffel := \Gamma^{i}_{j k} = (1/2)*g^{i l}*( \partial_{k}{g_{l j}} 
                                                 + \partial_{j}{g_{l k}} 
                                                 - \partial_{l}{g_{j k}} );

evaluate(Christoffel, metric, rhsonly=True);
# substitute(Christoffel, $V = @(V)$);
# map_sympy(Christoffel, 'simplify');

p_{#}::Depends(px, py); 

# Hamiltonian flow template
h_flow := \nabla{v?_{h}} -> g^{i \kappa}*p_{\kappa}*\partial_{i}{v?_{h}} + g^{\kappa \lambda} * \Gamma^{\mu}_{\lambda \nu}*p_{\kappa}*p_{\mu}*\partial_{\nu}{v?_{h}};

D1 := \nabla{q^{a}};

substitute(D1, h_flow);
substitute(D1, Christoffel);
evaluate(D1, metric, rhsonly=True);
in General questions by (220 points)
edited by

Latex display. Replacing the underscore symbol by by the ascii code for it fixed the display issue.

Or use backslash-underscore \_ instead of just a plain underscore _.

Got it. Thank you.

1 Answer

+1 vote

Because the a index is not in the index set of which h forms part, the pattern in h_flow will not match when you try to substitute h_flow into D1. You can do

D1 := \nabla{q_{h}};

and then the substitutes work. To make evaluate work, you need to declare that q_{h} depends on x and y,

q_{h}::Depends(x,y);

After that things seem to work as intended (have not verified the result).

by (76.4k points)

Thank you for your explanation and corrections. Evaluation now generates a result.

After that things seem to work as intended (have not verified the result).

Almost.

$X_k = g^{i j} p_j \frac{\partial } {\partial x^i} + g^{k j} \Gamma^l_{j i} p_k p_l \frac{\partial } {\partial p_i}$

The first term evaluates as expected, however, the second term of

h_flow := \nabla{v?^{h}} -> ( g^{i \kappa}*p_{\kappa}*\partial_{i}{v?^{h}} + 
                              g^{\kappa \lambda} * \Gamma^{\mu}_{\lambda \nu}*p_{\kappa}*p_{\mu}*\partial_{\nu}{v?^{h}});

with

q^{h}::Depends(x, y);
D1 := \nabla{q^{h}};

yields terms of the form

$p_x^2 \partial_x V \partial_x q^h$, $p_y^2 \partial_y V \partial_y q^h$, etc.

The partial derivatives are with respect $(x, y)$ rather then $(p_x, p_y)$.

I had guessed that using separate sets of indices for $(x, y)$ and $(p_x, p_y)$ would suffice and that metric terms such as $g^{\kappa \lambda}$ would be evaluated according the numerical values of the indices.

Is it possible to define and distinctly label two partial derivatives

\partial_{i, j, k, l, m, n, h#}::PartialDerivative;

and

\partial_{\iota, \kappa, \lambda, \mu, \nu, \xi, \eta#}::PartialDerivative;

so that the first would be with respect to $(x, y)$ and the second with respect to $(p_x, p_y)$?

...